Visualizacion del Ajuste de Parámetros
En la entrada anterior, vimos como ajustar los parámetros de gbm probando todas las combinaciones de valores de parámetros dentro de cierto rango. Esto se hacía con una búsqueda en rejilla (grid), haciendo en cada punto de la reja una validación cruzada (es decir, división en entrenamiento / validación múltiples veces (5 en este caso). Una vez que se habían determinado los mejores parámetros, se vuelve a construir el modelo con todos los datos y se evalúa con un conjunto de test adicional.
Pero en ocasiones, conviene visualizar como responde el algoritmo frente al cambio de parámetros, para ver si merecería comprobar otras combinaciones de parámetros o no. Esto es lo que vamos a hacer en esta entrada, usando los paquetes lattice y latticeExtra para visualizar resultados.
En primer lugar, cargamos los resultados para todos los parámetros que teníamos guardados de la entrada anterior y volvemos a calcular las medias para los 5 folds. Cuidado, porque parametros es un data.table y no un data.frame (esto último sería lo más habitual).
library(data.table)
library(lattice)
library(latticeExtra)
library(reshape2)
options(width=500)
semilla=1
nf = paste0("parametrosGBM-",semilla,".txt")
parametros = fread(nf)
nombres_parametros = c("interaction.depth", "shrinkage", "n.minobsinnode", "nFold")
claves = c(nombres_parametros, "n.trees")
parametros = unique(parametros)
# Ahora calculamos la media de los 5 folds
parametros = parametros[, list(Train=mean(Train), Validacion=mean(Validacion), Time=mean(Time)), by=eval(claves[claves!="nFold"])]Recordemos que la mejor combinación de parámetros era:
(p = parametros[which.min(Validacion)])## interaction.depth shrinkage n.minobsinnode n.trees Train Validacion Time
## 1: 4 0.05 10 100 7.227 10.51 1.645p = p[,1:(ncol(p)-3), with=FALSE] # Quitamos los valores de Train, Validación y Time que no son parámetros realmenteEs interesante también ver cuales son las combinaciones que ocupan los primeros puestos y si hay muchas diferencias entre ellas en términos de validación. Parece claro que los mejores resultados se obtienen con n.minobsinnode=10 y n.trees=100. Sin embargo, interaction.depth no parece ser tan determinante.
parametros[order(Validacion)][1:20]## interaction.depth shrinkage n.minobsinnode n.trees Train Validacion Time
## 1: 4 5e-02 10 100 7.227 10.51 1.6453
## 2: 6 5e-02 10 100 6.965 10.52 2.3046
## 3: 8 5e-02 10 100 6.965 10.52 2.3611
## 4: 10 5e-02 10 100 6.965 10.52 2.2805
## 5: 12 5e-02 10 100 6.965 10.52 2.2801
## 6: 14 5e-02 10 100 6.965 10.52 2.2790
## 7: 1 5e-04 10 9500 9.274 10.59 0.4785
## 8: 1 5e-04 10 9000 9.323 10.60 0.4785
## 9: 1 5e-04 10 8500 9.372 10.60 0.4785
## 10: 2 5e-02 10 100 8.257 10.60 0.8832
## 11: 1 5e-04 10 10000 9.234 10.60 0.4785
## 12: 1 1e-02 10 500 9.230 10.60 0.4735
## 13: 1 5e-04 10 10500 9.196 10.61 0.4785
## 14: 1 5e-04 10 8000 9.427 10.61 0.4785
## 15: 1 1e-03 5 4500 9.276 10.61 0.5173
## 16: 1 5e-04 5 9000 9.277 10.61 0.4784
## 17: 1 5e-04 5 9500 9.224 10.61 0.4784
## 18: 1 5e-04 5 8500 9.330 10.61 0.4784
## 19: 1 1e-03 10 4500 9.322 10.62 0.4810
## 20: 1 5e-04 10 11000 9.161 10.62 0.4785A continuación mostramos como cambian los resultados en Test para cada valor del parámetro. Como hay 4 parámetros (n.trees, shrinkage, interaction.depth, y n.minobsinnode), realmente necesitaríamos un plot en 5 dimensiones (4 para los parámetros y uno para el valor en Validación). Como esto no es posible se ha optado hacer plots bidimensionales, donde para cada valor de cada parámetro se visualizan los errores de Validación para los valores de los otros parámetros. La línea azul representa la tendencia media y la roja los valores mínimos. Podemos ver que no siempre la tendencia media no siempre es igual a la tendencia mínima, siendo esta última la que nos interesa, puesto que queremos seleccionar los parámetros que tengan el mínimo error en validación. Aparte de la gran variabilidad para cada valor de cada parámetro, no destacaremos nada más en estas gráficas.
xyplot(Validacion~value|variable, data=melt(parametros, id.vars=tail(names(parametros),3)), scales="free", type=c("p", "smooth"))+
as.layer(xyplot(Validacion~value|variable, data=melt(parametros, id.vars=tail(names(parametros),3))[,list(Validacion=min(Validacion)), by=c("variable","value")], scales=list(x="free"), type=c("b"), col="red"))
Debido a las grandes variabilidades en validación (eje y), en las figuras anteriores no se aprecia demasiado las curvas medias y mínimas. Vamos a representar a continuación sólo la evolución de la media y el mínimo. Podemos ver por ejemplo que aunque en media parece conveniente tener muchos árboles (n.trees), los mínimos aconsejan usar pocos árboles (100 en este caso). De igual manera, para el parámetro n.minobsinnode en media parece recomendable usar valores lo más pequeño posibles, en mínimo hay un óptimo en 10. En shrinkage también hay un comportamiento dispar entre media y mínimo para valores pequeños de este parámetro. En cualquier caso, insisto en que los valores que interesan son los mínimos, que se verán más claramente en las gráficas siguientes.
dt1 =melt(parametros, id.vars=tail(names(parametros),3))[,list(Validacion=mean(Validacion)), by=c("variable","value")]
dt2 = melt(parametros, id.vars=tail(names(parametros),3))[,list(Validacion=min(Validacion)), by=c("variable","value")]
dt1 = cbind(dt1, tipo="medias")
dt2 = cbind(dt2, tipo="minimos")
dtn = rbind(dt1,dt2)
dtn = dcast(data=dtn, formula = variable+value~tipo, value.var = "Validacion")
xyplot(medias+minimos~value|variable, data=dtn, scales=list(x="free"), type=c("b"), auto.key=TRUE)
Para ver más claro el comportamiento de cada parámetro, vamos a representar para cada valor de cada parámetro, los valores mínimos en validación para ese valor de ese parámetro (es lo que representaría la recta roja anterior).
(q=parametros[shrinkage==0.1])## interaction.depth shrinkage n.minobsinnode n.trees Train Validacion Time
## 1: 1 0.1 2 100 8.440 11.09 0.4853
## 2: 1 0.1 2 500 6.065 12.35 0.4853
## 3: 1 0.1 2 1000 4.511 13.02 0.4853
## 4: 1 0.1 2 1500 3.560 13.41 0.4853
## 5: 1 0.1 2 2000 2.861 13.71 0.4853
## ---
## 2228: 14 0.1 30 13000 4.369 17.47 0.7496
## 2229: 14 0.1 30 13500 4.346 17.57 0.7496
## 2230: 14 0.1 30 14000 4.288 17.58 0.7496
## 2231: 14 0.1 30 14500 4.206 17.56 0.7496
## 2232: 14 0.1 30 15000 4.133 17.72 0.7496Por ejemplo, para el parámetro 0.1 de shrinkage hay varios resultados en validación, uno para cada combinación de 0.1 con el resto de valores para los otros tres parámetros. Es decir, hay 2232 combinaciones que incluyen shrinkage=0.1 y para visualizar se escogerá aquella que tenga el mínimo. Vamos a ver un ejemplo para shrinkage=0.1.
q[which.min(Validacion)]## interaction.depth shrinkage n.minobsinnode n.trees Train Validacion Time
## 1: 1 0.1 10 100 8.762 10.8 0.4852Así pues, las siguientes figuras recogen el mejor comportamiento de los parámetros. La línea azul muestra los mínimos en validación. El cuadrado y la línea roja indican donde se produce el mínimo en validación, que lógicamente coincide con los valores óptimos que ya habíamos visto anteriormente. La línea rosa representa los valores de entrenamiento.
resultados_claves = list()
for(clave in claves[claves!="nFold"]){
resultados_claves[[clave]] = parametros[,.SD[which.min(Validacion)], by=clave, .SDcols=c("Train","Validacion", "Time")]
resultados_claves[[clave]][,variable:=clave]
setnames(resultados_claves[[clave]], 1, "valor")
}
resultados = rbindlist(resultados_claves)
resultadosMin = resultados[,.SD[which.min(Validacion)],by="variable"]
resultadosMinLinea = copy(resultados)[,Validacion:=min(Validacion),by="variable"]
xyplot(Validacion+Train~valor|variable, resultados, scales = "free", type=c("b"), auto.key=TRUE)+
as.layer(xyplot(Validacion~valor|variable, data=resultadosMin, pch=0, cex=2))+
as.layer(xyplot(Validacion~valor|variable, data=resultadosMinLinea, type="l", col="red"))
Es interesante ver que a medida que se incrementa la interaction.depth (profundidad del árbol, aproximadamente), el error de entrenamiento baja mucho, pero no ocurre lo mismo en validacion, el cual alcanza el mínimo en 4. La figura anterior no muestra demasiado bien las curvas de mínimos en validación porque también estabamos visualizando las curvas de entrenamiento y el rango de valores en el eje y es grande. Así que ya por fin, vamos a mostrar ahora sólo las de mínimos en validación que es la que realmente nos interesa:
xyplot(Validacion~valor|variable, resultados, scales = "free", type=c("b"), auto.key=TRUE)+
as.layer(xyplot(Validacion~valor|variable, data=resultadosMin, pch=0, cex=2))+
as.layer(xyplot(Validacion~valor|variable, data=resultadosMinLinea, type="l", col="red"))
El objetivo de visualizar estas curvas de mínimos es determinar si merece la pena probar otras combinaciones de parámetros, si detectamos alguna región donde podría haber mejores mínimos. ¿Es así?. El comportamiento de n.trees es peculiar porque inicialmente incrementar el número de árboles empeora el comportamiento en validación y a partir de cierto valor, empieza a disminuir. Así y todo, el mejor valor se alcanza en 100 árboles, con lo que podría ser interesante probar con valores menores a 100. En el caso de shrinkage, hay una tendencia descendiente hasta que se alcanza el mínimo en 0.05 y comienza a subir para valores menores a ese. Pero podría ser interesante ver si en la estrecha región ligeramente menor que 0.05 hay valores un poco mejores. También podrían encontrarse mejores valores alrededor de n.minobsinnode=10 y interaction.depth=4.
Las figuras anteriores tienen distintas escalas en el eje y (el error de validación). Para poder comparar de una manera absoluta la influencia de unos y otros parámetros, vamos a dejar libre sólo la escala del eje x, pero la del y será la misma para todas los gráficas.
xyplot(Validacion~valor|variable, resultados, scales = list(x="free"), type=c("b"), auto.key=TRUE)+
as.layer(xyplot(Validacion~valor|variable, data=resultadosMin, pch=0, cex=2))+
as.layer(xyplot(Validacion~valor|variable, data=resultadosMinLinea, type="l", col="red"))
Podemos ver que el parámetro que más influencia tiene es shrinkage seguido de n.minobsinnode y n.trees. Puede ser interesante por tanto ver la eficacia conjunta de shrinkage con n.trees y shrinkage con n.minobsinnode. Las siguientes cuatro figuras ilustran esto. La primera indica donde está el óptimo. La segunda y tercera dan una idea un poco mas clara con sólo tres colores simples, mostrando los mejores valores en rosa, los siguientes en blanco y los malos en azul. Esta segunda figura muestra que cuando shrinkage es grande, se requieren pocos árboles, pero cuando es más pequeño, hacen falta una gran cantidad de árboles para obtener buenos resultados. También se puede ver que cuando shrinkage es alto, sólo se pueden obtener buenos resultados con pocos árboles. Pero cuando shrinkage es bajo, se pueden obtener buenos resultados para grandes rangos de valores de n.trees. Por ejemplo, para shrinkage=0.001, un número de árboles entre 4500 y 11000 árboles va a dar buenos resultados. La última figura muestra la gradación completa de buenos a malos valores (azul, blanco, rosa).
#library(manipulate)
library(gridExtra)## Loading required package: gridplotea = function(nivel, v1, v2, salida){
param = parametros[,.SD[which.min(Validacion)], by=c(v1, v2), .SDcol="Validacion"]
rango = param[,range(Validacion)]
param = param[,v:=scale(Validacion,rango[1], diff(rango))]
if(nivel>=0){
param[,disc:=ifelse(v<=nivel/2,0, ifelse(v<=nivel, 0.5, 1.0))];
levelplot(disc ~ factor(get(v1))+factor(get(v2)), param, xlab=v1, ylab=v2)
} else {
levelplot(v ~ factor(get(v1))+factor(get(v2)), param, xlab=v1, ylab=v2)
}
}
#manipulate(plotea(nivel, v1, v2), nivel=slider(0,1))
v1 = "shrinkage"
v2 = "n.trees"
plotea(0,v1,v2)
plotea(0.044,v1,v2)
plotea(0.044*2,v1,v2)
plotea(-1, v1, v2)
A continuación se muestran las figuras de shrinkage vs. n.minobsinnode.
v1 = "shrinkage"
v2 = "n.minobsinnode";
plotea(0,v1,v2)
plotea(0.044,v1,v2)
plotea(0.044*2,v1,v2)
plotea(-1, v1, v2)
A continuación se muestran las figuras de shrinkage vs. interaction.depth.
v1 = "shrinkage"
v2 = "interaction.depth"
plotea(0,v1,v2)
plotea(0.044,v1,v2)
plotea(0.044*2,v1,v2)
plotea(-1, v1, v2)
Por último se muestran las figuras de n.trees vs. interaction.depth.
v1 = "n.trees"
v2 = "interaction.depth"
plotea(0,v1,v2)
plotea(0.044,v1,v2)
plotea(0.044*2,v1,v2)
plotea(-1, v1, v2)
Por último, y ya que hemos almacenado los tiempos que cuesta construir los modelos, vamos a ver que parámetro tiene más coste computacional. Obviaremos n.trees porque aunque está claro que cuantos más árboles tiene un modelo, más se tarda en construir el modelo, como siempre hemos construido el modelo con el máximo número de árboles, no vamos a poder saber cuanto se hubiera tardado en construir modelos con menor número de árboles. Vemos que el tiempo crece con interaction.depth y decrece con n.minobsinode. shrinkage no tiene mucha influencia en el tiempo, pero esto es debido a que siempre construimos modelos con el máximo número de árboles. Recordemos que pequeños valores de shrinkage necesitan de mayor número de árboles.
xyplot(Time ~ value|variable, melt(parametros[,!"n.trees", with=FALSE], id.vars=tail(names(parametros),3)), type=c("p"), scales=list(x="free"))
bwplot(Time~ factor(value)|variable , melt(parametros[,!"n.trees", with=FALSE], id.vars=tail(names(parametros),3)), scales=list(x="free"))