Tidymodels y machine learning en R

class: center, middle, inverse, title-slide

# Tidymodels y machine learning en R
### Rocío Joo
### Abril 2021

---

# Plan para este módulo

* Introducción a tidymodels

* Ajuste de un modelo de aprendizaje automático con tidymodels

* Calibración de hiperparámetros del modelo

* Gráficos para la interpretación del modelo

---

# ¿Qué es tidymodels?

Metapaquete para modelar de manera ordenada.

tidy data --> datos ordenados

Artwork by @allison_horst

---

# ¿Qué es tidymodels?

Imagen de https://github.com/data-datum/curso-r-analisis-datos

---

# ¿Por qué tidymodels?

Las funciones de aprendizaje automático en R varían respecto al formato de los
datos de entrada, dependiendo del paquete a utilizar.

Ejm:

```r
  randomForest::randomForest(x = x, y = y, data = data)
```

```r
  ranger::ranger(y ~ x, data = data)
```

---

# ¿Por qué tidymodels?

.pull-left[

**Flujo de trabajo** consistente para:

* Dividir en sets de entrenamiento y test

* Definir una "receta" para el pre-procesamiento de datos

* Especificar el modelo

* Calibrar los hiperparámetros

* Evaluar y comparar modelos

* Generar predicciones con el modelo final

]

.pull-right[

]

---

# ¿Por qué tidymodels?

.pull-left[

**Flujo de trabajo** consistente para:

* Dividir en sets de **entrenamiento** y **test**

* Definir una "receta" para el pre-procesamiento de datos

* Especificar el modelo

* Calibrar los hiperparámetros

* Evaluar y comparar modelos

* Generar predicciones con el modelo final

]

.pull-right[

]

---

# ¿Por qué tidymodels?

.pull-left[

**Flujo de trabajo** consistente para:

* Dividir en sets de **entrenamiento** y **test**

* Definir una "receta" para el **pre-procesamiento** de datos

* Especificar el modelo

* Calibrar los hiperparámetros

* Evaluar y comparar modelos

* Generar predicciones con el modelo final

]

.pull-right[

]

---

# ¿Por qué tidymodels?

.pull-left[

**Flujo de trabajo** consistente para:

* Dividir en sets de **entrenamiento** y **test**

* Definir una "receta" para el **pre-procesamiento** de datos

* Especificar el **modelo**

* Calibrar los hiperparámetros

* Evaluar y comparar modelos

* Generar predicciones con el modelo final

]

.pull-right[

]

---

# ¿Por qué tidymodels?

.pull-left[

**Flujo de trabajo** consistente para:

* Dividir en sets de **entrenamiento** y **test**

* Definir una "receta" para el **pre-procesamiento** de datos

* Especificar el **modelo**

* **Calibrar** los hiperparámetros

* Evaluar y comparar modelos

* Generar predicciones con el modelo final

]

.pull-right[

]

---

# ¿Por qué tidymodels?

.pull-left[

**Flujo de trabajo** consistente para:

* Dividir en sets de **entrenamiento** y **test**

* Definir una "receta" para el **pre-procesamiento** de datos

* Especificar el **modelo**

* **Calibrar** los hiperparámetros

* **Evaluar** y comparar modelos

* Generar predicciones con el modelo final

]

.pull-right[

]

---

# ¿Por qué tidymodels?

.pull-left[

**Flujo de trabajo** consistente para:

* Dividir en sets de **entrenamiento** y **test**

* Definir una "receta" para el **pre-procesamiento** de datos

* Especificar el **modelo**

* **Calibrar** los hiperparámetros

* **Evaluar** y comparar modelos

* Generar **predicciones** con el modelo final

]

.pull-right[

]

---

# Primer caso:

# Ajustando un modelo de aprendizaje automático y predecir con `tidymodels`.

* Ejercicio de clasificación
  * Bosques aleatorios
  * Entrenamiento y test
  * Parametrización por defecto
  * Flujo de trabajo simple
  * Matriz de confusión

---

# Los datos

Imaginemos que queremos crear un modelo para predecir **especie** de pingüinos en
la Antártida en base a su morfología u otras variables disponibles.

```r
library(datos)
pinguinos
```

```
## # A tibble: 344 x 8
## especie isla largo_pico_mm alto_pico_mm largo_aleta_mm masa_corporal_g sexo 
## <fct> <fct> <dbl> <dbl> <int> <int> <fct>
## 1 Adelia Torg… 39.1 18.7 181 3750 macho
## 2 Adelia Torg… 39.5 17.4 186 3800 hemb…
## 3 Adelia Torg… 40.3 18 195 3250 hemb…
## 4 Adelia Torg… NA NA NA NA <NA> 
## 5 Adelia Torg… 36.7 19.3 193 3450 hemb…
## 6 Adelia Torg… 39.3 20.6 190 3650 macho
## 7 Adelia Torg… 38.9 17.8 181 3625 hemb…
## 8 Adelia Torg… 39.2 19.6 195 4675 macho
## 9 Adelia Torg… 34.1 18.1 193 3475 <NA> 
## 10 Adelia Torg… 42 20.2 190 4250 <NA> 
## # … with 334 more rows, and 1 more variable: anio <int>
```

---

# Explorando los datos

```r
summary(pinguinos)
```

```
##     especie           isla     largo_pico_mm    alto_pico_mm   largo_aleta_mm 
##  Adelia :152   Biscoe   :168   Min.   :32.10   Min.   :13.10   Min.   :172.0  
##  Barbijo: 68   Dream    :124   1st Qu.:39.23   1st Qu.:15.60   1st Qu.:190.0  
##  Papúa  :124   Torgersen: 52   Median :44.45   Median :17.30   Median :197.0  
##                                Mean   :43.92   Mean   :17.15   Mean   :200.9  
##                                3rd Qu.:48.50   3rd Qu.:18.70   3rd Qu.:213.0  
##                                Max.   :59.60   Max.   :21.50   Max.   :231.0  
##                                NA's   :2       NA's   :2       NA's   :2      
##  masa_corporal_g     sexo          anio     
##  Min.   :2700    hembra:165   Min.   :2007  
##  1st Qu.:3550    macho :168   1st Qu.:2007  
##  Median :4050    NA's  : 11   Median :2008  
##  Mean   :4202                 Mean   :2008  
##  3rd Qu.:4750                 3rd Qu.:2009  
##  Max.   :6300                 Max.   :2009  
##  NA's   :2
```

---

# Explorando los datos

```r
library(magrittr) #pipes
library(tidyr) #drop_na
pinguinos <- pinguinos %>% 
 drop_na()

summary(pinguinos)
```

```
##     especie           isla     largo_pico_mm    alto_pico_mm   largo_aleta_mm
##  Adelia :146   Biscoe   :163   Min.   :32.10   Min.   :13.10   Min.   :172   
##  Barbijo: 68   Dream    :123   1st Qu.:39.50   1st Qu.:15.60   1st Qu.:190   
##  Papúa  :119   Torgersen: 47   Median :44.50   Median :17.30   Median :197   
##                                Mean   :43.99   Mean   :17.16   Mean   :201   
##                                3rd Qu.:48.60   3rd Qu.:18.70   3rd Qu.:213   
##                                Max.   :59.60   Max.   :21.50   Max.   :231   
##  masa_corporal_g     sexo          anio     
##  Min.   :2700    hembra:165   Min.   :2007  
##  1st Qu.:3550    macho :168   1st Qu.:2007  
##  Median :4050                 Median :2008  
##  Mean   :4207                 Mean   :2008  
##  3rd Qu.:4775                 3rd Qu.:2009  
##  Max.   :6300                 Max.   :2009
```

---

# Explorando los datos

```r
library(ggplot2)
pinguinos %>% 
  ggplot(aes(x = largo_pico_mm, y = largo_aleta_mm, color = especie)) +
  geom_point() +
  theme_classic()
```

![](tidymodels_files/figure-html/unnamed-chunk-14-1.png)

---

# Explorando los datos

```r
library(dplyr)
pinguinos %>% 
  select(-especie, -isla, -sexo, -anio) %>% 
  cor(use = "pairwise.complete.obs") %>% 
  round(digits = 2)
```

```
##                 largo_pico_mm alto_pico_mm largo_aleta_mm masa_corporal_g
## largo_pico_mm            1.00        -0.23           0.65            0.59
## alto_pico_mm            -0.23         1.00          -0.58           -0.47
## largo_aleta_mm           0.65        -0.58           1.00            0.87
## masa_corporal_g          0.59        -0.47           0.87            1.00
```
---

# Explorando los datos

```r
library(GGally)
ggpairs(pinguinos)
```

![](tidymodels_files/figure-html/unnamed-chunk-16-1.png)

---

# Dividiendo entre datos de entrenamiento y test

Usamos el paquete `rsample`

```r
library(rsample)

set.seed(4563) # for traceability
pinguinos_div <- initial_split(data = pinguinos, prop = 0.75)
pinguinos_div
```

```
## <Analysis/Assess/Total>
## <250/83/333>
```

```r
pinguinos_entren <- training(pinguinos_div)
pinguinos_test <- testing(pinguinos_div)
```

---

# Especificando el modelo

En general, hay que especificar:

* El tipo de modelo
* El modo (clasificación o regresión)
* El motor ("engine") o paquete

El [motor de búsqueda de tidymodels](https://www.tidymodels.org/find/parsnip/) 
ofrece muchas opciones.

Digamos que usaremos bosques aleatorios.

* El tipo de modelo es `rand_forest`

* El modo es `classification`

* El motor que usaremos es `ranger`

---

# Especificando el modelo

```r
library(parsnip)
rf_espec <- 
 # tipo de modelo
 rand_forest() %>% 
 # modo
 set_mode("classification") %>% 
 # motor
 set_engine("ranger")

rf_espec
```

```
## Random Forest Model Specification (classification)
## 
## Computational engine: ranger
```

---

# Creando una receta de pre-procesamiento

Artwork by @allison_horst

---

# Creando una receta de pre-procesamiento

```r
library(recipes)
receta <- recipe(especie ~ ., # modelar especie usando todas las otras variables
 data = pinguinos) %>% 
 # cambié de opinión respecto a anio 
 update_role(anio, new_role = "no usar") %>% 
 # remover variables con correlaciones > 0.8
 step_corr(all_predictors() & all_numeric(), threshold = 0.8)

summary(receta)
```

```
## # A tibble: 8 x 4
## variable type role source 
## <chr> <chr> <chr> <chr> 
## 1 isla nominal predictor original
## 2 largo_pico_mm numeric predictor original
## 3 alto_pico_mm numeric predictor original
## 4 largo_aleta_mm numeric predictor original
## 5 masa_corporal_g numeric predictor original
## 6 sexo nominal predictor original
## 7 anio numeric no usar original
## 8 especie nominal outcome original
```

---

# Armando el flujo de trabajo

```r
library(workflows)
# generando el flujo de trabajo con workflow
rf_flujo <- workflow() %>% 
 # añadir la receta
 add_recipe(receta) %>% 
 # añadir las especificaciones del modelo
 add_model(rf_espec)
```

---

# Armando el flujo de trabajo

```r
rf_flujo 
```

```
## ══ Workflow ════════════════════════════════════════════════════════════════════
## Preprocessor: Recipe
## Model: rand_forest()
## 
## ── Preprocessor ────────────────────────────────────────────────────────────────
## 1 Recipe Step
## 
## ● step_corr()
## 
## ── Model ───────────────────────────────────────────────────────────────────────
## Random Forest Model Specification (classification)
## 
## Computational engine: ranger
```

---

# Aplicando la receta y ajustando el modelo a los datos

```r
rf_ajuste <- fit(object = rf_flujo, 
 data = pinguinos_entren)
```

---

```
## ══ Workflow [trained] ══════════════════════════════════════════════════════════
## Preprocessor: Recipe
## Model: rand_forest()
## 
## ── Preprocessor ────────────────────────────────────────────────────────────────
## 1 Recipe Step
## 
## ● step_corr()
## 
## ── Model ───────────────────────────────────────────────────────────────────────
## Ranger result
## 
## Call:
##  ranger::ranger(x = maybe_data_frame(x), y = y, num.threads = 1,      verbose = FALSE, seed = sample.int(10^5, 1), probability = TRUE) 
## 
## Type:                             Probability estimation 
## Number of trees:                  500 
## Sample size:                      250 
## Number of independent variables:  5 
## Mtry:                             2 
## Target node size:                 10 
## Variable importance mode:         none 
## Splitrule:                        gini 
## OOB prediction error (Brier s.):  0.01904099
```

---

# Prediciendo sobre los datos de test

```r
predict(object = rf_ajuste, 
        new_data = pinguinos_test)
```

```
## # A tibble: 83 x 1
## .pred_class
## <fct> 
## 1 Adelia 
## 2 Adelia 
## 3 Adelia 
## 4 Adelia 
## 5 Adelia 
## 6 Adelia 
## 7 Adelia 
## 8 Adelia 
## 9 Adelia 
## 10 Adelia 
## # … with 73 more rows
```

---

# Prediciendo sobre los datos de test

```r
(predict(object = rf_ajuste, 
        new_data = pinguinos_test) %>% 
  transmute(pred = .pred_class,
            real = pinguinos_test$especie) -> pinguinos_test_pred)
```

```
## # A tibble: 83 x 2
## pred real 
## <fct> <fct> 
## 1 Adelia Adelia
## 2 Adelia Adelia
## 3 Adelia Adelia
## 4 Adelia Adelia
## 5 Adelia Adelia
## 6 Adelia Adelia
## 7 Adelia Adelia
## 8 Adelia Adelia
## 9 Adelia Adelia
## 10 Adelia Adelia
## # … with 73 more rows
```

---

# Calculando la matriz de confusión

```r
library(yardstick)
pinguinos_test_pred %>% 
  conf_mat(truth = real, estimate = pred,
           dnn = c("Predicción", "Real"))
```

```
##           Real
## Predicción Adelia Barbijo Papúa
##    Adelia      43       0     0
##    Barbijo      0      17     0
##    Papúa        0       0    23
```
Resultado perfecto.

---

# Segundo caso:

# Calibrando los hiperparámetros del modelo.

* Ejercicio de clasificación
  * Bosques aleatorios
  * Remuestreo en el set de entrenamiento
  * Calibración de hiperparámetros
  * Matriz de confusión e indicadores

---

Ahora buscaremos predecir el sexo de los pingüinos a partir morfológicos y otros.

# Cambiando la receta

```r
receta2 <- recipe(sexo ~ .,
 data = pinguinos) %>% 
 update_role(anio, new_role = "no usar") %>% 
 update_role(isla, new_role = "no usar") %>% 
 step_corr(all_predictors() & all_numeric(), threshold = 0.8)
```

---

# Especificando el modelo

Esta vez, especificando los hiperparámetros a calibrar para el tipo de modelo.

Los hiperparámetros a especificar se pueden ver con

```r
show_model_info("rand_forest")
```

```
## Information for `rand_forest`
##  modes: unknown, classification, regression 
## 
##  engines: 
##    classification: randomForest, ranger, spark
##    regression:     randomForest, ranger, spark
## 
##  arguments: 
##    ranger:       
##       mtry  --> mtry
##       trees --> num.trees
##       min_n --> min.node.size
##    randomForest: 
##       mtry  --> mtry
##       trees --> ntree
##       min_n --> nodesize
##    spark:        
##       mtry  --> feature_subset_strategy
##       trees --> num_trees
##       min_n --> min_instances_per_node
## 
##  fit modules:
##          engine           mode
##          ranger classification
##          ranger     regression
##    randomForest classification
##    randomForest     regression
##           spark classification
##           spark     regression
## 
##  prediction modules:
##              mode       engine                    methods
##    classification randomForest           class, prob, raw
##    classification       ranger class, conf_int, prob, raw
##    classification        spark                class, prob
##        regression randomForest               numeric, raw
##        regression       ranger     conf_int, numeric, raw
##        regression        spark                    numeric
```

---

# Especificando el modelo

```r
library(tune)
rf_espec2 <- rand_forest(
 mtry = tune(), # número de predictores a usar en cada división
 min_n = tune() # número de observaciones necesarias para seguir dividiendo en nodos
) %>% 
 set_mode("classification") %>% 
 set_engine("ranger")
```

---

# Entrenamiento y test

Usamos las particiones anteriores `pinguinos_entren` y `pinguinos_test` a partir
de `pinguinos_div`.

```r
pinguinos_div
```

```
## <Analysis/Assess/Total>
## <250/83/333>
```

* Esta vez no nos conformaremos con una sola división en entrenamiento y test
para probar modelos con parametrizaciones distintas.

* Haremos bootstrapping de los datos de entrenamiento para calibrar hiperparámetros

* Luego haremos la evaluación final del modelo sobre los datos de test.

---

# Entrenamiento y test

```r
set.seed(8723)
pinguinos_boot <- bootstraps(pinguinos_entren, strata = sexo)
pinguinos_boot
```

```
## # Bootstrap sampling using stratification 
## # A tibble: 25 x 2
## splits id 
## <list> <chr> 
## 1 <split [250/89]> Bootstrap01
## 2 <split [250/86]> Bootstrap02
## 3 <split [250/89]> Bootstrap03
## 4 <split [250/96]> Bootstrap04
## 5 <split [250/92]> Bootstrap05
## 6 <split [250/88]> Bootstrap06
## 7 <split [250/93]> Bootstrap07
## 8 <split [250/91]> Bootstrap08
## 9 <split [250/89]> Bootstrap09
## 10 <split [250/85]> Bootstrap10
## # … with 15 more rows
```

---

# Armando el flujo de trabajo

```r
rf_flujo2 <- workflow() %>% 
 add_recipe(receta2) %>% 
 add_model(rf_espec2)
```

---

# Armando el flujo de trabajo

```r
rf_flujo2
```

---

# Calibrando hiperparámetros

* Haremos cálculo en paralelo para optimizar procesos (opcional).
  
  * Si no tenemos idea de los valores a probar, fijamos un número de valores y 
  el algoritmo internamente los fija

```r
library(parallel)
all_cores <- parallel::detectCores(logical = FALSE) # paralelo
cl <- makePSOCKcluster(all_cores-3) # paralelo
doParallel::registerDoParallel(cl) # paralelo

set.seed(345)
calibrar_res <- tune_grid(
 rf_flujo2,
 resamples = pinguinos_boot,
 grid = 10 # 10 combinaciones de valores de hiperparámetros
)

stopCluster(cl) # paralelo
```

---

# Calibrando hiperparámetros

Veamos los valores probados y sus resultados respecto a la exactitud (*accuracy*)

```r
calibrar_res %>% 
  collect_metrics() 
```

```
## # A tibble: 20 x 8
## mtry min_n .metric .estimator mean n std_err .config 
## <int> <int> <chr> <chr> <dbl> <int> <dbl> <chr> 
## 1 2 30 accuracy binary 0.878 25 0.00667 Preprocessor1_Model01
## 2 2 30 roc_auc binary 0.956 25 0.00358 Preprocessor1_Model01
## 3 3 12 accuracy binary 0.883 25 0.00725 Preprocessor1_Model02
## 4 3 12 roc_auc binary 0.955 25 0.00351 Preprocessor1_Model02
## 5 2 27 accuracy binary 0.880 25 0.00655 Preprocessor1_Model03
## 6 2 27 roc_auc binary 0.957 25 0.00354 Preprocessor1_Model03
## 7 3 24 accuracy binary 0.874 25 0.00688 Preprocessor1_Model04
## 8 3 24 roc_auc binary 0.953 25 0.00372 Preprocessor1_Model04
## 9 3 39 accuracy binary 0.870 25 0.00680 Preprocessor1_Model05
## 10 3 39 roc_auc binary 0.951 25 0.00377 Preprocessor1_Model05
## 11 2 6 accuracy binary 0.890 25 0.00590 Preprocessor1_Model06
## 12 2 6 roc_auc binary 0.960 25 0.00308 Preprocessor1_Model06
## 13 3 16 accuracy binary 0.881 25 0.00723 Preprocessor1_Model07
## 14 3 16 roc_auc binary 0.954 25 0.00357 Preprocessor1_Model07
## 15 4 17 accuracy binary 0.870 25 0.00764 Preprocessor1_Model08
## 16 4 17 roc_auc binary 0.948 25 0.00460 Preprocessor1_Model08
## 17 2 7 accuracy binary 0.895 25 0.00632 Preprocessor1_Model09
## 18 2 7 roc_auc binary 0.960 25 0.00330 Preprocessor1_Model09
## 19 1 36 accuracy binary 0.894 25 0.00677 Preprocessor1_Model10
## 20 1 36 roc_auc binary 0.962 25 0.00346 Preprocessor1_Model10
```

---

# Calibrando hiperparámetros

Veamos los valores probados y sus resultados respecto a la exactitud (*accuracy*)

```r
calibrar_res %>% 
  collect_metrics() %>% 
  filter(.metric == "accuracy") 
```

```
## # A tibble: 10 x 8
## mtry min_n .metric .estimator mean n std_err .config 
## <int> <int> <chr> <chr> <dbl> <int> <dbl> <chr> 
## 1 2 30 accuracy binary 0.878 25 0.00667 Preprocessor1_Model01
## 2 3 12 accuracy binary 0.883 25 0.00725 Preprocessor1_Model02
## 3 2 27 accuracy binary 0.880 25 0.00655 Preprocessor1_Model03
## 4 3 24 accuracy binary 0.874 25 0.00688 Preprocessor1_Model04
## 5 3 39 accuracy binary 0.870 25 0.00680 Preprocessor1_Model05
## 6 2 6 accuracy binary 0.890 25 0.00590 Preprocessor1_Model06
## 7 3 16 accuracy binary 0.881 25 0.00723 Preprocessor1_Model07
## 8 4 17 accuracy binary 0.870 25 0.00764 Preprocessor1_Model08
## 9 2 7 accuracy binary 0.895 25 0.00632 Preprocessor1_Model09
## 10 1 36 accuracy binary 0.894 25 0.00677 Preprocessor1_Model10
```

---

# Calibrando hiperparámetros

Veamos los valores probados y sus resultados respecto a la exactitud (*accuracy*)

```r
calibrar_res %>% 
  collect_metrics() %>% 
  filter(.metric == "accuracy") %>% 
  select(mtry, min_n, mean)
```

```
## # A tibble: 10 x 3
## mtry min_n mean
## <int> <int> <dbl>
## 1 2 30 0.878
## 2 3 12 0.883
## 3 2 27 0.880
## 4 3 24 0.874
## 5 3 39 0.870
## 6 2 6 0.890
## 7 3 16 0.881
## 8 4 17 0.870
## 9 2 7 0.895
## 10 1 36 0.894
```

---

# Calibrando hiperparámetros

Veamos los valores probados y sus resultados respecto a la exactitud (*accuracy*)

```
## # A tibble: 20 x 3
## mean hiperparámetro valor
## <dbl> <chr> <int>
## 1 0.878 mtry 2
## 2 0.878 min_n 30
## 3 0.883 mtry 3
## 4 0.883 min_n 12
## 5 0.880 mtry 2
## 6 0.880 min_n 27
## 7 0.874 mtry 3
## 8 0.874 min_n 24
## 9 0.870 mtry 3
## 10 0.870 min_n 39
## 11 0.890 mtry 2
## 12 0.890 min_n 6
## 13 0.881 mtry 3
## 14 0.881 min_n 16
## 15 0.870 mtry 4
## 16 0.870 min_n 17
## 17 0.895 mtry 2
## 18 0.895 min_n 7
## 19 0.894 mtry 1
## 20 0.894 min_n 36
```

---

# Calibrando hiperparámetros

Veamos los valores probados y sus resultados respecto a la exactitud (*accuracy*)

```r
calibrar_res %>% 
  collect_metrics() %>% 
  filter(.metric == "accuracy") %>% 
  select(mtry, min_n, mean) %>% 
  pivot_longer(mtry:min_n,
               values_to = "valor",
               names_to = "hiperparametro") %>% 
  ggplot(aes(x = valor, y = mean)) +
  geom_point(show.legend = FALSE) + 
  ylab("exactitud") +
  facet_wrap(~hiperparametro, scales = "free_x") +
  theme_light()
```
---

# Calibrando hiperparámetros

![](tidymodels_files/figure-html/unnamed-chunk-42-1.png)

---

# Calibrando hiperparámetros

A partir de estos resultados podemos hacer una búsqueda de combinaciones de 
valores más sistemática.

```r
grilla <- expand.grid(
 mtry = c(1,2),
 min_n = seq(from = 10, to = 40, by=5)
)

grilla
```

```
##    mtry min_n
## 1     1    10
## 2     2    10
## 3     1    15
## 4     2    15
## 5     1    20
## 6     2    20
## 7     1    25
## 8     2    25
## 9     1    30
## 10    2    30
## 11    1    35
## 12    2    35
## 13    1    40
## 14    2    40
```

---

# Calibrando hiperparámetros

Y calibramos otra vez

```r
cl <- makePSOCKcluster(all_cores-3) # paralelo
doParallel::registerDoParallel(cl) # paralelo

set.seed(543)
calibrar2_res <- tune_grid(
 rf_flujo2,
 resamples = pinguinos_boot,
 grid = grilla
)

stopCluster(cl) # paralelo
```

---

# Calibrando hiperparámetros

Veamos los resultados

```r
calibrar2_res %>% 
  collect_metrics %>% 
  filter(.metric == "accuracy") 
```

```
## # A tibble: 14 x 8
## mtry min_n .metric .estimator mean n std_err .config 
## <dbl> <dbl> <chr> <chr> <dbl> <int> <dbl> <chr> 
## 1 1 10 accuracy binary 0.903 25 0.00670 Preprocessor1_Model01
## 2 2 10 accuracy binary 0.891 25 0.00685 Preprocessor1_Model02
## 3 1 15 accuracy binary 0.900 25 0.00615 Preprocessor1_Model03
## 4 2 15 accuracy binary 0.890 25 0.00679 Preprocessor1_Model04
## 5 1 20 accuracy binary 0.898 25 0.00637 Preprocessor1_Model05
## 6 2 20 accuracy binary 0.887 25 0.00636 Preprocessor1_Model06
## 7 1 25 accuracy binary 0.900 25 0.00671 Preprocessor1_Model07
## 8 2 25 accuracy binary 0.884 25 0.00642 Preprocessor1_Model08
## 9 1 30 accuracy binary 0.894 25 0.00612 Preprocessor1_Model09
## 10 2 30 accuracy binary 0.878 25 0.00652 Preprocessor1_Model10
## 11 1 35 accuracy binary 0.894 25 0.00638 Preprocessor1_Model11
## 12 2 35 accuracy binary 0.878 25 0.00630 Preprocessor1_Model12
## 13 1 40 accuracy binary 0.891 25 0.00654 Preprocessor1_Model13
## 14 2 40 accuracy binary 0.874 25 0.00706 Preprocessor1_Model14
```

---

# Calibrando hiperparámetros

Veamos los resultados

```r
calibrar2_res %>% 
  collect_metrics %>% 
  filter(.metric == "accuracy") %>% 
  select(mtry, min_n, mean, std_err)
```

```
## # A tibble: 14 x 4
## mtry min_n mean std_err
## <dbl> <dbl> <dbl> <dbl>
## 1 1 10 0.903 0.00670
## 2 2 10 0.891 0.00685
## 3 1 15 0.900 0.00615
## 4 2 15 0.890 0.00679
## 5 1 20 0.898 0.00637
## 6 2 20 0.887 0.00636
## 7 1 25 0.900 0.00671
## 8 2 25 0.884 0.00642
## 9 1 30 0.894 0.00612
## 10 2 30 0.878 0.00652
## 11 1 35 0.894 0.00638
## 12 2 35 0.878 0.00630
## 13 1 40 0.891 0.00654
## 14 2 40 0.874 0.00706
```

---

# Calibrando hiperparámetros

Veamos los resultados

```r
calibrar2_res %>% 
  collect_metrics %>% 
  filter(.metric == "accuracy") %>% 
  select(mtry, min_n, mean, std_err) %>% 
  mutate(min_n = factor(min_n)) %>% 
  ggplot(aes(x = mtry, y = mean, color = min_n)) +
  geom_line() +
  geom_point() +
  ylab("exactitud") +
  scale_color_brewer(type = "qual", palette = "Set1" ) +
  theme_bw()
```

---

# Calibrando hiperparámetros

Veamos los resultados

.pull-left[
![](tidymodels_files/figure-html/unnamed-chunk-48-1.png)
]

.pull-right[
Es obvio que mtry=1 da mejores resultados. No hay mucha diferencia en min_n;
nos quedamos con min_n=20.
]

---

# Actualizando las especificaciones del modelo

Lo podemos hacer directamente en nuestro flujo de trabajo

```r
rf_flujo2_cal <- finalize_workflow(
 x = rf_flujo2,
 parameters = list(mtry = 1, min_n = 20) # select_best(calibrar_res, "accuracy")
)

rf_flujo2_cal
```
---

# Actualizando las especificaciones del modelo

---

# Ajustando a los datos y prediciendo

```r
rf_ajuste2_cal <- rf_flujo2_cal %>% 
 fit(data = pinguinos_entren)

predict(object = rf_ajuste2_cal, new_data =  pinguinos_test) %>% 
  transmute(pred = .pred_class,
            real = pinguinos_test$sexo) -> pinguinos_sex_pred
```

---

# Calculando la matriz de confusión e indicadores

```r
pinguinos_sex_pred %>% 
conf_mat(truth = real, estimate = pred, 
         dnn = c("Predicción", "Real"))
```

```
##           Real
## Predicción hembra macho
##     hembra     33     4
##     macho       3    43
```

```r
multi_metric <- metric_set(accuracy, recall, precision)
multi_metric(data = pinguinos_sex_pred, truth = real, estimate = pred)
```

```
## # A tibble: 3 x 3
## .metric .estimator .estimate
## <chr> <chr> <dbl>
## 1 accuracy binary 0.916
## 2 recall binary 0.917
## 3 precision binary 0.892
```

---

# Tercer caso:

# Interpretando el modelo

* Ejercicio de clasificación
  * Bosques aleatorios
  * Hiperparametrización anterior
  * Importancia de variables
  * Dependencia parcial

---

# Especificando el modelo

Cuando nos interesa la interpretación, es mejor pedir el cálculo de 
la importancia al especificar el modelo.

Ahora usaremos los valores de hiperparámetros del caso anterior.

```r
rf_espec3 <- rand_forest() %>% 
 set_args(mtry = 1,
 min_n = 20) %>% 
 set_mode("classification") %>% 
 set_engine("ranger", importance = "permutation")

rf_flujo3 <- workflow() %>% 
 add_recipe(receta2) %>% 
 add_model(rf_espec3)
```

Esta vez ajustamos a todos los datos

```r
rf_ajuste3 <- fit(object = rf_flujo3,
 data = pinguinos)
```
---

# Graficando la importancia

Extraemos el modelo

```r
modelo_ajustado <- rf_ajuste3 %>% 
 pull_workflow_fit()
```

Y de él, la importancia de variables

```r
library(vip)
(importancia <- vi(modelo_ajustado))
```

```
## # A tibble: 4 x 2
## Variable Importance
## <chr> <dbl>
## 1 masa_corporal_g 0.129 
## 2 alto_pico_mm 0.110 
## 3 especie 0.0924
## 4 largo_pico_mm 0.0829
```

---

# Graficando la importancia

```r
importancia %>%
  mutate(Variable = factor(Variable, 
                           levels = Variable[order(Importance)])) %>% 
  ggplot(aes(y = Variable, x = Importance)) +
  geom_point() +
  theme_bw()
```

![](tidymodels_files/figure-html/unnamed-chunk-58-1.png)

---

# Graficando la dependencia parcial

Primero creamos una función para especificar 
* qué objeto se usa para predecir
* sobre qué datos
* qué clase de predicción
* y cómo se llama la columna de predicción

Es necesario para la función `partial`

```r
library(pdp)

pdp_pred_fun <- function(object, newdata) {
 predict(object, newdata, type = "class")$.pred_class
}
```

---

# Graficando la dependencia parcial

```r
partial(rf_ajuste3,
          pred.var = "masa_corporal_g",
          pred.fun = pdp_pred_fun,
          train = pinguinos, center = FALSE) %>% 
plotPartial(rug = FALSE, center = FALSE, levelplot = FALSE, 
            alpha = 0, pdp.col = "#1f78b4")
```

![](tidymodels_files/figure-html/unnamed-chunk-60-1.png)

---

# Graficando la dependencia parcial

```r
partial(rf_ajuste3,
          pred.var = "alto_pico_mm",
          pred.fun = pdp_pred_fun,
          train = pinguinos, center = FALSE) %>% 
plotPartial(rug = FALSE, center = FALSE, levelplot = FALSE, 
            alpha = 0, pdp.col = "#1f78b4")
```

![](tidymodels_files/figure-html/unnamed-chunk-61-1.png)

---

# Graficando la dependencia parcial

```r
partial(rf_ajuste3,
          pred.var = "largo_pico_mm",
          pred.fun = pdp_pred_fun,
          train = pinguinos, center = FALSE) %>% 
plotPartial(rug = FALSE, center = FALSE, levelplot = FALSE, 
            alpha = 0, pdp.col = "#1f78b4")
```

![](tidymodels_files/figure-html/unnamed-chunk-62-1.png)

---

# Bibliografía

Para preparar esta unidad se utilizó:

* [Tidy Modeling with R](https://www.tmwr.org/). Max Kuhn and Julia Silge. Book.

* [Tidymodels page](https://www.tidymodels.org/)

* [Un tour por tidymodels](https://www.youtube.com/watch?v=o38PVdS3k_Q). Ana Laura Diedrichs. RLadies Buenos Aires. Video.

* [Tidymodels R Ladies talk](https://github.com/rlbarter/tidymodels_talk). Rebecca Barter. RLadies Baltimore. Notebook.

* [Get started with tidymodels and #TidyTuesday Palmer penguins](https://juliasilge.com/blog/palmer-penguins/). Julia Silge. Blog post and video.

---

# Bibliografía

Para preparar esta unidad se utilizó:

* [Dials, Tune, and Parsnip: Tidymodels’ Way to Create and Tune Model Parameters](https://towardsdatascience.com/dials-tune-and-parsnip-tidymodels-way-to-create-and-tune-model-parameters-c97ba31d6173). Yu En Hsu. Blog post.

* [Tuning random forest hyperparameters with tidymodels](https://www.youtube.com/watch?v=ts5bRZ7pRKQ). Julia Silge. Video.

* [Tidymodels: tidy machine learning in R](http://www.rebeccabarter.com/blog/2020-03-25_machine_learning/#tune-the-parameters). Rebecca Barter. Blog post.

* [How to use `pdp:: partial()` in `tidymodels`](https://community.rstudio.com/t/how-to-use-pdp-partial-in-tidymodels/69735). Forum.

---

# Bibliografía

Más recomendaciones:

* [Different random forest packages in R](https://www.linkedin.com/pulse/different-random-forest-packages-r-madhur-modi/). LinkedIn post.

* [tidymodels, rpart, and pengnuins](https://koderkow.rbind.io/post/tidymodels-rpart-and-pengnuins/). Kyle Harris. Blog post.

---

# Ejercicio de aplicación

* Usen otro modelo para predecir `sexo` u otra variable

* Comparen al menos dos configuraciones de modelos

* Hagan su receta, división en entrenamiento y test, 
calculen los indicadores. Justifiquen sus elecciones.

* Interpreten

* Pueden compartir su archivo en .Rmd, o R.