Introducción

Vamos a trabajar con dos set de datos. El primero, es un dataset (.wb_bank_data_2019.csv) con información construida por el Banco Mundial acerca de la distribución de la población ocupada según grandes sectores de actividad y categoría ocupacional:

  • SL.AGR.EMPL.ZS: % total el empleo en agricultura: consistente en las actividades de agricultura, silvicultura, caza y pesca
  • SL.IND.EMPL.ZS: % total de empleo en industria: en actividades de minas y canteras, manufacturas, construcción y energía, gas y agua.
  • SL.SERV.EMPL.ZS % total de empleo el servicios comprende comercio al por mayor y menor y restaurantes y hoteles, transporte, almacenamiento y comunicaciones, finanzas, seguros, servicios inmobiliarios y a las empresas y servicios personales, sociales y comunales.
  • SL.FAM.WORK.ZS: % total de trabajadores familiares
  • SL.EMP.MPYR.ZS: % total de empleadores
  • SL.EMP.SELF.ZS: % total de independientes (cuenta propias)
  • SL.EMP.WORK.ZS: % total de trabajadores asalariados

El segundo (HDI_HDR2020_040722.csv) son las componentes del Human Development Index (HDI) elaborado por UNDP para el período 1992-2019. El HDI es un promedio a nivel país de cuatro variables fundamentales:

En el dataset va a encontrarlas llamadas de la siguiente forma: - le_2019: expectativa de vida al nacer de la población - eys_2019: años “esperados” de escolaridad de la población - mys_2019: años promedios de escolaridad - gnipc_2019: producto bruto interno per cápita - hdi_2019: HDI para el país

En ambos datasets van a encontrar los códigos de país (iso3c) y las etiquetas de país (country)

Preprocesamiento

Consigna 1.

Cargue los dos datasets. ¿Cuál es la estructura de los datos? ¿Qué formato tienen? ¿Cuál es la unidad de análisis y cuáles las variables? ¿Cuál es el alcance temporal y geográfico?

> df <- read_csv('./data/wb_bank_data_2019.csv')
> hdi <- read_csv('./data/HDI_HDR2020_040722.csv')
> 
> df <- df %>%
+      pivot_wider(id_cols = c(iso3c, iso2c, country), names_from = indicatorID, values_from = value)

Consigna 2.

Seleccione del dataset del HDI solamente las columnas que corresponden al año 2019.

> hdi <- hdi %>%
+         select(iso3c, country, contains("_2019"))
>

Consigna 3.

Con el dataset recortado, realice un join para tener todas las variables en una sola tabla. ¿Qué variable debería usar para vincular ambas tablas?

> final <- hdi %>%
+             left_join(
+               df %>% select(-country)
+               )
> 
> final <- final %>% drop_na()
>

Análisis descriptivo

Consigna 4.

’¿Qué relacion existe (si es que existe) entre el HDI y la proporción de población asalariada? ¿Y entre el HDI y la proporción de población en el sector agropecuario?

> final %>%
+   ggplot() + 
+     geom_point(aes(x=SL.EMP.WORK.ZS, y=hdi_2019)) + 
+     ylim(0,1) + 
+     xlim(0,100) +
+     theme_minimal()

> final %>%
+   ggplot() + 
+     geom_point(aes(x=SL.AGR.EMPL.ZS, y=hdi_2019)) + 
+     ylim(0,1) + 
+     xlim(0,100) +
+     theme_minimal()

Correlación

> final %>%
+   select(SL.AGR.EMPL.ZS,SL.EMP.WORK.ZS, hdi_2019) %>%
+   cor(., method='pearson')
##                SL.AGR.EMPL.ZS SL.EMP.WORK.ZS   hdi_2019
## SL.AGR.EMPL.ZS      1.0000000     -0.9044272 -0.8401293
## SL.EMP.WORK.ZS     -0.9044272      1.0000000  0.8735794
## hdi_2019           -0.8401293      0.8735794  1.0000000

Modelo de regresión lineal múltiple

Consigna 5

Queremos construir información sobre la relación entre el PBI per cápita (gnipc_2019) y variables relativas a la educación, a la salud y al grado de desarrollo capitalista del país. Entrenar una regresión lineal múltiple con las variables que considere relevantes.

> final <- final %>%
+                 mutate(gnipc_2019_log = log10(gnipc_2019))
> 
> lm_1 <- final %>%
+           lm(gnipc_2019_log ~ le_2019 + mys_2019 + SL.AGR.EMPL.ZS + SL.SRV.EMPL.ZS + SL.EMP.WORK.ZS -gnipc_2019_log, data = .)
> 
> summary(lm_1)
## 
## Call:
## lm(formula = gnipc_2019_log ~ le_2019 + mys_2019 + SL.AGR.EMPL.ZS + 
##     SL.SRV.EMPL.ZS + SL.EMP.WORK.ZS - gnipc_2019_log, data = .)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.52606 -0.09054  0.00492  0.11069  0.80886 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     1.860987   0.295885   6.290 2.57e-09 ***
## le_2019         0.020951   0.003577   5.857 2.36e-08 ***
## mys_2019        0.039773   0.008888   4.475 1.39e-05 ***
## SL.AGR.EMPL.ZS -0.002743   0.002492  -1.101 0.272569    
## SL.SRV.EMPL.ZS  0.001323   0.002404   0.550 0.582827    
## SL.EMP.WORK.ZS  0.005265   0.001518   3.468 0.000664 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2003 on 171 degrees of freedom
## Multiple R-squared:  0.8477, Adjusted R-squared:  0.8432 
## F-statistic: 190.3 on 5 and 171 DF,  p-value: < 2.2e-16

Análisis de supuestos del modelo

Consigna 6.

Chequeen los residuos del modelo. ¿Qué herramienta gráfica pueden usar para eso? ¿Son aleatorios?

> resids <- tibble(
+   iso3c = final$iso3c,
+   gnipc = final$gnipc_2019_log,
+   gnipc_log_pred = predict(lm_1, final),
+   resid =  gnipc - gnipc_log_pred
+ ) 
> 
> resids %>%
+   ggplot(aes(x=gnipc_log_pred, y=resid, color=resid)) + 
+     geom_point() +
+     scale_color_viridis_c() +
+     geom_hline(yintercept = 0, linetype='dashed') + 
+     theme_minimal()

Consigna 7.

Realicen un histograma de los residuos. ¿Qué forma tienen?

> resids %>%
+   ggplot() + 
+     geom_histogram(aes(x=resid)) +
+     theme_minimal()

> plot(lm_1)

Consigna 8.

¿Parece que la condición de variabilidad constante del modelo es respetada? ¿Por qué?

> ###

Consigna 9. Colinealdad

¿Hay multicolinealidad en el modelo? ¿De qué grado? ¿Es un problema?

> vif(lm_1)
##        le_2019       mys_2019 SL.AGR.EMPL.ZS SL.SRV.EMPL.ZS SL.EMP.WORK.ZS 
##       3.166152       3.442897      12.643221       7.636922       7.335243

Análisis de resultados

Consigna 10.

¿Cuál es el ajuste del modelo?

> summary(lm_1)$adj.r.squared
## [1] 0.8432293

Consigna 11.

¿Cómo interpreta los resultados del modelo? ¿Cuál es, en su opinion, la variable más importante?

> lm_scaled <- final %>%
+         select(gnipc_2019_log, le_2019, mys_2019, SL.AGR.EMPL.ZS, SL.SRV.EMPL.ZS, SL.EMP.WORK.ZS) %>%
+         mutate(across(everything(), scale)) %>%
+         lm(gnipc_2019_log ~ ., data=.)
> tidy(lm_scaled)
## # A tibble: 6 × 5
##   term            estimate std.error statistic      p.value
##   <chr>              <dbl>     <dbl>     <dbl>        <dbl>
## 1 (Intercept)    -4.73e-16    0.0298 -1.59e-14 1.000       
## 2 le_2019         3.11e- 1    0.0531  5.86e+ 0 0.0000000236
## 3 mys_2019        2.48e- 1    0.0554  4.47e+ 0 0.0000139   
## 4 SL.AGR.EMPL.ZS -1.17e- 1    0.106  -1.10e+ 0 0.273       
## 5 SL.SRV.EMPL.ZS  4.54e- 2    0.0825  5.50e- 1 0.583       
## 6 SL.EMP.WORK.ZS  2.80e- 1    0.0808  3.47e+ 0 0.000664