Tarea 5

Ilan Ben-Dov, Florencia Olbertz, Edwyn Turner

7 Junio, 2021

Pregunta #1

Cargue los paquetes “data.table”, “ggplot2” y “caret”, junto con la base de datos

library(data.table)
library(ggplot2)
library(moderndive)
library(knitr)
library(lattice)
library(caret)
seremi <- fread("datos_tarea_5.csv")
seremi <- as.data.table(seremi)
seremi <- seremi[,ss_t:=as.integer(seremi$ss_t)]
## Warning in eval(jsub, SDenv, parent.frame()): NAs introduced by coercion

Borramos todos los NA.

library(tidyr)
seremi <- seremi[!is.na(seremi$ss_t)]
seremi <- seremi[!is.na(seremi$sexo)]
seremi <- seremi[!is.na(seremi$nacionalidad)]
seremi <- seremi[!is.na(seremi$edad)]
seremi <- seremi[!is.na(seremi$r_p_c)]
seremi <- seremi[!is.na(seremi$region)]
seremi <- seremi[!is.na(seremi$cise)]

Pregunta 2

Muestre un histograma de ingreso por region. Limite el histograma a las observaciones con ingresos (ss_t) mayores a 0 y menores a 2 millones de pesos. Recuerde ver la clase de las variables

class(seremi$region)
## [1] "integer"
ingreso_por_region <- seremi[ss_t < 2000000 & ss_t > 0, ss_t, by = "region"]
#ippr <- ingreso_por_region[, (mean(ss_t)), by = "region"] ## areglar
ggplot(ingreso_por_region, aes(x=ss_t)) +
  geom_histogram(bins = 10) +  
  facet_wrap(facets="region")

Pregunta 3

Arregle los gráficos para que cada uno tenga un eje x legible, además arregle el eje y. Agregue color a cada gráfico: Hint: utilice la función scales = ‘free_y’ dentro de facet_wrap para dejar libre el eje y. También recuerde chequear la clase de la variable region para que se asigne correctamente un color a cada región.

ggplot(ingreso_por_region, aes(x=ss_t, fill=region)) +
  geom_histogram() +
  theme(axis.text.x = element_text(angle=75, vjust=0.5)) +
  facet_wrap(facets= "region", scales="free_y")
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Pregunta 4

  1. Realice un gráfico de dispersión que muestre la relación entre el ingreso y la edad.

Grafico con filtro en los ingresos (solicitado en preguntas anteriores)

ggplot(seremi[ss_t > 0 & ss_t < 2000000 & !is.na(ss_t)], aes(edad, ss_t)) +
  geom_point(color = "Blue", size = 0.8)

Grafico sin filtro, y mostrando la completa disperción de la edad~ingreso

ggplot(seremi[!is.na(ss_t)], aes(edad, ss_t)) +
  geom_point(color = "Blue", size = 0.8) +
  labs(x = "Edad", y = "Ingreso Promedio")

Pregunta 4b

Un miembro de su equipo propone un modelo que calcule la incidencia del sexo y la educación en el ingreso, sin considerar variables adicionales. Usted quiere mostrarle que aquel modelo está incompleto. Para lograr esto, haga el modelo de regresión anteriormente mencionado.

modelo1 <- lm(ss_t ~ sexo + nivel, data = seremi)
get_regression_table(modelo1) %>% kable()
term estimate std_error statistic p_value lower_ci upper_ci
intercept 169925.514 3498.895 48.565 0 163067.701 176783.326
sexo -54460.931 2162.967 -25.179 0 -58700.333 -50221.529
nivel 298.479 36.080 8.273 0 227.762 369.197

Pregunta 4.1

Calcule el salario predicho para cada observación.

seremi <- seremi[,predicciones_salario:=predict(modelo1)]

Pregunta 4.2

Dado el modelo anterior, calcule la predicción del salario para una mujer y para un hombre con educación Universitario Hombre

regressionHombre <- lm(ss_t ~ (sexo==1)+ (nivel==9) ,data = seremi)
get_regression_table(regressionHombre) %>% kable()
term estimate std_error statistic p_value lower_ci upper_ci
intercept 36690.96 1487.700 24.663 0 33775.08 39606.84
sexo == 1TRUE 53098.04 2107.393 25.196 0 48967.56 57228.52
nivel == 9TRUE 205527.95 3133.791 65.584 0 199385.73 211670.16 
seremi <- seremi[,predicciones_salariohombre_universitario:=predict(regressionHombre)]

Mujer

regressionMujer <- lm(ss_t ~  (sexo==2) + (nivel==9) ,data = seremi)
get_regression_table(regressionMujer) %>% kable()
term estimate std_error statistic p_value lower_ci upper_ci
intercept 89789.00 1598.816 56.160 0 86655.33 92922.67
sexo == 2TRUE -53098.04 2107.393 -25.196 0 -57228.52 -48967.56
nivel == 9TRUE 205527.95 3133.791 65.584 0 199385.73 211670.16 
seremi[,predicciones_salario_mujer_universitario:=predict(regressionMujer)]

Pregunta 5

Teniendo ya el modelo propuesto por su colega, contrarreste usted con un modelo de predicción en base a los datos con los que ya cuenta. Utilice el argumento na.action por los datos nulos en la función de la regresión.

propuesta <- lm(ss_t ~ sexo + 
                  nivel + 
                  nacionalidad + 
                  edad + 
                  r_p_c+region + 
                  cise  ,data = seremi)

get_regression_table(propuesta) %>% kable()
term estimate std_error statistic p_value lower_ci upper_ci
intercept 36646.039 4252.314 8.618 0 28311.529 44980.549
sexo -23329.579 1901.089 -12.272 0 -27055.701 -19603.456
nivel 125.410 31.477 3.984 0 63.715 187.105
nacionalidad -1.014 0.280 -3.615 0 -1.564 -0.464
edad -254.345 39.041 -6.515 0 -330.864 -177.825
r_p_c -35.683 5.782 -6.171 0 -47.016 -24.350
region 36757.988 5738.331 6.406 0 25510.893 48005.082
cise 107120.533 684.223 156.558 0 105779.460 108461.606

Pregunta 5.1

Calcule el salario predicho para cada observación

seremi[,Predic_5:=predict(propuesta)]
seremi[1:10,.(Predic_5)] %>% kable()
Predic_5
-18307.858
340877.452
-23676.547
6387.883
330798.095
307844.747
323160.700
-21021.787
302167.665
10487.218

Pregunta 6

Contando ya con su modelo, usted puede realizar una comparación directa entre ambos. Calcule los errores de predicción dentro de muestra para ambos modelos.

p1 <- predict(modelo1)
pr1 <-data.table(RMSE=RMSE(p1,seremi$ss_t), MAE=MAE(p1,seremi$ss_t))
pr1
##        RMSE      MAE
## 1: 303300.2 147359.3

Lo anterior los explica la regression en relación a sexo~y estudios. ahora volvemos a predecir con lo anterior y lo propuesto anteriormente

p2 <- predict(propuesta)
pr2 <-data.table(RMSE=RMSE(p2,seremi$ss_t),
                        MAE=MAE(p1,seremi$ss_t))
pr2
##        RMSE      MAE
## 1: 264053.7 147359.3

Pregunta 7

Interprete la diferencia en los errores de predicción entre el Modelo 1 y el Modelo 2. ¿Qué modelo hace una mejor predicción dentro de muestra?

Dentro de los dos modelos anteriores, podemos ver que el Error Absoluto Medio de la primera prediccion era mayor a el de la segunda. Esto quiere decir que la prediccion de el segundo modelo es “mejor” que el primero, dado los datos y el nivel predictivo obtenido.

Pregunta 8

Realice validación cruzada (CV) a los modelos de la pregunta anterior por el método K-folds con 5 folds. ¿Se mantienen las conclusiones obtenidas en el análisis dentro de muestra?

set.seed(12345)
setupKCV <- trainControl(method = "cv" , number = 5)
predkfolds1 <-train(ss_t ~ sexo + 
                      nivel ,data = seremi, method="lm", trControl= setupKCV)

predkfolds2<-train(ss_t ~ sexo + 
                     nivel + 
                     nacionalidad + 
                     edad + 
                     r_p_c + 
                     region + 
                     cise  ,data = seremi, method="lm", trControl= setupKCV)

predkfolds1
## Linear Regression 
## 
## 79088 samples
##     2 predictor
## 
## No pre-processing
## Resampling: Cross-Validated (5 fold) 
## Summary of sample sizes: 63270, 63270, 63271, 63271, 63270 
## Resampling results:
## 
##   RMSE      Rsquared     MAE     
##   302679.1  0.008774494  147373.4
## 
## Tuning parameter 'intercept' was held constant at a value of TRUE

Dado la siguente prediccion, y con el cross validation, determinamos que dejamos todo tal como fue dicho en la pregunta 7, ya que el cross validation nos da un valor de Error Cuadratico Medio menor a los anteriores.

Con esto mismo, nos damos cuetna que nuestro RMSE en ambos casos sigue un comportamiento similar, por lo que nuevamente podemos concluir que las concluciones anteriores siguen iguales.