[1.] Introducción

Los árboles de clasificación y regresión son modelos predictivos formados por reglas binarias (sí/no) o con las que se consigue repartir las observaciones en función de sus atributos y predecir así el valor de la variable respuesta.

Tomado de: https://i.pinimg.com

Muchos métodos predictivos generan modelos globales en los que una única ecuación se aplica a todo el espacio muestral. Cuando el caso de uso implica múltiples predictores, que interaccionan entre ellos de forma compleja y no lineal, es muy difícil encontrar un único modelo global que sea capaz de reflejar la relación entre las variables. Los métodos de machine learning basados en árboles engloban a un conjunto de técnicas supervisadas no paramétricas que consiguen segmentar el espacio de los predictores en regiones simples, dentro de las cuales es más sencillo manejar las interacciones. Es esta característica la que les proporciona gran parte de su potencial.

1.1 Consideraciones

  • Es sencillo, intuitivo y fácil de entender.

  • Los resultados pueden ser iguales a los de una regresión con muchas interacciones, el problema es que en una regresión tiene que definirse exante las interacciones. En arboles las interacciones las definen los datos.

  • Los arboles tienden a hacer overfit.

  • Son costosos computacionalmente.

  • Si la relación entre la variable dependiente y los predictores es lineal, los arboles no pueden ser la mejor opción.

  • No es necesario escalar las variables (a diferencia de Lasso)

1.2 Ventajas

  • Los árboles son fáciles de interpretar aun cuando las relaciones entre predictores son complejas.

  • Los modelos basados en un solo árbol se pueden representar gráficamente aun cuando el número de predictores es mayor de 3.

  • Los árboles pueden, en teoría, manejar tanto predictores numéricos como categóricos sin tener que crear variables dummy.

  • Al tratarse de métodos no paramétricos, no es necesario que se cumpla ningún tipo de distribución específica.

  • Por lo general, requieren mucha menos limpieza y preprocesado de los datos en comparación con otros métodos de aprendizaje estadístico (por ejemplo, no requieren estandarización).

  • No se ven muy influenciados por outliers en comparación con las regresiones lineales.

  • Si para alguna observación, el valor de un predictor no está disponible, a pesar de no poder llegar a ningún nodo terminal, se puede conseguir una predicción empleando todas las observaciones que pertenecen al último nodo alcanzado. La precisión de la predicción se verá reducida pero al menos podrá obtenerse.

  • Son muy útiles en la exploración de datos, permiten identificar de forma rápida y eficiente las variables (predictores) más importantes.

  • Son capaces de seleccionar predictores de forma automática.

  • Pueden aplicarse a problemas de regresión y clasificación.

1.3 Desventajas

  • La capacidad predictiva de los modelos basados en un único árbol es bastante inferior a la conseguida con otros modelos. Esto es debido a su tendencia al overfitting y alta varianza. Sin embargo, existen técnicas más complejas que, haciendo uso de la combinación de múltiples árboles (bagging, random forest, boosting), consiguen mejorar en gran medida este problema.

  • Son sensibles a datos de entrenamiento desbalanceados (una de las clases domina sobre las demás).

  • Cuando tratan con predictores continuos, pierden parte de su información al categorizarlos en el momento de la división de los nodos.

  • La creación de las ramificaciones de los árboles se consigue mediante el algoritmo de recursive binary splitting. Este algoritmo identifica y evalúa las posibles divisiones de cada predictor acorde a una determinada medida (RSS, Gini, entropía…). Los predictores continuos tienen mayor probabilidad de contener, solo por azar, algún punto de corte óptimo, por lo que suelen verse favorecidos en la creación de los árboles.

  • No son capaces de extrapolar fuera del rango de los predictores observado en los datos de entrenamiento.

[2.] Aplicación

2.0 Configuración inicial

Para replicar las secciones de esta clase, primero debe descargar el siguiente proyecto de R y abrir el archivo clase-06.Rproj.

2.0.0 Instalar/llamar las librerías de la clase

require(pacman) 
p_load(tidyverse,rio,janitor,
       caret, glmnet, elasticnet,
       yardstick, ## metric_set
       tune, ## tune_grid, select_best
       rsample, ## vfold_cv
       dials, ## grid_regular
       rattle, ## fancyRpartPlot
       parsnip,tree,rpart,rpart.plot) ## trees

2.0.1 fijar semilla

set.seed(12345)

2.0.2 leer datos

db <- import("input/data_regresiones.rds") %>% na.omit() %>% select(-id) %>% clean_names()

2.0.3 prepare db

### tibble to matrix
df <- model.matrix(~.+ I(surface_total^2) + I(dist_cbd^2) -1 -price , db) %>% as_tibble() %>% clean_names()

df$price <- db$price/1000000 %>% unlist()

### partir la muestra
sample_test <- sample(x = nrow(df) , size = nrow(df)*.3)

test <- df[sample_test,]

train <- df[-sample_test,]

2.1 Lasso

## predictores
x <- train %>% select(-price) %>% as.matrix()
y <- train$price

## estimations
lasso <- glmnet(x=x , y=y , alpha=1) ## (alpha = 1 para lasso)
dim(coef(lasso)) ## Número de variables y modelos
## [1] 22 67
plot(x=lasso, xvar="lambda" , label=TRUE) ## plot

## mejor lambda (CV)
lasso_cv <- cv.glmnet(x = x , y = y , alpha=1 , nfolds = 10)
min_lambda = lasso_cv$lambda.min

## hacer predicción con mejor lambda
coef(lasso)[,which(lasso$lambda==min_lambda)]
##                   (Intercept)                         rooms 
##                  5.288184e+01                  0.000000e+00 
##                     bathrooms                 surface_total 
##                  2.434626e+02                  1.334808e+00 
##     property_type_apartamento            property_type_casa 
##                  0.000000e+00                 -5.401139e+01 
##        property_type_deposito           property_type_finca 
##                  0.000000e+00                  0.000000e+00 
## property_type_local_comercial            property_type_lote 
##                  3.850758e+02                  8.256745e+02 
##         property_type_oficina            property_type_otro 
##                  2.356701e+02                  3.761056e+01 
##     property_type_parqueadero              property_type_ph 
##                  0.000000e+00                  0.000000e+00 
##                      dist_cbd                     dist_cole 
##                  0.000000e+00                 -3.912405e-02 
##                     dist_park                    luces_2019 
##                  0.000000e+00                  0.000000e+00 
##                    luces_2021                  cambio_activ 
##                  0.000000e+00                  0.000000e+00 
##             i_surface_total_2                  i_dist_cbd_2 
##                  7.981677e-07                 -9.927805e-07
test$price_lasso_cv <- predict(object=lasso , s=min_lambda , newx=select(test,-price) %>% as.matrix()) %>% unlist()

rmse_lasso = sqrt(with(test,mean((price-price_lasso_cv)^2)))

rmse_lasso
## [1] 2097.494

2.2 Arboles

## definir modelo
set_1 <- decision_tree() %>% set_engine("rpart") %>% set_mode("regression")

## ajustar primer modelo: 
model_1 <- fit(object = set_1 , formula = price ~ . , data = train)

model_1
## parsnip model object
## 
## n= 44355 
## 
## node), split, n, deviance, yval
##       * denotes terminal node
## 
##  1) root 44355 101837800000   829.6309  
##    2) surface_total< 18793.5 44346  69317590000   819.3645  
##      4) surface_total< 232.195 37365  20758580000   610.8041  
##        8) bathrooms< 3.5 30597   7714766000   490.0571 *
##        9) bathrooms>=3.5 6768  10580970000  1156.6810  
##         18) bathrooms< 8.5 6676   4361903000  1108.6200 *
##         19) bathrooms>=8.5 92   5084615000  4644.2930 *
##      5) surface_total>=232.195 6981  38234600000  1935.6600  
##       10) surface_total< 2837.5 6882  22686580000  1855.5930  
##         20) surface_total< 461 5291   7896822000  1577.1450 *
##         21) surface_total>=461 1591  13015290000  2781.5930 *
##       11) surface_total>=2837.5 99  12436970000  7501.5400  
##         22) dist_park>=61.66615 78   2709264000  4941.8690 *
##         23) dist_park< 61.66615 21   7318474000 17008.8900  
##           46) surface_total< 3889 9    665694900  6953.3850 *
##           47) surface_total>=3889 12   5060247000 24550.5100 *
##    3) surface_total>=18793.5 9   9485293000 51415.3300 *
fancyRpartPlot(model = model_1$fit)

## Importancia de las variables
import_1 <- varImp(model_1$fit) 
import_1
##                                  Overall
## bathrooms                     0.63936949
## dist_cbd                      0.40602043
## dist_cole                     0.05081865
## dist_park                     0.19371537
## i_dist_cbd_2                  0.36408096
## i_surface_total_2             0.88173502
## property_type_casa            0.06255768
## property_type_lote            0.03061077
## surface_total                 0.88173502
## rooms                         0.00000000
## property_type_apartamento     0.00000000
## property_type_deposito        0.00000000
## property_type_finca           0.00000000
## property_type_local_comercial 0.00000000
## property_type_oficina         0.00000000
## property_type_otro            0.00000000
## property_type_parqueadero     0.00000000
## property_type_ph              0.00000000
## luces_2019                    0.00000000
## luces_2021                    0.00000000
## cambio_activ                  0.00000000
## reescalar los pesos
import_1 <- import_1 %>% as.data.frame() %>%
            rownames_to_column(var = "var") %>%
            mutate(peso = Overall/sum(Overall)) %>%
            filter(peso > 0)

ggplot(import_1, aes(x=peso, y=reorder(var, peso))) +
geom_bar(stat="identity") + theme_test()

## predicción
test$price_tree <- predict(object=model_1 , new_data=test) %>% unlist()

rmse_tree = sqrt(with(test,mean((price-price_tree)^2)))

rmse_tree
## [1] 1600.939

2.3 Arboles: Tune

## definir modelo
set_2 <- decision_tree(cost_complexity = tune(), ## 
                       tree_depth = tune() , ## profundidad
                       min_n = tune() ## minima cantidad de observaciones por hoja
                       ) %>%  
         set_engine("rpart") %>% set_mode("regression")

## crear grilla con hiperparametos a probar
grid_regular(cost_complexity(), 
             min_n(),
             tree_depth(), 
             levels = 4)
## # A tibble: 64 × 3
##    cost_complexity min_n tree_depth
##              <dbl> <int>      <int>
##  1    0.0000000001     2          1
##  2    0.0000001        2          1
##  3    0.0001           2          1
##  4    0.1              2          1
##  5    0.0000000001    14          1
##  6    0.0000001       14          1
##  7    0.0001          14          1
##  8    0.1             14          1
##  9    0.0000000001    27          1
## 10    0.0000001       27          1
## # … with 54 more rows
tree_grid <- crossing(cost_complexity = c(0.0001), 
                      min_n = c(2,14,27),
                      tree_depth = c(2,4,8))

## Definimos CV
folds <- vfold_cv(data=train , v=5) ## resamples

## ajustar segundo modelo: 
model_2 <- tune_grid(set_2, price ~ .,
                     resamples = folds,
                     grid = tree_grid,
                     metrics = metric_set(rmse)
)

## ver los hiperparámetros
collect_metrics(model_2)
## # A tibble: 9 × 9
##   cost_complexity tree_depth min_n .metric .estimator  mean     n std_err
##             <dbl>      <dbl> <dbl> <chr>   <chr>      <dbl> <int>   <dbl>
## 1          0.0001          2     2 rmse    standard   1426.     5    162.
## 2          0.0001          4     2 rmse    standard   1409.     5    174.
## 3          0.0001          8     2 rmse    standard   1376.     5    157.
## 4          0.0001          2    14 rmse    standard   1377.     5    146.
## 5          0.0001          4    14 rmse    standard   1339.     5    158.
## 6          0.0001          8    14 rmse    standard   1284.     5    147.
## 7          0.0001          2    27 rmse    standard   1367.     5    123.
## 8          0.0001          4    27 rmse    standard   1328.     5    139.
## 9          0.0001          8    27 rmse    standard   1272.     5    144.
## # … with 1 more variable: .config <chr>
## visualizar los modelos
autoplot(model_2) + theme_test()

## best model
select_best(model_2)
## # A tibble: 1 × 4
##   cost_complexity tree_depth min_n .config             
##             <dbl>      <dbl> <dbl> <chr>               
## 1          0.0001          8    27 Preprocessor1_Model9
best_model_2 <- finalize_model(set_2, select_best(model_2))

## Entrenamos el mejor modelo
model_2_fit <- fit(object = best_model_2 , formula = price ~ . , data = train)

## cambia la importancia
varImp(model_2_fit$fit) %>% as.data.frame() %>%
rownames_to_column(var = "var") %>%
mutate(peso = Overall/sum(Overall)) %>%
filter(peso > 0) %>%
ggplot(., aes(x=peso, y=reorder(var, peso))) +
geom_bar(stat="identity") + theme_test()

## predicción
test$price_tree_tn <- predict(object=model_2_fit , new_data=test) %>% unlist()

rmse_tree_tn = sqrt(with(test,mean((price-price_tree_tn)^2)))

rmse_lasso
## [1] 2097.494
rmse_tree
## [1] 1600.939
rmse_tree_tn
## [1] 1541.632