library(paqueteMETODOS)
library(summarytools)
library(corrplot)
library(PerformanceAnalytics)
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 cientifico de datos en la construcción de un modelo que lo oriente sobre los precios de inmuebles.
data(vivienda4)
df = subset(vivienda4, vivienda4$tipo=="Apartamento" & vivienda4$estrato==4 & vivienda4$areaconst<=200)
# Antes de iniciar el analisis verificamos que el dataset contenga solo la informacion requerida (Apartamentos con menos de 200 m2 y estrato 4)
head(df,10)
## # A tibble: 10 × 5
## zona estrato preciom areaconst tipo
## <fct> <fct> <dbl> <dbl> <fct>
## 1 Zona Norte 4 220 52 Apartamento
## 2 Zona Norte 4 320 108 Apartamento
## 3 Zona Sur 4 290 96 Apartamento
## 4 Zona Norte 4 220 82 Apartamento
## 5 Zona Norte 4 220 75 Apartamento
## 6 Zona Norte 4 162 60 Apartamento
## 7 Zona Norte 4 225 84 Apartamento
## 8 Zona Norte 4 370 117 Apartamento
## 9 Zona Norte 4 155 60 Apartamento
## 10 Zona Norte 4 240 75 Apartamento
summary(df$estrato)
## 3 4 5 6
## 0 1363 0 0
summary(df$areaconst)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 40.00 60.00 70.00 75.48 84.00 200.00
summary(df$tipo)
## Apartamento Casa
## 1363 0
precio_zona = aggregate(preciom ~ zona, data=df, FUN = mean)
print(precio_zona)
## zona preciom
## 1 Zona Centro 236.0000
## 2 Zona Norte 204.5907
## 3 Zona Oeste 200.7692
## 4 Zona Oriente 262.5000
## 5 Zona Sur 201.7061
# Preciom
hist(df$preciom, freq = F, ylim = c(0,0.01))
lines(density(df$preciom))
shapiro.test(df$preciom) # Prueba de normalidad
##
## Shapiro-Wilk normality test
##
## data: df$preciom
## W = 0.89881, p-value < 2.2e-16
boxplot(df$preciom, main="Distribucion Precio m2", horizontal = T)
summarytools::descr(df$preciom) # Resumen estadisticos Precio
## Descriptive Statistics
## df$preciom
## N: 1363
##
## preciom
## ----------------- ---------
## Mean 202.44
## Std.Dev 65.29
## Min 78.00
## Q1 153.00
## Median 185.00
## Q3 240.00
## Max 645.00
## MAD 59.30
## IQR 86.50
## CV 0.32
## Skewness 1.44
## SE.Skewness 0.07
## Kurtosis 3.83
## N.Valid 1363.00
## Pct.Valid 100.00
Observando el histograma se puede notar una leve tendencia a la normalidad pero debe ser confirmado con las pruebas pertinentes, al realizar la prueba de normalidad tenemos un resultado de p-value < 2.2e-16 por lo que debemos rechazar la hipotesis nula para normalidad de la variable. Analizando la distribucion de quartiles se infiere una mayor concentracion de apartamentos hacia precios por encima de la media (el rango IQR es de solo 87 mientras que la distancia entre Q3 y el valor maximo es de 405), esto es confirmado al observar el boxplot
# Area Construida
hist(df$areaconst, freq = F, ylim = c(0,0.035))
lines(density(df$areaconst))
shapiro.test(df$areaconst) # Prueba de normalidad
##
## Shapiro-Wilk normality test
##
## data: df$areaconst
## W = 0.82347, p-value < 2.2e-16
boxplot(df$preciom, main="Distribucion Area Construida", horizontal = T)
summarytools::descr(df$areaconst) # Resumen estadisticos Area
## Descriptive Statistics
## df$areaconst
## N: 1363
##
## areaconst
## ----------------- -----------
## Mean 75.48
## Std.Dev 22.56
## Min 40.00
## Q1 60.00
## Median 70.00
## Q3 84.00
## Max 200.00
## MAD 14.83
## IQR 24.00
## CV 0.30
## Skewness 2.08
## SE.Skewness 0.07
## Kurtosis 6.32
## N.Valid 1363.00
## Pct.Valid 100.00
table(df$zona) # Distribucion Zona
##
## Zona Centro Zona Norte Zona Oeste Zona Oriente Zona Sur
## 7 237 52 2 1065
round(prop.table(table(df$zona)) * 100, 2) # Porcentajes Zona
##
## Zona Centro Zona Norte Zona Oeste Zona Oriente Zona Sur
## 0.51 17.39 3.82 0.15 78.14
Al igual que con el precio, para el area construida es dificil confirmar si tiene una distribucion normal o no solo con el histograma, por lo que debe ser sometida a pruebas de normalidad, en este caso obtenemos un resultado de p-value < 2.2e-16 por lo que rechazamos la hipotesis nula y la variable area construida no tiene normalidad. Con el area construida tambien hay una distribucion similar al precio, lo que sugiere una posible correlacion (el IQR es de 24, con una diferencia Max-Q3 = 116)
Ambas variables tienen un coeficiente de variacion similar (al rededor del 30%)
x = df$areaconst # Variable predictora
y = df$preciom # Variable respuesta
plot(x, y, main="Grafico de dispersion entre X y Y")
cor.test(x,y) # Correlacion de Pearson
##
## Pearson's product-moment correlation
##
## data: x and y
## t = 41.595, df = 1361, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.7237938 0.7706237
## sample estimates:
## cor
## 0.7481389
El resultado del coeficiente de relacion de Pearson confirma lo que se insinua en el grafico de dispersion, hay una correlacion fuerte (0.75)
ml = lm(y~x)
summary(ml)
##
## Call:
## lm(formula = y ~ x)
##
## 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 ***
## x 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
De acuerdo a este modelo se tiene \(Y=β0+β1X+ε\) donde:
confint(ml, level = 0.95)
## 2.5 % 97.5 %
## (Intercept) 31.00423 47.089340
## x 2.06264 2.266826
Como podemos ver aqui, el valor de β1 puede estar dentro del rango de 2.06 y 2.26 millones por cada metro cuadrado
En ninguno de los casos el coeficiente fue igual a 0, por lo que podemos confirmar que el precio y el area construida estan correlacionados.
cat("El R2 es igual a ", summary(ml)$r.squared)
## El R2 es igual a 0.5597117
El coeficiente de determinacion o \(R^2\), es una medida estadística que representa la proporción de la variabilidad de la variable dependiente que es explicada por las variables independientes en un modelo de regresión. Es una medida crucial para evaluar la bondad de ajuste de un modelo de regresión. en este case se tiene un \(R^2\) de 0.5597, lo que indica que la variable area explica un 55.97% de la variabilidad en el precio por metro cuadrado para los apartamentos estudiados
predict( ml, data.frame(x=110), interval="confidence", level = 0.95 )
## fit lwr upr
## 1 277.1674 272.9573 281.3775
predict( ml, data.frame(x=110), interval="prediction", level = 0.95 )
## fit lwr upr
## 1 277.1674 192.0449 362.2899
Como se puede observar en los resultados, el precio para un apartamento de 110 metros oscila entre 272.9 y 281.37 con un 95% de confianza, y el precio esperado con el modelo actual es de 277.1 millones de pesos por lo que una oferta de 200 millones es extremadamente atractiva.
El factor mas influyente a tener en cuenta es que como se vio en el analisis exploratorio, la zona tiene un gran impacto en el precio de los apartamentos, generando variaciones promedio de hasta 60 millones de pesos.
par(mfrow=c(2,2))
plot(ml)
Analizando los graficos se determina que los residuales tienen:
Linealidad: en el grafico de residuales vs valores ajustados se observa una dispersion homogenea en los valores, no siguen una funcion.
Normalidad: En el Q-Q plot se observa un ajuste a la linea de normalidad en la mayoria de los casos, pero para al ser una prueba subjetiva, se realiza prueba de normalidad sobre los residuales
shapiro.test(residuals(ml)) # Prueba de normalidad sobre residuos del modelo de regresion
##
## Shapiro-Wilk normality test
##
## data: residuals(ml)
## W = 0.96486, p-value < 2.2e-16
Se obtienen un p-value < 2.2e-16 por lo que se puede determinar que no hay normalidad, el supuesto no se cumple.
lmtest::bptest(ml)
##
## studentized Breusch-Pagan test
##
## data: ml
## BP = 292.99, df = 1, p-value < 2.2e-16
m1 = ml # Lin - Lin
m2 = lm(formula = y ~ log(x)) # Lin - Log
m3 = lm(formula = log(y) ~ x) # Log - Lin
m4 = lm(formula = log10(y) ~ log10(x)) # Log - Log
De ser necesario compare el ajuste y supuestos del modelo inicial y el transformado.
Estime varios modelos y compare los resultados obtenidos. En el mejor de los modelos, ¿se cumplen los supuestos sobre los errores?
# Se aplica Criterio de Informacion de Akaike para evaluar ajuste de cada modelo (menor es mejor)
aic_m1 = AIC(m1)
aic_m2 = AIC(m2)
aic_m3 = AIC(m3)
aic_m4 = AIC(m4)
AIC_values <- c(aic_m1, aic_m2, aic_m3, aic_m4)
nombres_modelos <- c("Lin - Lin", "Lin - Log", "Log - Lin", "Log - Log")
df_AIC <- data.frame(Modelo = nombres_modelos[order(AIC_values)], AIC = sort(AIC_values))
df_AIC
## Modelo AIC
## 1 Log - Log -2909.1990
## 2 Log - Lin -446.9303
## 3 Lin - Log 14059.7466
## 4 Lin - Lin 14146.4576
Por AIC, la transformacion Log-Log parece tener un mejor ajuste, se procede a evaluar los supuestos
lista_modelos <- list("Lin - Lin" = m1, "Lin - Log" = m2, "Log - Lin" = m3 , "Log - Log" = m4)
lista_residuos <- list()
#Se extraen los residuos para cada modelo
for (i in 1:length(lista_modelos)) {
modelo <- lista_modelos[[i]]
residuos <- residuals(modelo)
lista_residuos[[i]] <- residuos
}
# Realiza el análisis de residuos para cada modelo
for (i in 1:length(lista_residuos)) {
residuos <- lista_residuos[[i]]
# Gráfico de residuos vs. valores predichos
plot(predict(lista_modelos[[i]]), residuos, xlab = "Valores Predichos", ylab = "Residuos",
main = paste("Análisis de Residuos - Modelo", names(lista_modelos)[i]))
# Prueba de Breusch-Pagan para homocedasticidad
breusch_pagan_test <- lmtest::bptest(lista_modelos[[i]])
print(paste("Prueba de Breusch-Pagan - Modelo", names(lista_modelos)[i]))
print(breusch_pagan_test)
# Prueba de Shapiro-Wilk para normalidad
shapiro_test <- shapiro.test(residuos)
print(paste("Prueba de Shapiro-Wilk - Modelo", names(lista_modelos)[i]))
print(shapiro_test)
}
## [1] "Prueba de Breusch-Pagan - Modelo Lin - Lin"
##
## studentized Breusch-Pagan test
##
## data: lista_modelos[[i]]
## BP = 292.99, df = 1, p-value < 2.2e-16
##
## [1] "Prueba de Shapiro-Wilk - Modelo Lin - Lin"
##
## Shapiro-Wilk normality test
##
## data: residuos
## W = 0.96486, p-value < 2.2e-16
## [1] "Prueba de Breusch-Pagan - Modelo Lin - Log"
##
## studentized Breusch-Pagan test
##
## data: lista_modelos[[i]]
## BP = 214.66, df = 1, p-value < 2.2e-16
##
## [1] "Prueba de Shapiro-Wilk - Modelo Lin - Log"
##
## Shapiro-Wilk normality test
##
## data: residuos
## W = 0.95826, p-value < 2.2e-16
## [1] "Prueba de Breusch-Pagan - Modelo Log - Lin"
##
## studentized Breusch-Pagan test
##
## data: lista_modelos[[i]]
## BP = 150.38, df = 1, p-value < 2.2e-16
##
## [1] "Prueba de Shapiro-Wilk - Modelo Log - Lin"
##
## Shapiro-Wilk normality test
##
## data: residuos
## W = 0.99051, p-value = 1.035e-07
## [1] "Prueba de Breusch-Pagan - Modelo Log - Log"
##
## studentized Breusch-Pagan test
##
## data: lista_modelos[[i]]
## BP = 92.877, df = 1, p-value < 2.2e-16
##
## [1] "Prueba de Shapiro-Wilk - Modelo Log - Log"
##
## Shapiro-Wilk normality test
##
## data: residuos
## W = 0.98958, p-value = 2.857e-08
A pesar de que parece visualmente parece cumplirse la homoscedasticidad, la prueba de Breusch-Pagan arroja resultados muy inferiores a 0.05, lo mismo para Shapiro-wilk, por lo que ninguno de los suspuestos se cumple en ninguna transformacion
Como se puede observar en la comparacion de residuales y supuestos realizada en el punto inmediatamente anterior, ninguno de los modelos o transformaciones son validos para predecir valores de viviendas al no cumplir los supuestos de normalidad y homoscedasticidad.
Basado en los resultados de las pruebas realizadas a lo largo de este reporte, se recomienda utilizar metodos de analisis descriptivo tradicionales para estimar precios de viviendas, o utilizar modelos mas complejos que no dependan de normalidad de valores. El modelo utilizado en este ejercicio, y por ende sus transformaciones son dependientes en varios supuestos que simplemente no se cumplen en el dataset que se tiene disponible, por lo que los resultados que arroje, aunque a simple vista, cercanos, son arbitrarios y no deberian usarse como parte de un analisis o estrategia de mercado.