El ingreso per cápita varía entre municipios de Puerto Rico. Factores como la educación, desempleo, acceso a internet, pobreza y tamaño del hogar pueden ayudar a explicar estas diferencias. El objetivo es aplicar Regresión Lineal Múltiple (RLM) y técnicas de selección de variables para identificar los predictores más relevantes.

Las variables utilizadas son las siguientes:

ing_pc <- Ingreso per cápita.

bach_pct <- % con bachillerato o más.

desem_pct <- Tasa de desempleo.

banda_pct <- % de hogares con internet de banda ancha.

pobreza_pct <- % de población bajo pobreza.

tam_hogar <- Tamaño promedio del hogar.

Para llevar a cabo este ejercicio, trabajaremos con datos descargados desde la página del Censo de PR.

library(readr)
data_rlm <- read.csv("datos/data_rlm.csv")

Primero, analizaremos la correlación que poseen las variables entre sí. Esto lo podemos lograr ejecutando una gráfica de correlaciónes.

#correlaciones
library(corrplot)
R     <- cor(datos, use = "pairwise.complete.obs")

corrplot(R,
         method = "color",
         type = "lower",
         addCoef.col = "black",
         tl.col = "black",
         tl.srt = 45,
         diag = FALSE)

Podemos claramente observar, como el % de personas con alto nivel de educación y el % de hogares con internet de banda ancha, tienen una correlación positiva significativa al ingreso per cápita en el municipio. Esto nos sugiere que en Puerto Rico, las personas con alto nivel de educación, logran alcanzar mayor ingreso, lo que influye en que estas personas logren adquirir internet de banda ancha en sus hogares.

Por el contrario, notamos que el % de población bajo pobreza en los municipios tiene una correlacion negativa significativa, lo que indica que a mayor nivel de pobreza, menor ingreso per cápita en dichos municipios.

Una vez evaluada la correlacion entre las variables, podemos crear un modelo saturado para ver que tan significativas son las variabes respecto al ingreso per cápita.

saturado <- ing_pc ~ bach_pct + desem_pct + banda_pct + pobreza_pct + tam_hogar 

reg  <- lm(saturado, data=datos) # todos los predictores
summary(reg)
## 
## Call:
## lm(formula = saturado, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4585.2 -1168.6  -125.4  1208.2  5224.9 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  3757.23    3338.95   1.125 0.264211    
## bach_pct      264.14      68.73   3.843 0.000259 ***
## desem_pct      20.05      74.85   0.268 0.789533    
## banda_pct     125.32      45.30   2.766 0.007198 ** 
## pobreza_pct  -210.64      47.68  -4.417 3.45e-05 ***
## tam_hogar     816.67     622.90   1.311 0.193997    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1982 on 72 degrees of freedom
## Multiple R-squared:  0.8457, Adjusted R-squared:  0.835 
## F-statistic: 78.93 on 5 and 72 DF,  p-value: < 2.2e-16

Con ayuda de leaps, buscamos la mejor combinacion de variables para el modelo.

library(leaps)

ajuste <- regsubsets(
  saturado, data = datos,
  nvmax = 6,                 
  method = "exhaustive")    

saturado <- summary(ajuste)

names(saturado)
## [1] "which"  "rsq"    "rss"    "adjr2"  "cp"     "bic"    "outmat" "obj"

Luego, analizaremos y escogeremos el mejor número y combinación de variables. Para esto hay dos formas, utilizando el ajuste de modelos (hacia adelante, atrás o híbrido) usando BIC o ajustando modelos comparando el BIC y R^2 - ajustado. Comenzaremos con estos últimos.

par(mfrow = c(1,2))

plot(saturado$bic, type="b", col="black", pch=19,
     xlab="Número de predictores", ylab="BIC",
     main="Criterio BIC")
points(best_bic, saturado$bic[best_bic], pch=19, cex=1.5, col="red")

plot(saturado$adjr2, type="b", col="black", pch=19,
     xlab="Número de predictores", ylab="R^2 - ajustado",
     main="Criterio R^2 - ajustado")
points(best_r2, saturado$adjr2[best_r2], pch=19, cex=1.5, col="blue")

Notamos que hay una discrepancia entre ambos modelos. Esto se debe a que BIC suele premiar la simplicidad y penalizar la complejidad. Al contrario R^2 - ajustado premia al modelo más grande. Por eso BIC nos indica 3 y R^2 - Ajustado, 4.

Por consiguiente, podemos hacer el modelo ajustado final, extrayendo las variables de cada método.

BIC:

#Extracción de las variables del BIC
coef_best <- coef(ajuste, best_bic)
vars_best <- names(coef_best)[-1] # sin intercepto

# Modelo final ajustado
f_final <- as.formula(paste("ing_pc~", paste(vars_best, collapse = " + ")))

modelo_final <- lm(f_final, data = datos)
summary(modelo_final)
## 
## Call:
## lm(formula = f_final, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4644.3  -995.6   -35.5  1360.0  5031.6 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  6481.30    2522.94   2.569 0.012215 *  
## bach_pct      270.45      67.62   3.999 0.000149 ***
## banda_pct     119.54      44.32   2.697 0.008657 ** 
## pobreza_pct  -209.43      47.36  -4.422  3.3e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1979 on 74 degrees of freedom
## Multiple R-squared:  0.842,  Adjusted R-squared:  0.8356 
## F-statistic: 131.5 on 3 and 74 DF,  p-value: < 2.2e-16

Observando este modelo final del BIC, entendemos que las variables tienen una relación importante. Observamos que el ingreso per cápita en los municipios es de 6481 dólares, cuando todas las demás variables son 0. También reafirmamos la hipótesis de que el % de pobreza en la población y el % de personas con alta educacion son las variables más significativas en este modelo. Por cada aumento en el % de personas con nivel de alta educacion, el ingreso aumenta 270 dólares. Por el contrario por cada aumento en el % de personas en nivel de pobreza, el ingreso per cápita en los municipios disminuye 209 dólares.

R^2 - Ajustado:

coef_best <- coef(ajuste, best_r2)
vars_best <- names(coef_best)[-1] 

# modelo final ajustado
f_final <- as.formula(paste("ing_pc ~", paste(vars_best, collapse = " + ")))

modelo_final <- lm(f_final, data = datos)
summary(modelo_final)
## 
## Call:
## lm(formula = f_final, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4630.0 -1129.5   -83.4  1182.5  5207.2 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  4047.29    3138.40   1.290 0.201259    
## bach_pct      267.16      67.36   3.966 0.000169 ***
## banda_pct     123.03      44.21   2.783 0.006848 ** 
## pobreza_pct  -209.38      47.15  -4.441 3.12e-05 ***
## tam_hogar     790.86     611.48   1.293 0.199967    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1970 on 73 degrees of freedom
## Multiple R-squared:  0.8456, Adjusted R-squared:  0.8371 
## F-statistic: 99.92 on 4 and 73 DF,  p-value: < 2.2e-16

Observando este modelo final del R^2 - Ajustado, entendemos que las variables tienen una relación considerable. Observamos que el ingreso per cápita en los municipios es de 4047 dólares, cuando todas las demás variables son 0, menor que en el anterior modelo final. Volvemos a reafirmar la hipótesis de que el % de pobreza en la población y el % de personas con alta educacion son las variables más significativas en los modelos. Por cada aumento en el % de personas con nivel de alta educacion, el ingreso aumenta 267 dólares. Por el contrario por cada aumento en el % de personas en nivel de pobreza, el ingreso per cápita en los municipios tambien disminuye 209 dólares. Interesantemente, aprendemos que el tamaño del hogar no es una variable significativa en el modelo.

Como dijimos anteriormente, también podemos escoger el mejor número y combinación de variables con los paso a paso que se describen a continuacion:

Paso a paso hacia adelante:

# Hacia adelante 
m1_step_log <- step(m0_log,
                    scope     = list(lower = ~1, upper = formula(mF_log)),
                    direction = "forward",
                    k         = log(nrow(datos)),
                    trace     = TRUE)

Paso a paso hacia atrás:

# Hacia atras 
m2_step_log <- step(mF_log,
                    direction = "backward",
                    k         = log(nrow(datos)),
                    trace     = TRUE)

Paso a paso híbrido:

# híbrido 
m3_step_log <- step(m0_log,
                    scope     = list(lower = ~1, upper = formula(mF_log)),
                    direction = "both",
                    k         = log(nrow(datos)),
                    trace     = TRUE)

Debido a que los tres paso a paso llegaron a la misma conclusión, utilizaremos el método de selección híbrido, porque se supone que es el más complejo y adecuado para nuestro ejercicio:

summary(m3_step_log)
## 
## Call:
## lm(formula = ing_pc ~ bach_pct + pobreza_pct + banda_pct, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4644.3  -995.6   -35.5  1360.0  5031.6 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  6481.30    2522.94   2.569 0.012215 *  
## bach_pct      270.45      67.62   3.999 0.000149 ***
## pobreza_pct  -209.43      47.36  -4.422  3.3e-05 ***
## banda_pct     119.54      44.32   2.697 0.008657 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1979 on 74 degrees of freedom
## Multiple R-squared:  0.842,  Adjusted R-squared:  0.8356 
## F-statistic: 131.5 on 3 and 74 DF,  p-value: < 2.2e-16

Una vez más, recibimos un modelo con una relación significativa. Este método híbrido alcanzó unos números iguales a los del método BIC anterior, obteniendo las mismas conclusiones.

En definitiva, podemos concluir que el ingreso per cápita de los habitantes de los municipios está altamente relacionado al % de personas con bachillerato o más, de manera positiva, y al % de personas en nivel de pobreza, de manera negativa. También está medianamente relacionada al % de personas con internet de banda ancha en sus hogares, de manera positiva. Con esto, nos damos cuenta de la importancia de que los habitantes de Puerto Rico puedan recibir una buena educación, ya que esto aumenta el ingreso y la calidad de vida de los puertorriqueños.