TAREA 5
Luis Ochoa, Rafael Soqui, Cinthia Felix, Omar Cortez
28/2/2022
#Regresion lineal multiple
library(pacman)
p_load("prettydoc", "DT", "xfun", "dplyr", "psych", "GGally", "ggplot2","readr", "gridExtra", "corrplot")Un modelo de regresión lineal múltiple es un modelo estadístico versátil para evaluar las relaciones entre un destino continuo y los predictores. Los predictores pueden ser campos continuos, categóricos o derivados, de modo que las relaciones no lineales también estén soportadas.
Los modelos lineales múltiples siguen la siguiente ecuación:
\[ Y_{i}=(\beta_{0}+\beta_{1}X_{1i}+\beta_{2}X_{2i}+\cdots+\beta_{n}X_{ni})+e_{i} \] * \(\beta_{0}\) : es la ordenada en el origen, el valor de la variable dependiente Y cuando todos los predictores son cero.
\(\beta_{i}\) : es el efecto promedio que tiene el incremento en una unidad de la variable predictora Xi sobre la variable dependiente Y manteniéndose constantes el resto de variables. Se conocen como coeficientes parciales de regresión.
\(e_{i}\) : es el residuo o error, la diferencia entre el valor observado y el estimado por el modelo.
Objetivo de la asignación
En esta asignacion se nos señalo hacer un analisis de los datos de esperanza de vida en estados unidos y despues crear un modelo que pueda predecir la esperanza de vida promedio en funcion de las variables se dispone sobre:
- habitantes,
- analfabetismo,
- ingresos,
- esperanza de vida,
- asesinatos,
- universitarios,
- heladas,
- área
- y densidad poblacional.
Datos
datos <- as.data.frame(state.x77)
datatable(datos)Renombrar y modificar el conjunto de datos
datos <- rename(habitantes = Population, analfabetismo = Illiteracy,
ingresos = Income, esp_vida = `Life Exp`, asesinatos = Murder,universitarios = `HS Grad`, heladas = Frost, area = Area,.data = datos)
datos <- mutate(.data=datos, densidad_pobl = habitantes * 1000 / area)Analsis de datos
round(cor(x = datos, method = "pearson"), 3)## habitantes ingresos analfabetismo esp_vida asesinatos
## habitantes 1.000 0.208 0.108 -0.068 0.344
## ingresos 0.208 1.000 -0.437 0.340 -0.230
## analfabetismo 0.108 -0.437 1.000 -0.588 0.703
## esp_vida -0.068 0.340 -0.588 1.000 -0.781
## asesinatos 0.344 -0.230 0.703 -0.781 1.000
## universitarios -0.098 0.620 -0.657 0.582 -0.488
## heladas -0.332 0.226 -0.672 0.262 -0.539
## area 0.023 0.363 0.077 -0.107 0.228
## densidad_pobl 0.246 0.330 0.009 0.091 -0.185
## universitarios heladas area densidad_pobl
## habitantes -0.098 -0.332 0.023 0.246
## ingresos 0.620 0.226 0.363 0.330
## analfabetismo -0.657 -0.672 0.077 0.009
## esp_vida 0.582 0.262 -0.107 0.091
## asesinatos -0.488 -0.539 0.228 -0.185
## universitarios 1.000 0.367 0.334 -0.088
## heladas 0.367 1.000 0.059 0.002
## area 0.334 0.059 1.000 -0.341
## densidad_pobl -0.088 0.002 -0.341 1.000- Análisis con histogramas
multi.hist( x = datos, dcol = c("blue","red"), dlty = c("dotted", "solid"),
main = ""
)
ggpairs(datos, lower = list(continuous ="smooth"),
diag = list (continuos = "barDiag"), axisLabels = "none"
)
Entonces, de los análisis realizados hasta el momento, podemos obtener las siguientes conclusiones preliminares:
- Las variables que tienen una mayor relación lineal con la esperanza de vida son: asesinatos (r= -0.78), analfabetismo (r= -0.59) y universitarios (r= 0.58).
Se usara la variable de asesinatos para hacer el primer modelo
modelo <- lm(esp_vida ~ asesinatos, data=datos)
summary(modelo)##
## Call:
## lm(formula = esp_vida ~ asesinatos, data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.81690 -0.48139 0.09591 0.39769 2.38691
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 72.97356 0.26997 270.30 < 2e-16 ***
## asesinatos -0.28395 0.03279 -8.66 2.26e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.8473 on 48 degrees of freedom
## Multiple R-squared: 0.6097, Adjusted R-squared: 0.6016
## F-statistic: 74.99 on 1 and 48 DF, p-value: 2.26e-11El modelo creado tiene un coeficiente de determinacion ajustado de 0.6479 lo que nos indica que el 64.79% de los casos este modelo estaria en lo correcto.
asi mismo su p-value es bastante bajo lo que nos indicaria que no hay mucha probabilidad de datos extremos.
Tomando los resultados obtenidos el modelo quedaria asi:
\[ y = 72.97356 - 0.28395x \]
Evaluando el modelo
y = 72.97356 - (0.28395*11.6)
y## [1] 69.67974Conocer el residuo según este modelo
y - 67.96## [1] 1.71974El residuo de 1.72 nos indicaria que el modelo esta subestimando
Analisis de predictores
modelo2 <- lm(esp_vida ~ asesinatos + habitantes + ingresos + universitarios + area + densidad_pobl, data=datos)
summary(modelo2)##
## Call:
## lm(formula = esp_vida ~ asesinatos + habitantes + ingresos +
## universitarios + area + densidad_pobl, data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.64532 -0.45832 -0.00947 0.44113 2.15586
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 7.029e+01 1.260e+00 55.802 < 2e-16 ***
## asesinatos -2.707e-01 4.288e-02 -6.313 1.28e-07 ***
## habitantes 6.985e-05 2.888e-05 2.419 0.0199 *
## ingresos 1.180e-04 3.026e-04 0.390 0.6986
## universitarios 3.722e-02 2.280e-02 1.633 0.1098
## area -1.276e-06 1.752e-06 -0.728 0.4703
## densidad_pobl -7.872e-04 6.876e-04 -1.145 0.2586
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7702 on 43 degrees of freedom
## Multiple R-squared: 0.7111, Adjusted R-squared: 0.6708
## F-statistic: 17.64 on 6 and 43 DF, p-value: 3.519e-10Este modelo que incluye todos los predictores tiende a ser más exacto que el primero ya que su r cuadrada y su p-value son considerablemente mayor y menor respectivamente.
Selección de predictores
para la seleccion de los predictores se usara el metodo AIC que es una medida de la calidad relativa de un modelo estadístico, para un conjunto dado de datos, en este caso nos indicara cuales indicadores no son significativos para nuestro modelo.
step(object = modelo2, direction = "both", trace = 1)## Start: AIC=-19.66
## esp_vida ~ asesinatos + habitantes + ingresos + universitarios +
## area + densidad_pobl
##
## Df Sum of Sq RSS AIC
## - ingresos 1 0.0901 25.596 -21.480
## - area 1 0.3147 25.820 -21.043
## - densidad_pobl 1 0.7775 26.283 -20.155
## <none> 25.506 -19.656
## - universitarios 1 1.5814 27.087 -18.648
## - habitantes 1 3.4709 28.977 -15.277
## - asesinatos 1 23.6408 49.146 11.139
##
## Step: AIC=-21.48
## esp_vida ~ asesinatos + habitantes + universitarios + area +
## densidad_pobl
##
## Df Sum of Sq RSS AIC
## - area 1 0.2312 25.827 -23.0301
## - densidad_pobl 1 0.7427 26.338 -22.0496
## <none> 25.596 -21.4798
## + ingresos 1 0.0901 25.506 -19.6561
## - universitarios 1 3.0763 28.672 -17.8050
## - habitantes 1 3.8588 29.455 -16.4588
## - asesinatos 1 23.6039 49.200 9.1932
##
## Step: AIC=-23.03
## esp_vida ~ asesinatos + habitantes + universitarios + densidad_pobl
##
## Df Sum of Sq RSS AIC
## - densidad_pobl 1 0.6028 26.430 -23.877
## <none> 25.827 -23.030
## + area 1 0.2312 25.596 -21.480
## + ingresos 1 0.0067 25.820 -21.043
## - universitarios 1 3.0394 28.866 -19.467
## - habitantes 1 3.9396 29.767 -17.932
## - asesinatos 1 30.5433 56.370 13.996
##
## Step: AIC=-23.88
## esp_vida ~ asesinatos + habitantes + universitarios
##
## Df Sum of Sq RSS AIC
## <none> 26.430 -23.877
## + densidad_pobl 1 0.603 25.827 -23.030
## + area 1 0.091 26.338 -22.050
## + ingresos 1 0.091 26.339 -22.048
## - habitantes 1 3.341 29.770 -19.925
## - universitarios 1 4.015 30.445 -18.805
## - asesinatos 1 31.928 58.358 13.728##
## Call:
## lm(formula = esp_vida ~ asesinatos + habitantes + universitarios,
## data = datos)
##
## Coefficients:
## (Intercept) asesinatos habitantes universitarios
## 70.4146754 -0.2664148 0.0000625 0.0407496El metodo AIC considera a los indicadores analfabetismo, area, densidad poblacional e ingresos como no significativos para nuestro modelo
realizando el modelo excluyendo los indicadores indicados
modelo3 <- (lm(formula = esp_vida ~ habitantes + asesinatos + universitarios + heladas, data = datos))
summary(modelo3)##
## Call:
## lm(formula = esp_vida ~ habitantes + asesinatos + universitarios +
## heladas, data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.47095 -0.53464 -0.03701 0.57621 1.50683
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 7.103e+01 9.529e-01 74.542 < 2e-16 ***
## habitantes 5.014e-05 2.512e-05 1.996 0.05201 .
## asesinatos -3.001e-01 3.661e-02 -8.199 1.77e-10 ***
## universitarios 4.658e-02 1.483e-02 3.142 0.00297 **
## heladas -5.943e-03 2.421e-03 -2.455 0.01802 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7197 on 45 degrees of freedom
## Multiple R-squared: 0.736, Adjusted R-squared: 0.7126
## F-statistic: 31.37 on 4 and 45 DF, p-value: 1.696e-12Excluyendo las variables no significativas se logro una mejora en la r cuadrada ajustada de poco más del 1% lo cual no es mucho pero el p_value disminuyo bastante lo cual indica que se eliminaron datos atipicos del modelo haciendolo mucho mas exacto.
ahora, el intervalo de confianza
confint(lm(formula = esp_vida ~ habitantes + asesinatos + universitarios +
heladas, data = datos))## 2.5 % 97.5 %
## (Intercept) 6.910798e+01 72.9462729104
## habitantes -4.543308e-07 0.0001007343
## asesinatos -3.738840e-01 -0.2264135705
## universitarios 1.671901e-02 0.0764454870
## heladas -1.081918e-02 -0.0010673977plot1 <- ggplot(data = datos, aes(habitantes, modelo$residuals)) +
geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
theme_bw()
plot2 <- ggplot(data = datos, aes(asesinatos, modelo$residuals)) +
geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
theme_bw()
plot3 <- ggplot(data = datos, aes(universitarios, modelo$residuals)) +
geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
theme_bw()
plot4 <- ggplot(data = datos, aes(heladas, modelo$residuals)) +
geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
theme_bw()
grid.arrange(plot1, plot2, plot3, plot4)## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
El grafico de habitantes nos indica que la mayoria de las ciudades tenian menos de 6000 habitantes
El grafico de universitarios es el mas estable indicando un crecimiento ordenado, esto se debe al boom de educacion media- superior que hubo en estados unidos entre 1900 y 1970.
El grafico de asesinatos tiende a ser el mas disperso, esto es probablemente porque los asesinatos son algo poco comun
Las heladas sigen un patron de ascenso y descenso con respecto al tiempo.
Distribución normal de los residuos:
qqnorm(modelo3$residuals)
qqline(modelo3$residuals)
Prueba de normalidad de Shapiro-wilk
shapiro.test(modelo3$residuals)##
## Shapiro-Wilk normality test
##
## data: modelo3$residuals
## W = 0.97935, p-value = 0.525determinacion de homocedasticidad
ggplot(data = datos, aes(modelo$fitted.values, modelo$residuals)) +
geom_point() +
geom_smooth(color = "firebrick", se = FALSE) +
geom_hline(yintercept = 0) +
theme_bw()## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
La grafica prueba tener una distribucion aceptable lo que indicaria que hay una varianza relativamente constante en el modelo haciendolo fiable
Conclusión
Se puede concluir que el modelo para predecir tomara como predictores a la poblacion, la tasa de asesinatos, el porcentaje de graduados de preparatoria y las heladas
Se excluyeron los factores de area, analfabetismo e ingresos ya que estan bastante relacionados con la tasa de universitarios
el modelo final es el siguiente:
\[ E_vida = (7.103 + 5.014e^-5 * X1 -0.33 * X2 + 0.04658 * X3 -0.005.943 * X4)+0.7197 \]