1. El paquete caret

2. Visualización

Cargar conjunto de datos después de imputación KNN. Una vez cargada la tabla con los datos, podemos analizar los datos mediante gráficos de dispersión. Para realizar estos gráficos hay que utilizar la función featurePlot:

hepatitis.KnnImp = read.csv2(file="hepatitis.kNNImpute.csv", header=TRUE, sep = ",", dec = ".")
hepatitis.KnnImp = hepatitis.KnnImp[1:20]
#str(hepatitis)
#featurePlot(x = hepatitis.KnnImp[, 1:4], y = hepatitis.KnnImp$PRONOSTICO, %$ plot = "pairs",auto.key = list(columns = 2))
featurePlot(x = hepatitis.KnnImp[, 1:4], y = hepatitis.KnnImp$PRONOSTICO, plot = "pairs",auto.key = list(columns = 2))

Ejercicio 1.

Genera un gráfico de dispersión para las variables numericas FOSFATOalc, SGOT, ALBUMINA Y PROTIME.

  • 1.a) Prueba todas las opciones del parámetro pairs, box, strip, density y ellipse.
vars_1a = c("FOSFATOalc", "SGOT", "ALBUMINA", "PROTIME")
featurePlot(x=hepatitis.KnnImp[, vars_1a], y=hepatitis.KnnImp$PRONOSTICO, plot="pairs", scales=list(x=list(relation="free"), y=list(relation="free")), auto.key=list(columns=3))

#featurePlot(x = hepatitis.KnnImp[, vars_1a], y = hepatitis.KnnImp$PRONOSTICO, plot = "box",auto.key = list(columns = 2))
#Presentar en el eje Y la escala en función de la varaible presentada. "scales"
featurePlot(x=hepatitis.KnnImp[, vars_1a], y=hepatitis.KnnImp$PRONOSTICO, plot="box", scales=list(x=list(relation="free"), y=list(relation="free")), auto.key=list(columns=3))

featurePlot(x=hepatitis.KnnImp[, vars_1a], y=hepatitis.KnnImp$PRONOSTICO, plot="density", scales=list(x=list(relation="free"), y=list(relation="free")), auto.key=list(columns=3))

featurePlot(x=hepatitis.KnnImp[, vars_1a], y=hepatitis.KnnImp$PRONOSTICO, plot="ellipse", scales=list(x=list(relation="free"), y=list(relation="free")), auto.key=list(columns=3))

featurePlot(x=hepatitis.KnnImp[, vars_1a], y=hepatitis.KnnImp$PRONOSTICO, plot="strip", scales=list(x=list(relation="free"), y=list(relation="free")), auto.key=list(columns=3))

  • 1.b) >Se pueden extraer alguna conclusión de las gráficas anteriormente generadas?
    Conclusiones. Parece que la variable PROTIME es muy determinista para concluir si vive o fallece una persona, con un nivel de PROTIME sobre 30 hay mas probabilidad de que viva, por otroa lado con un bajo nivel más probabilidad de que fallesca. Un comportamiento similar se aprecia con la variable ALBUMINA, pero en esta con el poder determinista en los extremos, es decir, para valores muy altos de ALBUMINA (superior a 4.5) y fallece con valores bajos (inferior a ~ 2.5). SGOT no es una variable determinista, dado que se presentan las mismas curvas de densidad entre las personas que VIVEN y las que FALLECEN. Los mejores tipos gráficas para estas variables en mi opinión son los de box y density que permiten ver la distribución para ambas clases.

3. Preprocesamiento

#regresa los indices que corresponden a entrenamiento. (muestreo estratificado, manetniendo la proporción de las clases)
trainIndex <- createDataPartition(hepatitis.KnnImp$PRONOSTICO, p = 0.66, list = FALSE, times = 1)
hepatitisTrain <- hepatitis.KnnImp[trainIndex,]
hepatitisTest <- hepatitis.KnnImp[-trainIndex,]
#componentes principales
# no preprocesa sino genera la función que hará la transformación
hepatitisPCA <- preProcess(hepatitis.KnnImp[,1:ncol(hepatitis.KnnImp)-1],
method = "pca",thresh = 0.95)
#print(prostata.preprocess)
print(hepatitisPCA)
## Created from 155 samples and 19 variables
## 
## Pre-processing:
##   - centered (6)
##   - ignored (13)
##   - principal component signal extraction (6)
##   - scaled (6)
## 
## PCA needed 6 components to capture 95 percent of the variance
#escala y centra las 6 variables princiapales ignora 13 (por tipo de dato)
#se necesitan 6 transformaciones para capturar el 95% de la varianza. No son los datos originales, sino, combinaciones lineales de los originales. 
#Busca una orientación de los ejes, generando dimensiones ortogonales, donde las distancias se agrupen de una forma más sencilla
#con predict genera la transformación a partir de la funcion de transformacion y los datos que se quieren transformar #tanto para test como para training, excluyendo a la calse en caso de clasificaión.
PCATrain <- predict(hepatitisPCA,hepatitisTrain[,1:ncol(hepatitisTrain)-1])
PCATest <- predict(hepatitisPCA,hepatitisTest[,1:ncol(hepatitisTest)-1])
#se generan dataframe de test y train 
PCATrain$PRONOSTICO = hepatitisTrain$PRONOSTICO
PCATest$PRONOSTICO = hepatitisTest$PRONOSTICO
#PCATrain <- data.frame(PCATrain,hepatitisTrain$PRONOSTICO)
#PCATest <- data.frame(PCATest,hepatitisTest$PRONOSTICO)
#PC1-PC6 son una combinacion linea de las variables originales y se pierde el nombre original de las variables
str(PCATrain)
## 'data.frame':    104 obs. of  20 variables:
##  $ SEXO       : Factor w/ 2 levels "FEMENINO","MASCULINO": 2 1 1 1 1 1 1 1 1 1 ...
##  $ ESTEROIDES : logi  FALSE FALSE TRUE TRUE TRUE FALSE ...
##  $ ANTIVIRALES: logi  FALSE FALSE TRUE FALSE FALSE FALSE ...
##  $ FATIGA     : logi  FALSE TRUE FALSE FALSE FALSE TRUE ...
##  $ MALAISE    : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ ANOREXIA   : logi  FALSE FALSE FALSE FALSE FALSE TRUE ...
##  $ HIGgrande  : logi  FALSE FALSE TRUE TRUE TRUE TRUE ...
##  $ HIGfirme   : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ BAZOpalpa  : logi  FALSE FALSE FALSE FALSE FALSE TRUE ...
##  $ ARANIASvac : logi  FALSE FALSE FALSE FALSE FALSE TRUE ...
##  $ ASCITIS    : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ VARICES    : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ HISTIOLOGIA: logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ PC1        : num  -1.044 0.414 -1.565 -0.221 -0.993 ...
##  $ PC2        : num  -0.02431 0.65711 0.09062 -0.97581 0.00983 ...
##  $ PC3        : num  0.6393 0.5692 -0.0427 -0.8245 0.5449 ...
##  $ PC4        : num  -0.472 0.775 -0.509 0.208 -0.11 ...
##  $ PC5        : num  0.27704 -0.24859 0.37027 0.72329 0.00669 ...
##  $ PC6        : num  0.0341 -0.0356 0.624 0.2423 0.0588 ...
##  $ PRONOSTICO : Factor w/ 2 levels "FALLECE","VIVE": 2 2 2 2 2 1 2 2 2 2 ...

Ejercicio 2.

2.a) Utiliza la función preProcess para imputar, mediante el método kNN el conjunto de datos original antes de haberlo imputado. Observas alguna diferencia con el que generamos en la sesión anterior? A qué crees que es debido?

hepatitis = read.csv2(file="hepatitis.csv", header=FALSE, sep = ",", dec = ".",na.string="*")
atributos <- c("EDAD", "SEXO", "ESTEROIDES", "ANTIVIRALES",
               "FATIGA", "MALAISE", "ANOREXIA", "HIGgrande",
               "HIGfirme", "BAZOpalpa", "ARANIASvac", "ASCITIS",
               "VARICES", "BILIRRUBINA", "FOSFATOalc", "SGOT",
               "ALBUMINA", "PROTIME", "HISTIOLOGIA", "PRONOSTICO")
colnames(hepatitis) <- atributos
hepatitis$ANTIVIRALES <- as.logical(hepatitis$ANTIVIRALES)
hepatitis$ESTEROIDES <- as.logical(hepatitis$ESTEROIDES)
hepatitis$FATIGA <- as.logical(hepatitis$FATIGA)
hepatitis$MALAISE <- as.logical(hepatitis$MALAISE)
hepatitis$ANOREXIA <- as.logical(hepatitis$ANOREXIA)
hepatitis$HIGgrande <- as.logical(hepatitis$HIGgrande)
hepatitis$HIGfirme <- as.logical(hepatitis$HIGfirme)
hepatitis$BAZOpalpa <- as.logical(hepatitis$BAZOpalpa)
hepatitis$ARANIASvac <- as.logical(hepatitis$ARANIASvac)
hepatitis$ASCITIS <- as.logical(hepatitis$ASCITIS)
hepatitis$VARICES <- as.logical(hepatitis$VARICES)
hepatitis$HISTIOLOGIA <- as.logical(hepatitis$HISTIOLOGIA)
hepatitis$EDAD = as.integer(hepatitis$EDAD)
hepatitis$SGOT = as.integer(hepatitis$SGOT)
hepatitis$FOSFATOalc = as.integer(hepatitis$FOSFATOalc)
hepatitis$PROTIME = as.integer(hepatitis$PROTIME)
hepatitis$ALBUMINA = as.numeric(hepatitis$ALBUMINA)
hepatitis$BILIRRUBINA = as.numeric(hepatitis$BILIRRUBINA)
library(car)
hepatitis$PRONOSTICO <- as.factor(recode(hepatitis$PRONOSTICO,"0 = 'FALLECE'; 1 = 'VIVE'"))
hepatitis$SEXO <- as.factor(recode(hepatitis$SEXO,"0 = 'FEMENINO'; 1 = 'MASCULINO'"))
#hepatitisknn_function <- preProcess(hepatitis[,1:ncol(hepatitis)-1], method = "knnImpute",k=5)
#method = c("center", "scale", "knnImpute"), 
hepatitisknn_function <- preProcess(hepatitis, method = "knnImpute", k=5 )
str(hepatitisknn_function)
## List of 19
##  $ dim       : int [1:2] 80 20
##  $ bc        : NULL
##  $ yj        : NULL
##  $ et        : NULL
##  $ mean      : Named num [1:6] 41.2 1.43 105.33 85.89 3.82 ...
##   ..- attr(*, "names")= chr [1:6] "EDAD" "BILIRRUBINA" "FOSFATOalc" "SGOT" ...
##  $ std       : Named num [1:6] 12.566 1.212 51.508 89.651 0.652 ...
##   ..- attr(*, "names")= chr [1:6] "EDAD" "BILIRRUBINA" "FOSFATOalc" "SGOT" ...
##  $ ranges    : NULL
##  $ rotation  : NULL
##  $ method    :List of 4
##   ..$ knnImpute: chr [1:6] "EDAD" "BILIRRUBINA" "FOSFATOalc" "SGOT" ...
##   ..$ ignore   : chr [1:14] "SEXO" "ESTEROIDES" "ANTIVIRALES" "FATIGA" ...
##   ..$ center   : chr [1:6] "EDAD" "BILIRRUBINA" "FOSFATOalc" "SGOT" ...
##   ..$ scale    : chr [1:6] "EDAD" "BILIRRUBINA" "FOSFATOalc" "SGOT" ...
##  $ thresh    : num 0.95
##  $ pcaComp   : NULL
##  $ numComp   : NULL
##  $ ica       : NULL
##  $ wildcards :List of 2
##   ..$ PCA: chr(0) 
##   ..$ ICA: chr(0) 
##  $ k         : num 5
##  $ knnSummary:function (x, ...)  
##  $ bagImp    : NULL
##  $ median    : NULL
##  $ data      :'data.frame':  80 obs. of  6 variables:
##   ..$ EDAD       : num [1:80] -0.573 -0.1751 -0.7321 -0.0159 -0.8913 ...
##   ..$ BILIRRUBINA: num [1:80] -0.435 -0.105 -0.353 -0.435 0.637 ...
##   ..$ FOSFATOalc : num [1:80] -0.2 -0.531 -0.899 -0.472 -0.938 ...
##   ..$ SGOT       : num [1:80] -0.646 -0.623 1.819 -0.289 0.648 ...
##   ..$ ALBUMINA   : num [1:80] 0.28 0.894 -0.18 0.127 1.662 ...
##   ..$ PROTIME    : num [1:80] 0.575 1.012 -0.343 -0.431 0.706 ...
##  - attr(*, "class")= chr "preProcess"
#aplicar imputación con predict
#imputación con caret, solo considera las numericas
#hepatitisknn <- predict(hepatitisknn_function,hepatitis[,1:ncol(hepatitis)-1])
#####hepatitisknn <- predict(hepatitisknn_function,hepatitis, )

Nota: A diferencia de la función de imputación knn() del paquete VIM, el metodo de imputación knn de la función preProcess del paquete caret la imputación solo considera a las varaibles numericas; caret no es paquete que este orientado a la imputación, por tanto, solo cubre algunos aspectos básicos.

  • 2.b) Como podríamos generar un conjunto de datos qué solo contenga las componentes principales seleccionadas y el atributo que identifica la clase? Realiza el proceso para los conjuntos de entrenamiento y prueba.
#head(PCATrain[,14:19])
#levels(long$.imp) <- paste("Imputation",1:5)
pcaTrain = PCATrain[,14:19]
pcaTrain$PRONOSTICO = PCATrain[,20]
pcaTest = PCATest[,14:19]
pcaTest$PRONOSTICO = PCATest[,20]
head(pcaTrain)
##          PC1          PC2         PC3        PC4          PC5         PC6
## 1 -1.0442888 -0.024312322  0.63933295 -0.4724346  0.277036683  0.03405830
## 2  0.4142042  0.657114771  0.56919369  0.7747795 -0.248590807 -0.03561095
## 4 -1.5653359  0.090622445 -0.04265012 -0.5094899  0.370274377  0.62396707
## 5 -0.2208923 -0.975814208 -0.82452342  0.2078301  0.723287280  0.24228015
## 6 -0.9931569  0.009825509  0.54486786 -0.1100026  0.006694062  0.05877158
## 7  1.6672995  1.074564817  0.64572201  0.4996187  0.724300088  0.27155828
##   PRONOSTICO
## 1       VIVE
## 2       VIVE
## 4       VIVE
## 5       VIVE
## 6       VIVE
## 7    FALLECE

2.c) Crea un gráfico en el qué se muestren las cuatro primeras componentes principales del conjunto de entrenamiento.

featurePlot(x=pcaTrain[, 1:4], y=pcaTrain$PRONOSTICO, plot="density", scales=list(x=list(relation="free"), y=list(relation="free")), auto.key=list(columns=3))

featurePlot(x=pcaTrain[, 1:4], y=pcaTrain$PRONOSTICO, plot="ellipse", scales=list(x=list(relation="free"), y=list(relation="free")), auto.key=list(columns=3))

2.d) Repite el proceso seguido para el ACP para determinar las componentes independientes, tanto para el conjunto de entrenamiento como el de prueba.

#componentes independientes
hepatitisICA <- preProcess(hepatitis.KnnImp[,1:ncol(hepatitis.KnnImp)-1],
                           method = "ica",thresh = 0.95, n.com=5) #numero de componentes que se quieren elegir
print(hepatitisICA)
## Created from 155 samples and 19 variables
## 
## Pre-processing:
##   - centered (6)
##   - independent component signal extraction (6)
##   - ignored (13)
##   - scaled (6)
## 
## ICA used 5 components
#escala y centra las 6 variables para generar 5 independientes ignora 13
#se necesitan 6 (transformaciones) no son los datos originales sino combinaciones lineales de los originales. componentes para capturar el 95% de la varianza
ICATrain <- predict(hepatitisICA,hepatitisTrain[,1:ncol(hepatitisTrain)-1])
ICATest <- predict(hepatitisICA,hepatitisTest[,1:ncol(hepatitisTest)-1])
#se generan ddf de test y entrenamiento 
ICATrain$PRONOSTICO = hepatitisTrain$PRONOSTICO
ICATest$PRONOSTICO = hepatitisTest$PRONOSTICO
str(ICATrain)
## 'data.frame':    104 obs. of  19 variables:
##  $ SEXO       : Factor w/ 2 levels "FEMENINO","MASCULINO": 2 1 1 1 1 1 1 1 1 1 ...
##  $ ESTEROIDES : logi  FALSE FALSE TRUE TRUE TRUE FALSE ...
##  $ ANTIVIRALES: logi  FALSE FALSE TRUE FALSE FALSE FALSE ...
##  $ FATIGA     : logi  FALSE TRUE FALSE FALSE FALSE TRUE ...
##  $ MALAISE    : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ ANOREXIA   : logi  FALSE FALSE FALSE FALSE FALSE TRUE ...
##  $ HIGgrande  : logi  FALSE FALSE TRUE TRUE TRUE TRUE ...
##  $ HIGfirme   : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ BAZOpalpa  : logi  FALSE FALSE FALSE FALSE FALSE TRUE ...
##  $ ARANIASvac : logi  FALSE FALSE FALSE FALSE FALSE TRUE ...
##  $ ASCITIS    : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ VARICES    : logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ HISTIOLOGIA: logi  FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ ICA1       : num  0.0635 0.6412 0.1987 0.3294 0.2248 ...
##  $ ICA2       : num  0.83 0.434 0.447 -1.262 0.669 ...
##  $ ICA3       : num  0.809 -0.984 1.015 0.923 0.462 ...
##  $ ICA4       : num  0.0216 0.4082 -0.6341 -0.0393 0.1333 ...
##  $ ICA5       : num  0.1306 -0.0975 -0.0449 0.3706 -0.2283 ...
##  $ PRONOSTICO : Factor w/ 2 levels "FALLECE","VIVE": 2 2 2 2 2 1 2 2 2 2 ...

Tanto en ICA como PCA las clases estan muy solapadas, si se encontraran dos variables que genererarn separarción de clases se elegirian, pero en este caso ni PCA ni ICA generan variables que separen muy bien a las clases.
Generar dataFrame de solo las componentes independientes y de su clase:

icaTrain = ICATrain[,14:18]
icaTrain$PRONOSTICO = ICATrain[,19]
icaTest = ICATest[,14:18]
icaTest$PRONOSTICO = ICATest[,19]
str(icaTrain)
## 'data.frame':    104 obs. of  6 variables:
##  $ ICA1      : num  0.0635 0.6412 0.1987 0.3294 0.2248 ...
##  $ ICA2      : num  0.83 0.434 0.447 -1.262 0.669 ...
##  $ ICA3      : num  0.809 -0.984 1.015 0.923 0.462 ...
##  $ ICA4      : num  0.0216 0.4082 -0.6341 -0.0393 0.1333 ...
##  $ ICA5      : num  0.1306 -0.0975 -0.0449 0.3706 -0.2283 ...
##  $ PRONOSTICO: Factor w/ 2 levels "FALLECE","VIVE": 2 2 2 2 2 1 2 2 2 2 ...
  • 2.e) Crea un gráfico en el que se muestren las cuatro primeras componentes principales del conjunto de entrenamiento. ICA
featurePlot(x=icaTrain[, 1:4], y=icaTrain$PRONOSTICO, plot="density", scales=list(x=list(relation="free"), y=list(relation="free")), auto.key=list(columns=3))

featurePlot(x=icaTrain[, 1:4], y=icaTrain$PRONOSTICO, plot="ellipse", scales=list(x=list(relation="free"), y=list(relation="free")), auto.key=list(columns=3))

  • 2.f) Al realizar el análisis ACP o ICA ¿Qué otros métodos de preprocesamiento se han aplicado? Centrado y escalado, evidentemente sin considerar a las variables no númericas.

  • 2.g) Existe alguna variable que presente una varianza cercana a 0?

nearZeroVar(hepatitis.KnnImp)
## integer(0)

No existe alguna variable con varianza cercana a cero.

4. Selección de variables

4.1. Eliminación recursiva de variables

ctrl.rfe <- rfeControl(functions=treebagFuncs, method = "cv", number = 5, returnResamp="final", verbose = TRUE)
subsets <- c(3:19)
#rf.rfe4_e3 <- rfe(PRONOSTICO~., data=hepatitis.KnnImp,sizes=subsets, rfeControl=ctrl.rfe)
#saveRDS(rf.rfe4_e3, "rf.rfe4_e3.rds")
rf.rfe4_e3 = readRDS("rf.rfe4_e3.rds")
rf.rfe4_e3
## 
## Recursive feature selection
## 
## Outer resampling method: Cross-Validated (5 fold) 
## 
## Resampling performance over subset size:
## 
##  Variables Accuracy  Kappa AccuracySD KappaSD Selected
##          3   0.8769 0.6191    0.05442 0.17064         
##          4   0.8963 0.6747    0.04964 0.15846         
##          5   0.9094 0.7163    0.02844 0.08634         
##          6   0.9094 0.7221    0.02844 0.08692         
##          7   0.9159 0.7450    0.03791 0.11960         
##          8   0.9221 0.7620    0.03850 0.11960         
##          9   0.9026 0.6939    0.04790 0.15843         
##         10   0.9221 0.7587    0.03850 0.12190         
##         11   0.9090 0.7249    0.05014 0.14390         
##         12   0.9155 0.7390    0.04592 0.13046         
##         13   0.9155 0.7473    0.05127 0.14817         
##         14   0.9026 0.7037    0.04790 0.14530         
##         15   0.9155 0.7440    0.05127 0.14963         
##         16   0.9219 0.7643    0.05138 0.14747         
##         17   0.9155 0.7380    0.05127 0.14774         
##         18   0.9155 0.7361    0.05127 0.15297         
##         19   0.9286 0.7803    0.03722 0.12034        *
## 
## The top 5 variables (out of 19):
##    PROTIME, ALBUMINA, BILIRRUBINA, ASCITISTRUE, FOSFATOalc
rf.rfe4_e3$optVariable
##  [1] "PROTIME"         "ALBUMINA"        "BILIRRUBINA"    
##  [4] "ASCITISTRUE"     "FOSFATOalc"      "EDAD"           
##  [7] "VARICESTRUE"     "ARANIASvacTRUE"  "FATIGATRUE"     
## [10] "HISTIOLOGIATRUE" "MALAISETRUE"     "SEXOMASCULINO"  
## [13] "SGOT"            "HIGgrandeTRUE"   "BAZOpalpaTRUE"  
## [16] "ANOREXIATRUE"    "ANTIVIRALESTRUE" "HIGfirmeTRUE"   
## [19] "ESTEROIDESTRUE"
rf.rfe4_e3$fit
## 
## Bagging classification trees with 25 bootstrap replications

rf.profile$optVariables, lo que nos permite volver a generar los conjuntos de entrenamiento y prueba, pero sólo con las variables seleccionadas.

sel.cols = c(rf.rfe4_e3$optVariable, "PRONOSTICO")
#hepatitisTrain.sel = hepatitisTrain[,sel.cols]
#hepatitisTest.sel = hepatitisTest[,sel.cols]

Ejercicio 3.

¿Por qué no funcionan las operaciones anteriores? ¿Que hay que hacer para resolverlo? No funciona del todo bien, porque genera variables dummy para aquellas variables que son del tipo lógicas, esto depende del modelo clasificador. Por tanto, al seleccionar una variable dummy como puede ser ESTEROIDESTrue se genera error debido a que esta variable no existe en el conjunto de datos original. Para garantizar que funcione en todos los modelos clasificadores, se recomineda poner a 1 los valores que corresponden a true y a 0 los false. Duda ¿Pasa lo mismo para las variables categoricas?

#Dado que algunos metodos no funcionan correctamente con variables de tipo lógicas se hace el ajuste a entero
test = hepatitis.KnnImp
test$ANTIVIRALES <- as.integer(test$ANTIVIRALES)
test$ESTEROIDES <- as.integer(test$ESTEROIDES)
test$FATIGA <- as.integer(test$FATIGA)
test$MALAISE <- as.integer(test$MALAISE)
test$ANOREXIA <- as.integer(test$ANOREXIA)
test$HIGgrande <- as.integer(test$HIGgrande)
test$HIGfirme <- as.integer(test$HIGfirme)
test$BAZOpalpa <- as.integer(test$BAZOpalpa)
test$ARANIASvac <- as.integer(test$ARANIASvac)
test$ASCITIS <- as.integer(test$ASCITIS)
test$VARICES <- as.integer(test$VARICES)
test$HISTIOLOGIA <- as.integer(test$HISTIOLOGIA)
test$SEXO  = as.integer(as.character(recode(test$SEXO,"'FEMENINO' = 0; 'MASCULINO' = 1")))

Ejercicio 4

Realiza una selección de variables utilizando la función rfe con las funciones treebagFuncs y nbFuncs.

  • 4.a) Qué técnicas de clasificación utilizan? ¿Cuántas variables seleccionan cada técnica?
ctrl.rfe <- rfeControl(functions=treebagFuncs, method = "cv", number = 5, returnResamp="final", verbose = TRUE)
subsets <- c(3:19)
#rf.rfe1 <- rfe(PRONOSTICO~., data=test,sizes=subsets, rfeControl=ctrl.rfe)
#saveRDS(rf.rfe1, "rf.rfe4_treebagFuncs.rds")
rf.rfe4_treebagFuncs = readRDS("rf.rfe4_treebagFuncs.rds")
rf.rfe4_treebagFuncs$fit
## 
## Bagging classification trees with 25 bootstrap replications
length(rf.rfe4_treebagFuncs$optVariable)
## [1] 3
ctrl.rfe <- rfeControl(functions=nbFuncs, method = "cv", number = 5, returnResamp="final", verbose = TRUE)
#rf.rfe_nbFuncs <- rfe(PRONOSTICO~., data=test,sizes=subsets, rfeControl=ctrl.rfe)
#saveRDS(rf.rfe_nbFuncs, "rf.rfe_nbFuncs.rds")
rf.rfe_nbFuncs = readRDS("rf.rfe_nbFuncs.rds")
rf.rfe_nbFuncs$fit$call
## NaiveBayes.default(x = x, grouping = y, usekernel = TRUE, fL = 2)
length(rf.rfe_nbFuncs$optVariable)
## [1] 7
  • 4.b) Compara estos resultados aplicando la misma técnica pero con las técnicas de clasificación svmLinear y rpart (recuerda que esto se debe hacer a través de la opción functions=caretFuncs y el parámetro method de la función rfe.
ctrl.rfe <- rfeControl(functions=caretFuncs, method = "repeatedcv", number = 5, returnResamp="final", verbose = TRUE, allowParallel = TRUE)
#rf.rfe3 <- rfe(PRONOSTICO~., method = "svmLinear", data=test,sizes=subsets, rfeControl=ctrl.rfe)
#saveRDS(rf.rfe3, "rf.rfe_svmLinear.rds")
rf.rfe_svmLinear = readRDS("rf.rfe_svmLinear.rds")
rf.rfe_svmLinear$fit
## Support Vector Machines with Linear Kernel
## Loading required package: kernlab
## Warning: package 'kernlab' was built under R version 3.2.3
## 
## Attaching package: 'kernlab'
## 
## The following object is masked from 'package:ggplot2':
## 
##     alpha
## 155 samples
##  18 predictor
##   2 classes: 'FALLECE', 'VIVE' 
## 
## No pre-processing
## Resampling: Bootstrapped (25 reps) 
## Summary of sample sizes: 155, 155, 155, 155, 155, 155, ... 
## Resampling results
## 
##   Accuracy   Kappa      Accuracy SD  Kappa SD 
##   0.8547249  0.5312478  0.0405601    0.1110223
## 
## Tuning parameter 'C' was held constant at a value of 1
## 
rf.rfe_svmLinear$fit$call
## train.default(x = x, y = y, method = "svmLinear")
length(rf.rfe_svmLinear$optVariable)
## [1] 18

El objeto rf.rfe_svmLinear, nos indica el mejor clasificador y el mejor subconjunto de variables predcitoras que deben ser incorporadas. Dichas variables quedan almacenadas en el objeto rf.rfe_svmLinear$optVariable, lo que nos permite volver a generar los conjuntos de entrenamiento y prueba, pero sólo con las variables seleccionadas.

set.seed(342)
trainIndex <- createDataPartition(test$PRONOSTICO, p = 0.66, list = FALSE, times = 1)
hepatitisTrain <- test[trainIndex,]
hepatitisTest <- test[-trainIndex,]

sel.cols = c(rf.rfe_svmLinear$optVariable, "PRONOSTICO")
hepatitisTrain.sel = hepatitisTrain[,sel.cols]
hepatitisTest.sel = hepatitisTest[,sel.cols]
ctrl.rfe <- rfeControl(functions=caretFuncs, method = "repeatedcv", number = 5, returnResamp="final", verbose = TRUE, allowParallel = TRUE)
#rf.rfe_rpart <- rfe(PRONOSTICO~., method = "rpart", data=test,sizes=subsets, rfeControl=ctrl.rfe)
#saveRDS(rf.rfe_rpart, "rf.rfe_rpart.rds")
rf.rfe_rpart = readRDS("rf.rfe_rpart.rds")
rf.rfe_rpart$fit
## CART 
## 
## 155 samples
##   6 predictor
##   2 classes: 'FALLECE', 'VIVE' 
## 
## No pre-processing
## Resampling: Bootstrapped (25 reps) 
## Summary of sample sizes: 155, 155, 155, 155, 155, 155, ... 
## Resampling results across tuning parameters:
## 
##   cp       Accuracy   Kappa      Accuracy SD  Kappa SD 
##   0.00000  0.8666842  0.6032024  0.05028665   0.1449158
##   0.28125  0.8717041  0.6442335  0.05013984   0.1216461
##   0.56250  0.8346184  0.3598951  0.05412284   0.3320983
## 
## Accuracy was used to select the optimal model using  the largest value.
## The final value used for the model was cp = 0.28125.
rf.rfe_rpart$results[rf.rfe_rpart$bestSubset,]
##   Variables  Accuracy     Kappa AccuracySD    KappaSD
## 6         8 0.8965323 0.6902584  0.0275784 0.08756343

En cuanto a la eliminación recursiva de variables con caretFuncs y en concreto con el método svmLinear de las 19 variables predictoras se requieren 18, para generar un Accuracy de 0.85. Por otra parte, el método rpart sólo requiere de 6 variables predictoras para generar un Accuracy de 0.89, con lo cual, para este problema en concreto hepatitis se simplifica significativamente con el modelo final propuesto con selección de variables con el metodo rpart, si se toma como criterio de slección Accuracy.

4.2. Eliminación de variables por filtros (caret, sbf)

En definitiva, de lo que se trata es de encontrar aquellas variables que presentan difenrecias significativamente estadísticas entre las distintas clases o el resultado.

Ejercicio 5.

Aplica la selección por filtros utilizando las funciones de evaluación disponibles e indica cuantas variables se seleccionan en cada caso. Las funciones que están disponibles son: ldaSBF, rfSBF, nbSBF y nbSBF

ctrl.ranker <- sbfControl(functions = ldaSBF, method = "cv", number = 5)
#rf.ranker_ldaSBF <- sbf(PRONOSTICO~., data=test, sbfControl = ctrl.ranker)
#saveRDS(rf.ranker_ldaSBF, "rf.ranker_ldaSBF.rds")
rf.ranker_ldaSBF = readRDS("rf.ranker_ldaSBF.rds")
#names(rf.ranker_ldaSBF)
rf.ranker_ldaSBF$fit$call
## lda(x = x, y)
length(rf.ranker_ldaSBF$optVariables)
## [1] 13
rf.ranker_ldaSBF$optVariables
##  [1] "EDAD"        "SEXO"        "FATIGA"      "MALAISE"     "BAZOpalpa"  
##  [6] "ARANIASvac"  "ASCITIS"     "VARICES"     "BILIRRUBINA" "FOSFATOalc" 
## [11] "ALBUMINA"    "PROTIME"     "HISTIOLOGIA"
ctrl.ranker <- sbfControl(functions = rfSBF, method = "cv", number = 5)
#rf.ranker_rfSBF <- sbf(PRONOSTICO~., data=test, sbfControl = ctrl.ranker)
#saveRDS(rf.ranker_rfSBF, "rf.ranker_rfSBF.rds")
rf.ranker_rfSBF = readRDS("rf.ranker_rfSBF.rds")
#names(rf.ranker_rfSBF)
rf.ranker_rfSBF$fit$call
## randomForest(x = x, y = y)
length(rf.ranker_rfSBF$optVariables)
## [1] 13
ctrl.ranker <- sbfControl(functions = nbSBF, method = "cv", number = 5)
#rf.ranker_nbSBF <- sbf(PRONOSTICO~., data=test, sbfControl = ctrl.ranker)
#saveRDS(rf.ranker_nbSBF, "rf.ranker_nbSBF.rds")
rf.ranker_nbSBF = readRDS("rf.ranker_nbSBF.rds")
#names(rf.ranker_nbSBF)
rf.ranker_nbSBF$fit$call
## NaiveBayes.default(x = x, grouping = y, usekernel = TRUE, fL = 2)
length(rf.ranker_nbSBF$optVariables)
## [1] 13
ctrl.ranker <- sbfControl(functions = treebagSBF, method = "cv", number = 5)
#rf.ranker_treebagSBF <- sbf(PRONOSTICO~., data=test, sbfControl = ctrl.ranker)
#saveRDS(rf.ranker_treebagSBF, "rf.ranker_treebagSBF.rds")
rf.ranker_treebagSBF = readRDS("rf.ranker_treebagSBF.rds")
#names(rf.ranker_treebagSBF)
rf.ranker_treebagSBF$fit
## 
## Bagging classification trees with 25 bootstrap replications
length(rf.ranker_treebagSBF$optVariables)
## [1] 13

4.3. Eliminación de variables con búsqueda aleatoria.

Selección de variables de tipo wrapper con búsqueda aelatoria: Algoritmos genéticos, a través de la función gafs() y Enfriamiento simulado safs().

#"Muy costoso computacionalmente Disponible en AV para análisis"
rf.ga = readRDS("AV/rf.ga")
#str(rf.ga$internal)
#str(rf.ga$external)
plot(rf.ga)

#enfriamineto simulado
rf.sa = readRDS("AV/rf.sa")
plot(rf.sa)

5. Entrenamiento de clasificadores y ajuste de parámetros.

sel.cols = c(rf.rfe_svmLinear$optVariable, "PRONOSTICO")
hepatitisTrain.sel = hepatitisTrain[,sel.cols]
hepatitisTest.sel = hepatitisTest[,sel.cols]
str(hepatitisTrain.sel)
## 'data.frame':    104 obs. of  19 variables:
##  $ PROTIME    : int  63 74 63 75 84 54 52 78 85 64 ...
##  $ ALBUMINA   : num  3.5 4 4 4 4.4 3.7 3.9 4.9 4.1 4.2 ...
##  $ BILIRRUBINA: num  0.9 0.7 1 0.9 0.7 1 0.9 2.2 0.7 0.9 ...
##  $ ARANIASvac : int  0 0 0 0 0 1 0 0 0 0 ...
##  $ HISTIOLOGIA: int  0 0 0 0 0 0 0 0 1 0 ...
##  $ MALAISE    : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ FOSFATOalc : int  135 96 89 95 85 59 81 57 53 48 ...
##  $ EDAD       : int  50 78 34 34 39 32 41 30 38 22 ...
##  $ ESTEROIDES : int  0 1 1 1 1 1 1 1 1 1 ...
##  $ SGOT       : int  42 32 200 28 48 249 60 144 42 20 ...
##  $ SEXO       : int  0 0 0 0 0 0 0 0 0 1 ...
##  $ ANTIVIRALES: int  0 0 0 0 0 1 1 0 0 1 ...
##  $ HIGfirme   : int  0 0 0 0 1 1 1 1 0 0 ...
##  $ HIGgrande  : int  0 1 1 1 1 1 1 1 1 1 ...
##  $ ANOREXIA   : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ ASCITIS    : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ BAZOpalpa  : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ VARICES    : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ PRONOSTICO : Factor w/ 2 levels "FALLECE","VIVE": 2 2 2 2 2 2 2 2 2 2 ...
str(hepatitisTest.sel)
## 'data.frame':    51 obs. of  19 variables:
##  $ PROTIME    : int  70 80 35 90 100 85 72 46 74 63 ...
##  $ ALBUMINA   : num  4 4 2.8 4.3 3.9 4.4 4 2.9 4.3 4 ...
##  $ BILIRRUBINA: num  1 0.7 1 1 1 1.3 0.8 2 1.2 0.6 ...
##  $ ARANIASvac : int  0 0 1 0 0 0 0 0 0 0 ...
##  $ HISTIOLOGIA: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ MALAISE    : int  0 0 0 0 0 0 0 1 0 0 ...
##  $ FOSFATOalc : int  85 46 127 89 85 78 92 72 102 62 ...
##  $ EDAD       : int  30 31 51 23 30 39 47 38 66 40 ...
##  $ ESTEROIDES : int  0 1 0 1 1 0 0 0 1 0 ...
##  $ SGOT       : int  18 52 60 59 120 30 60 89 53 166 ...
##  $ SEXO       : int  1 0 0 0 0 0 0 0 0 0 ...
##  $ ANTIVIRALES: int  0 1 0 0 0 1 1 0 0 0 ...
##  $ HIGfirme   : int  0 0 0 0 0 1 0 0 0 1 ...
##  $ HIGgrande  : int  0 1 1 1 1 0 1 1 1 1 ...
##  $ ANOREXIA   : int  0 0 1 0 0 0 0 1 0 0 ...
##  $ ASCITIS    : int  0 0 0 0 0 0 0 1 0 0 ...
##  $ BAZOpalpa  : int  0 0 1 0 0 0 0 0 0 0 ...
##  $ VARICES    : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ PRONOSTICO : Factor w/ 2 levels "FALLECE","VIVE": 2 2 1 2 2 2 2 2 2 2 ...

Ejercicio 6.

Utiliza el conjunto de entrenamiento para generar los siguientes clasicadores: C5.0, svmLinear, knn y lda

  • 6.a) Los nombres de cada uno de los modelos deben ser C50Fit, svmFit, knnFit y ldaFit respectivamente
####Bloque preentrenamiento para evitar entrenar en tiempo de Rmarkdown
#fitcontrol <- trainControl(method = "LGOCV",p=.75,number=10,
#  returnResamp = "final",verboseIter=FALSE)
#set.seed(342)
#C50Fit <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                method="C5.0",
#                tuneLength=10,
#                trControl=fitcontrol)
#saveRDS(C50Fit, "C50Fit.rds")
C50Fit = readRDS("C50Fit.rds")
plot(C50Fit)

#summary(C50Fit)
####Bloque preentrenamiento para evitar entrenar en tiempo de Rmarkdown
#fitcontrol <- trainControl(method = "LGOCV",p=.75,number=10,
#  returnResamp = "final",verboseIter=FALSE)
#set.seed(342)
#svmFit <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                method="svmLinear",
#                tuneLength=10,
#                trControl=fitcontrol)
#saveRDS(svmFit, "svmFit.rds")
svmFit = readRDS("svmFit.rds")
svmFit
## Support Vector Machines with Linear Kernel 
## 
## 104 samples
##  18 predictor
##   2 classes: 'FALLECE', 'VIVE' 
## 
## No pre-processing
## Resampling: Repeated Train/Test Splits Estimated (10 reps, 0.75%) 
## Summary of sample sizes: 79, 79, 79, 79, 79, 79, ... 
## Resampling results
## 
##   Accuracy  Kappa     Accuracy SD  Kappa SD 
##   0.932     0.786623  0.06811755   0.2148844
## 
## Tuning parameter 'C' was held constant at a value of 1
## 
####Bloque preentrenamiento para evitar entrenar en tiempo de Rmarkdown
#fitcontrol <- trainControl(method = "LGOCV",p=.75,number=10,
#  returnResamp = "final",verboseIter=FALSE)
#set.seed(342)
#knnFit <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                method="knn",
#                tuneLength=10,
#                trControl=fitcontrol)
#saveRDS(knnFit, "knnFit.rds")
knnFit = readRDS("knnFit.rds")
plot(knnFit)

####Bloque preentrenamiento para evitar entrenar en tiempo de Rmarkdown
#fitcontrol <- trainControl(method = "LGOCV",p=.75,number=10,
#  returnResamp = "final",verboseIter=FALSE)
#set.seed(342)
#ldaFit <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                method="lda",
#                tuneLength=10,
#                trControl=fitcontrol)
#saveRDS(ldaFit, "ldaFit.rds")
ldaFit = readRDS("ldaFit.rds")
ldaFit
## Linear Discriminant Analysis 
## 
## 104 samples
##  18 predictor
##   2 classes: 'FALLECE', 'VIVE' 
## 
## No pre-processing
## Resampling: Repeated Train/Test Splits Estimated (10 reps, 0.75%) 
## Summary of sample sizes: 79, 79, 79, 79, 79, 79, ... 
## Resampling results
## 
##   Accuracy  Kappa      Accuracy SD  Kappa SD 
##   0.92      0.7650469  0.06798693   0.1934717
## 
##
  • 6.b) Crea una tabla en que contenga los resultados de la precisión y el índice Kappa para cada uno de los modelos. Las filas llevarán el nombre de cada uno de los modelos y las columnas el nombre del índice correspondiente.
cols = c("Accuracy","Kappa","AccuracySD","KappaSD")
C50 = C50Fit$results[rownames(C50Fit$bestTune),][,cols]
svm = svmFit$results[rownames(svmFit$bestTune),][,cols]
knn = knnFit$results[rownames(knnFit$bestTune),][,cols]
lda = ldaFit$results[rownames(ldaFit$bestTune),][,cols]
models = rbind.data.frame(C50, svm, knn, lda)
rownames(models) = c("C50","svm","knn","lda")
models
##     Accuracy      Kappa AccuracySD    KappaSD
## C50    0.936 0.79726587 0.03373096 0.11951821
## svm    0.932 0.78662301 0.06811755 0.21488443
## knn    0.804 0.02857143 0.01264911 0.09035079
## lda    0.920 0.76504687 0.06798693 0.19347166

En resumen, parece que el problema de clasificación hepatitis se ajusta mejor a los modelos lineales: svm y lda con lo cual, no hace falta agregar complejidad al clasificador, lo cual es conveniente para el tiempo de computo (comparando sólo contra knn).

Ejercicio 7.

  • 7.a) Vuelve a generar un clasificador C50 pero preprocesando primero aplicando al mismo tiempo un centrado y un escalado (center y scale). Repite la operación aplicando un cambio de rango (range) ¿Encuentras diferencias? ¿Qué opción da mejor resultado?
####Bloque preentrenamiento para evitar entrenar en tiempo de Rmarkdown
#si queremos hacer preprocesamiento de los atributos de nuestro conjunto
#de datos para que estén en el intervalo [-1, 1]: Range
fitcontrol <- trainControl(method = "cv",number=10,
                           returnResamp = "final",verboseIter=FALSE)
#set.seed(342)
#C50FitRange <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                method="C5.0",
#                tuneLength=10,
#                preProc= "range",
#                trControl=fitcontrol)
#saveRDS(C50FitRange, "C50FitRange.rds")
C50FitRange = readRDS("C50FitRange.rds")
plot(C50FitRange)

#preOpc = c("center","scale")
#set.seed(342)
#con preprocesamiento incluido
#C50Fit_CyS <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                method="C5.0",
#                tuneLength=10,
#                preProc= preOpc,
#                trControl=fitcontrol)
#saveRDS(C50Fit_CyS, "C50Fit_CyS.rds")
C50Fit_CyS = readRDS("C50Fit_CyS.rds")
plot(C50Fit_CyS)

En principio no se aprecian diferencias cuando se usa un preprocesamiento de centrado y escalado o range, dado que, en ambos casos la precisión se comporta de forma similar según se varian los valores de los hiperparámetros.

  • 7.a) Vuelve a generar una máquina de soporte de vectores con kernel lineal escalando y centrando los datos ¿Encuentras diferencias respecto a la calculada en el ejercicio anterior? ¿A qué crees que es debido?
#fitcontrol <- trainControl(method = "LGOCV",p=.75,number=10,
#  returnResamp = "final",verboseIter=FALSE)
#preOpc = c("center","scale")
#set.seed(342)
#svmFit_pre <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                method="svmLinear",
#                tuneLength=10,
#                preProc= preOpc,
#                trControl=fitcontrol)
#saveRDS(svmFit_pre, "svmFit_pre.rds")
svmFit_pre = readRDS("svmFit_pre.rds")
svmFit$results
##   C Accuracy    Kappa AccuracySD   KappaSD
## 1 1    0.932 0.786623 0.06811755 0.2148844
svmFit_pre$results
##   C Accuracy    Kappa AccuracySD   KappaSD
## 1 1    0.932 0.786623 0.06811755 0.2148844

No parece haber diferencias, considero que es porque algunos modelos ya incluyen algun preprocesamiento preconfigurado aunque no se le ordene explicitamente.

Ejercicio 8.

  • 8.a) Utiliza el código anterior y extiende el parámetro size para que busque la mejor red neuronales variando el numero de neuronas de la capa intermedia entre 10 y 20 neuronas.
#sizes = 10:20
#nnetGrid_8a <- expand.grid(.size=sizes,
#                        .decay=c(0.5,0.1,0.001,0.0001))
#set.seed(342)
#nnetFit_8a <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                 method="nnet",
#                 tuneGrid=nnetGrid_8a,
#                 trControl=fitcontrol)
#saveRDS(nnetFit_8a, "nnetFit_8a.rds")
nnetFit_8a = readRDS("nnetFit_8a.rds")
plot(nnetFit_8a)

  • 8.b) Repite el proceso anterior pero realizando un centrado y escalado de los datos.
#preOpc = c("center","scale")
#nnetGrid_8b <- expand.grid(.size=sizes,
#                           .decay=c(0.5,0.1,0.001,0.0001))
#set.seed(342)
#nnetFit_8b <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                    method="nnet",
#                    preProc= preOpc,
#                    tuneGrid=nnetGrid_8b,
#                    trControl=fitcontrol)
#saveRDS(nnetFit_8b, "nnetFit_8b.rds")
nnetFit_8b = readRDS("nnetFit_8b.rds")
plot(nnetFit_8b)

Para el caso de las redes neuronales parece que el clasificador en si, no incluye un preprocesamiento por defecto, con lo cual, parece que al aplicar un preprocesamiento de centrado y escalado da mejores valores de precisión para el problema de hepatitis.

  • 8.c) Elimina el grid definido, e introduce el parámetro ¿Sobre qué valores realiza la búsqueda?
#preOpc = c("center","scale")
#set.seed(342)
#nnetFit_8c <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                    method="nnet",
#                    preProc= preOpc,
#                    tuneLength=10,
#                    trControl=fitcontrol)
#saveRDS(nnetFit_8c, "nnetFit_8c.rds")
nnetFit_8c = readRDS("nnetFit_8c.rds")
plot(nnetFit_8c)

Realiza el entrenamiento con 1, 3, 5, 7, 9, 11, 13, 15, 17, 19. (10 valores diferentes de size), combinado con diez valores distintos de decay entre 0 y 0.1

  • 8.d) Entrena una máquina de soporte de vectores lineal haciendo que la función train() varíe el parámetro C entre los siguientes valores: 0.25, 0.5, 1, 2, 4, 8, 16, 32 y 64 ¿Con qué valor se obtiene el mejor resultado?.
preOpc = c("center","scale")
svmGrid_8d <- expand.grid(.C=c(0.1,0.25,0.5,1,2,4,8,16,32,64))
#set.seed(342)
#svmfit_8d_pre <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                     method="svmLinear",
#                     tuneGrid=svmGrid_8d,
#                     preProc= preOpc,
#                     trControl=fitcontrol)
#saveRDS(svmfit_8d_pre, "svmfit_8d_pre.rds")
svmfit_8d_pre = readRDS("svmfit_8d_pre.rds")
plot(svmfit_8d_pre)

svmfit_8d_pre$bestTune
##   C
## 4 1

Ejercicio 9.

  • 9.a) Genera predicciones de forma conjunta para los modelos que presenten mejor precisión.
#Los modelos dque se confrontan deberán estar entrenados bajo
#las mismas configuraciones de trainControl
#, C50=C50Fit_CyS esta definido diferente el train control
#fitcontrol <- trainControl(method = "LGOCV",p=.75,number=10,
#fitcontrol <- trainControl(method = "cv",number=10,
models <- list(SVM=svmfit_8d_pre,NNET=nnetFit_8b)
hepatitis.predict2 <- predict(models,newdata = hepatitisTest.sel)
## Loading required package: nnet
hepatitis.predict3 <- extractPrediction(models,
                                        testX = hepatitisTest.sel[,-ncol(hepatitisTest.sel)],
                                        testY = hepatitisTest.sel[,ncol(hepatitisTest.sel)])
#solo test
hepatitis.predict4 <- subset(hepatitis.predict3, dataType == "Test")
  • 9.a) Los objetos devueltos por la función pueden visualizarse a través de la función plotObsVsPred(). Visualiza los resultados para el objeto obtenido en el apartado anterior y explica qué significa dicha gráfica
plotObsVsPred(hepatitis.predict3)

Permite hacer una comparación de los modelos en cuanto a su precisión tanto con los datos de entrenamiento, como con los datos de test, para este problema de clasificación los dos modelos presentaron mejores valores de precisión con los datos de entrenamiento. En cuanto a su poder predictivo con datos de test, nnet es mejor que SVM Linear.

6. Área bajo la curva: ROC

En el caso de problemas de clasificiación en los que los objetos sólo pueden pertenecer a dos clases, se suele utilizar como medida de eficiencia del clasificador el área bajo la curva (ROC).

Ejercicio 10.

Calcula las probabiliadades de clase y la curva ROC para 2 de los tres modelos de clasificación generados en los ejercicios anteriores.

fitcontrolROC <- trainControl(method = "cv",number=10,
                              returnResamp = "final",
                              summaryFunction = twoClassSummary,
                              classProbs=TRUE,
                              verboseIter=TRUE)
#set.seed(342)
#svmFitROC <- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                   method="svmLinear",
#                   tuneLength=10,
#                   preProcess = c("center", "scale"),
#                   trControl=fitcontrolROC,
#                   metric = "ROC")
#saveRDS(svmFitROC, "svmFitROC.rds")
svmFitROC = readRDS("svmFitROC.rds")
svmProb <- predict(svmFitROC,newdata = hepatitisTest.sel, type = "prob")
head(svmProb)
##       FALLECE      VIVE
## 1 0.004177398 0.9958226
## 2 0.051582792 0.9484172
## 3 0.320459232 0.6795408
## 4 0.007569663 0.9924303
## 5 0.010336652 0.9896633
## 6 0.019473672 0.9805263
library(pROC)
## Warning: package 'pROC' was built under R version 3.2.3
## Type 'citation("pROC")' for a citation.
## 
## Attaching package: 'pROC'
## 
## The following objects are masked from 'package:stats':
## 
##     cov, smooth, var
svmROC <- roc(hepatitisTest$PRONOSTICO, svmProb$VIVE, dataGrid = TRUE,
  gridLength = 100)
svmROC
## 
## Call:
## roc.default(response = hepatitisTest$PRONOSTICO, predictor = svmProb$VIVE,     dataGrid = TRUE, gridLength = 100)
## 
## Data: svmProb$VIVE in 10 controls (hepatitisTest$PRONOSTICO FALLECE) < 41 cases (hepatitisTest$PRONOSTICO VIVE).
## Area under the curve: 0.7634
plot(svmROC)

## 
## Call:
## roc.default(response = hepatitisTest$PRONOSTICO, predictor = svmProb$VIVE,     dataGrid = TRUE, gridLength = 100)
## 
## Data: svmProb$VIVE in 10 controls (hepatitisTest$PRONOSTICO FALLECE) < 41 cases (hepatitisTest$PRONOSTICO VIVE).
## Area under the curve: 0.7634
#set.seed(342)
#nnetFitROC<- train(PRONOSTICO~.,data=hepatitisTrain.sel,
#                   method="nnet",
#                   tuneLength=10,
#                   preProcess = c("center", "scale"),
#                   trControl=fitcontrolROC,
#                   metric = "ROC")
#saveRDS(nnetFitROC, "nnetFitROC.rds")
nnetFitROC = readRDS("nnetFitROC.rds")
nnetProb <- predict(nnetFitROC,newdata = hepatitisTest.sel, type = "prob")
head(nnetProb)
##         FALLECE      VIVE
## 1  0.0001419806 0.9998580
## 4  0.0004257661 0.9995742
## 7  0.6636062335 0.3363938
## 8  0.0001407746 0.9998592
## 10 0.0001499132 0.9998501
## 11 0.0001033157 0.9998967
nnetROC <- roc(hepatitisTest$PRONOSTICO, nnetProb$VIVE, dataGrid = TRUE,
  gridLength = 100)
nnetROC
## 
## Call:
## roc.default(response = hepatitisTest$PRONOSTICO, predictor = nnetProb$VIVE,     dataGrid = TRUE, gridLength = 100)
## 
## Data: nnetProb$VIVE in 10 controls (hepatitisTest$PRONOSTICO FALLECE) < 41 cases (hepatitisTest$PRONOSTICO VIVE).
## Area under the curve: 0.8341
plot(nnetROC)

## 
## Call:
## roc.default(response = hepatitisTest$PRONOSTICO, predictor = nnetProb$VIVE,     dataGrid = TRUE, gridLength = 100)
## 
## Data: nnetProb$VIVE in 10 controls (hepatitisTest$PRONOSTICO FALLECE) < 41 cases (hepatitisTest$PRONOSTICO VIVE).
## Area under the curve: 0.8341

7. Evaluación de los conjuntos de prueba

Con la función confusionMatrix podemos crear un objeto que contiene todas las medidas incluidas en una matriz de confusión.

#matriz de confusión de las pariciones de entrenamiento con LOGCV, no aplica para boot o LOOCV
hepatitis.conf <- confusionMatrix(nnetFit_8b)
hepatitis.conf
## Repeated Train/Test Splits Estimated (10 reps, 0.75%) Confusion Matrix 
## 
## (entries are percentages of table totals)
##  
##           Reference
## Prediction FALLECE VIVE
##    FALLECE    17.2  4.8
##    VIVE        2.8 75.2

Ejercicio 11.

Elige dos modelos de los anteriormente generados (ambos tienen que haber sido evaluados con la misma técnica)

  • 11.a) Calcula la matriz de confusión y los principales índices de eficiencia de ambos modelos.
#otra forma de calcular una matriz de confusión es la de utilizar directamente 
#las clases predichas para el conjunto de prueba
hepatitis.pred.nnet <- predict(nnetFit_8b,hepatitisTest.sel)
hepatitis.conf.nnet <- confusionMatrix(hepatitis.pred.nnet,
hepatitisTest.sel[,ncol(hepatitisTest.sel)])
hepatitis.conf.nnet
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction FALLECE VIVE
##    FALLECE       7    4
##    VIVE          3   37
##                                          
##                Accuracy : 0.8627         
##                  95% CI : (0.7374, 0.943)
##     No Information Rate : 0.8039         
##     P-Value [Acc > NIR] : 0.1911         
##                                          
##                   Kappa : 0.5805         
##  Mcnemar's Test P-Value : 1.0000         
##                                          
##             Sensitivity : 0.7000         
##             Specificity : 0.9024         
##          Pos Pred Value : 0.6364         
##          Neg Pred Value : 0.9250         
##              Prevalence : 0.1961         
##          Detection Rate : 0.1373         
##    Detection Prevalence : 0.2157         
##       Balanced Accuracy : 0.8012         
##                                          
##        'Positive' Class : FALLECE        
## 
hepatitis.pred.svm <- predict(svmfit_8d_pre,hepatitisTest.sel)
hepatitis.conf.svm <- confusionMatrix(hepatitis.pred.svm,
hepatitisTest.sel[,ncol(hepatitisTest.sel)])
hepatitis.conf.svm
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction FALLECE VIVE
##    FALLECE       5    4
##    VIVE          5   37
##                                          
##                Accuracy : 0.8235         
##                  95% CI : (0.6913, 0.916)
##     No Information Rate : 0.8039         
##     P-Value [Acc > NIR] : 0.4441         
##                                          
##                   Kappa : 0.4183         
##  Mcnemar's Test P-Value : 1.0000         
##                                          
##             Sensitivity : 0.50000        
##             Specificity : 0.90244        
##          Pos Pred Value : 0.55556        
##          Neg Pred Value : 0.88095        
##              Prevalence : 0.19608        
##          Detection Rate : 0.09804        
##    Detection Prevalence : 0.17647        
##       Balanced Accuracy : 0.70122        
##                                          
##        'Positive' Class : FALLECE        
##

Comparando las matrices de confusión para cada modelo, considero que la que arroja nnet es mejor, dado que va a generar menos casos en los que predice que VIVE cuando en realidad FALLECE. Y evidentemente, con una mejor tasa de acierto.

  • 11.b) ¿Cuáles son los atributos del objeto generado por la función confusionMatrix()?
names(hepatitis.conf.nnet)
## [1] "positive" "table"    "overall"  "byClass"  "dots"
  • 11.c) ¿Cómo podríamos generar una tabla con los nombres de los modelos en las filas y el de los indicadores en las columnas?
#what: data to conver to matrix. Either "xtabs", "overall" or "classes"
overall.nnet = t(as.matrix(hepatitis.conf.nnet,what = "overall"))
overall.svm = t(as.matrix(hepatitis.conf.svm,what = "overall"))
tablesCM = rbind.data.frame(overall.nnet, overall.svm)
rownames(tablesCM) = c("NNET","SVM")
tablesCM
##       Accuracy     Kappa AccuracyLower AccuracyUpper AccuracyNull
## NNET 0.8627451 0.5804935     0.7374485     0.9429882    0.8039216
## SVM  0.8235294 0.4182510     0.6912740     0.9159913    0.8039216
##      AccuracyPValue McnemarPValue
## NNET      0.1910541             1
## SVM       0.4440568             1

8. Comparación de diferentes modelos

Ejercicio 12

  • 12.a) Remuestrea los dos modelos y muestra un resumen con las medidas utilizadas para su evaluación.
#agrupar todos los resultados de remuestreo utilizando la funcion resamples()
hepatitis.resample <- resamples(models)
summary(hepatitis.resample)
## 
## Call:
## summary.resamples(object = hepatitis.resample)
## 
## Models: SVM, NNET 
## Number of resamples: 10 
## 
## Accuracy 
##      Min. 1st Qu. Median  Mean 3rd Qu. Max. NA's
## SVM  0.84    0.88   0.94 0.932       1    1    0
## NNET 0.80    0.88   0.92 0.924       1    1    0
## 
## Kappa 
##        Min. 1st Qu. Median   Mean 3rd Qu. Max. NA's
## SVM  0.5000  0.5867 0.8169 0.7866       1    1    0
## NNET 0.4898  0.6087 0.7674 0.7770       1    1    0
  • 12.b) Genera algunas gráficas que muetren el resultado de la operación anterior.
bwplot(hepatitis.resample, main="bwplot")

densityplot(hepatitis.resample, metric = "Accuracy",auto.key=TRUE)

  • 12.c) Aplica el test de estudent e interpreta el resultado.

Dado que los modelos obtenidos han sido generados a partir de los mismos datos de entrenamiento, podemos realizar algún tipo de inferencia estadística con ellos. Por ejemplo, podemos calcular la diferencia entre modelos y después realizar la prueba t de Student para evaluar la hipótesis nula de que no hay diferencia entre los distintos modelos.

difValues <- diff(hepatitis.resample)
summary(difValues)
## 
## Call:
## summary.diff.resamples(object = difValues)
## 
## p-value adjustment: bonferroni 
## Upper diagonal: estimates of the difference
## Lower diagonal: p-value for H0: difference = 0
## 
## Accuracy 
##      SVM    NNET 
## SVM         0.008
## NNET 0.1679      
## 
## Kappa 
##      SVM    NNET    
## SVM         0.009579
## NNET 0.5433

Estadisticamente se comportan de forma similar ambos modelos, por tanto se desacarta la hipotesis nula de que los modelos oferecen diferencias significativas con base en la precisión.

#Los resultados sobre las diferencias entre métodos se pueden mostrar de forma
#gráfica a través de diferentes tipos de gráficas.
densityplot(difValues, metric = "Accuracy",auto.key=TRUE, pch = "|")

levelplot(difValues, what="diferences")

La diferencia de las medias sería cercana a cero en caso de que los modelos sean estadísticamente similares. En este caso de comparación entorno (0.008 y 0.009579) estan coincidiendo mas,diagonal superior de la matriz. t-test se acepta la hipotesis nula aunque hay una tendencia a mejor comportamiento en la SVM. T-test se comprueba que la hipotesis nula es verdadera (que no hay diferencia entre los distintos modelos) en este caso dice que si porque en menor a 0.5.

  • 12.d) La función compare_models() también realiza el mismo test. Aplícala y analiza los resultados.
#compare_models is a shorthand function to compare two models using a single metric. It returns the results of t.test on the differences.
compare_models(a = nnetFit_8b, b = svmfit_8d_pre)
## 
##  One Sample t-test
## 
## data:  x
## t = -1.5, df = 9, p-value = 0.1679
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.020064838  0.004064838
## sample estimates:
## mean of x 
##    -0.008

Es un método “rápido” para comparar dos modelos. En lugar de hacer las llamdas manuales a resamples y difValues. Indica que no hay diferencias significativas estadísticamente en cuanto a diferencias precisión en estos dos modelos clasificadores.