Maestría en Explotación de datos y Descubrimiento de conocimiento


Materia: Enfoque Estadistico del Aprendizaje
Trabajo práctico 1: Regresion Lineal (Parte 2)

2. Modelo inicial

2.1. Se plantea la siguiente primera alternativa para modelar el peso:

\(E(peso) = \beta_0 + \beta_1 * altura + \beta_2 * edad + \beta_3 * genero + \beta4 * diasActividadFisicaSemanal + \beta5 * consumoDiarioAlcohol\)

Primero se cargan las librerías necesarias:

options(warn=-1)
rm(list=ls())
gc()

##          used (Mb) gc trigger (Mb) max used (Mb)
## Ncells 469536 25.1    1008755 53.9   665551 35.6
## Vcells 884367  6.8    8388608 64.0  1813089 13.9

options(warn=-2)

# install.packages("pacman") -- Descomentar par instalar pacman
library(pacman)
p_load_gh('adrianmarino/commons')
import('../src/dataset.R')

## [1] "-> '../src/dataset.R' script loadded successfuly!"

import('../src/preprocessing.R')

## [1] "-> '../src/preprocessing.R' script loadded successfuly!"

import('../src/model.R')

## [1] "-> '../src/model.R' script loadded successfuly!"

import('../src/plot.R')

## [1] "-> '../src/plot.R' script loadded successfuly!"

A continuación se carga los conjuntos de entrenamiento y test. También se resumen los valores de las variables categóricas y transforman todo los valores faltantes a NA’s ya que el modelo lineal los excluye por defecto. Fuera de esto tenemos muy pocos valores faltantes como ya analizamos en la notebook eea-tp1-eda.

train_set <- load_train_set() %>% 
  preprocess() %>% 
  shorten_values() %>%
  process_missings()

test_set <- load_test_set() %>% 
  preprocess() %>% 
  shorten_values() %>% 
  process_missings()

Validamos la estrutura de los datos:

glimpse(train_set)

## Rows: 7,024
## Columns: 15
## $ edad                          <int> 17, 15, 15, 16, 17, 15, 13, 17, 17, 16, …
## $ genero                        <fct> Femenino, Masculino, Masculino, Masculin…
## $ nivel_educativo               <ord> 2, 1, 2, 1, 2, 1, 9, 9, 1, 3, 3, 8, 9, 3…
## $ altura                        <int> 165, 178, 172, 170, 170, 178, 156, 163, …
## $ peso                          <int> 62, 62, 62, 65, 75, 88, 46, 60, 57, 51, …
## $ frecuencia_hambre_mensual     <ord> Rara vez, Rara vez, Nunca, Nunca, Rara v…
## $ dias_consumo_comida_rapida    <int> 0, 0, 3, 1, 1, 2, 0, 0, 0, 3, 4, 2, 1, 1…
## $ edad_consumo_alcohol          <ord> 14-15, <=7, 0, 14-15, 16-17, 8-9, 10-11,…
## $ consumo_diario_alcohol        <dbl> 5.0, 4.0, 0.0, 0.0, 0.0, 5.0, 1.0, 0.5, …
## $ dias_actividad_fisica_semanal <int> 7, 7, 7, 7, 0, 7, 0, 2, 7, 3, 2, 2, 7, 1…
## $ consumo_semanal_frutas        <ord> 0, 0, 0, 4-6, 14, 7, 14, 21, 0, 14, <=3,…
## $ consumo_semanal_verdura       <ord> 4-6, 4-6, 7, >=28, <=3, 14, 4-6, 7, 0, 4…
## $ consumo_semanal_gaseosas      <ord> <=3, <=3, 4-6, <=3, 7, 4-6, 0, 7, <=3, 4…
## $ consumo_semanal_snacks        <ord> <=3, 0, 4-6, <=3, 0, 4-6, 0, <=3, 0, <=3…
## $ consumo_semanal_comida_grasa  <ord> 0, 4-6, 0, 0, <=3, 4-6, <=3, 7, 0, <=3, …

Las variables categóricas se transforman a factores ordinales o no ordinales según sea el caso.

A contunuación se fija la semilla y se validan las proporciones de los conjuntos de entrenamiento y test:

set.seed(25)
show_train_test_props(train_set, test_set)

## [1] "Train: 70%, Test: 30%"

Modelo 1

Se plantea el primer modelo lineal:

model_1 <- lm(
  peso ~ altura + edad + genero + dias_actividad_fisica_semanal + consumo_diario_alcohol, 
  data = train_set
)

2.2. ¿Cuál es la interpretación de cada uno de los coeficientes estimados?

Veamos a continuación un resumen de los coeficiente del modelo 1:

coefficients_summary(model_1)

## _________________________________________________________________________________________________________ 
## term                               estimate  std.error   statistic       p.value    conf.low    conf.high 
## ========================================================================================================= 
## (Intercept)                   -68.922688070 2.33805445 -29.4786497 3.614866e-180 -73.5059810 -64.33939510 
## altura                          0.650606544 0.01437975  45.2446353  0.000000e+00   0.6224179   0.67879520 
## edad                            1.406727060 0.09385081  14.9889709  5.121599e-50   1.2227511   1.59070300 
## generoMasculino                 1.262643558 0.27282821   4.6279802  3.758831e-06   0.7278179   1.79746926 
## dias_actividad_fisica_semanal  -0.087391031 0.04992917  -1.7503000  8.011025e-02  -0.1852673   0.01048523 
## consumo_diario_alcohol          0.007271379 0.06138558   0.1184542  9.057112e-01  -0.1130629   0.12760566 
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ 
## ____________________________________________________________________________ 
## Termino                         Coeficiente Signiticativo IC incluye al cero 
## ============================================================================ 
## (Intercept)                   -68.922688070            Si                 No 
## altura                          0.650606544            Si                 No 
## edad                            1.406727060            Si                 No 
## generoMasculino                 1.262643558            Si                 No 
## dias_actividad_fisica_semanal  -0.087391031            No                 Si 
## consumo_diario_alcohol          0.007271379            No                 Si 
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

Al analizar cada coeficiente se encuentra que:

• \(\hat{\beta_0}\) (Ordenada al origen) de valor -68.92 Kg, es el peso esperado o promedio de un individuo de genero femenino que tiene cero altura, edad, actividad física y consumo diario de alcohol. Esto no es interpretable, ya que una persona tiene que tener una altura superior a cero y no puede tener un peso negativo, pero si podría no realizar actividad física ni consumir alcohol.

• El coeficiente \(\hat{\beta_1}\) de valor 650 gramos, corresponde a la altura del individuo. Este coeficiente indica que dada una edad, genero, consumo de alcohol diario y días de actividad física semanal fijos, cada incremento en 1 cm adicional en la altura del individuo implica un aumento de su peso esperado o promedio de 650 gramos.

• El coeficiente \(\hat{\beta_2}\) de valor 1.4 kg, corresponde a la edad del individuo. Este coeficiente indica que dada una altura, genero, días de actividad física y consumo de alcohol diario fijos, cada vez que el individuo cumple un año su peso esperado o promedio aumenta en 1.4 kg.

• El coeficiente \(\hat{\beta_3}\) de valor 1.26 kg, corresponde a los individuos de genero masculinos. Este coeficiente indica que dada una altura, edad, consumo de alcohol diario y días de actividad física semanal fijos, el peso promedio o esperado para el genero masculino es 1.26 kg mayor al peso femenino (categoría basal). Por otro lado, el coeficientes nos indica cuanto mas alto es el peso del genero masculino respecto del femenino al fijar los demás coeficientes.

• El coeficiente \(\hat{\beta_4}\) de valor 87.1 gramos, corresponde a los días de actividad física semanal que realiza el individuo. Este coeficiente indica que dada una altura, edad, genero y consumo de alcohol diario fijos, cada vez que un individuo realiza un día mas de actividad física semanal su peso esperado o promedio disminuye en 87.1 gramos.

• El coeficiente \(\hat{\beta_5}\) de valor 7 gramos, corresponde al nivel de consumo diario de alcohol del individuo. Este coeficiente indica que dada una altura, edad, genero y días de actividad física semanal fijos, cada vez que el individuo consume un trago de alcohol su peso esperado o promedio aumenta en 8 gramos.

2.3. ¿Son significativos los coeficientes?

Para determina si los coeficientes son aptos para explicar el peso de un individuo se realiza un \({T}\) test para cada coeficiente en el cual se evalúan las siguientes hipótesis:

• \({H_0: \beta_i = 0}\)
• \({H_1: \beta_i \neq 0}\)

Si \({\beta_i \neq 0}\) podemos decir que existe una diferencia estadisticamente significativas del cero para el coeficiente \({\beta_i}\), y por lo tanto el coeficiente \({\beta_i}\) explicar la variable \({y}\) (Peso en nuestro caso).

Luego analizando la salida de coefficients_summary concluimos que:

• Los coeficientes correspondientes a la categoria basal \({\beta_1}\)(Genero femenino), altura(\({\beta_1}\)), edad(\({\beta_2}\)) y genero masculino (\({\beta_3}\)) tienen \(p-valor < 0.05\). Por lo tanto, se rechaza la hipótesis nula y resultan estadistitamente significativos para explicar el peso. Por otro lado, se puede apreciar que los intervalos de confianza del 95% no incluyen al cero.
• Lo contrario sucede con días de actividad física semanal(\({\beta_4}\)) y consumo de alcohol diario (\({\beta_5}\)), dado que ambos no rechazar la hipótesis de nulidad (\(p-valor > 0.05\)) y por lo tanto no existe una diferencia significativa del cero. Finalmente, no hay evidencia estadistitamente significativas pensar que estos coeficientes expliquen el peso.

2.4. ¿El modelo resulta significativo para explicar el peso?

Para determinar si es modelo es significativo para explicar el peso de un individuo se realiza un \(F\) test con las siguientes hipótesis:

• \(H_0: {\beta_1} = {\beta_2} = · · · = {\beta_{p−1}} = 0\)
• \(H_1:\) Por lo menos un \({\beta_k}\) (\(k = 1, 2,..., p−1\)) es distinto de 0.

Donde:

• \(H_0\) afirma que no hay vinculo entre la variable \({y}\)(Peso) y las variables regresoras.
• \(H_1\) afirma que al menos una de las variables regresoras sirve para predecir la variable \({y}\)(Peso).

Veamos los resultados del \(F\) test:

glance(model_1)

Podemos apreciar que el \(p-valor < 0.05\) y igual a 0. Con mucha certeza podemos decir que al menos una de las variables regresoras permite explicar el peso. Esto concuerda con los resultados de los \(T\) test para los coeficientes correspondientes a altura, edad y genero femenino(basa) y masculino.

2.5. ¿Qué porcentaje de la variabilidad explica el modelo?

Según el valor de \(R^2\) ajustado (adj.r.squared), este modelo llega a explica el 35.39% de la variabilidad del dataset de entrenamiento, lo cuales un valor bajo pero tampoco es despreciable.

2.6. ¿Que sucede si ponemos al genero masculino como variable basal?

train_set_genero <- data.frame(train_set) 
train_set_genero$genero <- factor(
  train_set_genero$genero,
  levels=c('Masculino', 'Femenino'), 
  ordered=FALSE
)
table(train_set_genero$genero)

## 
## Masculino  Femenino 
##      3260      3764

model_genero <- lm(
  peso ~ altura + edad + genero + dias_actividad_fisica_semanal + consumo_diario_alcohol, 
  data = train_set_genero
)

coefficients_summary(model_genero)

## _________________________________________________________________________________________________________ 
## term                               estimate  std.error   statistic       p.value    conf.low    conf.high 
## ========================================================================================================= 
## (Intercept)                   -67.660044511 2.44965480 -27.6202364 1.682519e-159 -72.4621079 -62.85798114 
## altura                          0.650606544 0.01437975  45.2446353  0.000000e+00   0.6224179   0.67879520 
## edad                            1.406727060 0.09385081  14.9889709  5.121599e-50   1.2227511   1.59070300 
## generoFemenino                 -1.262643558 0.27282821  -4.6279802  3.758831e-06  -1.7974693  -0.72781785 
## dias_actividad_fisica_semanal  -0.087391031 0.04992917  -1.7503000  8.011025e-02  -0.1852673   0.01048523 
## consumo_diario_alcohol          0.007271379 0.06138558   0.1184542  9.057112e-01  -0.1130629   0.12760566 
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ 
## ____________________________________________________________________________ 
## Termino                         Coeficiente Signiticativo IC incluye al cero 
## ============================================================================ 
## (Intercept)                   -67.660044511            Si                 No 
## altura                          0.650606544            Si                 No 
## edad                            1.406727060            Si                 No 
## generoFemenino                 -1.262643558            Si                 No 
## dias_actividad_fisica_semanal  -0.087391031            No                 Si 
## consumo_diario_alcohol          0.007271379            No                 Si 
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

glance(model_genero)

Observaciones:

• En este caso, cambia el valor del coeficiente \(\hat\beta_0\) para la categoría basal, manteniéndose los demás coeficientes.
• Siguen siendo significativos los mismo coeficientes. La única diferencia es que \(\hat\beta_0\)(basal) representa al genero masculino y \(\hat\beta_3\) al femenino.
• Se mantiene la significatividad global del modelo y se obtiene el mismo \(R^2\) ajustado.
• Como ultimo, cabe aclarar que al cambiar la categoría basar cambian las interpretaciones del modelo.

3. Modelo categóricas

3.1. Se sugiere probar un modelo que incorpore el consumo semanal de snacks y una interacción entre el género y la edad, en lugar de actividad física y consumo de alcohol. Además se pide explicitamente que la categoría “No comí comida salada o snacks en los últimos 7 días” de la variable consumo_semanal_snacks se encuentre como nivel/categoría basal.

\(E(peso) = \beta_0 + \beta_1 * altura + \beta_2 * edad + \beta_3 * genero + \beta4 * consumoSemanalSnacks + \beta_5 * genero * edad\)

Primero validamos que las primeras categorías en cada variable de tipo factor sean las correctas, ya que esta sera la que el modelo defina como categoría basal:

table(train_set$consumo_semanal_snacks)

## 
##    0  <=3  4-6    7   14   21 >=28 
## 2162 3144  623  604  231  100  134

table(train_set$genero)

## 
##  Femenino Masculino 
##      3764      3260

Se puede apreciar que la primeras categorías corresponden a 0 consumo de snacks semanal y genero femenino. Por otro lado la categoría genero se encuentra balanceada en gran medida.

A continuación se cambiar la variable factor a no ordenada para que aparezcan correctamente valores de las categorías en los coeficientes:

train_set_model_1 <- train_set 

train_set_model_1$consumo_semanal_snacks <- factor(
  train_set_model_1$consumo_semanal_snacks,
  levels=c("0", "<=3", "4-6", "7", "14", "21", ">=28"), 
  ordered=FALSE
)

Modelo 2

Definimos el nuevo modelo:

model_2 <- lm(
  peso ~ altura + edad + genero + consumo_semanal_snacks +  genero * edad, 
  data = train_set_model_1
)

3.2. ¿Cuál es la interpretación de los coeficientes estimados para las categorías de consumo_semanal_snacks y genero * edad? ¿Son significativas?

coefficients_summary(model_2, show_tidy_sumamry = FALSE)

## _______________________________________________________________________ 
## Termino                    Coeficiente Signiticativo IC incluye al cero 
## ======================================================================= 
## (Intercept)                -64.2548500            Si                 No 
## altura                       0.6431229            Si                 No 
## edad                         1.2253900            Si                 No 
## generoMasculino             -4.6046463            No                 Si 
## consumo_semanal_snacks<=3   -1.3515386            Si                 No 
## consumo_semanal_snacks4-6   -2.2695264            Si                 No 
## consumo_semanal_snacks7     -0.6080931            No                 Si 
## consumo_semanal_snacks14    -1.0943997            No                 Si 
## consumo_semanal_snacks21    -1.2761227            No                 Si 
## consumo_semanal_snacks>=28  -2.5682972            Si                 No 
## edad:generoMasculino         0.3892757            Si                 No 
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

Si interpretamos los coeficientes que son significativos para el \(T\) test:

1. Coeficiente correspondiente al consumo_semanal_snacks<=3: Si fijamos los coeficientes correspondientes a la altura, edad, generoMasculino y edad:generoMasculino; el peso promedio o esperado de un individuo de consume snacks hasta 3 veces por semana es 1.33 kg menor que aquellos que no consumen snacks.

2. Sucede algo similar con las categorías consumo_semanal_snacks4-6 y consumo_semanal_snacks>=28 donde:

• El peso para un individuo que consume de 4 a 6 veces por semana es 2.26 kg menor que aquellos que no consumen snacks.
• El peso para un individuo que consume 28 veces a mas por semana es 2.56 kg menor que aquellos que no consumen snacks.

Puede parecer contradictorio que os que mas consumen snacks menos pesan, pero pensemos que no tenemos desglose de esta categorías. Hoy en día se consumen muchos snacks saludables, los cuales no tienen grasas, pero también existen los tradicionales ricos en grasas como papas fritas, palitos con queso, 3D, etc. Según los coeficientes que obtuvimos con este dataset, se podría pensar que los individuos están consumiendo en mayor medida snacks saludables, pero solo podemos asegurarlo si incluimos la variables tipo de snack.

3. Coeficiente correspondiente al edad * genero:

Dado el modelo original:

\(E(peso) = \beta_0 + \beta_1 * altura + \beta_2 * edad + \beta_3 * genero + \ \beta4 * consumoSemanalSnacks + \beta_{2,3} * genero * edad\)

y sabiendo que el genero femenino toma el valor 0 y masculino 1. Si reemplazamos estos valores en el modelo original encontramos que:

\(E_f(peso) = \beta_0 + \beta_1 * altura + \beta_2 * edad + \beta_4 * consumoSemanalSnacks\)

El genero femenino tiene la ordenada \(\beta0\) y las pendientes determinada por \(\beta_1\), \(\beta_2\) y \(\beta_4\).

\(E_m(peso) = (\beta_0 + \beta_3) + \beta_1 * altura + (\beta_2 + \beta_{2,3}) * edad + \ \beta_4 * consumoSemanalSnacks\)

El genero masculino tiene una ordenada que es la suma de la ordenada del genero femenino \(\beta_0\) mas \(\beta_3\). Tambien cambia la pendiente de \(\beta_2\) de la edad, a la cual se le suma \(\beta_{2,3}\).

Ahora, si solo cambian los coeficientes correspondientes al genero y edad, es decir si mantenemos contantes los demás coeficientes obtenemos:

• \(E_f(peso) = \beta_0 + \beta_2 * edad + cte\)
• \(E_m(peso) = (\beta_0 + \beta_3) + (\beta_2 + \beta_{2,3}) * edad\)

Reemplazamos por los coeficientes por lo valores que encontró el modelo:

• Femenino:

  • \(E_f(peso) = -64.2548500 + 1.2253900 * edad + cte\)

• Masculino:

  1. \(E_m(peso) = (-64.2548500 + -4.6046463) + (1.2253900 + 0.3892757) * edad + cte\)
  2. \(E_m(peso) = -68.8595 + 1.614666 * edad + cte\)

Finalmente tenemos:

• \(E_f(peso) = -64.2548500 + 1.2253900 * edad + cte\)
• \(E_m(peso) = -68.8595 + 1.614666 * edad + cte\)

Finalmente, graficamos ambas rectas definiendo la \(cte\) con un valor que de pesos positivos para tener una gráfica consistente:

cte = 100

train_set %>% 
  mutate(
    peso = ifelse(
      genero=='Femenino', 
      (-64.2548500 + 1.2253900 * edad) + cte,
      (-68.8595 + 1.614666 * edad) + cte
    ) 
  ) %>%
  ggplot(aes(x = edad, y = peso, colour=genero)) +
  geom_line() +
  ylab('Peso') +
  xlab('Edad')

Finalmente, se puede apreciar que las ordenadas de ambos géneros son distintas, donde el genero femenino inicia desde un peso 4kg menor al masculino. Luego si variamos únicamente la edad, se aprecia que el peso del genero masculino es mayor al femenino para la misma edad en todo los casos. Esto se debe a que la recta correspondiente al genero masculina esta por arriba de la resta correspondiente al genero femenino.

3.3. ¿Qué porcentaje de la variabilidad explica el modelo? En caso de detectar que existen categorías no significativas de la variable consumo_semanal_snacks evaluar si la variable es significativa en su conjunto y, en caso afirmativo, proponer una redefinición de las mismas que permita obtener una mayor proporción de categorías significativas individualmente. Luego, analizar si existen cambios en la variabilidad explicada por el modelo.

Viendo el resultado de coefficients_summary se aprecia que las siguientes categorías de consumo_semanal_snacks no son significativas:

• 2 veces al día (14 veces/semana).
• 4 a 6 veces durante los últimos 7 (4 a 6 veces/semana).
• 3 veces al día (21 veces/semana).

Pero si son significativos los extremos:

• De 1 a 3 veces/semana.
• De 28 o mas veces/semana.

A continuación se realiza un \(F\) test para evaluar la significatividad conjunta de las categóricas de la variable consumo_semanal_snacks para explicar el peso.

El \(F\) test también llamando ANOVA (Análisis de la varianza) se realiza para probar la significatividad conjunta de todos los valores de una variable categórica.

Las hipótesis son las siguientes:

• \(H_0: β_q = β_{q+1} = · · · = β_{p−1} = 0\)
• \(H_1:\) por lo menos uno de los \(β_k\) (con \(k\) entre \(q\) y \(p−1\)) es tal que \(β_k \neq 0\).

Luego si todos los coeficientes asociados a los valores de variable categórica son cero, se rechaza la hipótesis nula y por lo tanto la variable no es significartiva para explicar el peso en nuestro caso.

A continuación veremos el p-valor resultado de aplicar \(F\) test para cada variable del modelo:

anova_summary(model_2)

Podemos apreciar que el \(p-value < 0.005\) para la variable consumo_semanal_snacks. Por lo tanto se rechaza la hipótesis nula y podemos decir que en su conjunto resulta estadísticamente significativa para explicar el peso. Luego, como la variable consumo_semanal_snacks es significativa vale la pena re-definirla. Por otro lado, la combinación de variables genero-edad no es estadísticamente significativa para explicar el peso, pero si lo es el genero en forma separada. Finalmente, como ya vimos en pasos anteriores, edad y altura son significativas.

Modelo 2: Redefinición 1

Dado que no todas las categorías de la variable consumo_semanal_snacks son significativas, a continuación se propone una re-definición, en la que se reagrupan las categorías en nuevas categorías que si resulten significativas para el modelo 2:

train_set_snack_1 <- train_set %>% mutate(consumo_semanal_snacks = case_when(
  consumo_semanal_snacks %in% c('<=3', '4-6' , '7')  ~ '<=7',
  consumo_semanal_snacks %in% c('14', '21', '>=28')  ~ '>=14',
  TRUE ~ as.character(consumo_semanal_snacks)
))
train_set_snack_1$consumo_semanal_snacks <- factor(
  train_set_snack_1$consumo_semanal_snacks,
  levels=c('0', '<=7', '>=14'),
  ordered=FALSE
)

test_set_snack_1 <- test_set %>% mutate(consumo_semanal_snacks = case_when(
  consumo_semanal_snacks %in% c('<=3', '4-6' , '7')  ~ '<=7',
  consumo_semanal_snacks %in% c('14', '21', '>=28')  ~ '>=14',
  TRUE ~ as.character(consumo_semanal_snacks)
))

table(train_set_snack_1$consumo_semanal_snacks)

## 
##    0  <=7 >=14 
## 2162 4371  465

table(test_set_snack_1$consumo_semanal_snacks) 

## 
##  <=7 >=14    0 
## 1934  192  875

Como se puede ver, se reagruparan las categorías en 3 grupos:

• 0 veces/semana.
• Hasta 7 veces/semana.
• Mas de 14 veces/semana.

Luego entrenamos el Modelo 2 usando esta nueva re-definición:

model_2_redefinicion_1 <- lm(
  peso ~ altura + edad + genero + consumo_semanal_snacks +  genero * edad, 
  data = train_set_snack_1
)

coefficients_summary(model_2_redefinicion_1, show_tidy_sumamry = FALSE)

## _______________________________________________________________________ 
## Termino                    Coeficiente Signiticativo IC incluye al cero 
## ======================================================================= 
## (Intercept)                -64.0389167            Si                 No 
## altura                       0.6419791            Si                 No 
## edad                         1.2234673            Si                 No 
## generoMasculino             -4.6552035            No                 Si 
## consumo_semanal_snacks<=7   -1.3800004            Si                 No 
## consumo_semanal_snacks>=14  -1.5602043            Si                 No 
## edad:generoMasculino         0.3928142            Si                 No 
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

anova_summary(model_2_redefinicion_1)

glance(model_2_redefinicion_1)

Modelo 2: Redefinición 2

En este caso se propone calcular la media del ratio altura/edad para cada categoría de la variables consumo_semanal_snacks. Luego calculamos los cuantiles de esta nueva distribución y los utilizamos para crear una nueva categorización: Los individuos que tenga un ratio menor al cuantíl 2 (Mediana) tendran el valor Bajo y Alto en caso contrario. Se intento llevar a mas niveles pero el test \(T\) no daba significativo para todos los coeficientes.

1. Definimos la districión altura/edad por categoria de consumo_semanal_snacks:

train_set_snack_2 <- train_set %>% 
  mutate(alt_edad_ratio = round(altura/edad, 0))

avg_train_set_snack_2 <- train_set_snack_2 %>% 
   group_by(consumo_semanal_snacks) %>% 
   summarise(avg_alt_edad_ratio = mean(alt_edad_ratio))

ggplot(data = avg_train_set_snack_2, aes(x = avg_alt_edad_ratio)) + 
  geom_boxplot(alpha = 0.75, fill="blue") +
  theme_bw()

DE los siguientes cuantíles utilizaremos el cuantíl 2(Mediana o 50%):

quantiles_avg_alt_edad_ratio <- quantile(avg_train_set_snack_2$avg_alt_edad_ratio)
quantiles_avg_alt_edad_ratio

##       0%      25%      50%      75%     100% 
## 10.73077 10.89892 10.97694 11.03058 11.09952

2. Creamos un dataset intermedio donde están mapeadas las categorias originales vs las nuevas:

q2 <- quantiles_avg_alt_edad_ratio[3]

snack_level_mapping <- avg_train_set_snack_2 %>% 
  mutate(level = case_when(
    avg_alt_edad_ratio < q2  ~ 'Bajo',
    avg_alt_edad_ratio >=  q2 ~ 'Alto'
  )) %>% select(consumo_semanal_snacks, level)

snack_level_mapping %>%
  arrange(consumo_semanal_snacks)

3. Creamos nuevos dataset de train/test con la re-definición de la variable consumo_semanal_snacks:

train_set_snack_2 <- train_set %>%
  inner_join(snack_level_mapping, by = 'consumo_semanal_snacks') %>%
  mutate(consumo_semanal_snacks = level) %>% 
  select(-level)

test_set_snack_2 <- test_set %>%
  inner_join(snack_level_mapping, by = 'consumo_semanal_snacks') %>%
  mutate(consumo_semanal_snacks = level) %>% 
  select(-level)

train_set_snack_2 %>%  
  group_by(consumo_semanal_snacks) %>% 
  tally()

test_set_snack_2 %>%  
  group_by(consumo_semanal_snacks) %>% 
  tally()

4. Evaluamos el nuevo train_set con el modelo 2:

model_2_redefinicion_2 <- lm(
  peso ~ altura + edad + genero + consumo_semanal_snacks +  genero * edad, 
  data = train_set_snack_2
)

coefficients_summary(model_2_redefinicion_2, show_tidy_sumamry = FALSE)

## _______________________________________________________________________ 
## Termino                    Coeficiente Signiticativo IC incluye al cero 
## ======================================================================= 
## (Intercept)                -65.6507415            Si                 No 
## altura                       0.6435785            Si                 No 
## edad                         1.2217889            Si                 No 
## generoMasculino             -4.6387828            No                 Si 
## consumo_semanal_snacksBajo   1.1231822            Si                 No 
## edad:generoMasculino         0.3929346            Si                 No 
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

anova_summary(model_2_redefinicion_2)

glance(model_2_redefinicion_2)

Conclusión

models <- list(
  'Modelo 1' = model_1, 
  'Modelo 2' = model_2, 
  'Modelo 2 - Re-definición 1' = model_2_redefinicion_1, 
  'Modelo 2 - Re-definición 2' = model_2_redefinicion_2
)

models %>% 
  map_df(glance, .id = "model") %>%
  arrange(desc(adj.r.squared))

Ambos modelos son significativos para explicar el peso. El modelo Modelo 2 - Re-definición 1 es mas explicativo, ya que \(R^2\) ajustado es mayor par este dataset de entrenamiento. Por otro lado, ambos modelos son menos explicativos que el Modelo 2 pero mas explicativos que el Modelo 1 (inicial). Finalmente es importante apreciar que las diferencian del \(R^2\) ajustado son muy bajas entre los modelos, deberiamos evaluarlo en test para entnder cual modelo podria performar menor.

4. Modelos propios y evaluación

4.1. Realizar 2 modelos lineales múltiples adicionales y explicar breve-mente la lógica detrás de los mismos (se valorará la creación y/o inclusión de variables nuevas).

Al continuación se define 2 modelos.

Modelo 3

En base a los resultados de importancia de variables realizados en la notebook eea-tp1-eda, \(R^2\) ajustado de los resultados y la significación de los coeficientes se seleccionaron las siguientes variables:

• Altura: Ya que la altura influye en el peso del individuo directamente.
• Genero: No tiene un relación directa, pero lo individuos masculinos tienden a tener mas altura y por ende mayor peso.
• consumo_semanal_snacks y consumo_semanal_comida_grasa: Puede considerarse que estas variables afectan en gran media al aumento de peso.
• altura * edad: Ya que existe una relación importante entre la edad y la altura en el rango de edades del dataset (12 a 18 años). Esta relación se entiende que afecta al peso al igual que la altura.
• altura * genero: Parece haber una relación definida entre la altura y el genero, ya que al agregar esta variable aumenta ligeramente el \(R^2\) Ajustado.

model_3 <- lm(
  peso~ 
    altura + 
    genero + 
    consumo_semanal_snacks +
    consumo_semanal_comida_grasa + 
    altura * edad + 
    altura * genero,
  data = train_set
)

coefficients_summary(model_3, show_tidy_sumamry = FALSE)

## ____________________________________________________________________________ 
## Termino                         Coeficiente Signiticativo IC incluye al cero 
## ============================================================================ 
## (Intercept)                    -15.74275494            No                 Si 
## altura                           0.32562490            Si                 No 
## generoMasculino                -15.67600313            Si                 No 
## consumo_semanal_snacks.L        -1.35103936            No                 Si 
## consumo_semanal_snacks.Q        -0.19162670            No                 Si 
## consumo_semanal_snacks.C        -1.66725662            Si                 No 
## consumo_semanal_snacks^4         0.07791653            No                 Si 
## consumo_semanal_snacks^5         0.21636253            No                 Si 
## consumo_semanal_snacks^6        -0.69763266            No                 Si 
## consumo_semanal_comida_grasa.L  -0.11654538            No                 Si 
## consumo_semanal_comida_grasa.Q   0.35094641            No                 Si 
## consumo_semanal_comida_grasa.C   1.06360279            No                 Si 
## consumo_semanal_comida_grasa^4   1.01675663            No                 Si 
## consumo_semanal_comida_grasa^5   0.09779008            No                 Si 
## consumo_semanal_comida_grasa^6  -1.51128497            Si                 No 
## edad                            -1.56720061            No                 Si 
## altura:edad                      0.01777853            No                 Si 
## altura:generoMasculino           0.10259715            Si                 No 
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

anova_summary(model_3)

glance(model_3)

Modelo 4

En este modelo principalmente se busco sumar explicabilidad al modelo, es decir el mayor \(R^2\) ajustado posible sin prestar atención si las variables son explicativas o sus categorìas en caso de ser categòricas. Para esto se agregaron todas las variables del train_set y luego quellas combinaciòn que sumaron \(R^2\) en el modelo 3.

model_4 <- lm(
  peso ~
  edad +
  genero +
  nivel_educativo + 
  altura + 
  frecuencia_hambre_mensual + 
  dias_consumo_comida_rapida +
  edad_consumo_alcohol + 
  consumo_diario_alcohol + 
  dias_actividad_fisica_semanal +
  consumo_semanal_frutas +
  consumo_semanal_verdura +
  consumo_semanal_gaseosas +
  consumo_semanal_snacks +
  consumo_semanal_comida_grasa +
  altura * edad + 
  altura * genero + 
  altura * nivel_educativo,
  data = train_set
)

coefficients_summary(model_4, show_tidy_sumamry = FALSE)

## _____________________________________________________________________________ 
## Termino                          Coeficiente Signiticativo IC incluye al cero 
## ============================================================================= 
## (Intercept)                    -81.322024010            Si                 No 
## edad                             2.866682357            No                 Si 
## generoMasculino                -15.340016623            Si                 No 
## nivel_educativo.L              -21.003289330            Si                 No 
## nivel_educativo.Q                2.892897973            No                 Si 
## nivel_educativo.C                3.678132249            No                 Si 
## nivel_educativo^4               -9.809417974            Si                 No 
## altura                           0.725953429            Si                 No 
## frecuencia_hambre_mensual.L      0.419188819            No                 Si 
## frecuencia_hambre_mensual.Q      0.931304812            No                 Si 
## frecuencia_hambre_mensual.C      0.005646566            No                 Si 
## frecuencia_hambre_mensual^4      0.005481189            No                 Si 
## dias_consumo_comida_rapida      -0.089875750            No                 Si 
## edad_consumo_alcohol.L          -1.651306167            No                 Si 
## edad_consumo_alcohol.Q          -2.183122338            No                 Si 
## edad_consumo_alcohol.C          -1.339149569            No                 Si 
## edad_consumo_alcohol^4          -1.215373662            No                 Si 
## edad_consumo_alcohol^5          -0.281883308            No                 Si 
## edad_consumo_alcohol^6          -0.219142126            No                 Si 
## edad_consumo_alcohol^7          -0.574866030            No                 Si 
## consumo_diario_alcohol           0.027146287            No                 Si 
## dias_actividad_fisica_semanal   -0.116645647            Si                 No 
## consumo_semanal_frutas.L        -0.631395889            No                 Si 
## consumo_semanal_frutas.Q        -0.735855107            No                 Si 
## consumo_semanal_frutas.C        -0.091526430            No                 Si 
## consumo_semanal_frutas^4         0.428976864            No                 Si 
## consumo_semanal_frutas^5         0.189968669            No                 Si 
## consumo_semanal_frutas^6         0.651034410            No                 Si 
## consumo_semanal_verdura.L        0.060649375            No                 Si 
## consumo_semanal_verdura.Q        0.436999347            No                 Si 
## consumo_semanal_verdura.C        0.770944414            No                 Si 
## consumo_semanal_verdura^4       -0.084173295            No                 Si 
## consumo_semanal_verdura^5        0.260346330            No                 Si 
## consumo_semanal_verdura^6        0.182684068            No                 Si 
## consumo_semanal_gaseosas.L      -1.328574893            Si                 No 
## consumo_semanal_gaseosas.Q       0.136219732            No                 Si 
## consumo_semanal_gaseosas.C      -0.296205585            No                 Si 
## consumo_semanal_gaseosas^4      -0.527578792            No                 Si 
## consumo_semanal_gaseosas^5      -0.866338586            Si                 No 
## consumo_semanal_gaseosas^6      -0.617701274            No                 Si 
## consumo_semanal_snacks.L        -0.739112310            No                 Si 
## consumo_semanal_snacks.Q        -0.263445557            No                 Si 
## consumo_semanal_snacks.C        -1.800008228            Si                 No 
## consumo_semanal_snacks^4         0.114405964            No                 Si 
## consumo_semanal_snacks^5         0.241512402            No                 Si 
## consumo_semanal_snacks^6        -0.688788449            No                 Si 
## consumo_semanal_comida_grasa.L   0.609882988            No                 Si 
## consumo_semanal_comida_grasa.Q   0.527005294            No                 Si 
## consumo_semanal_comida_grasa.C   0.811837477            No                 Si 
## consumo_semanal_comida_grasa^4   0.817166249            No                 Si 
## consumo_semanal_comida_grasa^5  -0.117789611            No                 Si 
## consumo_semanal_comida_grasa^6  -1.573730114            Si                 No 
## edad:altura                     -0.009310081            No                 Si 
## generoMasculino:altura           0.100919188            Si                 No 
## nivel_educativo.L:altura         0.128387205            Si                 No 
## nivel_educativo.Q:altura        -0.020600451            No                 Si 
## nivel_educativo.C:altura        -0.020533157            No                 Si 
## nivel_educativo^4:altura         0.062739320            Si                 No 
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

anova_summary(model_4)

glance(model_4)

Finalmente, si comparamos los modelos por \(R^2\) Ajustado, se puede apreciar que el Modelo 4 llega a captar la mayor varianza explicada sobre el conjunto de entrenamiento. Por supuesto esto no dice nada acerca de la performance del modelo en test, pero si que tiene la mejor capacidad para extraer información de los dato de entrenamiento.

4.2. Evaluar la performance del modelo inicial, el modelo con las categorías redefinidas de la variable consumo_semanal_snacks y los modelos desarrollados en este punto en el dataset de entrenamiento y evaluación (usar dataset “encuesta_salud_test.csv”). La evaluación de performance consiste en comparar en ambos sets la performance en términos del R cuadrado ajustado, RMSE y MAE ¿Cuál es el mejor modelo para nuestro objetivo de predecir el peso? ¿Por qué?

Ahora comparamos la performance de todo los modelos evaliando el error de los mismo al predecir el peso en el conjunto de train y test tanto para RMSE como MAE:

RMSE

eval_summary <- function(metric_fn = rmse) {
  summary_part_1 <- list(
    'Modelo 1' = model_1, 
    'Modelo 2' = model_2, 
    'Modelo 3' = model_3,
    'Modelo 4' = model_4
  ) %>% train_test_eval_metric_summary(test_set = test_set, metric_fn = metric_fn)
  
  summary_part_2 <- list('Modelo 2 - Re-definición 1' = model_2_redefinicion_1) %>% 
    train_test_eval_metric_summary(test_set = test_set_snack_1, metric_fn = metric_fn)
  
  summary_part_3 <- list('Modelo 2 - Re-definición 2' = model_2_redefinicion_2) %>% 
    train_test_eval_metric_summary(test_set = test_set_snack_2, metric_fn = metric_fn)

  summary <- summary_part_1 %>% 
    union_all(summary_part_2) %>% 
    union_all(summary_part_3)
  
  min_test_error <- summary %>% 
    pull(test_error) %>% 
    min()
  
  max_r_2_adjusted <- summary %>% 
    pull(r_2_adjusted) %>% 
    max()

  summary %>% mutate(
    min_test_error_diff = abs(test_error - min_test_error),
    max_r2_diff         = abs(r_2_adjusted - max_r_2_adjusted)
  ) %>%
    select(model, r_2_adjusted, max_r2_diff, train_error , test_error, error_diff, min_test_error_diff) %>%
    arrange(test_error)
}

rmse_eval_summary_result <- eval_summary(rmse)
rmse_eval_summary_result

Si utilizamos la métrica RMSE podemos ver que el modelo Modelo 2 - Re-definición 1 tiene el menor error al evaluarlo con el conjunto de test. Por otro lado, es el que tiene menor diferencia de error entre test y entrenamiento. Esto nos dice que tiene el menor grado de sobre-ajuste (overfitting) al conjunto de entrenamiento. Otro punto importante, es que es uno de los modelos con mayor diferencian en \(R^2\) ajustado contra el Modelo 4 el cual tiene el mayor \(R^2\) ajustado, siendo este es otro indicador de que el Modelo 4 se ajusta mas al conjunto de entrenamiento, cayendo en el mayor sobre ajuste posible. Sucede algo muy similar con el modelo Modelo 3 pero en menor medida es cual es el segundo sobre ajustando y el que tiene mayor error al evaluar en test seguido de los modelos 1 y 4.

MAE

mae_eval_summary_result <- eval_summary(mae)
mae_eval_summary_result

Usan la métrica MAE el Modelo 4 es el que tiene el menor error en test seguido de Modelo 2 - Re-definición 1. El Modelo 2 - Re-definición 1 parece ser el modelo que menos sobre-ajuste y sorprendente-mente el Modelo 4 a pesar de ser el que menos error tiene el test también es el que mas sobre ajusta al conjunto de test dado. Por esta ultima cuestión la mejor opción es el modelo Modelo 2 - Re-definición 1 ya que el modelo Modelo 4 es el que peor generaliza.

Finalmente, según ambas métricas Modelo 2 - Re-definición 1 es el modelo mas adecuado para ambas metricas (RMSE y MAE).

5. Diagnóstico del modelo

Analizar en profundidad el cumplimiento de los supuestos del modelo lineal para el modelo inicial.

plot(model_1)

Homocedasticidad

Al visualizar el primer gráfico (Residuos vs. Valores ajustados) se puede apreciar que hay presencia de homocedasticidad, ya que a medida que aumentan los valores predichos la variabilidad o amplitud de los residuos parece mantenerse en los mismo niveles. Dadas esta condiciones podemos decir que se cumple el supuesto de varianza constante.

Normalidad

Al visualizar el diagrama QQ-Plot podemos observas que en el extremo derecho, el modelo sobre-estima el peso del los individuos ya que hay una gran diferencia entre los valores predichos y los valores esperados teóricos. Lo contrario sucede a izquierda, donde el modelo subestima el valor de peso en comparación al valor esperado teórico, aunque los valores de los residuos son menores en este caso. Finalmente el QQ-Plot muestra un grado de alejamiento pronunciado de una distribución normal teórica y por lo tanto no se cumple el supuesto de normalidad del modelo.

Apalancamiento (Leverage)

Si observamos el gráfico de Residuos vs Apalacamiento vemos varias observaciones o individuos que se alejan a derecha del cumulo principal. Estos ejercen un alejamiento de las prediciones del modelo vs los valores reales a partir de un apalancamiento(leverage) 0.0025 y es mas pronunciado desde 0.0035. Finalmente, vemos un grado importante de desvió de las predicciones vs valores reales y por ente un grado importante de apalancamiento(leverage).

A continuación se pueden ver lo individuos que producen mayor apalancamiento(leverage) y por ende sesgo en las predicciones del modelo:

augment(model_1) %>%
  filter(.hat>0.0025) %>%
  arrange(.hat)

6. Modelo Robusto

6.1. Leer el archivo “encuesta_salud_modelo6.csv”. Este último consiste en el dataset original de train

con la incorporación de algunas observaciones adicionales que pueden incluir valores atípicos. En particular, observar la relación entre peso y altura ¿Qué ocurre con estos nuevos datos?

A continuación se carga el conjunto de entrenamiento en crudo, es decir sin pre-procesamiento. Luego se resumen los valores de las variables categóricas y se transforman lo misssings a NA’s:

original_train_set <- load_original_train_set() %>% 
  preprocess() %>% 
  shorten_values() %>% 
  process_missings()

original_train_set %>% missings_summary()

Se aprecia que hay muy pocos missing values. Por otro lado estos no se eliminan ya que el modelo lineal por defecto los filtra. Comparemos las distribuciones del peso vs. altura en ambos conjunto de entrenamiento:

box_plots(train_set %>% select(altura, peso), title = 'Peso vs. Altura en train')

box_plots(original_train_set %>% select(altura, peso), title = 'Peso vs. Altura en train(Original)')

outliers(original_train_set, 'altura')

## $inf
## [1] 125 130 132 133 134 135 136 138
## 
## $sup
##  [1] 189 190 191 192 193 194 195 197 198 200

outliers(train_set, 'altura')

## $inf
## [1] 125 130 132 133 134 135 136 138
## 
## $sup
##  [1] 189 190 191 192 193 194 195 197 198 200

outliers(original_train_set, 'peso')

## $inf
## [1] 27
## 
## $sup
##  [1]  90  91  92  93  94  95  96  97  98  99 100 101 102 103 105 106 107 108 109
## [20] 110 112 115 117 120 121 122 123 126 130 131 135 143 144 146 148 149 152 153
## [39] 154 155 157 159 160 161 163 164 165 166 169 172 175 180 187

outliers(train_set, 'peso')

## $inf
## [1] 27
## 
## $sup
##  [1]  90  91  92  93  94  95  96  97  98  99 100 101 102 103 105 106 107 108 109
## [20] 110 112 115 117 120 121 122 123 130 131 144

En el dataset de entrenamiento original la variable peso tiene outliers mas extremos que el dataset pre-procesado, los cuales se encuentran entre los 144 y 187 kg. por otro lado la altura tiene los mismos outliers en ambos dataset.

6.2. Entrenar el modelo inicial con estos nuevos datos y comentar qué se observa en los coeficientes

estimados y las métricas de evaluación (R cuadrado ajustado, RMSE y MAE) respecto al modelo entrenado con el set de entrenamiento original.

Modelo 5

Definimos un modelo igual al modelo 1 pero entrenando en el dataset de entrenamiento original.

model_5 <- lm(
  peso ~ altura + edad + genero + dias_actividad_fisica_semanal + consumo_diario_alcohol, 
  data = original_train_set
)

coefficients_summary(model_5, show_tidy_sumamry = FALSE)

## ___________________________________________________________________________ 
## Termino                        Coeficiente Signiticativo IC incluye al cero 
## =========================================================================== 
## (Intercept)                   -70.03018295            Si                 No 
## altura                          0.66195138            Si                 No 
## edad                            1.40007823            Si                 No 
## generoMasculino                 1.03107761            Si                 No 
## dias_actividad_fisica_semanal  -0.11320339            No                 Si 
## consumo_diario_alcohol          0.02083144            No                 Si 
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

anova_summary(model_5)

glance(model_5)

compare_r2(model_1, model_5)

## [1] "R2 - Modelo A:  0.354 Modelo B: 0.273 Diff:  22.799 %"

Hay una diferencia porcentual muy alta del \(R^2\) ajustado en ambos modelos. El hecho de que existan outliers en el conjunto de entrenamiento original produce que disminuya la explicabilidad del mismo.

models <- list('Modelo 1' = model_1, 'Modelo 5' = model_5)
models %>% eval_models_summary(test_set,  metric_fn = rmse)

models %>% eval_models_summary(test_set,  metric_fn = mae)

En este caso los resultados son muy consistentes con ambas métricas (RMSE y MAE). En ambos caso el Modelo 1 es el que logra un menor error en test, ademas de tener la menor diferencia entre train y test. Por otro lado, con ambas métricas se aprecia que el Modelo 5 tiene peor performance en train que en test. Esto podría deberse a que el modelo tiene ejemplos mas difíciles de ajustar en train que en test, ya sea por efecto de los outliers o por que existe algún ejemplo extra el cual es difícil de ajustar para el modelo. A pesar de esto, hay muy poca diferencia de error en test entre ambos modelos.

6.3. Entrenar un modelo robusto con la misma especificación que el modelo inicial

sobre los nuevos datos. Comparar los coeficientes y su performance (RMSE y MAE) respecto al modelo inicial no robusto entrenado en este punto. ¿Qué puede concluir al respecto?

Modelo 7

Definimos un modelo igual al modelo 1 entrenando en el dataset de entrenamiento original y usamos un modelo lineal robusto.

model_robusto <- rlm(
  peso ~ altura + edad + genero + dias_actividad_fisica_semanal + consumo_diario_alcohol, 
  data = original_train_set
)

coefficients_summary(model_robusto, show_tidy_sumamry = FALSE)

## [1] "WARN: p.value column is required to make model coefficients summary!\n"
## [1] "WARN: p.value column is required to plot tidy coefficients!\n"

## NULL

models <- list('Modelo 1' = model_1, 'Modelo Robusto' = model_robusto)
models %>% eval_models_summary(test_set,  metric_fn = rmse, include_r2 = FALSE)

models %>% eval_models_summary(test_set,  metric_fn = mae, include_r2  = FALSE)

Con RMSE los resultados son muy similares al punto anterior. En este caso el Modelo Robusto pierde con respecto al Modelo 1, inclusive por una diferencia mucho mas grande en test que el Modelo 5. El modelo robusto sufre del mismo problema que el Modelo 5 donde tenemos un error en train mayor que test. Ahora Midiendo con MAE los resultados se invierten y el Modelo Robusto gana por una diferencia poco significativa.

Como conclusión, Los modelos Modelo Robusto y Modelo 5 on son mejores que el Modelo 1 sobre otro conjunto este conjunto de test. Por otro lado, los outliers presentes en el dataset original dificultan el ajuste de los modelos en el entrenamiento. Por otro lado, no se aprecia una diferencia significativa en la performance del Modelo 5 y el Modelo Robusto.

Evaluación de todos los modelos

models <- list('Modelo 5' = model_5, 'Modelo Robusto' = model_robusto)

models %>% 
  eval_models_summary(test_set,  metric_fn = rmse) %>% 
  join_eval_summaries(rmse_eval_summary_result) %>%
  select(model, r_2_adjusted, train_error, test_error, error_diff)

models %>% 
  eval_models_summary(test_set,  metric_fn = mae) %>% 
  join_eval_summaries(mae_eval_summary_result) %>%
  select(model, r_2_adjusted, train_error, test_error, error_diff)

Si comparamos todos los modelos, se puede ver con claridad que el Modelo Robusto queda ultimo al medir su performance con RMSE, pero es el mejor modelo al medir con MAE ¿Por que sucede esto?. El Modelo Robusto fue entrenada con el dataset original de entrenamiento, este dataset tiene mas outliers que el dataset original. Luego esto repercute en el modelo, ya que vamos a tener residuos con un mayor desvío. Por otro lado, la métricas RMSE penaliza los errores, es decir que le da un peso mayor con respeto al MAE. Esto causa que aumente mucho el RMSE del Modelo Robusto.

Finalmente Modelo 2 - Re-definición 1 sigue teniendo el menor RMSE en test y ademas están entre los modelos que tiene menor overfitting.

 

---

Realizado por Adrian Marino

adrianmarino@gmail.com