En este análisis de regresión lineal múltiple se busca hallar un modelo confiable para predecir el valor medio de las casas en Boston.
knitr::opts_chunk$set(echo = TRUE)
library (MASS)
library(PerformanceAnalytics)
## Loading required package: xts
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
## Attaching package: 'PerformanceAnalytics'
## The following object is masked from 'package:graphics':
##
## legend
library(corrplot)
## corrplot 0.95 loaded
library(nortest)
datos <- (Boston)
y <- datos$medv
x1 <- datos$crim
x2 <- datos$zn
x3 <- datos$indus
x4 <- datos$chas
x5 <- datos$nox
x6 <- datos$rm
x7 <- datos$age
x8 <- datos$dis
x9 <- datos$rad
x10 <- datos$tax
x11 <- datos$ptratio
x12 <- datos$black
x13 <- datos$lstat
lillie.test(y)$p.value
## [1] 1.671245e-30
lillie.test(x1)$p.value
## [1] 1.655823e-165
lillie.test(x2)$p.value
## [1] 1.214282e-261
lillie.test(x3)$p.value
## [1] 1.56259e-67
lillie.test(x5)$p.value
## [1] 1.006836e-14
lillie.test(x6)$p.value
## [1] 1.252791e-08
lillie.test(x7)$p.value
## [1] 4.928761e-30
lillie.test(x8)$p.value
## [1] 1.043784e-23
lillie.test(x9)$p.value
## [1] 6.250103e-145
lillie.test(x10)$p.value
## [1] 8.538468e-61
lillie.test(x11)$p.value
## [1] 1.097664e-43
lillie.test(x12)$p.value
## [1] 1.020562e-173
lillie.test(x13)$p.value
## [1] 1.639952e-10
Como se puede observar, ninguna de las variables siguen una distribución porque el p valor es menor a 0.05, la variable chas no fui incluida porque esta es categórica.
chart.Correlation(datos, histogram = TRUE, method = "p")
Esta gráfica será usada más adelante para tomar las variables
explicativas con mayor corelación a la valiablre medv.
mat_cor <- cor(datos[-14], method = "p")
corrplot(mat_cor, method = 'circle')
Aquí se puede observar como algunas de las varibales explicativas tienen
una moderada correlación entre sí, entonces se tratarán de evitar usar
estas en el mismo modelo.
m1 <- lm(y ~ x1 + x2 +x3 +x4 +x5 +x6 +x7 +x8 +x9 +x10 +x11 +x12 +x13)
summary(m1)
##
## Call:
## lm(formula = y ~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 +
## x10 + x11 + x12 + x13)
##
## Residuals:
## Min 1Q Median 3Q Max
## -15.595 -2.730 -0.518 1.777 26.199
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.646e+01 5.103e+00 7.144 3.28e-12 ***
## x1 -1.080e-01 3.286e-02 -3.287 0.001087 **
## x2 4.642e-02 1.373e-02 3.382 0.000778 ***
## x3 2.056e-02 6.150e-02 0.334 0.738288
## x4 2.687e+00 8.616e-01 3.118 0.001925 **
## x5 -1.777e+01 3.820e+00 -4.651 4.25e-06 ***
## x6 3.810e+00 4.179e-01 9.116 < 2e-16 ***
## x7 6.922e-04 1.321e-02 0.052 0.958229
## x8 -1.476e+00 1.995e-01 -7.398 6.01e-13 ***
## x9 3.060e-01 6.635e-02 4.613 5.07e-06 ***
## x10 -1.233e-02 3.760e-03 -3.280 0.001112 **
## x11 -9.527e-01 1.308e-01 -7.283 1.31e-12 ***
## x12 9.312e-03 2.686e-03 3.467 0.000573 ***
## x13 -5.248e-01 5.072e-02 -10.347 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.745 on 492 degrees of freedom
## Multiple R-squared: 0.7406, Adjusted R-squared: 0.7338
## F-statistic: 108.1 on 13 and 492 DF, p-value: < 2.2e-16
Si se usan todas las variables de la base de datos se alcanza un R ajustado de 0.73, ahora veamos si aliminando algunas variables con altos p valor, se consigue mejorar el modelo.
m2 <- lm(y ~ x1 + x2 + x3 + x4 + x5 + x6 + x8 + x9 + x10 + x11 + x12 + x13)
summary(m2)
##
## Call:
## lm(formula = y ~ x1 + x2 + x3 + x4 + x5 + x6 + x8 + x9 + x10 +
## x11 + x12 + x13)
##
## Residuals:
## Min 1Q Median 3Q Max
## -15.6054 -2.7313 -0.5188 1.7601 26.2243
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 36.436927 5.080119 7.172 2.72e-12 ***
## x1 -0.108006 0.032832 -3.290 0.001075 **
## x2 0.046334 0.013613 3.404 0.000719 ***
## x3 0.020562 0.061433 0.335 0.737989
## x4 2.689026 0.859598 3.128 0.001863 **
## x5 -17.713540 3.679308 -4.814 1.97e-06 ***
## x6 3.814394 0.408480 9.338 < 2e-16 ***
## x8 -1.478612 0.190611 -7.757 5.03e-14 ***
## x9 0.305786 0.066089 4.627 4.75e-06 ***
## x10 -0.012329 0.003755 -3.283 0.001099 **
## x11 -0.952211 0.130294 -7.308 1.10e-12 ***
## x12 0.009321 0.002678 3.481 0.000544 ***
## x13 -0.523852 0.047625 -10.999 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.74 on 493 degrees of freedom
## Multiple R-squared: 0.7406, Adjusted R-squared: 0.7343
## F-statistic: 117.3 on 12 and 493 DF, p-value: < 2.2e-16
m3 <- lm(y ~ x1 + x2 + x4 + x5 + x6 + x8 + x9 + x10 + x11 + x12 + x13)
summary(m3)
##
## Call:
## lm(formula = y ~ x1 + x2 + x4 + x5 + x6 + x8 + x9 + x10 + x11 +
## x12 + x13)
##
## Residuals:
## Min 1Q Median 3Q Max
## -15.5984 -2.7386 -0.5046 1.7273 26.2373
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 36.341145 5.067492 7.171 2.73e-12 ***
## x1 -0.108413 0.032779 -3.307 0.001010 **
## x2 0.045845 0.013523 3.390 0.000754 ***
## x4 2.718716 0.854240 3.183 0.001551 **
## x5 -17.376023 3.535243 -4.915 1.21e-06 ***
## x6 3.801579 0.406316 9.356 < 2e-16 ***
## x8 -1.492711 0.185731 -8.037 6.84e-15 ***
## x9 0.299608 0.063402 4.726 3.00e-06 ***
## x10 -0.011778 0.003372 -3.493 0.000521 ***
## x11 -0.946525 0.129066 -7.334 9.24e-13 ***
## x12 0.009291 0.002674 3.475 0.000557 ***
## x13 -0.522553 0.047424 -11.019 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.736 on 494 degrees of freedom
## Multiple R-squared: 0.7406, Adjusted R-squared: 0.7348
## F-statistic: 128.2 on 11 and 494 DF, p-value: < 2.2e-16
Después de eliminar a x7 (age) y a x3(indus) se alcanza un R ajustado de 0.7348, esto quiere decir que el modelo predice el 73.48% de lavariabilidad del precio medio de las casas a partir de las variables explicativas (crim, zn,chas,nox,rm, etc) y el 25.52% se explica a partir de otras variables. Ahora se comprobarán supuestos relacionados con los residuos del modelo.
plot(m3)
## Análisis de supuestos Como se pude ver en la gráfica Residuals vs
Fitted, los puntos tienden a tomar la forma de una párabola que abre
hacia arriba, esto incumple el supuesto de linealidad de los
residuos.
En el Q-Q Residuals, se puede observar como muchos valores en la cola derecha sobresalen de la línea teórica, para asegurarnos de si hay o no normalidad, haremos un lillie.test.
lillie.test(m3$residuals)$p.value
## [1] 2.938098e-21
Con un p valor menor a 0.05 se conlcuye que los residuos no siguen una distribución normal.
En la gráfica Scale-Location se puede observar un comportamiento similar a la de Residuals vs Fitted, con una forma de parábola, entonces se conyucle que hay heterocedasticidad en los residuos.
Y en la última gráfica, Residuals vs Leverange, se puede ver como hay observaciones que se alejan mucho del resto, esto podría indicuar puntos con un alto Leverage.
mean(m3$residuals)
## [1] 3.028395e-17
t.test(m3$residuals)$p.value
## [1] 1
Como se puede ver, con un p valor mayor a 0.05, se puede conluir que los residuos sí tienen media igual a 0.
plot(m1)
plot(m2)
Como se observan en estas gráficas, se pueden llegar a las mismas concluiones que en el modelo que habías seleccionado, entonces estos modelos tampoco cumplen con los supuestos requeridos.
No se pudo encontrar un modelo que fue confiable para el análisis de la variabilidad del precio medio de las casas en Boston, esto porque los modelos que se tuvieron en cuenta no cumplían con los supuestos exigidos por el mismo modelo.