¿Como modelar con tidymodels?

Se presentan algunos ejemplos de modelos de clasificación no lineales para que se pueda observar la utilización del nuevo conjunto de herramientas para Ciencia de Datos.

Primer ejemplo

Para el ejemplo inicial se ha transcrito, con importantes ajustes, la explicación asociada a un comando de rsample1, la cual utiliza la base de datos attrition. Proviene del IBM Watson Analytics Lab. Se describe como data que devela los factores que llevan a los empleados a retirarse del empleo. Consta de 1470 registros.

data("attrition", package = "modeldata")
names(attrition)
##  [1] "Age"                      "Attrition"               
##  [3] "BusinessTravel"           "DailyRate"               
##  [5] "Department"               "DistanceFromHome"        
##  [7] "Education"                "EducationField"          
##  [9] "EnvironmentSatisfaction"  "Gender"                  
## [11] "HourlyRate"               "JobInvolvement"          
## [13] "JobLevel"                 "JobRole"                 
## [15] "JobSatisfaction"          "MaritalStatus"           
## [17] "MonthlyIncome"            "MonthlyRate"             
## [19] "NumCompaniesWorked"       "OverTime"                
## [21] "PercentSalaryHike"        "PerformanceRating"       
## [23] "RelationshipSatisfaction" "StockOptionLevel"        
## [25] "TotalWorkingYears"        "TrainingTimesLastYear"   
## [27] "WorkLifeBalance"          "YearsAtCompany"          
## [29] "YearsInCurrentRole"       "YearsSinceLastPromotion" 
## [31] "YearsWithCurrManager"
table(attrition$Attrition)
## 
##   No  Yes 
## 1233  237

El primer paso es separar la base de datos en dos subconjuntos, el que se utiliza para el entrenamiento (training) y el de validación (validate). Se realiza con el comando initial_split() del paquete rsample. Luego se podrá separar el conjunto de datos training, el que propiamente se utiliza para el entrenamiento del modelo y el que se utiliza para validar las métricas del modelo con el objeto de ayudar a seleccionar, de entre un conjunto de modelos aplicados sobre una misma data, cuál funciona mejor. Podría separarse en dos subconjuntos (training y test), pero es mejor práctica realizar una validación cruzada, así que se aplicará ésta.

att_split <- initial_split(attrition, prop = .8)
data_train <- training(att_split)
data_validate <- testing(att_split)

El comando initial_split() separa por defecto en train=75/validate=25. Aquí se especificó train=80/validate=20.

Por ser la variable de interés attrition una variable categórica, se elige un modelo de regresión logística como el adecuado para iniciar. Será la línea de base. Por simplicidad, se modela con sólo tres variables.

variables <- c('Attrition', 'Gender', 'MonthlyIncome', 'JobSatisfaction')
data_train <- data_train %>% dplyr::select(all_of(variables))
data_validate <- data_validate %>% dplyr::select(all_of(variables))
glm(Attrition ~ JobSatisfaction + Gender + MonthlyIncome, data = data_train, family = binomial)
## 
## Call:  glm(formula = Attrition ~ JobSatisfaction + Gender + MonthlyIncome, 
##     family = binomial, data = data_train)
## 
## Coefficients:
##       (Intercept)  JobSatisfaction.L  JobSatisfaction.Q  JobSatisfaction.C  
##         -1.106682          -0.676307          -0.055543          -0.152725  
##        GenderMale      MonthlyIncome  
##          0.236362          -0.000126  
## 
## Degrees of Freedom: 1176 Total (i.e. Null);  1171 Residual
## Null Deviance:       1017 
## Residual Deviance: 963.5     AIC: 975.5

Por conveniencia, se creará un objeto que guarda la fórmula para ser usada posteriormente:

model_formula <- as.formula(Attrition ~ JobSatisfaction + Gender + MonthlyIncome)

Para evaluar el modelo se usan 10 validaciones cruzadas, cada una de 10 folds, por tanto se obtendrán 100 muestras para evaluar la precisión del modelo.

set.seed(4622)
folds_with_repeats <- rsample::vfold_cv(data = data_train, v = 10, repeats = 10)
print(folds_with_repeats)
## #  10-fold cross-validation repeated 10 times 
## # A tibble: 100 x 3
##    splits             id       id2   
##    <list>             <chr>    <chr> 
##  1 <split [1.1K/118]> Repeat01 Fold01
##  2 <split [1.1K/118]> Repeat01 Fold02
##  3 <split [1.1K/118]> Repeat01 Fold03
##  4 <split [1.1K/118]> Repeat01 Fold04
##  5 <split [1.1K/118]> Repeat01 Fold05
##  6 <split [1.1K/118]> Repeat01 Fold06
##  7 <split [1.1K/118]> Repeat01 Fold07
##  8 <split [1.1K/117]> Repeat01 Fold08
##  9 <split [1.1K/117]> Repeat01 Fold09
## 10 <split [1.1K/117]> Repeat01 Fold10
## # ... with 90 more rows

El comando rsample::vfold_cv() genera los 100 subconjuntos de igual tamaño, variando la semilla de la subdivisión. El resultado genera un objeto con tres variables: la repetición (‘id’), el fold (‘id2’) y splits. splits es una variable que almacena una lista de tamaño (v*repeats)*3, donde v es el número de folds y repeat el número de repeticiones. Cada elemento de la lista es un objeto de la clase tibble con varias variables:

data: matriz de dimensión n*p, dónde n es el número de registros y p el número de variables.

in_id: índice lógico de la posición de los registros de la data que pertenecen a los registros que se analizan. Denomina analisys a este conjunto de registros, y al dejado aparte (hold out) lo denomina assessment2. Ver diagrama.

id: almacena el id del fold y de la repetición a la que corresponde la data.

Para accesar al subconjunto denominado analysis, se usa el comando analysis() y para el otro subconjunto, el comando assessment().

Imagen tomada de https://rviews.rstudio.com/2020/04/21/the-case-for-tidymodels/

Primero se define una función que aplique sobre cada muestra los siguientes pasos:

  1. La regresión logística al conjunto de datos de analisys.
  2. Un ajuste del modelo logístico.
  3. Una predicción sobre el conjunto de datos assessment, ayudado por medio del paquete broom.
  4. Asignación lógica asociada a si el registro fue predicho correctamente.
## splits es una de las muestras generadas en la cv
res_set_aparte <- function(splits, modelo) {
  # Ajusta el modelo
  model <- glm(modelo, data = analysis(splits), family = binomial)
  # Identifique el conjunto de datos dejado aparte.
  holdout <- assessment(splits)
  # El comando augment() realiza la prediccion sobre el conjunto de datos dejado aparte.
  res <- broom::augment(model, newdata = holdout)
  # lvls seran los niveles del factor con las predicciones
  lvls <- levels(holdout$Attrition)
  predictions <- factor(ifelse(res$.fitted > 0, lvls[2], lvls[1]), levels = lvls)
  # Calcule la variable lógica que establece para cada registro si la prediccion es correcta
  res$correcto <- predictions == holdout$Attrition
  # Retorna el conjunto de datos con las columnas adicionales
  res
}

Obsérvese que la predicción (res$.fitted) se transforma en una variable nominal dependiendo de si es mayor o menor a cero. Si mayor a cero, entonces No (se va del trabajo), de lo contrario, Yes.

Por ejemplo, para la primera de las 100 muestras:

ejemplo <- res_set_aparte(folds_with_repeats$splits[[1]], model_formula)

print(dim(ejemplo))
## [1] 118   7
dim(rsample::assessment(folds_with_repeats$splits[[1]]))
## [1] 118   4
## nombre de las columnas añadidas:
print(ejemplo[1:6, setdiff(names(ejemplo), names(attrition))])
## # A tibble: 6 x 3
##   .rownames .fitted correcto
##   <chr>       <dbl> <lgl>   
## 1 10         -0.863 TRUE    
## 2 19         -1.12  FALSE   
## 3 26         -1.60  TRUE    
## 4 39         -1.63  TRUE    
## 5 51         -0.972 TRUE    
## 6 68         -1.27  TRUE

Para este modelo, el valor .fitted es el predictor lineal del log-odds.

Para calcular el resultado sobre las 100 muestras, se usa el comando map() del paquete purrr. Es un comando muy útil cuando se tiene una lista de objetos. Aplica lo solicitado sobre cada lista.

# La función *map* transforma la entrada aplicando la funcion 
# a cada elemento de la lista 
folds_with_repeats$resultados <- purrr::map(folds_with_repeats$splits, res_set_aparte,
                                model_formula)
# Genera una lista por cada split
print(folds_with_repeats)
## #  10-fold cross-validation repeated 10 times 
## # A tibble: 100 x 4
##    splits             id       id2    resultados        
##    <list>             <chr>    <chr>  <list>            
##  1 <split [1.1K/118]> Repeat01 Fold01 <tibble [118 x 7]>
##  2 <split [1.1K/118]> Repeat01 Fold02 <tibble [118 x 7]>
##  3 <split [1.1K/118]> Repeat01 Fold03 <tibble [118 x 7]>
##  4 <split [1.1K/118]> Repeat01 Fold04 <tibble [118 x 7]>
##  5 <split [1.1K/118]> Repeat01 Fold05 <tibble [118 x 7]>
##  6 <split [1.1K/118]> Repeat01 Fold06 <tibble [118 x 7]>
##  7 <split [1.1K/118]> Repeat01 Fold07 <tibble [118 x 7]>
##  8 <split [1.1K/117]> Repeat01 Fold08 <tibble [117 x 7]>
##  9 <split [1.1K/117]> Repeat01 Fold09 <tibble [117 x 7]>
## 10 <split [1.1K/117]> Repeat01 Fold10 <tibble [117 x 7]>
## # ... with 90 more rows

Ahora podemos calcular la métrica para todos los conjuntos de datos de assessment. Se calcula el porcentaje de registros con una predicción correcta:

folds_with_repeats$resultados <- map_dbl(folds_with_repeats$resultados, function(x) {mean(x$correcto)})
boxplot(folds_with_repeats$resultados, main = 'Distribución de la predicción de acuerdo a los cien subconjuntos', col = 'wheat', ylab = 'Porcentaje', cex.main = 0.8, cex.lab = 0.8, cex.axis = 0.8)

La accuracy a sobrepasar es esta línea de base, 0.8445024. Queda ilustrado un modelamiento básico.

Segundo ejemplo

El segundo ejemplo trabaja sobre un conjunto de datos simulados, para comparar diferentes métodos de clasificación. Es una transcripción del ejemplo asociado a la librería discrim3, pero ajustada, por ejemplo, al uso de la librería ggformula en vez de ggplot2.

Este ejemplo utiliza más ampliamente el paradigma tidymodels.

library(AppliedPredictiveModeling)
set.seed(855)
parabolic <- AppliedPredictiveModeling::quadBoundaryFunc(1500)
print(head(parabolic))
##          X1         X2       prob  class
## 1 0.6131216  0.1272546 0.09371977 Class2
## 2 1.7617845 -0.7437384 0.01732629 Class2
## 3 2.4981652  2.2086611 0.92496990 Class2
## 4 3.0761547  3.0271284 0.99990694 Class1
## 5 1.6600466  1.0245654 0.05886856 Class2
## 6 0.6092829 -0.4227656 0.12615966 Class2
ggformula::gf_point(X2 ~ X1, data = parabolic, color = ~class, alpha = 0.3)

Se dividen los datos. Esta vez en train=66/validate=33.

set.seed(115)
data_split <- rsample::initial_split(parabolic, prop = 2/3)
data_train <- rsample::training(data_split)
data_validate <- rsample::testing(data_split)

Se aplica el modelo de análisis discriminante flexible:

fda_mod <- discrim::discrim_flexible() %>% parsnip::set_engine("earth")
fda_mod %>% parsnip::translate()
## Flexible Discriminant Model Specification (classification)
## 
## Computational engine: earth 
## 
## Model fit template:
## mda::fda(formula = missing_arg(), data = missing_arg(), method = earth::earth)
model_formula <- as.formula(class ~ X1 + X2)

fda_fit <- fda_mod %>% fit(model_formula, data = data_train)
fda_fit 
## parsnip model object
## 
## Fit time:  81ms 
## Call:
## mda::fda(formula = class ~ X1 + X2, data = data, method = earth::earth)
## 
## Dimension: 1 
## 
## Percent Between-Group Variance Explained:
##  v1 
## 100 
## 
## Training Misclassification Error: 0.125 ( N = 1000 )

Obsérvese que el comando de la librería parsnip set_engine() especifica que se use la librería earth. Es una librería que ejecuta un algoritmo MARS (Multivariate Adaptative Regression Splines).

El comando translate() presenta el esquema funcional, sin referencia a la data. translate() traduce la especificación del modelo en un objeto de código que es específico para el motor particular especificado en set_engine().

La última línea de la salida presenta que el error es 12.5%, o en otras palabras, que el acierto (accuracy) es del 87.5%.

Como no se especificaron parámetros de ajuste del modelo, el algoritmo MARS sigue su propio método interno para optimizar la cantidad de características que se incluyen en el modelo. El modelo MARS subyacente es:

summary(fda_fit$fit$fit)
## Call: earth(x=x, y=Theta, weights=weights)
## 
##                  coefficients
## (Intercept)        -1.8008471
## h(X1- -0.820471)    0.6364784
## h(X1-0.599259)     -0.3075588
## h(X2- -1.76861)     1.3797082
## h(X2- -0.774763)   -1.3781369
## h(X2-0.85974)      -1.3368848
## h(X2-2.55583)       1.4516810
## 
## Selected 7 of 12 terms, and 2 of 2 predictors
## Termination condition: RSq changed by less than 0.001 at 12 terms
## Importance: X2, X1
## Number of terms at each degree of interaction: 1 6 (additive model)
## GCV 0.4269824    RSS 415.953    GRSq 0.5738712    RSq 0.584047

El límite clasificatorio, superpuesto al conjunto de prueba, muestra una serie de líneas segmentadas:

pred_grid <- expand.grid(X1 = seq(-5, 5, length = 100), X2 = seq(-5, 5, length = 100))

pred_grid <- dplyr::bind_cols(pred_grid, 
    data.frame(predict(fda_fit, pred_grid, type = "prob")) %>% 
      dplyr::select(.pred_Class1) %>% setNames("fda_pred"))

ggformula::gf_point(X2 ~ X1, data = parabolic, color = ~class, alpha = 0.3) +
 geom_contour(data = pred_grid, aes(z = fda_pred), breaks = .5, col = "black", linetype = 2)

La librería parsnip devuelve nombres de columnas estándar, independientemente del modelo que aplique. Es una manera de estandarizar la salida. Si esta realizando una regresión, devuelve variables \(.pred\). En clasificación devuelve variables \(.Class_i\)

La frontera presentada en la figura parece bastante razonable. No obstante, parece demasiado ajustada a los datos. Con nuevos datos, tal vez no funcione igual de bien.

Para demostrar cómo realizar el afinamiento (tuning), se optimizarán los parámetros del modelo mediante una búsqueda en una grilla por medio del paquete tune. Se cambiará de especificación: de MARS a RDA, buscando una solución menos ajustada a los datos, es decir, evitando el overfitting.

Primero, se establecen los parámetros a afinar:

rda_mod <- discrim::discrim_regularized(frac_common_cov = tune('frac_covar_en_comun'), 
          frac_identity = tune('frac_independencia')) %>% parsnip::set_engine("klaR")
rda_mod
## Regularized Discriminant Model Specification (classification)
## 
## Main Arguments:
##   frac_common_cov = tune("frac_covar_en_comun")
##   frac_identity = tune("frac_independencia")
## 
## Computational engine: klaR

El comando discrim_regularized() es una forma de generar la especificación de un modelo de análisis discriminante regularizado (RDA) antes de ajustarlo. RDA es una generalización del análisis discriminante lineal (LDA) y el análisis discriminante cuadrático (QDA), y permite generar una especificación intermedia. Sea el parámetro frac_identity (gamma en la librería klaR) igual a cero. Si el parámetro frac_common_cov (lambda en la librería klaR) se establece en 1, se está especificando un modelo LDA, pero si el parámetro alfa se establece en 0, este operador especifica un modelo QDA.

En vez de especificar el valor de los parámetros, se hace un llamado al comando tune() etiquetando de este modo qué variables deben afinarse. Ambas en este caso. Usualmente el comando tune() no contiene parámetros, pero se puede especificar un label asociado, como se realizó en el código anterior.

El comando set_engine() especifica que se use la implementación de RDA de la librería klaR.

Para el ajuste del modelo, requerimos una especificación de validación cruzada con remuestreo.

set.seed(2014)
folds <- rsample::vfold_cv(data_train, repeats = 1)

Se podrían haber especificado 10 repeticiones, pero no mejora el cálculo, así que se dejó en una sola.

Y una cuadrícula de valores candidatos.

# Cuadra la grilla de candidatos
candidatos <- rda_mod %>% dials::parameters() %>% 
  dials::grid_max_entropy(size = 10)
candidatos
## # A tibble: 10 x 2
##    frac_covar_en_comun frac_independencia
##                  <dbl>              <dbl>
##  1              0.891              0.983 
##  2              0.936              0.0132
##  3              0.0708             0.0157
##  4              0.720              0.697 
##  5              0.374              0.0817
##  6              0.0117             0.873 
##  7              0.625              0.361 
##  8              0.988              0.530 
##  9              0.125              0.430 
## 10              0.355              0.869

El comando parameters() de la librería dials contiene información acerca de los posibles valores, rangos, tipos, y otros aspectos del modelo de entrada.

También se procede a configurar un objeto de criterio de medida, el área bajo la curva ROC, por medio del paquete yardstick:

roc_values <- yardstick::metric_set(roc_auc)

Procédase a afinar el modelo:

rda_res <- tune::tune_grid(object = rda_mod, 
           preprocessor = model_formula, resamples = folds, 
           grid = candidatos, metrics = roc_values)

La afinación toma su tiempo en sacar los resultados. Obsérvese que el paradigma tidymodels ha postergado el trabajo con los datos en tanto se procedió a realizar la afinación. Previamente sólo se trataba de especificación del modelo.

auc_values <- tune::show_best(rda_res, 'roc_auc') %>% arrange(desc(mean)) 
auc_values %>% dplyr::slice(1:5)
## # A tibble: 5 x 8
##   frac_covar_en_c~ frac_independen~ .metric .estimator  mean     n std_err
##              <dbl>            <dbl> <chr>   <chr>      <dbl> <int>   <dbl>
## 1           0.0708           0.0157 roc_auc binary     0.939    10 0.00717
## 2           0.125            0.430  roc_auc binary     0.920    10 0.00686
## 3           0.374            0.0817 roc_auc binary     0.904    10 0.00827
## 4           0.0117           0.873  roc_auc binary     0.904    10 0.00762
## 5           0.355            0.869  roc_auc binary     0.881    10 0.0101 
## # ... with 1 more variable: .config <chr>

Visualización:

autoplot(rda_res)

El comando de visualización por defecto está en proceso de mejora. Por ahora, grafica cada variable, pero en este caso es mejor la visualización bivariada, así que se realiza una gráfica ad hoc:

ggplot(auc_values, aes(x = frac_covar_en_comun, y = frac_independencia, size = mean)) + 
  geom_point(alpha = .5) + coord_equal()

Obsérvese que el nombre de las variables en el objeto de salida del tune() es el fijado como parámetro.

Hay un amplio rango de combinaciones de parámetros asociados a una buen desempeño. Los bajos desempeños ocurren cuando la matriz de covarianza trata de asemejarse a un clasificador lineal LDA (frac_covar_en_comun > 0.9) al tiempo que frac_independencia > 0.6. Este último parámetro induce la matriz de covarianza hacia uno, es decir, a considerar los parámetros casi como independientes.

El modelo generado se actualiza ahora con la mejor combinación de los parámetros mediante el comando tune::select_best() (0.0708109, 0.0156591). Estos valores son casi un modelo Quadratic Discriminat Analysis (QDA) puro.

Con el comando merge() se añaden los valores al objeto parsnip4:

final_param <- mejores %>% 
  dplyr::select(frac_covar_en_comun, frac_independencia)
rda_mod <- rda_mod %>% merge(final_param) %>% 
  dplyr::pull(x) %>% purrr::pluck(1)
rda_fit <- rda_mod %>% fit(model_formula, data = data_train)
rda_fit$spec$method$libs
## [1] "klaR"
rda_fit$spec$method$fit$args
## $formula
## missing_arg()
## 
## $data
## missing_arg()
## 
## $lambda
## <quosure>
## expr: ^0.0708109484985471
## env:  empty
## 
## $gamma
## <quosure>
## expr: ^0.0156591301783919
## env:  empty

Visualización de la nueva frontera:

pred_grid <- bind_cols(pred_grid, 
    predict(rda_fit, pred_grid, type = "prob") %>% 
      dplyr::select(.pred_Class1) %>% setNames("rda_pred"))

ggformula::gf_point(X2 ~ X1, data = parabolic, color = ~class, alpha = 0.3) +
  geom_contour(data = pred_grid, aes(z = fda_pred), breaks = .5, 
          col = "black", alpha = .5, lty = 2) + 
  geom_contour(data = pred_grid, aes(z = rda_pred), breaks = .5, col = "black")

Esta bastante cercana a la frontera simulada, que es de naturaleza parabólica.

Los resultados del test son:

probs_test <- predict(rda_fit, data_validate, type = "prob") %>% 
  bind_cols(data_validate %>% dplyr::select(class))
probs_test
## # A tibble: 500 x 3
##    .pred_Class1 .pred_Class2 class 
##           <dbl>        <dbl> <fct> 
##  1        0.172        0.828 Class2
##  2        0.488        0.512 Class2
##  3        0.164        0.836 Class2
##  4        0.255        0.745 Class2
##  5        0.175        0.825 Class2
##  6        0.320        0.680 Class1
##  7        0.125        0.875 Class2
##  8        0.480        0.520 Class1
##  9        0.246        0.754 Class2
## 10        0.440        0.560 Class2
## # ... with 490 more rows
probs_test %>% roc_auc(class, .pred_Class1)
## # A tibble: 1 x 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 roc_auc binary         0.924
autoplot(roc_curve(probs_test, class, .pred_Class1))

Mucho mejor!

El caso es que el conjunto de piezas del rompecabezas de tidymodels está en desarrollo, así que puede ir variando5. Este documento es un abrebocas para los entusiastas de R.

¿Qué otros paquetes están asociados a tidymodels?

modeltime6: diseñado para desarrollar y probar rápidamente modelos de series de tiempo utilizando modelos de aprendizaje automático, modelos clásicos y modelos automatizados. Trabaja a su vez con timetk7, una colección de herramientas para trabajar con series temporales en R.

embed8: es un paquete que contiene pasos adicionales para el paquete recipes para incrustar predictores en una o más columnas numéricas. Todos los métodos de pre procesamiento son supervisados.

themis9: Una colección de nuevos pasos para el paquete recipes que maneja datos no balanceados.

rules10: paquete adyacente a parsnip que permite tres modelos específicos:

baguette11: funciones básicas y envolturas para parsnip para bagging (también conocido como agregación bootstrap).

poissonreg12: Envoltura para parsnip que tiene métodos para modelos lineales, logísticos y multinomiales, en donde el resultado es un recuento. Hay motores para glm, rstanarm, glmnet, hurdle y zeroinfl. Los dos últimos permiten modelos de Poisson cero inflados desde el paquete pscl.

plsmod13: envolturas para parsnip que tiene modelos de regresión y clasificación de mínimos cuadrados parciales (PLS) basados en el trabajo del paquete mixOmics de Bioconductor. Este paquete facilita los modelos PLS comunes y las versiones sparse. Además, también se puede utilizar para modelos multivariados.

discrim, ya mencionado, que contiene enlaces de parsnip para modelos de clasificación.

vip14: Algunos algoritmos modernos, como Random Forest y los árboles de decisión impulsados por gradientes (gradient boosted decision trees), tienen una forma natural de cuantificar la importancia o la influencia relativa de cada variable (feature). Otros algoritmos, como los clasificadores ingenuos de Bayes (naive bayes) y las máquinas de soporte vectorial (SVM), no son capaces de hacerlo y, por lo general, se utilizan acercamientos independientes del modelo para medir la importancia de cada predictor. vip es un paquete para construir puntuaciones / plots de importancia de las variables (VI) para muchos tipos de algoritmos de aprendizaje supervisado utilizando enfoques agnósticos novedosos y específicos a cada modelo.

  1. https://cran.r-project.org/web/packages/rsample/vignettes/Working_with_rsets.html↩︎

  2. En inglés, algunas veces la llaman muestra “out-of-bag” (OOB).↩︎

  3. https://www.tidyverse.org/blog/2019/10/discrim-0-0-1/↩︎

  4. dplyr::pull() es similar al $, pero se ve mejor en los comandos pipe del paquete magrittr. purrr::pluck () es similar a doble corchete cuadrado el cual permite indexar de manera profunda y flexible las estructuras de datos↩︎

  5. De hecho la página fuente usaba un comando denominado estimate que no existe en la versión actual y debí indagar para encontrar el equivalente.↩︎

  6. https://www.r-bloggers.com/introducing-modeltime-tidy-time-series-forecasting-using-tidymodels/↩︎

  7. https://cran.r-project.org/web/packages/timetk/vignettes/TK03_Forecasting_Using_Time_Series_Signature.html↩︎

  8. https://www.tidyverse.org/blog/2020/06/embed-0-1-0/↩︎

  9. https://www.tidyverse.org/blog/2020/02/themis-0-1-0/↩︎

  10. https://www.tidyverse.org/blog/2020/05/rules-0-0-1/↩︎

  11. https://www.tidyverse.org/blog/2020/04/parsnip-adjacent/↩︎

  12. https://www.tidyverse.org/blog/2020/04/parsnip-adjacent/↩︎

  13. https://www.tidyverse.org/blog/2020/04/parsnip-adjacent/↩︎

  14. https://koalaverse.github.io/vip/articles/vip.html↩︎