Regresión lineal múltiple

setwd("~/ea9am")
library(pacman)
p_load("prettydoc", "DT", "xfun", "dplyr", "psych", "GGally", "ggplot2","readr", "gridExtra","leaflet","TSstudio", "dplyr", "corrplot", "lmtest","car")

La regresión lineal múltiple permite generar un modelo lineal en el que el valor de la variable dependiente o respuesta (Y) se determina a partir de un conjunto de variables independientes llamadas predictores X1X2X3.

Los modelos de regresión múltiple pueden emplearse para predecir el valor de la variable dependiente o para evaluar la influencia que tienen los predictores sobre ella (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} \]

  • $ _{0} $ : es la ordenada en el origen, el valor de la variable dependiente Y cuando todos los predictores son cero.

  • $ _{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.

Condiciones para la regresión lineal múltiple

Los modelos de correlación lineal múltiple requieren de las mismas condiciones que los modelos lineales simples más otras adicionales.

Parsimonia

Este término hace referencia a que el mejor modelo es aquel capaz de explicar con mayor precisión la variabilidad observada en la variable respuesta empleando el menor número de predictores, por lo tanto, con menos asunciones.

Tamaño de la muestra

No se trata de una condición de por sí pero, si no se dispone de suficientes observaciones, predictores que no son realmente influyentes podrían parecerlo. En el libro Hanbook of biological statistics recomiendan que el número de observaciones sea como mínimo entre 10 y 20 veces el número de predictores del modelo.

La gran mayoría de condiciones se verifican utilizando los residuos, por lo tanto, se suele generar primero el modelo y posteriormente validar las condiciones. De hecho, el ajuste de un modelo debe verse como un proceso iterativo en el que se ajusta el modelo, se evalúan sus residuos y se mejora. Así hasta llegar a un modelo óptimo.

Caso de estudio

La gestión de la contaminación del aire se ha vuelto importante en el último periodo debido a las siguientes afirmaciones: “A causa del desarrollo industrial y avance tecnológico se estima aproximadamente 1200 millones de personas están expuestas a niveles de dióxido de azufre (SO2), muy por encima de por directrices de la Organización Mundial de la Salud (OMS) y aproximadamente 1400 millones de personas están expuestas a niveles excesivos de humo y material articulado (PM)” (Rico, 2018).

La importancia de respirar aire limpio sin contaminantes es crucial para todo ser viviente. Por lo cual es importante que este sea aire limpio, pero ¿Cómo detectar si el aire que respiramos no tiene contaminantes? ¿Cómo nos afecta a nosotros?. La calidad del aire está directamente relacionada con la movilidad de carros, autobuses y hasta la propia movilidad de la ciudad, donde las personas realizan sus actividades diarias, por lo tanto, entre más movilidad haya en la ciudad, aumentan los contaminantes en el aire.

Los principales contaminantes en el aire han sido creados principalmente por la actividad económica humana, la concentración de estas sustancias es altamente nociva para la salud del ser humano. Este estudio tuvo en cuenta los siguientes contaminantes:

  • Material particulado \((PM10)\)
  • Ozono \((O3)\)
  • Dióxido de azufre \((SO2)\)

Contaminantes.

  1. PM10: son aquellas partículas sólidas o líquidas de polvo, cenizas, hollín, partículas metálicas, cemento o polen, dispersas en la atmósfera y cuyo diámetro varía entre 2,5 y 10 µm. La exposición repetitiva de este contaminante puede provocar defectos nocivos en el sistema respiratorio de la persona. a mayoría de éstas partículas precipitan en la tierra, provocando una capa de polvo en la superficie que puede afectar seriamente a la salud tanto de los organismos terrestres como los organismos acuáticos.

  2. O3 : es una molécula triatómica formada por átomos de oxígeno. El ozono es un componente natural de la atmósfera que se encuentra en bajas concentraciones y es vital para la vida. Al nivel de la superficie, el ozno en altas concentraciones es un contaminante del aire que provoca efectos nocivos en la salud humana, las plantas y los animales. Además, tiene una contribución al calentamiento global. A concentraciones elevadas de ozono se puede agravar el asma, debido a que el ozono favorece la susceptibilidad del organismo a los alérgenos (agentes que provocan los ataques).

Sin embargo, para este análisis nos enfocaremos en el siguiente contaminante \(SO2\):

  1. SO2: este es un gas incoloro, irritante, con un olor penetrante que se comienza a percibir on 0,3 a 1,4 ppm y es perfectamente distinguible a partir de 3 ppm -partes por millón. Este no es un gas inflamable, explosivo y tiene mucha estabilidad. La principal fuente de emisión de dióxido de azufre es la combustión de productos petrolíferos y la quema de carbón en centrales eléctricas y calefacciones centrales. Existen también algunas fuentes naturales, como es el caso de los volcanes.

¿Qué efectos tiene en la salud?

La contaminación del aire por SO2 puede provocar:

  • Dificultad para respirar.

  • Inflamación de las vías respiratorias.

  • irritación ocular por formación de ácido sulfuroso sobre las mucosas húmedas.

  • Alteraciones psíquicas.

  • Edema pulmonar.

  • Paro cardíaco.

  • Entre otras.

El dióxido de azufre también se ha asociado a problemas de asma y bronquitis crónica, aumentando la morbilidad y mortalidad en personas mayores y niños. os asmáticos y las personas con enfermedades pulmonares obstructivas crónicas (EPOC) y con problemas cardíacos son los más sensibles a los efectos del SO2.

Curiosamente, a pesar de ser un contaminante muy perjudicial, es también un aditivo alimenticio muy utilizado al tener características conservantes y antibacterianas. Se conoce como E220 y se emplea en vinos- los famosos sulfitos- cervezas, zumos, caramelos, yemas de huevo y otros productos con huevo.

En este caso en particular se analizará la relación de la variable de Dioxido de azufre con las demás variables de movilidad, para posteriormente realizar un modelo de regresion lineal múltiple con las variables necesarias.

Mapa de donde se obtuvieron los datos

content <- paste(sep = "<br/>",
  "<b><a href='https://www.ruoa.unam.mx/index.php?page=estaciones&id=6#datos'>ERNO</a></b>","Lng: -110.9706, Lat: 29.0814")


m <- leaflet() %>%
  addTiles() %>%  
  addMarkers(lng=-110.9706, lat= 29.0814, popup= content)

m

Datos

library(readxl)
datos <- read_excel("Concentracion.xlsx")
datatable(datos)

1.Analizar la relación entre variables

El primer paso a la hora de establecer un modelo lineal múltiple es estudiar la relación que existe entre variables. Esta información es crítica a la hora de identificar cuáles pueden ser los mejores predictores para el modelo, qué variables presentan relaciones de tipo no lineal (por lo que no pueden ser incluidas) y para identificar colinialidad entre predictores. A modo complementario, es recomendable representar la distribución de cada variable mediante histogramas.

round(cor(x = datos, method = "pearson"), 3)
##                            O3    SO2   PM10 Reactivacion_Comercial
## O3                      1.000 -0.220 -0.229                 -0.193
## SO2                    -0.220  1.000  0.612                  0.404
## PM10                   -0.229  0.612  1.000                  0.281
## Reactivacion_Comercial -0.193  0.404  0.281                  1.000
## Supermercado_Farmacia  -0.047  0.443  0.280                  0.875
## Parques_Centros        -0.182  0.019  0.065                  0.814
## Estaciones_Transito    -0.174  0.399  0.270                  0.944
## Lugares_Trabajo        -0.119  0.138  0.100                  0.581
## Residencia              0.130 -0.371 -0.238                 -0.842
##                        Supermercado_Farmacia Parques_Centros
## O3                                    -0.047          -0.182
## SO2                                    0.443           0.019
## PM10                                   0.280           0.065
## Reactivacion_Comercial                 0.875           0.814
## Supermercado_Farmacia                  1.000           0.655
## Parques_Centros                        0.655           1.000
## Estaciones_Transito                    0.872           0.770
## Lugares_Trabajo                        0.492           0.371
## Residencia                            -0.729          -0.617
##                        Estaciones_Transito Lugares_Trabajo Residencia
## O3                                  -0.174          -0.119      0.130
## SO2                                  0.399           0.138     -0.371
## PM10                                 0.270           0.100     -0.238
## Reactivacion_Comercial               0.944           0.581     -0.842
## Supermercado_Farmacia                0.872           0.492     -0.729
## Parques_Centros                      0.770           0.371     -0.617
## Estaciones_Transito                  1.000           0.474     -0.758
## Lugares_Trabajo                      0.474           1.000     -0.853
## Residencia                          -0.758          -0.853      1.000
  • Análisis con histogramas
multi.hist( x = datos, dcol = c("blue","red"), dlty = c("dotted", "solid"),
            main = "")

En el anterior gráfico, se observan 9 histogramas, cada uno, representa las variables antes vistas (03, SO2, PM10, Reactivación comercial, supermercados y farmacias, parques y centros, estaciones de transito, lugares de trabajo y residencia). De los histogramas podemos rescatar que la mayoría de las variables tienen una distribución normal, en especial: el O3 (ozono), PM10 (material particulado) y las residencias. Dicho esto, inferimos que nuestro modelo será bueno, ya que nuestros datos son “normales”, es decir, están alrededor de la media.

Ahora representando esto mismo, más la recta utilizando ggplot y ggally

ggpairs(datos, lower = list(continuous ="smooth"),
        diag = list (continuos = "barDiag"), axisLabels = "none")

Este gráfico nos permite conocer la correlación entre las variables continuas, así como los diagramas de disperción y la diagonal de densidad. Podemos observar como el SO2, se relaciona principalmente con el material particulado. Ambas variables suelen tener relación, ya que ambos, suelen compartir sus origenes, así como sus consecuencias, como por ejemplo, la exposición a estos contaminantes se relaciona con un deterioro de la salud respiratoria. Esto, explica de igual manera, la gran correlación de SO2 con supermercados y farmacias. Por último, se puede ver, como la correlación entre residencias y SO2 es negativa, ya que, al disminuir la estancia en nuestras casas, los niveles de SO2 aumentan, esto debido a los medios de transporte y nuestras actividades cotidianas.

Entonces, de los análisis realizados hasta el momento, podemos obtener las siguientes conclusiones preliminares:

  • Las variables que que tienen una mayor relación con el dioxido de azufre son: la reactivacion comercial (0.404), supermercados y farmacias (0.443).

  • La reactivación comercial claramente tiene un impacto debido a que aumenta la movilidad de automoviles y la producción en fábricas que generan este tipo de gas. Pues la combustion de carbono y petróleo generan \((SO2)\). Del mismo modo, los supermercados y farmacias son lugares concurridos en los que hay tambien mucho tráfico.

  • Por otro lado podemos ver que la variable residencia tiene una relación inversa con el \((SO2)\) pues cuando las personas se quedan en casa no generan este tipo de contaminante.

2.Generar el primer modelo de regresión lineal simple basado en la ecuación de la recta de mínimos cuadrados para una relación de X y Y

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 <- lm(SO2 ~ Supermercado_Farmacia, data=datos)
summary(modelo)
## 
## Call:
## lm(formula = SO2 ~ Supermercado_Farmacia, data = datos)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.77447 -0.16040 -0.00167  0.15632  0.98634 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)           0.285796   0.015300   18.68   <2e-16 ***
## Supermercado_Farmacia 0.011437   0.001098   10.42   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2544 on 445 degrees of freedom
## Multiple R-squared:  0.1961, Adjusted R-squared:  0.1943 
## F-statistic: 108.5 on 1 and 445 DF,  p-value: < 2.2e-16

Bajo este planteamiento el modelo tomaría la siguiente forma:

\[ y = 0.285796 - 0.011437x \]

Evaluando el modelo

y = 0.285796 - ( 0.011437*4)
y
## [1] 0.240048

Conocer el residuo según este modelo

y - (0.105416670)
## [1] 0.1346313
  • Aquí podemos notar que nos da un residuo algo grande para estos datos, lo que nos indica que este modelo no es del todo correcto.

Evaluando gráficamente el modelo

plot(datos$Supermercado_Farmacia, datos$SO2)
abline(modelo)

Esta gráfica intenta evaluar el modelo de manera gráfica, representando así, nuestro resultado esperado, y lo que realmente obtuvimos. En un caso ideal, cuando nuestra X esté en 10, el SO2 debería estar en -0.202083330. En el gráfico se observa que se encuentra debajo del 0.0, muy cerca del valor esperado, sin embargo el modelo sobreestima por 0.1346313 los valores.

Ahora, con múltiples predictores

modelo2 <- lm(SO2 ~ Supermercado_Farmacia + Reactivacion_Comercial + Residencia + Parques_Centros + Estaciones_Transito + Lugares_Trabajo, data=datos)
summary(modelo2) 
## 
## Call:
## lm(formula = SO2 ~ Supermercado_Farmacia + Reactivacion_Comercial + 
##     Residencia + Parques_Centros + Estaciones_Transito + Lugares_Trabajo, 
##     data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.7163 -0.1199 -0.0009  0.1168  0.7773 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             0.090469   0.027457   3.295  0.00106 ** 
## Supermercado_Farmacia   0.003807   0.001791   2.126  0.03405 *  
## Reactivacion_Comercial  0.015459   0.002329   6.638 9.36e-11 ***
## Residencia             -0.036921   0.004868  -7.584 2.00e-13 ***
## Parques_Centros        -0.020616   0.001226 -16.811  < 2e-16 ***
## Estaciones_Transito    -0.001840   0.001605  -1.146  0.25226    
## Lugares_Trabajo        -0.011906   0.001136 -10.482  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1879 on 440 degrees of freedom
## Multiple R-squared:  0.5661, Adjusted R-squared:  0.5602 
## F-statistic: 95.68 on 6 and 440 DF,  p-value: < 2.2e-16

El modelo con todas las variables introducidas como predictores tiene un R2 alta (0.5661), es capaz de explicar el 56.61% de la variabilidad observada en el SO2. El p-value del modelo es significativo (2.2e-16) 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.

Seleccionando los mejores predictores

En este caso 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).

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.

AIC maneja un trade-off entre la bondad de ajuste del modelo y la complejidad del modelo. Se basa en la entropía de información: se ofrece una estimación relativa de la información perdida cuando se utiliza un modelo determinado para representar el proceso que genera los datos.

AIC no proporciona una prueba de un modelo en el sentido de probar una hipótesis nula, es decir AIC no puede decir nada acerca de la calidad del modelo en un sentido absoluto. Si todos los modelos candidatos encajan mal, AIC no dará ningún aviso de ello.

step(object = modelo2, direction="both", trace = 1)
## Start:  AIC=-1487.46
## SO2 ~ Supermercado_Farmacia + Reactivacion_Comercial + Residencia + 
##     Parques_Centros + Estaciones_Transito + Lugares_Trabajo
## 
##                          Df Sum of Sq    RSS     AIC
## - Estaciones_Transito     1    0.0464 15.589 -1488.1
## <none>                                15.543 -1487.5
## - Supermercado_Farmacia   1    0.1597 15.703 -1484.9
## - Reactivacion_Comercial  1    1.5567 17.099 -1446.8
## - Residencia              1    2.0320 17.575 -1434.5
## - Lugares_Trabajo         1    3.8814 19.424 -1389.8
## - Parques_Centros         1    9.9827 25.526 -1267.7
## 
## Step:  AIC=-1488.12
## SO2 ~ Supermercado_Farmacia + Reactivacion_Comercial + Residencia + 
##     Parques_Centros + Lugares_Trabajo
## 
##                          Df Sum of Sq    RSS     AIC
## <none>                                15.589 -1488.1
## + Estaciones_Transito     1    0.0464 15.543 -1487.5
## - Supermercado_Farmacia   1    0.1245 15.714 -1486.6
## - Reactivacion_Comercial  1    1.9777 17.567 -1436.7
## - Residencia              1    2.0131 17.602 -1435.8
## - Lugares_Trabajo         1    3.8591 19.448 -1391.3
## - Parques_Centros         1    9.9995 25.589 -1268.6
## 
## Call:
## lm(formula = SO2 ~ Supermercado_Farmacia + Reactivacion_Comercial + 
##     Residencia + Parques_Centros + Lugares_Trabajo, data = datos)
## 
## Coefficients:
##            (Intercept)   Supermercado_Farmacia  Reactivacion_Comercial  
##               0.091438                0.003223                0.013836  
##             Residencia         Parques_Centros         Lugares_Trabajo  
##              -0.036727               -0.020632               -0.011663

Según la llamada final de este proceso (AIC) tenemos que el mejor modelo es:

modelo3 <- lm(formula = SO2 ~ Supermercado_Farmacia + Reactivacion_Comercial + 
    Residencia + Parques_Centros + Lugares_Trabajo, data = datos)
summary(modelo3)
## 
## Call:
## lm(formula = SO2 ~ Supermercado_Farmacia + Reactivacion_Comercial + 
##     Residencia + Parques_Centros + Lugares_Trabajo, data = datos)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.75212 -0.12115  0.00057  0.11759  0.75632 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             0.091438   0.027454   3.331 0.000939 ***
## Supermercado_Farmacia   0.003223   0.001717   1.877 0.061187 .  
## Reactivacion_Comercial  0.013836   0.001850   7.480 4.06e-13 ***
## Residencia             -0.036727   0.004867  -7.546 2.59e-13 ***
## Parques_Centros        -0.020632   0.001227 -16.819  < 2e-16 ***
## Lugares_Trabajo        -0.011663   0.001116 -10.448  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.188 on 441 degrees of freedom
## Multiple R-squared:  0.5648, Adjusted R-squared:  0.5599 
## F-statistic: 114.5 on 5 and 441 DF,  p-value: < 2.2e-16

Es recomendable mostrar el intervalo de confianza para cada uno de los coeficientes parciales de regresión:

confint(lm(formula = SO2 ~ Supermercado_Farmacia + Reactivacion_Comercial + 
    Residencia + Parques_Centros + Lugares_Trabajo, data = datos))
##                                2.5 %       97.5 %
## (Intercept)             0.0374815632  0.145394155
## Supermercado_Farmacia  -0.0001518597  0.006598654
## Reactivacion_Comercial  0.0102008433  0.017472070
## Residencia             -0.0462915208 -0.027161549
## Parques_Centros        -0.0230430192 -0.018221112
## Lugares_Trabajo        -0.0138566915 -0.009469078

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(Supermercado_Farmacia, modelo3$residuals)) +
    geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
    theme_bw()
plot2 <- ggplot(data = datos, aes(Residencia, modelo3$residuals)) +
    geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
    theme_bw()
plot3 <- ggplot(data = datos, aes(Reactivacion_Comercial, modelo3$residuals)) +
    geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
    theme_bw()
plot4 <- ggplot(data = datos, aes(Parques_Centros, modelo3$residuals)) +
    geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
    theme_bw()
plot5 <- ggplot(data = datos, aes(Lugares_Trabajo, modelo3$residuals)) +
    geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
    theme_bw()
grid.arrange(plot1, plot2,plot3,plot4,plot5)
## `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'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

Estas gráficas permiten validar mediante diagramas de dispersión entre la variable dependiente y cada uno de los predictores o con diagramas de dispersión entre cada uno de los predictores y los residuos del modelo. Se puede observar que, se cumple la linealidad para todos los predictores. Ya que la mayoría de ellos, se encuentran distribuidos sobre el eje de la x, o bien, muy cerca de él. Se puede destacar el caso de residencia, esto significaría que su relación es lineal.

Distribución normal de los residuos:

qqnorm(modelo3$residuals)
qqline(modelo3$residuals)

El gráfico Q-Q Plot, nos permite comprobar si dos variables provienen de la misma distribución, es decir, si nuestros datos son normales. La anterior gráfica, nos dice que nuestras variables (modelo3 y residuals) son datos normales, es decir, se encuentran alrededor de la media.

  • Prueba de shapiro wilk para determinar normalidad
shapiro.test(modelo3$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo3$residuals
## W = 0.99415, p-value = 0.084

Tanto el análisis gráfico como es test de hipótesis aprueban la normalidad.

  • 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.

ggplot(data = datos, aes(modelo3$fitted.values, modelo3$residuals)) +
geom_point() +
geom_smooth(color = "firebrick", se = FALSE) +
geom_hline(yintercept = 0) +
theme_bw()
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

El anterior gráfico representa los residuos frente a los valores ajustados por el modelo. Nuestra línea, comienza cerca de nuestro punto de origen, manteniendo su distribución la mayoría del tiempo sobre el eje de las X, además, no se observa ningún patrón específico. Esto quiere decir que, existe la homocedasticidad, es decir nuestros residuos de nuestro modelo son constantes para todos los valores.

Probando el modelo de regresión múltiple

bptest(modelo3)
## 
##  studentized Breusch-Pagan test
## 
## data:  modelo3
## BP = 9.1954, df = 5, p-value = 0.1015

Matriz de correlación entre las variables predictoras

corrplot(cor(dplyr::select(datos,Supermercado_Farmacia,Residencia,
                          Reactivacion_Comercial,Parques_Centros,Lugares_Trabajo)), 
        method = "number", tl.col = "black")

Aquí se puede notar que hay predictores que estan fuertemente relacionados, residencia tiene una grán correlación negativa con todos los demás predictotes. Así mismo, reactivación comercial influye positivamente en supermercado-farmacia y parques-centros.

vif(modelo3)
##  Supermercado_Farmacia Reactivacion_Comercial             Residencia 
##               4.479966              12.488933              10.974351 
##        Parques_Centros        Lugares_Trabajo 
##               3.270771               4.825273

Por concecuencia el factor de inflación de la varianza de reactivacion comercial y residencia es bastante alto. Esto podría deberse a que al estar en casa(residencia) no nos trasladamos a otros lugares como parques, supermercados. Por otro lado, la reactivación comercial aumenta la actividad esta otras zonas.

Autocorrelación

dwt(modelo, alternative = "two.sided")
##  lag Autocorrelation D-W Statistic p-value
##    1       0.9041928     0.1809021       0
##  Alternative hypothesis: rho != 0

Identificación de los posibles valores atípicos o influyentes

datos$studentized_residual <- rstudent(modelo3)
ggplot(data = datos, aes(x = predict(modelo3), y = abs(studentized_residual))) +
geom_hline(yintercept = 3, color = "grey", linetype = "dashed") +

geom_point(aes(color = ifelse(abs(studentized_residual) > 3, 'red', 'black'))) +
scale_color_identity() +
labs(title = "Distribución de los residuos studentized",
     x = "predicción modelo") + 
theme_bw() + theme(plot.title = element_text(hjust = 0.5))

Como lo marcan los puntos rojos podemos observar que si existen valores atípicos.

which(abs(datos$studentized_residual) > 3)
##  27  33 315 
##  27  33 315

Cuantificación de las influencias de los predictores

summary(influence.measures(modelo))
## Potentially influential observations of
##   lm(formula = SO2 ~ Supermercado_Farmacia, data = datos) :
## 
##     dfb.1_ dfb.Sp_F dffit   cov.r   cook.d hat    
## 9   -0.14  -0.08    -0.14    0.99_*  0.01   0.00  
## 23  -0.14  -0.08    -0.14    0.98_*  0.01   0.00  
## 24  -0.11  -0.04    -0.11    0.98_*  0.01   0.00  
## 25  -0.16  -0.11    -0.16    0.98_*  0.01   0.00  
## 26  -0.16  -0.11    -0.16    0.98_*  0.01   0.00  
## 27  -0.14  -0.07    -0.14    0.98_*  0.01   0.00  
## 28  -0.29  -0.29    -0.32_*  0.97_*  0.05   0.01  
## 29  -0.30  -0.30    -0.34_*  0.98_*  0.06   0.01  
## 30  -0.25  -0.24    -0.27_*  0.98_*  0.04   0.01  
## 31  -0.30  -0.31    -0.33_*  0.98_*  0.05   0.01_*
## 32  -0.21  -0.19    -0.22_*  0.99    0.02   0.01  
## 58   0.01   0.03    -0.04    1.01_*  0.00   0.01  
## 59   0.01   0.02    -0.03    1.02_*  0.00   0.02_*
## 60   0.01   0.04    -0.05    1.02_*  0.00   0.01  
## 64   0.00  -0.01     0.01    1.02_*  0.00   0.01  
## 65  -0.01  -0.04     0.04    1.02_*  0.00   0.02_*
## 66   0.00  -0.01     0.01    1.02_*  0.00   0.01  
## 67  -0.01  -0.02     0.02    1.02_*  0.00   0.01_*
## 68   0.00   0.00     0.00    1.02_*  0.00   0.01  
## 69   0.00   0.01    -0.01    1.01_*  0.00   0.01  
## 70   0.00   0.01    -0.01    1.01_*  0.00   0.01  
## 71   0.00   0.01    -0.01    1.01_*  0.00   0.01  
## 72   0.00  -0.01     0.01    1.02_*  0.00   0.02_*
## 73   0.00   0.02    -0.02    1.02_*  0.00   0.01  
## 74   0.00   0.00     0.00    1.01_*  0.00   0.01  
## 79   0.01   0.02    -0.02    1.02_*  0.00   0.01  
## 296  0.14   0.08     0.14    0.98_*  0.01   0.00  
## 312  0.06   0.07     0.07    1.02_*  0.00   0.01_*
## 313 -0.02  -0.02    -0.02    1.04_*  0.00   0.03_*
## 314  0.29   0.36     0.37_*  1.02_*  0.07   0.03_*
## 315 -0.14  -0.47     0.51_*  0.95_*  0.13   0.02_*
## 320 -0.06  -0.07    -0.07    1.02_*  0.00   0.02_*
## 321 -0.07  -0.08    -0.08    1.03_*  0.00   0.03_*
## 322 -0.09  -0.27     0.29_*  1.00    0.04   0.02_*
## 333  0.06  -0.06     0.13    0.98_*  0.01   0.00  
## 335  0.04  -0.10     0.17    0.97_*  0.01   0.00  
## 336  0.05  -0.07     0.13    0.98_*  0.01   0.00  
## 339  0.04  -0.07     0.13    0.98_*  0.01   0.00  
## 371  0.09   0.02     0.11    0.98_*  0.01   0.00  
## 411 -0.01  -0.01    -0.01    1.02_*  0.00   0.02_*
## 412  0.04   0.05     0.05    1.02_*  0.00   0.02_*
## 435 -0.01  -0.01    -0.01    1.01_*  0.00   0.01  
## 440 -0.06  -0.06    -0.07    1.02_*  0.00   0.01_*
## 441 -0.10  -0.12    -0.13    1.03_*  0.01   0.02_*
## 442 -0.02  -0.03    -0.03    1.02_*  0.00   0.01_*
## 444 -0.06  -0.06    -0.06    1.01_*  0.00   0.01  
## 446 -0.01  -0.01    -0.01    1.01_*  0.00   0.01  
## 447 -0.02  -0.02    -0.02    1.02_*  0.00   0.01

Descargas.

  • Código.
xfun::embed_file("A6U1.Rmd")

Download A6U1.Rmd

  • Datos.
xfun::embed_file("Concentracion.xlsx")

Download Concentracion.xlsx