---
title: "AjusteHiperRF"
author: "Ricardo Aler"
date: "1/8/2019"
output:
  pdf_document: default
  html_notebook: default
  word_document: default
---

## Ajuste hiper-parámetros Random Forests

```{r}
library(mlr)
library(mlbench)
library(mlrHyperopt) # Esta librería ayuda en la obtención de hiper-parámetros 
library(ggplot2)     # Librería para hacer plots
library(mlrMBO)      # Para hacer ajuste de hiper-parámetros con model-based optimization

# Esto es para que no haya demasiados mensajes de información en el documento
configureMlr(show.info = FALSE, on.learner.warning = "quiet")

data(BostonHousing)
# Primero, definimos una tarea de regresión, cuyos datos están contenidos en el data.frame BostonHousing y cuyo target es "mdev"
task_bh <- makeRegrTask(data= BostonHousing, target="medv")

# Segundo, definimos el nombre delgoritmo de aprendizaje ("learner"). Hay muchos de random forests, probemos ranger.
learner_name <- "regr.ranger"

# Tercero, vemos que hiper-parámetros ajustables tiene
filterParams(getParamSet(learner_name), tunable = TRUE)

# Nuestro learner inicial va a ser una SVM con kernel gausiano (radial). Lo definimos así:

learner_rf <- makeLearner(learner_name)


```

Vemos que tiene muchos hiper-parámetros, pero sabemos que los más importantes son el número de árboles **num.trees** y el **mtry**. También **min.nodesize**,  ue controla la profundidad, aunque de manera indirecta. Más información aquí [https://www.rdocumentation.org/packages/ranger/versions/0.10.1]. Más adelante tomaremos una decisión acerca de cuales ajustar.

Es interesante saber que podemos conseguir ayuda sobre el método y de sus parámetros así:

```{r, eval=FALSE}
helpLearner(learner_name)
helpLearnerParam(learner_name, "num.random.splits")
```


```{r}

# Ahora habría que definir un espacio de búsqueda para los hiper-parámetros
# Una opción es ir a la página web # http://mlrhyperopt.jakob-r.de/parconfigs

# Una segunda opción es descargar directamente el espacio de búsqueda desde la web así:

pc = downloadParConfigs(learner.name=learner_name) # http://mlrhyperopt.jakob-r.de/parconfigs
print(pc)
if(length(pc)>0) {
  ps = getParConfigParSet(pc[[1]], task=task_bh)
  print(ps)
}

# Desgraciadamente, en este caso no hay ninguno exactamente para ese learner
# Tendríamos que ir a la página web http://mlrhyperopt.jakob-r.de/parconfigs
# y buscarlo a mano, para randomForest.

```

En este caso, usaremos una tercera opción, más sencilla: la función **generateParConfig** nos da un espacio de búsqueda básico para un learner concreto.

```{r}
pc_rf <- generateParConfig(learner=learner_name, task=task_bh)
ps_rf <- getParConfigParSet(pc_rf, task=task_bh)
print(ps_rf)
```

Vemos que aparecen sólo 2 hiper-parámetro. La razón de que no aparezca **num.trees** es que su valor por defecto es de 500 y el autor ha considerado que son más que suficientes. También hay que recordar RandomForest no entra en sobreaprendizaje por más árboles que introduzcamos en el ensemble (a diferencia de Boosting). 

Ahora definimos cómo se búsca en el espacio de hiper-parámetros (con **randomsearch** en este caso y un **budget** de 50 evaluaciones), y cómo se evalúan los distintos hiper-parámetros (validación cruzada de tres folds).

```{r}
# Ojo, he puesto 1000 evaluaciones para hacer una buena visualización después, pero puede tardar mucho. 
# Posiblemente con budget = 50 sea suficiente

control_grid <- makeTuneControlRandom(budget=1000)
inner_desc <- makeResampleDesc("CV", iter=3)

# Aquí construimos una secuencia de ajuste seguida de construcción de modelo
learner_tune_rf <- makeTuneWrapper(learner_rf, resampling = inner_desc, par.set = ps_rf, control = control_grid, measures = list(rmse), show.info = FALSE)

# También tenemos que definir cómo se va a evaluar el modelo final. En este caso con train/test (holdout)
outer_desc <-  makeResampleDesc("Holdout")
set.seed(0)
outer_inst <- makeResampleInstance(outer_desc, task_bh)

# Por último, usamos resample para entrenar y evaluar learner_tune_rf
set.seed(0)
error_tune_rf <- resample(learner_tune_rf, task_bh, outer_inst, measures = list(rmse), extract = getTuneResult)
```

Vemos que los hiper-parámetros que se eligieron en la partición de entrenamiento son los siguientes (también se muestra el error que producen dichos hiper-parámetros en la validación cruzada de 3 folds que se usó para seleccionarlos).

```{r}
error_tune_rf$extract
```

y el error del modelo final es:

```{r}
error_tune_rf$aggr
```


```{r}

data <- generateHyperParsEffectData(error_tune_rf, include.diagnostics = FALSE, trafo=TRUE)

# Las siguientes dos líneas son para solventar un bug de MLR
names(data$data)[names(data$data)=="rmse.test.rmse"] <- "rmse.test.mean"
data$measures[data$measures=="rmse.test.rmse"] <- "rmse.test.mean"
```


Aquí vemos el espacio de búsqueda que hemos explorado:

```{r}
plt = plotHyperParsEffect(data, x = "mtry", y = "min.node.size", z = "rmse.test.mean",
  plot.type = "contour", interpolate = "regr.earth", show.experiments = TRUE)

plot(plt + geom_point(x=error_tune_rf$extract[[1]]$x$mtry, y=error_tune_rf$extract[[1]]$x$min.node.size, size=5))
```


A continuación podemos ver como se ha ido reduciendo el error con las iteraciones.

```{r}
plt = plotHyperParsEffect(data, x = "iteration", y = "rmse.test.mean",  plot.type = "line")
plt
```


## Vamos a hacer el ajuste de hiper-parámetros con Model Based Optimization, del páquete **mbo**

```{r}

# Dedica 80 iteraciones al ajuste
control = makeMBOControl()
control = setMBOControlTermination(control, iters = 80)
control = setMBOControlInfill(control, crit = makeMBOInfillCritEI())

control_grid <-makeTuneControlMBO(mbo.control = control)

inner_desc <- makeResampleDesc("CV", iter=3)

# Aquí construimos una secuencia de ajuste seguida de construcción de modelo
learner_tune_rf <- makeTuneWrapper(learner_rf, resampling = inner_desc, par.set = ps_rf, control = control_grid, measures = list(rmse))

# También tenemos que definir cómo se va a evaluar el modelo final. En este caso con train/test (holdout)
outer_desc <-  makeResampleDesc("Holdout")
set.seed(0)
outer_inst <- makeResampleInstance(outer_desc, task_bh)

# Por último, usamos resample para entrenar y evaluar learner_tune_rf
set.seed(0)
error_tune_rf <- resample(learner_tune_rf, task_bh, outer_inst, measures = list(rmse), extract = getTuneResult )
```

Vemos que los hiper-parámetros que se eligieron en la partición de entrenamiento son los siguientes (también se muestra el error que producen dichos hiper-parámetros en la validación cruzada de 3 folds que se usó para seleccionarlos).

```{r}
error_tune_rf$extract
```

y el error del modelo final es muy parecido también al que obtuvimos con Random Search:

```{r}
error_tune_rf$aggr
```

Lo que si parece claro que este tipo de búsqueda se centra mucho mejor en los valores adecuados de n.trees y 
Por último, podemos ver un plot que nos muestra cómo cambia el error con los distintos valores de los hiper-parámetros. 

```{r}
data <- generateHyperParsEffectData(error_tune_rf, include.diagnostics = FALSE)

# Las siguientes dos líneas son para solventar un bug de MLR
names(data$data)[names(data$data)=="rmse.test.rmse"] <- "rmse.test.mean"
data$measures[data$measures=="rmse.test.rmse"] <- "rmse.test.mean"
```


Aquí vemos el espacio de búsqueda que hemos explorado:

```{r}
plt = plotHyperParsEffect(data, x = "mtry", y = "min.node.size", z = "rmse.test.mean",
  plot.type = "contour", interpolate = "regr.earth", show.experiments = TRUE)

plot(plt + geom_point(x=error_tune_rf$extract[[1]]$x$mtry, y=error_tune_rf$extract[[1]]$x$min.node.size, size=5))
```


A continuación podemos ver como se ha ido reduciendo el error con las iteraciones. Parece que con pocas iteraciones se alcanzan buenos mínimos.

```{r}
plt = plotHyperParsEffect(data, x = "iteration", y = "rmse.test.mean",
  plot.type = "line")
plt
```


Por último, puede ser interesante saber qué hubiera pasado de haber utilizado un número de árboles diferente a 500. Definimos un learner con los mejores hiper-parámetros encontrados hasta el momento y un nuevo espacio de búsqueda donde sólo cambie **num.trees**. Ahora utilizaremos la búsqueda sistemática de grid-search para hacerlo.


```{r}

learner_rf <- setHyperPars(makeLearner(learner_name), par.vals = error_tune_rf$extract[[1]]$x)

ps_rf <- makeParamSet(
    makeDiscreteParam("num.trees", values = c(10, seq(50, 2000, by=100)))
)

control_grid <- makeTuneControlGrid()
inner_desc <- makeResampleDesc("CV", iter=3)

# Aquí construimos una secuencia de ajuste seguida de construcción de modelo
# Nótese que inner_Desc, control_grid, outer_inst, etc no cambién, con lo que usamos los mismos que antes.
learner_tune_rf <- makeTuneWrapper(learner_rf, resampling = inner_desc, par.set = ps_rf, control = control_grid, measures = list(rmse))

# Por último, usamos resample para entrenar y evaluar learner_tune_rf
set.seed(0)
error_tune_rf <- resample(learner_tune_rf, task_bh, outer_inst, measures = list(rmse), extract = getTuneResult )




```

```{r}
data <- generateHyperParsEffectData(error_tune_rf, include.diagnostics = FALSE)

# Las siguientes dos líneas son para solventar un bug de MLR
names(data$data)[names(data$data)=="rmse.test.rmse"] <- "rmse.test.mean"
data$measures[data$measures=="rmse.test.rmse"] <- "rmse.test.mean"
```


Aquí vemos el espacio de búsqueda que hemos explorado. Aunque la dependencia del error respecto al número de árboles es algo ruidosa, parece que 500 árboles era un valor razonable, aunque tal vez nos hubieramos podido beneficiar con algo mas de 500 árboles.

```{r}
plt = plotHyperParsEffect(data, x = "num.trees", y = "rmse.test.mean",
  plot.type = "line", show.experiments = TRUE)

plot(plt + geom_point(x=error_tune_rf$extract[[1]]$x$num.trees, y=error_tune_rf$extract[[1]]$y, size=5))
```

```{r}
plt = plotHyperParsEffect(data, x = "iteration", y = "rmse.test.mean",
  plot.type = "line")
plt
```
```{r}

```