Modelo de Regresion Lineal Multiple. Datos Esperanza de Vida

Descripción

Un estudio quiere generar un modelo que permita predecir la esperanza de vida media de los habitantes de una ciudad en función de diferentes variables.

Se dispone de información sobre: habitantes, analfabetismo, ingresos, esperanza de vida, asesinatos, universitarios, heladas, área y densidad poblacional.

Objetivo

Realizar un análisis de datos identificando coeficientes de correlación entre variables así como un modelo de regresión lineal múltiple para determianr predicciones con nuevos datos.

Objetivos Específicos

Cargar los datos
Explorar los datos
Limpiar los datos
Determinar e Interpretar Coeficientes de Correlación
Determinar el modelo de Regresión Lineal Múltiple (RLM)
Determinar y generar el mejor modelo de RLM
Interpretar estadísticos del mejor modelo
Realizar predicciones con el mejor modelo
Interpretar predicciones con el mejor modelo

Las Librerías

library(dplyr)   # selet,filter,mutate,rename,group by,summarize,...
library(readr)   # Leer datos
library(knitr)   # tablas amigables
library(ggplot2) # gráficos más amigables y potentes
library(knitr)   # Para tablas amigables

Los datos

Se cargan datos llamados state.77 que contiene algunas características de las ciudades de USA. Estos datos provienen de la librería dplyr

datos <- as.data.frame(state.x77)

Explorar datos

str(datos)

## 'data.frame':    50 obs. of  8 variables:
##  $ Population: num  3615 365 2212 2110 21198 ...
##  $ Income    : num  3624 6315 4530 3378 5114 ...
##  $ Illiteracy: num  2.1 1.5 1.8 1.9 1.1 0.7 1.1 0.9 1.3 2 ...
##  $ Life Exp  : num  69 69.3 70.5 70.7 71.7 ...
##  $ Murder    : num  15.1 11.3 7.8 10.1 10.3 6.8 3.1 6.2 10.7 13.9 ...
##  $ HS Grad   : num  41.3 66.7 58.1 39.9 62.6 63.9 56 54.6 52.6 40.6 ...
##  $ Frost     : num  20 152 15 65 20 166 139 103 11 60 ...
##  $ Area      : num  50708 566432 113417 51945 156361 ...

summary(datos)

##    Population        Income       Illiteracy       Life Exp    
##  Min.   :  365   Min.   :3098   Min.   :0.500   Min.   :67.96  
##  1st Qu.: 1080   1st Qu.:3993   1st Qu.:0.625   1st Qu.:70.12  
##  Median : 2838   Median :4519   Median :0.950   Median :70.67  
##  Mean   : 4246   Mean   :4436   Mean   :1.170   Mean   :70.88  
##  3rd Qu.: 4968   3rd Qu.:4814   3rd Qu.:1.575   3rd Qu.:71.89  
##  Max.   :21198   Max.   :6315   Max.   :2.800   Max.   :73.60  
##      Murder          HS Grad          Frost             Area       
##  Min.   : 1.400   Min.   :37.80   Min.   :  0.00   Min.   :  1049  
##  1st Qu.: 4.350   1st Qu.:48.05   1st Qu.: 66.25   1st Qu.: 36985  
##  Median : 6.850   Median :53.25   Median :114.50   Median : 54277  
##  Mean   : 7.378   Mean   :53.11   Mean   :104.46   Mean   : 70736  
##  3rd Qu.:10.675   3rd Qu.:59.15   3rd Qu.:139.75   3rd Qu.: 81163  
##  Max.   :15.100   Max.   :67.30   Max.   :188.00   Max.   :566432

Limpiar datos

Dentro de un proceso de limpieza de datos, además de encontrar anomalías, valores nulos, vacíos, datos atípicos, transformar valores, o quitar errores, entre otras cosas, también se puede modificar nombres de variables o columnas a gusto del analista, así como quitar o agregar columnas.

La recomendación es dejar en el idioma original para efecgos de brcar más en la profesión de científico de atos internacional.

Renombrar algunas los nombres de las variables con la funcipon rename() de la librería ggplot().
Agregar densidad de poblacion mediante \(densidad_pobl = habitantes * 1000 / area\) mediante la función mutate de la librería ggplot()

datos <- rename(habitantes = Population, 
                analfabetismo = Illiteracy,
                ingresos = Income, 
                esperanza_vida = `Life Exp`, 
                asesinatos = Murder,
                universitarios = `HS Grad`, 
                heladas = Frost, 
                area = Area,
                .data = datos)

datos <- mutate(.data = datos, densidad_pob = habitantes * 1000 / area)


kable(head(datos))

habitantes	ingresos	analfabetismo	esperanza_vida	asesinatos	universitarios	heladas	area	densidad_pob
3615	3624	2.1	69.05	15.1	41.3	20	50708	71.2905261
365	6315	1.5	69.31	11.3	66.7	152	566432	0.6443845
2212	4530	1.8	70.55	7.8	58.1	15	113417	19.5032491
2110	3378	1.9	70.66	10.1	39.9	65	51945	40.6198864
21198	5114	1.1	71.71	10.3	62.6	20	156361	135.5708904
2541	4884	0.7	72.06	6.8	63.9	166	103766	24.4877898

kable(tail(datos))

	habitantes	ingresos	analfabetismo	esperanza_vida	asesinatos	universitarios	heladas	area	densidad_pob
45	472	3907	0.6	71.64	5.5	57.1	168	9267	50.933420
46	4981	4701	1.4	70.08	9.5	47.8	85	39780	125.213675
47	3559	4864	0.6	71.72	4.3	63.5	32	66570	53.462521
48	1799	3617	1.4	69.48	6.7	41.6	100	24070	74.740341
49	4589	4468	0.7	72.48	3.0	54.5	149	54464	84.257491
50	376	4566	0.6	70.29	6.9	62.9	173	97203	3.868193

Coeficiente de Correlación

Se determina el Coeficiente de Correlación de todas las variables numéricas del conjunto de datos.

Usando la funcióm cor() de todos los datos y clasificando conforme a su clasificación:

−1.00 = correlación negativa perfecta. (“A mayor X, menor Y”, de manera proporcional. Es decir, cada vez que X aumenta una unidad, Y disminuye siempre una cantidad constante).Esto también se aplica “a menor X, mayor Y”.
−0.90 = Correlación negativa muy fuerte.
−0.75 = Correlación negativa considerable.
−0.50 = Correlación negativa media.
−0.25 = Correlación negativa débil.
−0.10 = Correlación negativa muy débil.
0.00 = No existe correlación alguna entre las variables.
+0.10 = Correlación positiva muy débil.
+0.25 = Correlación positiva débil.
+0.50 = Correlación positiva media.
+0.75 = Correlación positiva considerable.
+0.90 = Correlación positiva muy fuerte.
+1.00 = Correlación positiva perfecta (“A mayor X, mayor Y” o “a menor X, menor Y”, de manera proporcional. Cada vez que X aumenta, Y aumenta siempre una cantidad constante).

r <- cor(datos)

r <- as.data.frame(r)

kable(r)

	habitantes	ingresos	analfabetismo	esperanza_vida	asesinatos	universitarios	heladas	area	densidad_pob
habitantes	1.0000000	0.2082276	0.1076224	-0.0680520	0.3436428	-0.0984897	-0.3321525	0.0225438	0.2462279
ingresos	0.2082276	1.0000000	-0.4370752	0.3402553	-0.2300776	0.6199323	0.2262822	0.3633154	0.3299683
analfabetismo	0.1076224	-0.4370752	1.0000000	-0.5884779	0.7029752	-0.6571886	-0.6719470	0.0772611	0.0092743
esperanza_vida	-0.0680520	0.3402553	-0.5884779	1.0000000	-0.7808458	0.5822162	0.2620680	-0.1073319	0.0910618
asesinatos	0.3436428	-0.2300776	0.7029752	-0.7808458	1.0000000	-0.4879710	-0.5388834	0.2283902	-0.1850352
universitarios	-0.0984897	0.6199323	-0.6571886	0.5822162	-0.4879710	1.0000000	0.3667797	0.3335419	-0.0883672
heladas	-0.3321525	0.2262822	-0.6719470	0.2620680	-0.5388834	0.3667797	1.0000000	0.0592291	0.0022767
area	0.0225438	0.3633154	0.0772611	-0.1073319	0.2283902	0.3335419	0.0592291	1.0000000	-0.3413885
densidad_pob	0.2462279	0.3299683	0.0092743	0.0910618	-0.1850352	-0.0883672	0.0022767	-0.3413885	1.0000000

Interpretación del Coeficiente de Correlación

Con respecto a la variable **esperanza de vida*:
La variable de analfebitismo tiene un CR de -0.5884779 que significa negativa media.
La variable de asesinatos tiene un CR de -0.5884779 que significa negativa considerable.
La variable de universitarios tiene un CR de -0.5822162 que significa negativa considerable.

Estas variables pueder ser guía para un modelo de regresión sin embargo no son absolutamente las que más inciden toda vez elaborado el modelo.

Modelo de Regresión Lineal Múltiple (RLM)

Variable Depeneiente: Esperanza de Vida (esperanza_vida)
Variables Independientes: Todas las demás (habitantes + ingresos + analfabetismo + asesinatos + universitarios + heladas + area + densidad_pob)
La fórmula en Regresión Lineal Simple es : \[y = \beta0 + \beta1x_! + \varepsilon\]
\[y = a + bx + \varepsilon\]
La fórmula en la Regresión Lineal Múltiple para determinar \(y\) es \[ y = \beta0 + \beta1x_1 +\beta2x_2 + .... + \beta_nx_n + \varepsilon\]

\[y = \beta0 + \sum_{i=1}^n\beta_ix_i + \varepsilon\]

modelo <- lm(esperanza_vida ~ habitantes + ingresos + analfabetismo + asesinatos +
               universitarios + heladas + area + densidad_pob, data = datos )
summary(modelo)

## 
## Call:
## lm(formula = esperanza_vida ~ habitantes + ingresos + analfabetismo + 
##     asesinatos + universitarios + heladas + area + densidad_pob, 
##     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_pob   -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-10

El modelo con todas las variables introducidas como predictores tiene un \(R^2\) 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.
Muchos de ellos no son significativos, lo que es un indicativo de que podrían no contribuir al modelo.

Determinar mejores predictores

Se hace una recorrido de manera interna con la función step() para obtener como resultado el mejor modelo, se busca identifica cuales variabes representan variabilidad con respecto a la esperanza de vida, o lo que es lo mismo cuáles SI son importantes y cuáles NO con respecto a la variable dependiente esperanza de vida

Se van a emplear la estrategia de stepwise mixto. El valor matemático empleado para determinar la calidad del modelo va a ser Akaike(AIC).

step(object = modelo, direction = "both", trace = 1)

## Start:  AIC=-22.89
## esperanza_vida ~ habitantes + ingresos + analfabetismo + asesinatos + 
##     universitarios + heladas + area + densidad_pob
## 
##                  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_pob    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
## esperanza_vida ~ habitantes + ingresos + asesinatos + universitarios + 
##     heladas + area + densidad_pob
## 
##                  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_pob    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
## esperanza_vida ~ habitantes + ingresos + asesinatos + universitarios + 
##     heladas + densidad_pob
## 
##                  Df Sum of Sq    RSS     AIC
## - ingresos        1     0.132 22.648 -27.598
## - densidad_pob    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
## esperanza_vida ~ habitantes + asesinatos + universitarios + heladas + 
##     densidad_pob
## 
##                  Df Sum of Sq    RSS     AIC
## - densidad_pob    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
## esperanza_vida ~ habitantes + asesinatos + universitarios + heladas
## 
##                  Df Sum of Sq    RSS     AIC
## <none>                        23.308 -28.161
## + densidad_pob    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 = esperanza_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-03

Generando el mejor modelo basado en Akaike (AIK)

modelo <- (lm(formula = esperanza_vida ~ habitantes + asesinatos + universitarios + heladas, data = datos))
summary(modelo)

## 
## Call:
## lm(formula = esperanza_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-12

Interpretación del Modelo

\[Esperanza-de-vida = 7.103e01 + 5.014e05 * habitantes - 3.001e01 * aesinatos + 4.658e02 * universitarios -5.943e03 * heladas \]

En el modelo el valor de Multiple R-squared: 0.736 significa que las variables representan el 73.60% de la variabilida de los datos, es decir la lcntidad de habitantes, el índicador de asesinatos,, el indicador de universitaios y el indicador de heladas representan el 73.6% de la esperanza de vida en una ciudad.

El test F muestra que es significativo (p-value: 1.696e-12). Se satisfacen todas las condiciones para este tipo de regresión múltiple.

Predicciones

Cuál sería la esperanza de vida para tres ciudades que tienen los siguientes valores: \(habitantes = {500, 2000, 10000}\) \(asesnatos = {4, 12, 15}\) \(universitarios = {60,40,20}\) \(heladas = {150,100,50}\)

Los nuevos valores que participan en la predicción

nuevos.Valores <- data.frame(habitantes = c(500, 2000, 10000),
                             asesinatos = c(4, 12, 15),
                             universitarios = c(60,40,20),
                             heladas = c(150,100,50))
nuevos.Valores

##   habitantes asesinatos universitarios heladas
## 1        500          4             60     150
## 2       2000         12             40     100
## 3      10000         15             20      50

La ejecución de la predicción

prediccion <- predict(modelo, newdata = nuevos.Valores)
prediccion

##        1        2        3 
## 71.75504 68.79458 67.66078