Bookings Modeling
Andrés Expósito
Carga de librerías
library(tidyverse)
library(tidymodels)
library(rpart)
library(rpart.plot)
library(caret)
library(themis)
library(kknn)
library(randomForest)
library(pROC)Fase 3: Modelo
Vamos a utilizar el fichero que creamos en el análisis, con la limpieza y preprocesamiento ya hechos.
bookings <- read.csv("data/bookings_for_models.csv")Selección de variables
Para la selección de variables, vamos a escoger algunas de las que vimos que tenían más influencia en el análisis exploratorio, y también la que creamos en el preprocesamiento. De momento, escogeremos 12 variables predictoras.
bookings_for_model <- bookings %>%
mutate_if(is.character, as.factor) %>%
mutate(is_canceled = as.factor(is_canceled)) %>%
select(assigned_room_type, adr, days_in_waiting_list, lead_time, stays_in_week_nights, stays_in_weekend_nights, adults, booking_changes, deposit_type, market_segment, previous_cancellations, is_portuguese, is_canceled) Particionamiento entrenamiento vs test 70/30
Estableceremos una semilla antes de realizar el particionamiento. Después, utilizaremos sample_frac() para generar los datos de entrenamiento y setdiff() para crear los datos de test con los datos que no se incluyeron en los datos de entrenamiento.
set.seed(1234)
bookings_train <- sample_frac(bookings_for_model, .7)
bookings_test <- setdiff(bookings_for_model, bookings_train)Receta
Crearemos una receta con la función recipe() utilizando los datos de entrenamiento. Vamos a realizar downsampling a la variable is_canceled, y, como dijimos en la parte de preprocesamiento, debemos factorizar las variables categóricas y normalizar las variables numéricas por rango.
También vamos a utilizar la función bake(), esta función nos permite aplicar la receta que hemos creado a los datos que le indiquemos.
bookings_recipe <- recipe(is_canceled ~ ., data = bookings_train) %>%
step_downsample(is_canceled) %>%
step_dummy(all_nominal(), - all_outcomes()) %>%
step_zv(all_numeric()) %>%
prep()
bookings_recipe## ## ── Recipe ──────────────────────────────────────────────────────────────────────## ## ── Inputs## Number of variables by role## outcome: 1
## predictor: 12## ## ── Training information## Training data contained 83447 data points and no incomplete rows.## ## ── Operations## • Down-sampling based on: is_canceled | Trained## • Dummy variables from: assigned_room_type and deposit_type, ... | Trained## • Zero variance filter removed: assigned_room_type_L | Trainedtest_proc <- bake(bookings_recipe, new_data = bookings_test)Validación cruzada
Vamos a utilizar la validación cruzada Monte Carlo, con una proporción de 0.75, lo que indica que dividirá los datos en 75/25, y times=2, para que repita dos veces el sampling, lo que es especialmente importante para que no aparezcan NAs en las métricas de los modelos.
validation_splits <- mc_cv(juice(bookings_recipe), prop = 0.75, strata = is_canceled, times = 2)Clasificación por regresión logística
El primero modelo que vamos a crear va a ser un modelo de regresión logística, ya que el problema planteado es un problema de clasificación. Creamos un primer modelo con la función logistic_reg().
log_spec <- logistic_reg() %>%
set_engine("glm") %>%
set_mode("classification")Ahora vamos a generar las predicciones, aplicando al modelo la receta que hemos creado.
log_fit <- log_spec %>%
fit(is_canceled ~ .,
data = juice(bookings_recipe))
log_fit## parsnip model object
##
##
## Call: stats::glm(formula = is_canceled ~ ., family = stats::binomial,
## data = data)
##
## Coefficients:
## (Intercept) adr
## -2.509513 0.005470
## days_in_waiting_list lead_time
## -0.002869 0.004466
## stays_in_week_nights stays_in_weekend_nights
## 0.061301 0.097837
## adults booking_changes
## 0.187486 -0.408738
## previous_cancellations is_portuguese
## 0.903677 1.531922
## assigned_room_type_B assigned_room_type_C
## -0.389592 -0.902529
## assigned_room_type_D assigned_room_type_E
## -0.457857 -0.606012
## assigned_room_type_F assigned_room_type_G
## -0.644760 -0.620807
## assigned_room_type_H assigned_room_type_I
## -0.703269 -3.245616
## assigned_room_type_K deposit_type_Non.Refund
## -1.760004 5.337590
## deposit_type_Refundable market_segment_Complementary
## 0.039155 -0.846921
## market_segment_Corporate market_segment_Direct
## -0.840499 -0.611087
## market_segment_Groups market_segment_Offline.TA.TO
## -0.370090 -0.842693
## market_segment_Online.TA market_segment_Undefined
## 0.824148 10.145153
##
## Degrees of Freedom: 61793 Total (i.e. Null); 61766 Residual
## Null Deviance: 85660
## Residual Deviance: 58500 AIC: 58560Resultados de la primera clasificación por regresión logística
Para observar las métricas de nuestro modelo y poder realizar una primera valoración, vamos a utilizar la función fit_resamples().
log_res <- fit_resamples(
log_spec,
is_canceled ~ .,
validation_splits,
control = control_resamples(save_pred = TRUE)
)
log_res %>%
collect_metrics()## # A tibble: 3 × 6
## .metric .estimator mean n std_err .config
## <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 accuracy binary 0.731 2 0.00172 Preprocessor1_Model1
## 2 brier_class binary 0.163 2 0.000743 Preprocessor1_Model1
## 3 roc_auc binary 0.837 2 0.00139 Preprocessor1_Model1La precisión es de 0.73 y el área bajo la curva es de 0.83, mientras que el error estándar es muy bajo. Estos datos no están mal para ser un primer modelo, pero añadiendo alguna variable más hay mucho margen de mejora. Por lo tanto vamos a añadir más variables y a repetir el proceso, probando además otros modelos como árbol de decisión o random forest.
Creación de un segundo modelo
Vamos a tomar como base las variables que hemos escogido para el primero, pero además añadiremos total_of_special_requests y required_car_parking_spaces. En total tendremos un total de 14 variables predictoras.
bookings_for_model_imp <- bookings %>%
mutate_if(is.character, as.factor) %>%
mutate(is_canceled = as.factor(is_canceled)) %>%
select(assigned_room_type, adr, days_in_waiting_list, lead_time, stays_in_week_nights, stays_in_weekend_nights, adults, booking_changes, deposit_type, market_segment, previous_cancellations, total_of_special_requests, required_car_parking_spaces, is_portuguese, is_canceled)Particionamiento entrenamiento vs test 70/30
Volvemos a realizar el particionamiento, pero esta vez con los datos del nuevo dataframe.
set.seed(1234)
bookings_train <- sample_frac(bookings_for_model_imp, .7)
bookings_test <- setdiff(bookings_for_model_imp, bookings_train)Receta
Repetimos la creación de la receta con los nuevos datos de entrenamiento y volvemos a crear la variable test_proc.
bookings_recipe <- recipe(is_canceled ~ ., data = bookings_train) %>%
step_downsample(is_canceled) %>%
step_dummy(all_nominal(), - all_outcomes()) %>%
step_zv(all_numeric()) %>%
prep()
bookings_recipe## ## ── Recipe ──────────────────────────────────────────────────────────────────────## ## ── Inputs## Number of variables by role## outcome: 1
## predictor: 14## ## ── Training information## Training data contained 83447 data points and no incomplete rows.## ## ── Operations## • Down-sampling based on: is_canceled | Trained## • Dummy variables from: assigned_room_type and deposit_type, ... | Trained## • Zero variance filter removed: assigned_room_type_L | Trainedtest_proc <- bake(bookings_recipe, new_data = bookings_test)Validación cruzada
Establecemos la validación cruzada con los mismos parámetros que anteriormente.
validation_splits <- mc_cv(juice(bookings_recipe), prop = 0.75, strata = is_canceled, times = 2)Clasificación por regresión logística mejorada
Creamos de nuevo el modelo de regresión logística:
log_spec <- logistic_reg() %>%
set_engine("glm") %>%
set_mode("classification")Y generamos las predicciones en base a la nueva receta.
log_fit <- log_spec %>%
fit(is_canceled ~ .,
data = juice(bookings_recipe))## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurredlog_fit## parsnip model object
##
##
## Call: stats::glm(formula = is_canceled ~ ., family = stats::binomial,
## data = data)
##
## Coefficients:
## (Intercept) adr
## -2.717495 0.006662
## days_in_waiting_list lead_time
## -0.003896 0.004856
## stays_in_week_nights stays_in_weekend_nights
## 0.071203 0.116174
## adults booking_changes
## 0.268183 -0.347177
## previous_cancellations total_of_special_requests
## 0.998528 -0.797612
## required_car_parking_spaces is_portuguese
## -21.018089 1.751220
## assigned_room_type_B assigned_room_type_C
## -0.388061 -0.790203
## assigned_room_type_D assigned_room_type_E
## -0.445113 -0.425685
## assigned_room_type_F assigned_room_type_G
## -0.691895 -0.598387
## assigned_room_type_H assigned_room_type_I
## -0.710786 -3.218589
## assigned_room_type_K deposit_type_Non.Refund
## -1.681254 5.121188
## deposit_type_Refundable market_segment_Complementary
## 0.046776 -0.325603
## market_segment_Corporate market_segment_Direct
## -0.806856 -0.305556
## market_segment_Groups market_segment_Offline.TA.TO
## -0.546881 -0.892088
## market_segment_Online.TA market_segment_Undefined
## 1.360160 19.443906
##
## Degrees of Freedom: 61793 Total (i.e. Null); 61764 Residual
## Null Deviance: 85660
## Residual Deviance: 52420 AIC: 52480Clasificación por KNN
Vamos a crear un modelo KNN, para ello utilizaremos la función nearest_neighbor() y en set_engine le indicaremos ‘kknn’. El modo no va a variar, ya que es un problema de clasificación.
knn_spec <- nearest_neighbor() %>%
set_engine("kknn") %>%
set_mode("classification")Generamos las predicciones del mismo modo que para regresión logística:
knn_fit <- knn_spec %>%
fit(is_canceled ~ .,
data = juice(bookings_recipe))
knn_fit## parsnip model object
##
##
## Call:
## kknn::train.kknn(formula = is_canceled ~ ., data = data, ks = min_rows(5, data, 5))
##
## Type of response variable: nominal
## Minimal misclassification: 0.1911189
## Best kernel: optimal
## Best k: 5Clasificación por árbol de decisión
La clasificación por árbol de decisión también puede ofrecer muy buenos resultados en problemas de este tipo, por lo que también haremos un modelo de árbol de decisión, con la función decision_tree():
tree_spec <- decision_tree() %>%
set_engine("rpart") %>%
set_mode("classification")Creamos las predicciones:
tree_fit <- tree_spec %>%
fit(is_canceled ~ .,
data = juice(bookings_recipe))
tree_fit## parsnip model object
##
## n= 61794
##
## node), split, n, loss, yval, (yprob)
## * denotes terminal node
##
## 1) root 61794 30897 0 (0.500000000 0.500000000)
## 2) deposit_type_Non.Refund< 0.5 51615 20755 0 (0.597888211 0.402111789)
## 4) lead_time< 11.5 10427 1810 0 (0.826412199 0.173587801) *
## 5) lead_time>=11.5 41188 18945 0 (0.540035933 0.459964067)
## 10) is_portuguese< 0.5 28493 10758 0 (0.622433580 0.377566420)
## 20) market_segment_Online.TA< 0.5 9039 634 0 (0.929859498 0.070140502) *
## 21) market_segment_Online.TA>=0.5 19454 9330 1 (0.479592886 0.520407114)
## 42) total_of_special_requests>=0.5 11392 4213 0 (0.630179073 0.369820927) *
## 43) total_of_special_requests< 0.5 8062 2151 1 (0.266807244 0.733192756) *
## 11) is_portuguese>=0.5 12695 4508 1 (0.355100433 0.644899567)
## 22) required_car_parking_spaces>=0.5 621 0 0 (1.000000000 0.000000000) *
## 23) required_car_parking_spaces< 0.5 12074 3887 1 (0.321931423 0.678068577) *
## 3) deposit_type_Non.Refund>=0.5 10179 37 1 (0.003634935 0.996365065) *Clasificación por random forest
Por último, haremos un modelo utilizando random forest. La función a utilizar en este caso es rand_forest():
random_forest_spec <- rand_forest() %>%
set_engine("ranger", importance = "impurity") %>%
set_mode("classification")Generamos las predicciones como hemos hecho hasta ahora:
rf_fit <- random_forest_spec %>%
fit(is_canceled ~ .,
data = juice(bookings_recipe))
rf_fit## parsnip model object
##
## Ranger result
##
## Call:
## ranger::ranger(x = maybe_data_frame(x), y = y, importance = ~"impurity", num.threads = 1, verbose = FALSE, seed = sample.int(10^5, 1), probability = TRUE)
##
## Type: Probability estimation
## Number of trees: 500
## Sample size: 61794
## Number of independent variables: 29
## Mtry: 5
## Target node size: 10
## Variable importance mode: impurity
## Splitrule: gini
## OOB prediction error (Brier s.): 0.1106386Fase 4: Evaluación de resultados
Una vez creados los cuatro modelos, llega el turno de evaluar los resultados de cada uno de ellos.
- Resultados regresión logística
log_res <- fit_resamples(
log_spec,
is_canceled ~ .,
validation_splits,
control = control_resamples(save_pred = TRUE)
)
log_res %>%
collect_metrics()## # A tibble: 3 × 6
## .metric .estimator mean n std_err .config
## <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 accuracy binary 0.785 2 0.000227 Preprocessor1_Model1
## 2 brier_class binary 0.143 2 0.0000497 Preprocessor1_Model1
## 3 roc_auc binary 0.878 2 0.0000955 Preprocessor1_Model1- Resultados KNN
knn_res <- fit_resamples(
knn_spec,
is_canceled ~ .,
validation_splits,
control = control_resamples(save_pred = TRUE)
)
knn_res %>%
collect_metrics()## # A tibble: 3 × 6
## .metric .estimator mean n std_err .config
## <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 accuracy binary 0.799 2 0.00188 Preprocessor1_Model1
## 2 brier_class binary 0.142 2 0.00116 Preprocessor1_Model1
## 3 roc_auc binary 0.889 2 0.00146 Preprocessor1_Model1- Resultados árbol de decisión
tree_res <- fit_resamples(
tree_spec,
is_canceled ~ .,
validation_splits,
control = control_resamples(save_pred = TRUE)
)
tree_res %>%
collect_metrics()## # A tibble: 3 × 6
## .metric .estimator mean n std_err .config
## <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 accuracy binary 0.792 2 0.00191 Preprocessor1_Model1
## 2 brier_class binary 0.146 2 0.00116 Preprocessor1_Model1
## 3 roc_auc binary 0.863 2 0.00204 Preprocessor1_Model1- Resultados random forest
random_forest_res <- fit_resamples(
random_forest_spec,
is_canceled ~ .,
validation_splits,
control = control_resamples(save_pred = TRUE)
)## Warning: package 'ranger' was built under R version 4.2.3random_forest_res %>%
collect_metrics()## # A tibble: 3 × 6
## .metric .estimator mean n std_err .config
## <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 accuracy binary 0.839 2 0.000259 Preprocessor1_Model1
## 2 brier_class binary 0.112 2 0.000822 Preprocessor1_Model1
## 3 roc_auc binary 0.927 2 0.00129 Preprocessor1_Model1Como hemos podido ver, el modelo creado con random forest es el que mejores resultados ofrece, por lo tanto nos centraremos en él para sacar conclusiones.
Matriz de confusión - Datos de entrenamiento
Vamos a sacar la matriz de confusión para ver la cantidad de falsos negativos y falsos positivos.
En primer lugar vamos a almacenar en una variable las predicciones que sacamos de los resultados.
rf_pred <- random_forest_res %>%
collect_predictions()A partir de esa variable podremos obtener la matriz de confusión:
rf_pred %>%
conf_mat(truth = is_canceled, estimate = .pred_class) %>%
autoplot("heatmap") +
scale_fill_gradient(low = "salmon", high = "blue")## Scale for fill is already present.
## Adding another scale for fill, which will replace the existing scale.
Vemos que de la cantidad de falsos positivos/negativos es ampliamente inferior a la de verdaderos positivos/negativos.
Prueba del modelo con los datos de test
Lo último que haremos para evaluar nuestro modelo, será realizar pruebas con los datos de test. De modo que seguiremos realizando los pasos que hemos hecho hasta ahora.
Particionamiento entrenamiento vs test 70/30
set.seed(1234)
bookings_train <- sample_frac(bookings_for_model_imp, .7)
bookings_test <- setdiff(bookings_for_model_imp, bookings_train)Receta
bookings_recipe_test <- recipe(is_canceled ~ ., data = bookings_test) %>%
step_downsample(is_canceled) %>%
step_dummy(all_nominal(), - all_outcomes()) %>%
step_zv(all_numeric()) %>%
prep()
bookings_recipe_test## ## ── Recipe ──────────────────────────────────────────────────────────────────────## ## ── Inputs## Number of variables by role## outcome: 1
## predictor: 14## ## ── Training information## Training data contained 22868 data points and no incomplete rows.## ## ── Operations## • Down-sampling based on: is_canceled | Trained## • Dummy variables from: assigned_room_type and deposit_type, ... | Trained## • Zero variance filter removed: <none> | Trainedtest_proc <- bake(bookings_recipe_test, new_data = bookings_test)Validación cruzada
validation_splits <- mc_cv(juice(bookings_recipe_test), prop = 0.75, strata = is_canceled, times = 2)Predicciones con los datos de test
rf_fit_test <- random_forest_spec %>%
fit(is_canceled ~ .,
data = juice(bookings_recipe_test))
rf_fit_test## parsnip model object
##
## Ranger result
##
## Call:
## ranger::ranger(x = maybe_data_frame(x), y = y, importance = ~"impurity", num.threads = 1, verbose = FALSE, seed = sample.int(10^5, 1), probability = TRUE)
##
## Type: Probability estimation
## Number of trees: 500
## Sample size: 12544
## Number of independent variables: 30
## Mtry: 5
## Target node size: 10
## Variable importance mode: impurity
## Splitrule: gini
## OOB prediction error (Brier s.): 0.1587809Resultados del modelo aplicado a los datos de test
random_forest_res_test <- fit_resamples(
random_forest_spec,
is_canceled ~ .,
validation_splits,
control = control_resamples(save_pred = TRUE)
)
random_forest_res_test %>%
collect_metrics()## # A tibble: 3 × 6
## .metric .estimator mean n std_err .config
## <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 accuracy binary 0.759 2 0.000957 Preprocessor1_Model1
## 2 brier_class binary 0.158 2 0.00177 Preprocessor1_Model1
## 3 roc_auc binary 0.850 2 0.00491 Preprocessor1_Model1Matriz de confusión - Datos de test
rf_pred_test <- random_forest_res_test %>%
collect_predictions()rf_pred_test %>%
conf_mat(truth = is_canceled, estimate = .pred_class) %>%
autoplot("heatmap") +
scale_fill_gradient(low = "salmon", high = "blue")## Scale for fill is already present.
## Adding another scale for fill, which will replace the existing scale.
Aquí, al haber menos datos, se aprecia una menor diferencia que con los datos de entrenamiento pero siguen siendo verdaderos positivos y negativos la mayoría.
Fase 5: Elección del modelo con mejores resultados y conclusiones
Hemos utilizado 4 modelos diferentes de clasificación para realizar pruebas y determinar cuál es el que ofrece mejores resultados para el problema planteado, los 4 modelos utilizados han sido:
- Regresión logística
- KNN
- Árbol de decisión
- Random fores
En primer lugar, creamos un modelo con 12 variables predictoras, que nos daba un área bajo la curva alrededor de 0.84, lo cual no está nada mal para ser un primer modelo, pero todavía había margen de mejora.
Por tanto, después creamos un segundo modelo mejorado, añadiendo dos variables predictoras más, lo que hace un total de 14 varables predictoras. En este caso, el área bajo la curva para regresión logística ha aumentado hasta casi 0.88.
Los resultados para KNN y árbol de decisión han sido 0.888 y 0.86 respectivamente, por lo que KNN se situaría como el segundo mejor modelo mientras que el de árbol de decisión en este caso ha obtenido los peores resultados.
Finalmente, el que mejores resultados nos ofrece, es Random forest, con una precisión de casi 0.84, un área bajo la curva de 0.927 y un error estándar también algo inferior en comparación con los modelos anteriores.