Construcción de modelo de regresión lineal para análisis del mercado inmobiliario de la empresa A&C

Con base en los datos de ofertas de vivienda descargadas del portal Fincaraiz para apartamento de estrato 4 con área construida menor a 200 m2 (vivienda4.RDS) la inmobiliaria A&C requiere el apoyo de un científico de datos en la construcción de un modelo que lo oriente sobre los precios de inmuebles.

Con este propósito el equipo de asesores ha diseñado los siguientes pasos para obtener un modelo y así poder a futuro determinar los precios de los inmuebles a negociar

Primer punto

Realice un análisis exploratorio de las variables precio de vivienda (millones de pesos COP) y área de la vivienda (metros cuadrados) - incluir gráficos e indicadores apropiados interpretados.

Se cuenta con una base de datos que contiene 1706 observaciones y 5 variables relacionadas con el mercado inmobiliario de tipo cuantitativo (precio, area construida) y variables de tipo cualitativo (tipo de vivienda, zona de ubicación y estrato sociodemografico), la base no presenta valores perdidos por lo cual no se realiza ningun tipo de transformación de los datos.

En relacion al análisis exploratorio, se encuentra que predomina el tipo de vivienda apartamento 79.89% (n=1363).

data(vivienda4)
pie(table(vivienda4$tipo), col=c("blue", "grey"), main="Distribuciones por tipo de vivienda")

La totalidad de las viviendas se encuentran en el estrato 4, principalmente en la zona Sur 78.78% (n=1344) y Norte 16.88 (n=288) de la ciudad.

barplot(table(vivienda4$zona), col=c("orange","blue", "yellow", "red"), main="Distribucion de viviendas por zonas" )

En cuanto a las variables cuantitativas, se encuentra que el precio minimo de las viviendas es de 78 Millones, con una mediana de 210 Millones y un Valor máximo de 760 Millones.

Con relación al área construida, se encuentra que el valor minimo es de 40 metros cuadrados con una Mediana de 75 metros cuadrados y un valor maximo de 200 metros cuadrados.

summary(vivienda4)
##            zona      estrato     preciom        areaconst     
##  Zona Centro :   8   3:   0   Min.   : 78.0   Min.   : 40.00  
##  Zona Norte  : 288   4:1706   1st Qu.:160.0   1st Qu.: 60.00  
##  Zona Oeste  :  60   5:   0   Median :210.0   Median : 75.00  
##  Zona Oriente:   6   6:   0   Mean   :225.4   Mean   : 87.63  
##  Zona Sur    :1344            3rd Qu.:265.0   3rd Qu.: 98.00  
##                               Max.   :760.0   Max.   :200.00  
##           tipo     
##  Apartamento:1363  
##  Casa       : 343  
##                    
##                    
##                    
##

Al analizar la variable precio de manera general, se encuentra que presenta una media de 225 Millones, con una desviación estandar de +-85 millones. El 90% de las observaciones se encuentran por debajo de los 340 Millones de pesos.

vivienda4|> summarise(media_preciom = mean(preciom),
                varianza_preciom = var(preciom),
                desvi_preciom = sd(preciom),
                Q1 = quantile(preciom, probs=0.25),
                P90 = quantile(preciom, probs=0.90))
## # A tibble: 1 × 5
##   media_preciom varianza_preciom desvi_preciom    Q1   P90
##           <dbl>            <dbl>         <dbl> <dbl> <dbl>
## 1          225.            7376.          85.9   160   340

Por otro lado, la variable area construida presenta de manera general, una media de 87.6 metros cuadrados, con una desviación estandar de +-36 metros cuadrados. El 90% de las observaciones se encuentran por debajo de los 145 metros cuadrados.

vivienda4|> summarise(media_area = mean(areaconst),
                varianza_area = var(areaconst),
                desvi_area = sd(areaconst),
                Q1 = quantile(areaconst, probs=0.25),
                D4 = quantile(areaconst, probs=0.40),
                P90 = quantile(areaconst, probs=0.90))
## # A tibble: 1 × 6
##   media_area varianza_area desvi_area    Q1    D4   P90
##        <dbl>         <dbl>      <dbl> <dbl> <dbl> <dbl>
## 1       87.6         1321.       36.3    60    70  144.

Segundo punto

Realice un análisis exploratorio bivariado de datos, enfocado en la relación entre la variable respuesta (precio) en función de la variable predictora (área construida) - incluir gráficos e indicadores apropiados interpretados.

Al analizar la distribución del precio por el tipo de vivienda, se evidencia una asimetría hacia la izquierda dada principalmente por el precio de las viviendas tipo apartamento que es predominante; de igual manera, se observa un sesgo hacia la derecha de la gráfica dado por precios atípicos (mayores a 600 Millones).

ggplot(vivienda4, aes(x = preciom, fill = tipo)) + 
  geom_histogram (bins=30) +  labs (title =  "Distribucion de los precios (millones COP) segun tipo de vivienda")

Al comparar el area construida con el tipo de vivienda se encuentra que las casas presentan mayor área construida y que en los apartamentos predominan areas menores a 80 metros cuadrados.

ggplot(vivienda4, aes(x = areaconst, fill = tipo)) + 
  geom_histogram (bins=30)  + labs (title =  "Distribucion del area construida (metros cuadrados) segun tipo de vivienda")

En comparación de los precios por tipo de vivienda, se encuentra que los apartamentos son mas economicos que las casas, tienen un precio promedio de 202 Millones con una Desviación estándar de 65 millones y mediana de 185 millones, el 90% de los datos se ubica por debajo de los 285 millones. Por otro lado, las viviendas tipo casa presentan un precio promedio de 316 Millones con una Desviación estándar de 97 millones y mediana de 300 millones.

aptos=subset(vivienda4, vivienda4$tipo=="Apartamento")
casas=subset(vivienda4, vivienda4$tipo=="Casa")

aptos|> summarise(media_preciom = mean(preciom),
                  mediana_preciom = median(preciom),
                varianza_preciom = var(preciom),
                desvi_preciom = sd(preciom),
                Q1 = quantile(preciom, probs=0.25),
                P90 = quantile(preciom, probs=0.90))
## # A tibble: 1 × 6
##   media_preciom mediana_preciom varianza_preciom desvi_preciom    Q1   P90
##           <dbl>           <dbl>            <dbl>         <dbl> <dbl> <dbl>
## 1          202.             185            4263.          65.3  154.   285
casas|> summarise(media_preciom = mean(preciom),
                  mediana_preciom = median (preciom),
                varianza_preciom = var(preciom),
                desvi_preciom = sd(preciom),
                Q1 = quantile(preciom, probs=0.25),
                P90 = quantile(preciom, probs=0.90))
## # A tibble: 1 × 6
##   media_preciom mediana_preciom varianza_preciom desvi_preciom    Q1   P90
##           <dbl>           <dbl>            <dbl>         <dbl> <dbl> <dbl>
## 1          317.             300            9368.          96.8   255   440

Se observa una mayor variabilidad de los datos, en las casas, sin embargo en los apartamentos se encuentran mayor numero de observaciones atípicas por encima del tercer cuartil.

boxplot (vivienda4$preciom~vivienda4$tipo, main="Precio en millones (COP) segun tipo de vivienda")

Al analizar la dispersión de los precios según el area construida, se observa una una relacion entre el área construida y el precio de la vivienda la cual aparentemente es de tipo lineal, es decir que a mayor área mayor es el precio de la vivienda. Se observa una mayor concentración de apartamentos con área menor a 100 metros cuadrados con precio entre 100 y 300 Millones.

ggplot2::ggplot(aptos, aes(x=areaconst, y=preciom)) + geom_point(colour=6) + labs (title = "Distribucion de los precios de apartamentos (millones COP) segun area construida (m2)")

En relación a las casas, se observa la misma tendencia lineal entre precio y área contruida, no obstante al ser las casas el 20% de los datos se evidencia mas dispersión en los puntos. Se observa un mayor numero de datos atipicos respecto a los apartamentos.La mayor concentración se presenta en casas con area entre 100 y 200 metros cuadrados.

ggplot2::ggplot (casas, aes(x=areaconst, y=preciom)) + geom_point(colour=4) + labs (title = "Distribucion de los precios de casas (millones COP) segun area construida (m2)")

Para confirmar esta relación que se observa en los diagramas de disperión se procede a aplicar pruebas de correlación entre estas variables, encontrando una correlacion positiva debil de 0.763.

ggpairs(vivienda4[,3:4], title="correlacion") 

cor.test(x = vivienda4$areaconst, y = vivienda4$preciom, method = "pearson", digits = 3)
## 
##  Pearson's product-moment correlation
## 
## data:  vivienda4$areaconst and vivienda4$preciom
## t = 48.728, df = 1704, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.7424432 0.7821521
## sample estimates:
##       cor 
## 0.7630166

Tercer punto

Estime el modelo de regresión lineal simple entre precio=f(area)+ε. Interprete los coeficientes del modelo β0, β1 en caso de ser correcto.

Dado que el 80% de la poblacion se concentra en la vivienda tipo apartamento, se procede a realizar el modelo de regresión lineal, teniendo en cuenta unicamente los apartamentos.

modeloap=lm(preciom ~ areaconst, data=aptos)
summary(modeloap)
## 
## Call:
## lm(formula = preciom ~ areaconst, data = aptos)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -225.404  -23.902   -4.754   25.763  209.021 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 39.04679    4.09977   9.524   <2e-16 ***
## areaconst    2.16473    0.05204  41.595   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 43.34 on 1361 degrees of freedom
## Multiple R-squared:  0.5597, Adjusted R-squared:  0.5594 
## F-statistic:  1730 on 1 and 1361 DF,  p-value: < 2.2e-16

Con el modelo de regresión lineal aplicado a los apartamentos, se observa que β0 (intercepto) es de 39 millones, es decir que ese sería el precio sin área construida, para efectos reales corresponde al precio del terreno donde se construirá el apartamento. Con relacion a β1 el cual es de 2.16 corresponde a la pendiente, es decir que por cada metro cuadrado que se aumente en el apartamento el precio aumentaría en 2.16 Millones.

Cuarto punto

Construir un intervalo de confianza (95%) para el coeficiente β1, interpretar y concluir si el coeficiente es igual a cero o no. Compare este resultado con una prueba de hipótesis t.

confint(modeloap, level = 0.95)
##                2.5 %    97.5 %
## (Intercept) 31.00423 47.089340
## areaconst    2.06264  2.266826

El intervalo de confianza para la pendiente β1 en el modeloap se encuentra entre 2.06 y 2.26; es decir, por cada metro cuadrado adicional, el precio promedio del apartamento aumenta entre 2.06 y 2.26 millones con una confianza del 95%.

En las pruebas t realizadas en el modeloap se tiene como hipótesis nula que los coeficientes β0 y β1 son iguales a cero. Al obtener valores inferiores a 0.05 tanto para β0 como β1, decimos que los coeficientes son significativos y por tanto se rechaza hipotesis nula y se acepta la hipótesis alterna, donde ambos valores (estimadores) son diferentes de cero.

Quinto punto

Calcule e interprete el indicador de bondad R2.

Con el modeloap se encuentra que el indicador de ajuste R cuadrado es de 0.559, esto indica que el precio de los apartamentos es explicado por el área construida en el 56%, es decir que existen otras variables (covariables) que influyen en el precio y que no han sido tenidas en cuenta para este modelo.

Sexto punto

¿Cuál sería el precio promedio estimado para un apartamento de 110 metros cuadrados? Considera entonces con este resultado que un apartamento en la misma zona con 110 metros cuadrados en un precio de 200 millones sería una atractiva esta oferta? ¿Qué consideraciones adicionales se deben tener?.

nuevos_datos <- data.frame(areaconst= 110)
predict (modeloap, nuevos_datos, interval = "confidence", level = 0.95)
##        fit      lwr      upr
## 1 277.1674 272.9573 281.3775

Acorde con el modelo, el precio promedio estimado para los apartamentos de 110 metros cuadrados sería de 277 Millones de pesos, con un limite inferior 273 Millones y limite superior de 281 Millones (IC 95%). Por lo anterior, se determina que el precio de 200 Millones por un apartamento de 110 metros cuadrados en la misma zona es una oferta muy atractiva para cualquier comprador.

Septimo punto

Realice la validación de los supuestos del modelo por medio de gráficos apropiados, interpretarlos y sugerir posibles soluciones si se violan algunos de ellos. Utilice las pruebas de hipótesis para la validación de supuestos y compare los resultados con lo observado en los gráficos asociados.

NORMALIDAD

Dado que los datos son mayores a 50, se debe aplicar el test de kolmogorov-Smirnov para verificar la normalidad. No obstante, dado que no se conoce la media y la varianza, se debe utilizar la modificación test Lilliefors el cual asume que la media y varianza son desconocidas, permitiendo contrastar la normalidad.

#normalidad de los residuos
lillie.test(modeloap$residuals)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  modeloap$residuals
## D = 0.060955, p-value = 6.975e-13

Acorde con los resultados del Test Lilliefors, se encuentra que la hipotesis nula (los datos de los residuos tienen una distribución normal) se rechaza y se acepta la hipotesis alterna, los datos de los residuos no siguen una distribucion normal.

HOMOCEDASTICIDAD

En la grafica se encuentra que la varianza de los residuos del modeloap alrededor de la linea de regresión se distribuye de manera heterogenea se dispersan mas a medida que aumenta el valor de datos ajustados, esto es un signo revelador de que existe heterocedasticidad. Tambien se observa como la varianza aumenta a medida que lo hacen los valores ajustados, decimos así que el modelo no cumple el supuesto de homocedasticidad.

#Homocedasticidad
res.stud<-studres(modeloap)
mod.fit<-modeloap$fitted.values

par(mfrow=c(1,2))
plot(mod.fit,res.stud,ylab="Residuos estudentizados",xlab="Valores ajustados")
abline(h=0,lty=2)
lines(lowess(res.stud~mod.fit), col = 2)
plot(mod.fit,abs(res.stud),ylab="|Residuos estudentizados|",xlab="Valores ajustados")
lines(lowess(abs(res.stud)~mod.fit), col = 2)

Para confirmar este supuesto se aplica el Test de Breush-Pagan, en el cual la hipotesis nula contempla que la homocedasticidad está presente (los residuos se distribuyen con la misma varianza), encontrando que se rechaza hipotesis nula y se acepta la hipotesis alterna acorde con valor de p. Es decir, que está presente la Heteroscedasticidad (los residuos no se distribuyen con la misma varianza alrededor de la linea de regresion)

# Test de Breush-Pagan (homocedasticidad de los residuos)
library(lmtest)
bptest(modeloap)
## 
##  studentized Breusch-Pagan test
## 
## data:  modeloap
## BP = 292.99, df = 1, p-value < 2.2e-16

LINEALIDAD

Se evidencia que la distribución de los residuos no siguen una relacion lineal.La grafica Q-Q de los residuos muestra que no se alinean perfectamente alrededor de la linea intercuartílica y que esta situacion se presenta especialmente hacia los extremos. En la distancia Cook, se encuentra que existe gran diferencia entre los estimadores B0 y B1 respecto a cada observación.

#Linealidad
par(mfrow=c(1,2))
plot(modeloap)

Al analizar mas detenidamente las observaciones y la distancia Cook respecto a los estimadores B0 y B1, se encuentran varios valores atípicos (outliers) los cuales pueden influir en la estimación del modeloap y por lo tanto, requiere ser ajustado añadiendo otras covariables o dando algun tipo de tratamiento a estos datos atípicos, como se muestra en la siguiente gráfica:

ols_plot_cooksd_bar(modeloap)

NO AUTOCORRELACION DE ERRORES

Los residuos deben ser independientes entre sí, para probarlo se usa el Test de Durbin-Watson, en el cual la hipotesis nula plantea que No existe correlación entre los residuos y la alterna que si hay autocorrelación. En este caso, dado que el resultado obtenido es menor de 2, indica que los errores de los residuos presentan correlación serial positiva y es significante (p valor menor a 0.05). Para corregirlo se debe considerar agregar rezagos de la variable dependiente y/o independiente al modelo.

lmtest::dwtest(modeloap)
## 
##  Durbin-Watson test
## 
## data:  modeloap
## DW = 1.443, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0

OUTLIERS

El análisis de los residuos detecta observaciones atípicas que pueden sesgar a los estimadores los coeficientes del modelo, las cuales se encuentran por encima de 150 y -150. Es decir que se encuentran varios apartamentos con precios y areas atípicos, representados por la diferencia de los residuos y la estimación del modelo, como se muestra en la gráfica:

ggplot(data = aptos, aes(x = seq_along(modeloap$residuals), 
                         y = modeloap$residuals)) +
geom_point(aes(color = modeloap$residuals)) +
scale_color_gradient2(low = "blue3", mid = "grey", high = "red") +
geom_line(linewidth = 0.3) + 
labs(title = "Distribución de los residuos Modeloap", x = "observaciones aptos", y = "residuo")+ 
geom_hline(yintercept = 0) + 
theme(plot.title = element_text(hjust = 0.5), legend.position = "none")

El análisis de los residuos estudentizados detecta varias observaciones atípicas respecto a la linea de regresión, con apartamentos por encima de los 400 metros cuadrados con área construida mayor a 175 metros cuadrados. Así como apartamentos de 100 metros cuadrados con precio menor a cien millones de pesos.

ggplot(data = aptos, mapping = aes(x = areaconst, y = preciom)) +
geom_point(color = "grey50", size = 2) +
#recta de regresión con todas las observaciones
geom_smooth(method = "lm", se = FALSE, color = "black") +
geom_point(data = aptos[,3:4], color = "red", size = 2) +
geom_smooth(data = aptos[,3:4], method="lm", se =FALSE,color = "blue") +
labs(x = "area construida (M2)", y = "precio (Millones)") +
theme_bw() + theme(plot.title = element_text(hjust = 0.5), legend.position = "none")
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'

Octavo punto

De ser necesario realice una transformación apropiada para mejorar el ajuste y supuestos del modelo.

La transformación de la variable respuesta, en este caso precio en Millones es util para solucionar problemas de validación de supuestos, que como ya vimos no tienen criterios de validez, para de esta manera mejorar el nivel de ajuste del modelo final. Para ello se procede a encontrar el valor de lambda optimo. El valor maximo de este valor puede ayudar a orientar el tipo de transformación requerida en el modeloap.

boxcox(lm(aptos$preciom ~ aptos$areaconst, data=aptos), lambda = -3:3)

bc<-boxcox(lm(aptos$preciom ~ aptos$areaconst), lambda = -1:1)

lambda <- bc$x[which.max(bc$y)]
lambda
## [1] -0.1313131

El valor de lambda de -0.1313 indica que la variable dependiente (precio) debe estar expresada en escala logarítmica. La línea vertical punteada central representa el parámetro estimado λ mientras que las otras dos líneas representan su intervalo de confianza al 95%. Como el gráfico muestra que el 0 está dentro del intervalo de confiaza del λ óptimo y la estimación está realmente cerca del 0, en este caso la mejor opción es aplicar la transformación logarítmica:

nuevo_y = ((aptos$preciom ^ lambda) - 1 )/ lambda
modBox=lm(nuevo_y ~ areaconst, data=aptos)
summary(modBox)
## 
## Call:
## lm(formula = nuevo_y ~ areaconst, data = aptos)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.48740 -0.06527 -0.00535  0.07957  0.32740 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 3.4459708  0.0097192  354.55   <2e-16 ***
## areaconst   0.0046628  0.0001234   37.79   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1027 on 1361 degrees of freedom
## Multiple R-squared:  0.5121, Adjusted R-squared:  0.5117 
## F-statistic:  1428 on 1 and 1361 DF,  p-value: < 2.2e-16

Con este modelo transformado ModBox no se aumenta el valor de R cuadrado respecto al modelo inicial, de hecho disminuye a 51% y este resultado obtenido es significativamente estadístico (p valor menor 0.05).

Noveno punto

De ser necesario compare el ajuste y supuestos del modelo inicial y el transformado.

Se procede a validar supuestos del modelo transformado ModBox para determinar su validez.

NORMALIDAD

Acorde con los resultados del Test de Lilliefors, se encuentra que la hipotesis nula (los datos de los residuos tienen una distribución normal) se rechaza y se acepta la hipotesis alterna, los residuos no siguen una distribucion normal, igual que en el modeloap; es decir, que este supuesto no cambia con la transformación realizada.

lillie.test(modBox$residuals)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  modBox$residuals
## D = 0.03952, p-value = 3.526e-05

HOMOCEDASTICIDAD

En la grafica se encuentra que la varianza de los residuos del modelo inicial modeloap alrededor de la linea de regresión se distribuye de manera heterogenea. En el modelo transformado modBox se observa que los residuos se dispersan menos a valores más altos en la gráfica, mejorando el resultado respecto al modeloap inicial.Tambien se observa como la varianza se mantiene mas constante y aumenta a medida que lo hacen los valores ajustados.

res.stud.box<-studres(modBox)
mod.fit.box<-modBox$fitted.values

par(mfrow=c(1,2))
plot(mod.fit.box,res.stud.box,ylab="Residuos estudentizados",xlab="Valores ajustados")
abline(h=0,lty=2)
lines(lowess(res.stud.box~mod.fit.box), col = 2)
plot(mod.fit.box,abs(res.stud.box),ylab="|Residuos estudentizados|",xlab="Valores ajustados")
lines(lowess(abs(res.stud.box)~mod.fit.box), col = 2)

bptest(modBox)
## 
##  studentized Breusch-Pagan test
## 
## data:  modBox
## BP = 130.48, df = 1, p-value < 2.2e-16

Para confirmar este supuesto se aplica el Test de Breush-Pagan, en el cual la hipotesis nula contempla que la homocedasticidad está presente (los residuos se distribuyen con la misma varianza), encontrando que se rechaza hipotesis nula y se acepta la hipotesis alterna acorde con valor de p. Es decir, que está presente la Heteroscedasticidad (los residuos no se distribuyen con la misma varianza alrededor de la linea de regresion). Respecto al modelo inicial (modeloap) el valor de BP es mucho menor pasando de 292.99 a 130.48, no obstante los p valor en ambos casos es menor a 0.05 confirmando la Heteroscedasticidad en ambos modelos, por ende, este supuesto no cambia con la transformación realizada.

LINEALIDAD

se evidencia que los residuos del modelo transformado ModBox no siguen una relacion lineal aunque mejora ligeramente respecto al modelo inicial modeloap, segun la grafica Q-Q. En la distancia Cook, se encuentra que mejora la diferencia entre los B0 y B1 respecto a cada observación en el modelo ajustado ModBox respecto al modelo inicial modeloap, no obstante, persisten la presencia de valores atípicos (779, 376) afectando a los estimadores del modelo transformado.

par(mfrow=c(1,2))
plot(modBox)

En la distancia Cook, se encuentra que mejora la diferencia entre los B0 y B1 respecto a cada observación en el modelo ajustado ModBox respecto al modelo inicial modeloap, no obstante, persisten la presencia de valores atípicos afectando a los estimadores del modelo transformado y por lo tanto, a pesar de la transformación se requiere añadir otras covariables o dar algun tipo de tratamiento a estos datos atípicos, como se muestra en la siguiente gráfica:

ols_plot_cooksd_bar(modBox)

NO AUTOCORRELACION DE ERRORES

Los residuos del modelo transformado ModBox deben ser independientes entre sí, para probarlo se usa el Test de Durbin-Watson,dado que el resultado obtenido es menor de 2, indica que los residuos presentan correlación serial positiva y es significante (p valor menor a 0.05) resultado también obtenido para el modelo inicial modeloap. No mejora este supuesto con la transformación del modelo.

lmtest::dwtest(modBox)
## 
##  Durbin-Watson test
## 
## data:  modBox
## DW = 1.3101, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0

OUTLIERS

El análisis de los residuos del modelo transformado ModBox mejora respecto al modelo inicial modeloap, como se muestra en la gráfica:

ggplot(data = aptos, aes(x = seq_along(modBox$residuals), 
                         y = modBox$residuals)) +
geom_point(aes(color = modBox$residuals)) +
scale_color_gradient2(low = "blue3", mid = "grey", high = "red") +
geom_line(linewidth = 0.3) + 
labs(title = "Distribución de los residuos modBox", x = "observaciones aptos", y = "residuo")+ 
geom_hline(yintercept = 0) + 
theme(plot.title = element_text(hjust = 0.5), legend.position = "none")

Decimo punto

Estime varios modelos y compare los resultados obtenidos. En el mejor de los modelos, ¿se cumplen los supuestos sobre los errores?

modelo1=lm (preciom ~ areaconst, data=aptos)          # Lin - Lin
modelo2=lm(preciom ~ log(areaconst), data=aptos)      # Lin - Log
modelo3=lm(log(preciom) ~ areaconst, data=aptos)      # Log - Lin
modelo4=lm(log(preciom) ~ log(areaconst), data=aptos) # Log - Log
stargazer(modelo1, modelo2, modelo3, modelo4, type="text", df=FALSE)
## 
## =======================================================================
##                                     Dependent variable:                
##                     ---------------------------------------------------
##                              preciom                log(preciom)       
##                         (1)          (2)          (3)          (4)     
## -----------------------------------------------------------------------
## areaconst             2.165***                  0.009***               
##                       (0.052)                   (0.0002)               
##                                                                        
## log(areaconst)                    195.419***                 0.882***  
##                                    (4.445)                   (0.020)   
##                                                                        
## Constant             39.047***   -635.532***    4.551***     1.484***  
##                       (4.100)      (19.092)     (0.019)      (0.087)   
##                                                                        
## -----------------------------------------------------------------------
## Observations           1,363        1,363        1,363        1,363    
## R2                     0.560        0.587        0.520        0.582    
## Adjusted R2            0.559        0.587        0.519        0.582    
## Residual Std. Error    43.339       41.982       0.205        0.191    
## F Statistic         1,730.157*** 1,933.199*** 1,473.424*** 1,894.288***
## =======================================================================
## Note:                                       *p<0.1; **p<0.05; ***p<0.01

Al realizar transformaciones de las variables independiente (areaconst) y dependiente (Preciom) para el tipo de vivienda apartamentos y comparalos entre sí, se encuentra que el mejor modelo es el 4 donde se utiliza logaritmo para la variable dependiente (preciom) tal como se expresó al encontrar el valor de lambda de -0.1313 con la transformación BOxCoX. En el modelo 4, se aumenta ligeramente el valor de R2 que explica que el area construida influye en el precio de apartamentos en un 58%, es decir un 3% mas que en el modelo 1. Por otro lado, es el que presenta un error estandar de los residuos mas cercano de cero. En los cuatro modelos se cuenta con pruebas F con significancia estadística.

Se procede a probar los supuestos con cada uno de los modelos propuestos:

Modelo 1

NORMALIDAD: Los residuos del modelo1 no provienen de una distribución normal (p valor menor a 0.05).

lillie.test(modelo1$residuals)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  modelo1$residuals
## D = 0.060955, p-value = 6.975e-13

HOMOCEDASTICIDAD: En la grafica se encuentra que las varianzas de los residuos del modelo 1, se dispersan bastante alrededor de la linea.Al aplicar el Test de Breush-Pagan, se rechaza la hipotesis nula la cual contempla que la homocedasticidad está presente (los residuos se distribuyen con la misma varianza), dado que el valor de p es significativo (menor a 0.05)

#Homocedasticidad
res.stud<-studres(modelo1)
mod.fit<-modelo1$fitted.values

par(mfrow=c(1,2))
plot(mod.fit,res.stud,ylab="Residuos estudentizados Modelo 1",xlab="Valores ajustados Modelo 1")
abline(h=0,lty=2)
lines(lowess(res.stud~mod.fit), col = 2)
plot(mod.fit,abs(res.stud),ylab="|Residuos estudentizados Modelo 1|",xlab="Valores ajustados Modelo 1")
lines(lowess(abs(res.stud)~mod.fit), col = 2)

bptest(modelo1)
## 
##  studentized Breusch-Pagan test
## 
## data:  modelo1
## BP = 292.99, df = 1, p-value < 2.2e-16

LINEALIDAD: se evidencia que los residuos del modelo1 no siguen una relacion lineal..La grafica Q-Q demuestra que los residuos en los extremos no se aliean alrededorde la linea intercuartílica. En la distancia Cook, se encuentra que existe gran diferencia entre los B0 y B1, respecto a cada observación.

par(mfrow=c(1,2))
plot(modelo1)

NO AUTOCORRELACION DE ERRORES: el Test de Durbin-Watson demuestra que los errores de los residuos no son independientes entre sí, dado que el resultado obtenido es menor de 2, indica que los residuos presentan correlación serial positiva y es significante (p valor menor a 0.05).

lmtest::dwtest(modelo1)
## 
##  Durbin-Watson test
## 
## data:  modelo1
## DW = 1.443, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0

OUTLIERS: El análisis de los residuos detecta observaciones atípicas que pueden sesgar a los estimadores los coeficientes del modelo 1 (150 y -150).

ggplot(data = aptos, aes(x = seq_along(modelo1$residuals), 
                         y = modelo1$residuals)) +
geom_point(aes(color = modelo1$residuals)) +
scale_color_gradient2(low = "blue3", mid = "grey", high = "red") +
geom_line(linewidth = 0.3) + 
labs(title = "Distribución de los residuos Modelo1", x = "observaciones aptos", y = "residuo")+ 
geom_hline(yintercept = 0) + 
theme(plot.title = element_text(hjust = 0.5), legend.position = "none")

Modelo 2

NORMALIDAD: Los residuos en el modelo 2 no provienen de una distribución normal acorde con el valor de p (menor a 0.05) del test Lilliefors.

lillie.test(modelo2$residuals)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  modelo2$residuals
## D = 0.060407, p-value = 1.22e-12

HOMOCEDASTICIDAD: En la grafica se muestra que la varianza de los residuos del modelo 2 tiene una distribución mas uniforme que el modelo1.Al aplicar el Test de Breush-Pagan, se rechaza hipotesis nula la cual contempla que los residuos se distribuyen con la misma varianza(homocedasticidad) dado que el valor de p significativo.

#Homocedasticidad
res.stud<-studres(modelo2)
mod.fit<-modelo2$fitted.values

par(mfrow=c(1,2))
plot(mod.fit,res.stud,ylab="Residuos estudentizados Mod 2",xlab="Valores ajustados Mod 2")
abline(h=0,lty=2)
lines(lowess(res.stud~mod.fit), col = 2)
plot(mod.fit,abs(res.stud),ylab="|Residuos estudentizados Mod 2|",xlab="Valores ajustados Mod 2")
lines(lowess(abs(res.stud)~mod.fit), col = 2)

bptest(modelo2)
## 
##  studentized Breusch-Pagan test
## 
## data:  modelo2
## BP = 214.66, df = 1, p-value < 2.2e-16

LINEALIDAD: se evidencia que los residuos en el Modelo 2 no siguen una relacion lineal, acrode con la grafica Q-Q. En la distancia Cook, se identifican puntos potencialmente influyentes los cuales pueden generar sesgos en los estimadores del modelo2.

par(mfrow=c(1,2))
plot(modelo2)

NO AUTOCORRELACION DE ERRORES: Los errores de los residuos del modelo 2 deben ser independientes entre sí, acorde con el resultado del Test de Durbin-Watson, menor de 2, indica que los residuos presentan correlación serial positiva y es significante (p valor menor a 0.05).

lmtest::dwtest(modelo2)
## 
##  Durbin-Watson test
## 
## data:  modelo2
## DW = 1.4775, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0

OUTLIERS: El análisis de los residuos detecta observaciones atípicas que pueden sesgar a los estimadores los coeficientes del modelo2, los cuales se observan por encima de 100 y -100.

ggplot(data = aptos, aes(x = seq_along(modelo2$residuals), 
                         y = modelo2$residuals)) +
geom_point(aes(color = modelo2$residuals)) +
scale_color_gradient2(low = "blue3", mid = "grey", high = "red") +
geom_line(linewidth = 0.3) + 
labs(title = "Distribución de los residuos Modelo2", x = "observaciones aptos", y = "residuo")+ 
geom_hline(yintercept = 0) + 
theme(plot.title = element_text(hjust = 0.5), legend.position = "none")

Modelo 3

NORMALIDAD: Los residuos en el modelo 3 no provienen de una distribución Normal (se rechaza hipotesis nula acorde con valores de p menor a 0.05) acorde con el Test Lilliefors.

lillie.test(modelo3$residuals)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  modelo3$residuals
## D = 0.033644, p-value = 0.001022

HOMOCEDASTICIDAD: En la grafica se observa que la varianza de los residuos del modelo 3, se distribuyen de manera dispersa. Al aplicar el Test de Breush-Pagan, se obtiene un p valor que rechaza la hipotesis nula (homocedasticidad está presente en los residuos), y confirmando la Heteroscedasticidad en el modelo3.

#Homocedasticidad
res.stud<-studres(modelo3)
mod.fit<-modelo3$fitted.values

par(mfrow=c(1,2))
plot(mod.fit,res.stud,ylab="Residuos estudentizados Mod3",xlab="Valores ajustados Mod3")
abline(h=0,lty=2)
lines(lowess(res.stud~mod.fit), col = 2)
plot(mod.fit,abs(res.stud),ylab="|Residuos estudentizados Mod3|",xlab="Valores ajustados Mod3")
lines(lowess(abs(res.stud)~mod.fit), col = 2)

bptest(modelo3)
## 
##  studentized Breusch-Pagan test
## 
## data:  modelo3
## BP = 150.38, df = 1, p-value < 2.2e-16

LINEALIDAD: se evidencia que los residuos mejoran la distribución sobre la linea intercuartílica de la grafica Q-Q. En la distancia Cook, se encuentra que persisten valores atípicos los cuales influyen en los estimadores de este modelo.

par(mfrow=c(1,2))
plot(modelo3)

NO AUTOCORRELACION DE ERRORES: Los errores de los residuos del modelo 3 presentan correlación serial positiva y es significante (p valor menor a 0.05), acorde con el resultado del Test de Durbin-Watson, el cual es menor de 2, por lo anterior, se rechaza la hipotesis nula (No existe correlación entre los residuos).

lmtest::dwtest(modelo3)
## 
##  Durbin-Watson test
## 
## data:  modelo3
## DW = 1.3187, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0

OUTLIERS: El análisis de los residuos en el modelo3 mejoran la distribución de las observaciones atípicas respecto a los dos modelos anteriores.

ggplot(data = aptos, aes(x = seq_along(modelo3$residuals), 
                         y = modelo3$residuals)) +
geom_point(aes(color = modelo3$residuals)) +
scale_color_gradient2(low = "blue3", mid = "grey", high = "red") +
geom_line(linewidth = 0.3) + 
labs(title = "Distribución de los residuos Mod3", x = "observaciones aptos", y = "residuo")+ 
geom_hline(yintercept = 0) + 
theme(plot.title = element_text(hjust = 0.5), legend.position = "none")

Modelo 4

NORMALIDAD: Los residuos en modelo 4, no provienen de una distribución Normal (valor de p menor a 0.05) acorde con el test de Lilliefors.

lillie.test(modelo4$residuals)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  modelo4$residuals
## D = 0.035893, p-value = 0.0003038

HOMOCEDASTICIDAD: En la grafica se encuentra que la varianza de los residuos del modelo 4, se dispersan menos a valores más altos en la gráfica respecto al modelo ModBox y mejorando el resultado respecto al modeloap inicial.Tambien se observa como la varianza se mantiene mas constante aumenta a medida que lo hacen los valores ajustados.Al aplicar el Test de Breush-Pagan, en el cual la hipotesis nula contempla que la homocedasticidad está presente (los residuos se distribuyen con la misma varianza), ésta se rechaza y se acepta la hipótesis alterna, donde el valor de BP es mucho menor para el modelo 4 (92.87), que el modelo ModBox (130) y mucho menor al modelo inicial modeloap (292.99). No obstante los p valor en los tres casos es menor a 0.05 confirmando la Heteroscedasticidad en los tres modelos mencionados.

#Homocedasticidad
res.stud<-studres(modelo4)
mod.fit<-modelo4$fitted.values

par(mfrow=c(1,2))
plot(mod.fit,res.stud,ylab="Residuos estudentizados",xlab="Valores ajustados")
abline(h=0,lty=2)
lines(lowess(res.stud~mod.fit), col = 2)
plot(mod.fit,abs(res.stud),ylab="|Residuos estudentizados|",xlab="Valores ajustados")
lines(lowess(abs(res.stud)~mod.fit), col = 2)

bptest(modelo4)
## 
##  studentized Breusch-Pagan test
## 
## data:  modelo4
## BP = 92.877, df = 1, p-value < 2.2e-16

LINEALIDAD: se evidencia que los parámetros no siguen una relacion lineal hacia los extremos de la línea intercuartílica; no obstante, en la grafica Q-Q presenta mejor alineación de los residuos alrededor de la linea intercuartílica respecto a los modelos anteriores. En la distancia Cook, se encuentra que continuan presentandose valores influyentes.

par(mfrow=c(1,2))
plot(modelo4)

NO AUTOCORRELACION DE ERRORES: Los residuos deben ser independientes entre sí, para probarlo se usa el Test de Durbin-Watson, en el cual la hipotesis nula plantea que No existe correlación entre los residuos y la alterna que si hay autocorrelación. En este caso, dado que el resultado obtenido es menor de 2, indica que los residuos presentan correlación serial positiva y es significante (p valor menor a 0.05), resultados similares obtenido en los modelos anteriores.

lmtest::dwtest(modelo4)
## 
##  Durbin-Watson test
## 
## data:  modelo4
## DW = 1.3214, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0

OUTLIERS: El análisis de los residuos detecta observaciones atípicas que pueden sesgar a los estimadores los coeficientes del modelo 4, los cuales se encuentran por encima de 0.5 y -0.5, sin embargo presenta menor numero de outliers que en los modelos inicial (modeloap), transformado (ModBox), modelo1, y modelo2.

ggplot(data = aptos, aes(x = seq_along(modelo4$residuals), 
                         y = modelo4$residuals)) +
geom_point(aes(color = modelo4$residuals)) +
scale_color_gradient2(low = "blue3", mid = "grey", high = "red") +
geom_line(linewidth = 0.3) + 
labs(title = "Distribución de los residuos Modelo4", x = "observaciones aptos", y = "residuo")+ 
geom_hline(yintercept = 0) + 
theme(plot.title = element_text(hjust = 0.5), legend.position = "none")

Undécimo punto

Con los resultados obtenidos construya un informe para los directivos de la inmobiliaria, indicando el modelo apropiado y sus principales características. A este informe se deben añadir los anexos como evidencia de la realización de los pasos anteriores.

A nivel predictivo los modelos propuestos no serían un buenos modelos a ser utilizados, dado que no se cumple los requisitos para la mayoría de supuestos. A continuación se realiza un comparativo de los supuestos para los 4 modelos propuestos:

#Linealidad comparativo
par(mfrow=c(1,4))
plot(modelo1, 1, caption = "Modelo 1")
plot(modelo2, 1, caption = "Modelo 2")
plot(modelo3, 1, caption = "Modelo 3")
plot(modelo4, 1, caption = "Modelo 4")

Lo que debemos esperar de este gráfico es que no aparezca un patrón lineal entre los residuales y los valores pronosticados del precio de apartamentos (Y). En tal sentido el modelo 2 y 4, se comportan mejor en este supuesto que los modelos 1 y 3.

#Normalidad comparativo
par(mfrow=c(2,2))
plot(modelo1, 2, caption = "Modelo 1")
plot(modelo2, 2, caption = "Modelo 2")
plot(modelo3, 2, caption = "Modelo 3")
plot(modelo4, 2, caption = "Modelo 4")

plot(density(modelo1$residuals))
plot(density(modelo2$residuals))
plot(density(modelo3$residuals))
plot(density(modelo4$residuals))

Para los primeros cuatro gráficos Q-Q plot, en el eje vertical se muestran los residuos estandarizados y se espera que los puntos se acercen lo más posible a la línea oblicua del gráfico que representa una distribución normal. En tal sentido, el modelo 3 presenta una mejor normalidad. Los modelos 1 y 2 tienen una distribución de residuos asimétrica respecto a la línea intercuartílica.

Como puede verse en los gráficos de densidad, la distribución de los residuales del los modelos 1, 3 y 4 tienen cierta asimetría positiva, mientras que en el modelo 2 es mas simetrico.

#Homocedasticidad comparativo
par(mfrow=c(1,4))
plot(modelo1, 3, caption = "Modelo 1")
plot(modelo2, 3, caption = "Modelo 2")
plot(modelo3, 3, caption = "Modelo 3")
plot(modelo4, 3, caption = "Modelo 4")

Se espera que la dispersión de la raiz cuadrada de los residuos estandarizados no presente muchas variaciones según los valores esperados de Y (precio). En este caso, los cuatro modelos propuestos presentan problemas de heterocedasticidad.

# Observaciones Influyentes: “Leverage” (palanca) y Valores Atípicos
par(mfrow=c(1,4))
plot(modelo1, 5, caption = "Modelo 1")
plot(modelo2, 5, caption = "Modelo 2")
plot(modelo3, 5, caption = "Modelo 3")
plot(modelo4, 5, caption = "Modelo 4")

Al analizar las observaciones y la distancia Cook respecto a los estimadores de los cuatro modelos, se encuentra que en todos ellos existen valores atípicos (outliers) que pueden influir en la estimación de los modelos; en el modelo 2 mejoran las distancias cook de algunos de ellos, como se muestra en las siguientes gráficas:

par(mfrow=c(1,4))
ols_plot_cooksd_bar(modelo1)

ols_plot_cooksd_bar(modelo2)

ols_plot_cooksd_bar(modelo3)

ols_plot_cooksd_bar(modelo4)

Los residuos estandarizados nos indican qué tan extremos son los valores o casos en la variable dependiente precio de los apartamentos (Y), por lo general residuos que exceden más de 3 desviaciones estándar indican casos muy atípicos.

Por otro lado el indicador de “leverage”, nos indica qué tanta “palanca” pueden tener los valores de los casos en la variable independientes (area construida) en la recta de regresion. La influencia o “leverage” mide la proporción de la suma total de cuadrados de la variable dependiente precio de los apartamentos (Y) a la cual contibuye el valor en la variable independiente area construida (X). Casos que sean muy extremos en Y y mucha palanca en X son los casos influyentes y pueden causar problemas al modelo. Estos casos se ubican en los extremos inferior y superior derecho del gráfico como en los modelos 1 y 2. Es decir que en este caso se comportan mejor los modelos 3 y 4 para este supuesto.

Según los resultados obtenidos, se evidencia que el área construida pese a tener una correlacion demostrada (debil positiva Coeficiente Pearson) con el precio de los apartamentos, no es la unica variable explicativa (R cuadrado 51 modBox a 56% modeloap). Existen otro tipo de covariables, que influyen en el precio final de los apartamentos las cuales no han sido contempladas en los modelos, debido a que no fueron suministradas en la base de datos entregada por la compañía. Es importante tener en cuenta, que el precio puede verse afectado por el equipamiento de los apartamentos (numero de habitaciones, numero de baños, numero de parqueaderos, acabados, piso, estrato, entre otros). Seria valioso contar con esta información para generar nuevos modelos explicativos del precio en función de estas variables y complementar el análisis realizado.

Dentro de la base de datos se encuentran datos atípicos en el precio y area construida de los apartamentos, los cuales generan sesgo en los estimadores de los coeficientes para los modelos de regresión lineal (B0 y B1). Por lo anterior, es importante considerar no tenerlos en cuenta para generar proximos modelos predictivos o generar modelos por rangos de precios y areas para optimizar los modelos de regresión, Ej: generar un modelo para areas construidas por debajo de los 150 metros cuadrados y los 400 Millones de pesos (concentración de las observaciones).

Se realizaron transformaciones del modelo inicial (ModeloAp) con el fin de mejorar supuestos: normalidad, linealidad, homocedasticidad, autocorrelación de errores, outliers; sin embargo, solamente se obtuvo mejoria en la linealidad y los outliers, los demas supuestos no se cumplen. Por lo anterior, se sugiere aplicar otro tipo de modelo que no exija normalidad de las observaciones.

No se recomienda excluir observaciones influyentes ya que puede afectar de manera sustancial la estimación de los coeficientes de regresión en el modelo. Sería importante para mejorar el modelo aumentar el tamaño del dataset, ya que entre mayor numero de observaciones analizadas menor será la posibilidad de que una de ellas tenga demasiada influencia.