Regresión Lineal Múltiple
La regresión lineal múltiple nos permite generar un modelo lineal dónde el valor de la variable dependiente es determinado a partir de un conjunto de variables independientes llamadas predictores (esto último se debe que analizar con cautela para no malinterpretar causa-efecto).
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 las 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.
Caso de Estudio
Esperanza de Vida en EU En el presente documento se mostrará información estadística a partir de un análisis previo acerca de la esperanza de vida de los habitantes de Estados Unidos.
Esperanza de Vida
Definición
La esperanza de vida es una estimación del número promedio de años de vida adicionales que una persona podría esperar vivir si las tasas de mortalidad por edad específica para un año determinado permanecieran durante el resto de su vida. La esperanza de vida es una medida hipotética porque se basa en las tasas de mortalidad actuales, pero las tasas de mortalidad reales cambian durante el transcurso de la vida de una persona. La esperanza de vida de cada persona cambia según la misma envejece y a medida que cambian las tendencias de mortandad.
Debido a que la esperanza de vida varía de manera significante según el sexo, la edad actual y la raza, normalmente estas categorías se nombran por separado. La esperanza de vida al nacer es la medida más comúnmente citada. Es un buen indicador de las condiciones de salud actuales.
Factores
Para poder hacer un cálculo exacto se estudia dos factores específicos que explican el aumento o disminución de la esperanza de vida de un grupo humano.
Genética. Para empezar está el factor de los genes que influye en un tercio del estudio. La genética nos demostrará si el grupo humano que estamos analizando tiene propensión a vivir más o menos.
Ambiente. Estos factores corresponden a los dos tercios restantes del estudio. Dentro de los factores ambientales está un amplio abanico de causas que influyen en el aumento o disminución de la esperanza de vida.
Factores que influyen en la esperanza de vida:
Dieta
Salud pública
Educación
Actividad física
Factores de apoyo
Paquetes
library(pacman)
p_load("prettydoc", "DT", "xfun", "dplyr", "psych", "GGally", "ggplot2","readr","stringr","vembedr", "xfun", "gridExtra","corrplot")Datos
Los datos utilizados dentro de este documento fueron obtenidos desde las librerias de R.
A continuación se muestran los datos obtenidos de nuestra fuente mencionada anteriormente:
-Área: Territorio en millas cuadradas.
-Heladas: Media de días con temperaturas mínimas menores a 0°C en la Capital o Ciudad más Grande.
-Habitantes: Población estimada a 1 de Julio de 1975.
-Ingresos: Ingreso per capita.
-Asesinatos: Ratio de asesinatos por cada 100,000 habitantes.
-Analfabetismo: Analfabetismo.
-Esperanza de Vida: Esperanza de vida en los años.
-Universitarios: Porcentaje de universitarios.
datos <- as.data.frame(state.x77)
datatable(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)1. Análisis
Análisis de Correlación
El análisis de correlación consiste en un procedimiento estadístico para determinar si dos variables están relacionadas o no. El resultado del análisis es un coeficiente de correlación que puede tomar valores entre -1 y +1. El signo indica el tipo de correlación entre las dos variables. A continuación se realizará el análisis de correlación:
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.000Análisis con histogramas
multi.hist( x = datos, dcol = c("red","green"), dlty = c("dotted", "solid")
)
Análisis de Dispersión de los datos utilizando ggplot y ggally
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).
Asesinatos y analfabetismo están medianamente correlacionados (r = 0.7) por lo que posiblemente no sea útil introducir ambos predictores en el modelo.
Las variables habitantes, área y densidad poblacional muestran una distribución exponencial, una transformación logarítmica posiblemente haría más normal su distribución.
2.Generar modelo de regresión lineal.
Como se ha explicado en la introducción, hay diferentes formas de llegar al modelo final más adecuado. En este caso se va a emplear el método mixto iniciando el modelo con todas las variables como predictores y realizando la selección de los mejores predictores con la medición Akaike(AIC)
Modelo general de correalación
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-11Ecuación de la recta de los mínimos cuadrados:
\[ y=72.97356−0.28395x \]
Evaluando gráficamente el modelo.
plot(datos$asesinatos, datos$esp_vida)
abline(modelo)
Ahora substituyendo X en nuestro modelo
y = 72.97356 - (0.28395*11.5)
y## [1] 69.70814Para conocer el residuo según este modelo:
y - 69.03## [1] 0.678135Modelo de regresión lineal múltiple
Se evaluaran las variables predictoras (X).
modelo2 <- (lm(formula = esp_vida ~ habitantes +ingresos+analfabetismo+ asesinatos + universitarios + heladas+area+densidad_pobl, data = datos))
summary(modelo2) ##
## Call:
## lm(formula = esp_vida ~ habitantes + ingresos + analfabetismo +
## asesinatos + universitarios + heladas + area + densidad_pobl,
## data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.47514 -0.45887 -0.06352 0.59362 1.21823
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.995e+01 1.843e+00 37.956 < 2e-16 ***
## habitantes 6.480e-05 3.001e-05 2.159 0.0367 *
## ingresos 2.701e-04 3.087e-04 0.875 0.3867
## analfabetismo 3.029e-01 4.024e-01 0.753 0.4559
## asesinatos -3.286e-01 4.941e-02 -6.652 5.12e-08 ***
## universitarios 4.291e-02 2.332e-02 1.840 0.0730 .
## heladas -4.580e-03 3.189e-03 -1.436 0.1585
## area -1.558e-06 1.914e-06 -0.814 0.4205
## densidad_pobl -1.105e-03 7.312e-04 -1.511 0.1385
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7337 on 41 degrees of freedom
## Multiple R-squared: 0.7501, Adjusted R-squared: 0.7013
## F-statistic: 15.38 on 8 and 41 DF, p-value: 3.787e-10El modelo tiene un R2 alta (0.7501), es capaz de explicar el 75,01% de la variabilidad observada en la esperanza de vida. El p-value del modelo es significativo (3.787e-10) por lo que se puede aceptar que el modelo no es por azar, al menos uno de los coeficientes parciales de regresión es distinto de 0.
Evaluacion del modelo(AIC)
El criterio de información de Akaike (AIC) es una medida de la calidad relativa de un modelo estadístico, para un conjunto dado de datos. Como tal, el AIC proporciona un medio para la selección del modelo.
step(object = modelo2, direction = "both", trace = 1)## Start: AIC=-22.89
## esp_vida ~ habitantes + ingresos + analfabetismo + asesinatos +
## universitarios + heladas + area + densidad_pobl
##
## Df Sum of Sq RSS AIC
## - analfabetismo 1 0.3050 22.373 -24.208
## - area 1 0.3564 22.425 -24.093
## - ingresos 1 0.4120 22.480 -23.969
## <none> 22.068 -22.894
## - heladas 1 1.1102 23.178 -22.440
## - densidad_pobl 1 1.2288 23.297 -22.185
## - universitarios 1 1.8225 23.891 -20.926
## - habitantes 1 2.5095 24.578 -19.509
## - asesinatos 1 23.8173 45.886 11.707
##
## Step: AIC=-24.21
## esp_vida ~ habitantes + ingresos + asesinatos + universitarios +
## heladas + area + densidad_pobl
##
## Df Sum of Sq RSS AIC
## - area 1 0.1427 22.516 -25.890
## - ingresos 1 0.2316 22.605 -25.693
## <none> 22.373 -24.208
## - densidad_pobl 1 0.9286 23.302 -24.174
## - universitarios 1 1.5218 23.895 -22.918
## + analfabetismo 1 0.3050 22.068 -22.894
## - habitantes 1 2.2047 24.578 -21.509
## - heladas 1 3.1324 25.506 -19.656
## - asesinatos 1 26.7071 49.080 13.072
##
## Step: AIC=-25.89
## esp_vida ~ habitantes + ingresos + asesinatos + universitarios +
## heladas + densidad_pobl
##
## Df Sum of Sq RSS AIC
## - ingresos 1 0.132 22.648 -27.598
## - densidad_pobl 1 0.786 23.302 -26.174
## <none> 22.516 -25.890
## - universitarios 1 1.424 23.940 -24.824
## + area 1 0.143 22.373 -24.208
## + analfabetismo 1 0.091 22.425 -24.093
## - habitantes 1 2.332 24.848 -22.962
## - heladas 1 3.304 25.820 -21.043
## - asesinatos 1 32.779 55.295 17.033
##
## Step: AIC=-27.6
## esp_vida ~ habitantes + asesinatos + universitarios + heladas +
## densidad_pobl
##
## Df Sum of Sq RSS AIC
## - densidad_pobl 1 0.660 23.308 -28.161
## <none> 22.648 -27.598
## + ingresos 1 0.132 22.516 -25.890
## + analfabetismo 1 0.061 22.587 -25.732
## + area 1 0.043 22.605 -25.693
## - habitantes 1 2.659 25.307 -24.046
## - heladas 1 3.179 25.827 -23.030
## - universitarios 1 3.966 26.614 -21.529
## - asesinatos 1 33.626 56.274 15.910
##
## Step: AIC=-28.16
## esp_vida ~ habitantes + asesinatos + universitarios + heladas
##
## Df Sum of Sq RSS AIC
## <none> 23.308 -28.161
## + densidad_pobl 1 0.660 22.648 -27.598
## + ingresos 1 0.006 23.302 -26.174
## + analfabetismo 1 0.004 23.304 -26.170
## + area 1 0.001 23.307 -26.163
## - habitantes 1 2.064 25.372 -25.920
## - heladas 1 3.122 26.430 -23.877
## - universitarios 1 5.112 28.420 -20.246
## - asesinatos 1 34.816 58.124 15.528##
## Call:
## lm(formula = esp_vida ~ habitantes + asesinatos + universitarios +
## heladas, data = datos)
##
## Coefficients:
## (Intercept) habitantes asesinatos universitarios heladas
## 7.103e+01 5.014e-05 -3.001e-01 4.658e-02 -5.943e-03Despues del proceso AIC se muestra el modelo ajustado, donde se quita la variable Área, Analfabetismo e Ingreso.
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-12El modelo ajustado describe la esperanza de vida con un 73% de exactitud.
Intervalo de Confianza de los coeficientes
El intervalo de confianza para los coeficientes parciales de regresión se puede observar en la siguiente tabla:
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.0010673977Cada una de las pendientes de un modelo de regresión lineal múltiple (coeficientes parciales de regresión de los predictores) se define del siguiente modo: Si el resto de variables se mantienen constantes, por cada unidad que aumenta el predictor en cuestión, la variable (Y) varía en promedio tantas unidades como indica la pendiente. Para este ejemplo, por cada unidad que aumenta el predictor universitarios, la esperanza de vida aumenta en promedio 0.04658 unidades, manteniéndose constantes el resto de predictores
Validación de condiciones para la regresión múltiple lineal.
Relación lineal entre los predictores numéricos y la variable respuesta:
Esta condición se puede validar bien mediante diagramas de dispersión entre la variable dependiente y cada uno de los predictores (como se ha hecho en el análisis preliminar) o con diagramas de dispersión entre cada uno de los predictores y los residuos del modelo. Si la relación es lineal, los residuos deben de distribuirse aleatoriamente en torno a 0 con una variabilidad constante a lo largo del eje X. Esta última opción suele ser más indicada ya que permite identificar posibles datos atípicos.
plot1 <- 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'
En la grafica de habitantes los datos se agrupan a la izquierda entre 0 y 5000, los asesinatos se encuentran agrupados entre Y -1 y 1, los universitarios cuenta con una concentracion entre 50 y 60, las heladas estan dispersas, pero estas disminuyen con el tiempo.
Distribución normal de los residuos
qqnorm(modelo3$residuals)
qqline(modelo3$residuals)
* Prueba de shapiro wilk para determinar normalidad
shapiro.test(modelo$residuals)##
## Shapiro-Wilk normality test
##
## data: modelo$residuals
## W = 0.98431, p-value = 0.7416Tanto el análisis gráfico como es test de hipótesis confirman la normalidad.
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'
Variabilidad constante de los residuos (homocedasticidad)
Al representar los residuos frente a los valores ajustados por el modelo, los primeros se tienen que distribuir de forma aleatoria en torno a cero, manteniendo aproximadamente la misma variabilidad a lo largo del eje X. Si se observa algún patrón específico, por ejemplo forma cónica o mayor dispersión en los extremos, significa que la variabilidad es dependiente del valor ajustado y por lo tanto no hay homocedasticidad.
corrplot(cor(dplyr::select(datos, habitantes,asesinatos,universitarios,heladas,esp_vida)),
method = "number", tl.col = "black")
Esta grafica nos muestra la relación que existe entre las variables. Lo que se ve rojo nos dice que tienen una relacion mas fuerte.
Conclusion
La esperanza de vida es el tiempo de vida estimado de una persona, lo que puede ser influido por muchos factores, y como vimos en este caso, teniamos varios factores que afectaban a la esperanza de vida, los cuales algunos se fueron descartando como el ingreso monetario, ya que esta variable no tenia mucha importancia para determinar la esperanza de vida.
Este modelo esta hecho para determinar la esperanza de vida que había en estados unidos en 1970, si se quisiera determinar con algun otro año, o en la actualidad, se necesitarian mas datos, ya que con el tiempo mas factores influyen en la esperanza de vida.