#install.packages("printr", type = "source", repos = c("http://yihui.name/xran"))
library(printr)
#install.packages("devtools") # solo la primera vez
#devtools::install_github("dgonxalex80/paqueteMODELOS", force =TRUE)
library(paqueteMODELOS)
data("vivienda")
head(vivienda, 5)
| id | zona | piso | estrato | preciom | areaconst | parqueaderos | banios | habitaciones | tipo | barrio | longitud | latitud |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1147 | Zona Oriente | NA | 3 | 250 | 70 | 1 | 3 | 6 | Casa | 20 de julio | -76.51168 | 3.43382 |
| 1169 | Zona Oriente | NA | 3 | 320 | 120 | 1 | 2 | 3 | Casa | 20 de julio | -76.51237 | 3.43369 |
| 1350 | Zona Oriente | NA | 3 | 350 | 220 | 2 | 2 | 4 | Casa | 20 de julio | -76.51537 | 3.43566 |
| 5992 | Zona Sur | 02 | 4 | 400 | 280 | 3 | 5 | 3 | Casa | 3 de julio | -76.54000 | 3.43500 |
| 1212 | Zona Norte | 01 | 5 | 260 | 90 | 1 | 2 | 3 | Apartamento | acopi | -76.51350 | 3.45891 |
Ahora se describe cada una de las columnas de la tabla.
Realice un filtro a la base de datos e incluya solo las ofertas de : base1: casas, de la zona norte de la ciudad. Presente los primeros 3 registros de las bases y algunas tablas que comprueben la consulta.
# Filtrar data frame
base1 <- subset(vivienda, tipo == "Casa" & zona == "Zona Norte")
#Mostrar los primeros 3 registros de la base
head(base1 , 3)
| id | zona | piso | estrato | preciom | areaconst | parqueaderos | banios | habitaciones | tipo | barrio | longitud | latitud |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1209 | Zona Norte | 02 | 5 | 320 | 150 | 2 | 4 | 6 | Casa | acopi | -76.51341 | 3.47968 |
| 1592 | Zona Norte | 02 | 5 | 780 | 380 | 2 | 3 | 3 | Casa | acopi | -76.51674 | 3.48721 |
| 4057 | Zona Norte | 02 | 6 | 750 | 445 | NA | 7 | 6 | Casa | acopi | -76.52950 | 3.38527 |
El resumen de la base de datos:
#install.packages("skimr")
skimr::skim(base1)
| Name | base1 |
| Number of rows | 722 |
| Number of columns | 13 |
| _______________________ | |
| Column type frequency: | |
| character | 4 |
| numeric | 9 |
| ________________________ | |
| Group variables | None |
Variable type: character
| skim_variable | n_missing | complete_rate | min | max | empty | n_unique | whitespace |
|---|---|---|---|---|---|---|---|
| zona | 0 | 1.00 | 10 | 10 | 0 | 1 | 0 |
| piso | 372 | 0.48 | 2 | 2 | 0 | 5 | 0 |
| tipo | 0 | 1.00 | 4 | 4 | 0 | 1 | 0 |
| barrio | 0 | 1.00 | 4 | 25 | 0 | 103 | 0 |
Variable type: numeric
| skim_variable | n_missing | complete_rate | mean | sd | p0 | p25 | p50 | p75 | p100 | hist |
|---|---|---|---|---|---|---|---|---|---|---|
| id | 0 | 1.0 | 2574.65 | 1986.26 | 58.00 | 766.25 | 2257.00 | 4225.00 | 8319.00 | ▇▃▃▂▁ |
| estrato | 0 | 1.0 | 4.20 | 0.98 | 3.00 | 3.00 | 4.00 | 5.00 | 6.00 | ▇▅▁▇▂ |
| preciom | 0 | 1.0 | 445.91 | 268.36 | 89.00 | 261.25 | 390.00 | 550.00 | 1940.00 | ▇▃▁▁▁ |
| areaconst | 0 | 1.0 | 264.85 | 167.17 | 30.00 | 140.00 | 240.00 | 336.75 | 1440.00 | ▇▃▁▁▁ |
| parqueaderos | 287 | 0.6 | 2.18 | 1.40 | 1.00 | 1.00 | 2.00 | 3.00 | 10.00 | ▇▂▁▁▁ |
| banios | 0 | 1.0 | 3.56 | 1.52 | 0.00 | 2.00 | 3.00 | 4.00 | 10.00 | ▅▇▃▁▁ |
| habitaciones | 0 | 1.0 | 4.51 | 1.83 | 0.00 | 3.00 | 4.00 | 5.00 | 10.00 | ▁▇▅▂▁ |
| longitud | 0 | 1.0 | -76.52 | 0.02 | -76.59 | -76.53 | -76.52 | -76.50 | -76.47 | ▁▁▇▆▂ |
| latitud | 0 | 1.0 | 3.46 | 0.03 | 3.33 | 3.45 | 3.47 | 3.48 | 3.50 | ▁▁▁▂▇ |
Teniendo en cuenta el resumen y el head anterior se evidencia que en la variable zona y en la variable tipo existe solo un número unico (n_unique), lo que se concluye que el filtro que se pretendía quedo correcto.
Se realiza un grafico de mapa para validar si todas las viviendas según la latitud y longitud si se encuentran en la zona norte de la ciudad de Cali como lo especifica el data set.
#install.packages("leaflet")
library(leaflet)
coordenadas <- subset(base1,
select = c( latitud,longitud)
)
#crear mapa
map <- leaflet(coordenadas) %>%
addTiles() %>%
addCircleMarkers(lng = ~longitud, lat = ~latitud, radius = 5, color = "mediumseagreen")%>%
addProviderTiles("Stamen.Toner") %>%
addMiniMap(position = "bottomright")
map # Muestra el mapa
En el anterior mapa se muestra cómo se distribuyen los puntos elegidos en las bases de datos. Sin embargo, al revisar el mapa, se descubrió que ciertos lugares no pertenecen a la región del Norte, lo que plantea una dificultad importante. Esta novedad podría ser el resultado de errores en la clasificación inicial de las propiedades en la base de datos o problemas de geolocalización potenciales. Este hallazgo enfatiza la importancia de limpiar y verificar minuciosamente los datos antes de realizar análisis y tomar decisiones basadas en ellos. Para garantizar la precisión de futuros análisis, es fundamental considerar la necesidad de revisar y corregir la información proporcionada sobre la ubicación geográfica.
Realice un análisis exploratorio de datos enfocado en la correlación entre la variable respuesta (precio de la casa) en función del área construida, estrato, numero de baños, numero de habitaciones y zona donde se ubica la vivienda. Use gráficos interactivos con el paquete plotly e interprete los resultados.
Adicional, se observa que para todas las variables se encuentran datos vacíos pero para las variables piso y parquea hay una cantidad considerable de variables nulas.
#install.packages("heatmaply")
apply(X= is.na(base1),
MARGIN = 2,
FUN = sum)
## id zona piso estrato preciom areaconst
## 0 0 372 0 0 0
## parqueaderos banios habitaciones tipo barrio longitud
## 287 0 0 0 0 0
## latitud
## 0
apply(X= is.na(base1),
MARGIN = 2,
FUN = mean)
## id zona piso estrato preciom areaconst
## 0.0000000 0.0000000 0.5152355 0.0000000 0.0000000 0.0000000
## parqueaderos banios habitaciones tipo barrio longitud
## 0.3975069 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
## latitud
## 0.0000000
A partir de lo anteriormente expuesto, se evidencia que un 51.52% de los datos correspondientes a la variable “piso” se encuentran en estado de falta (Null), mientras que un 39.75% de los datos asociados a la variable “parquea” presentan la misma condición de ausencia.
Se sugiere que el responsable de la recolección de datos valide y describa la naturaleza de los valores faltantes (NA), así como proponga posibles métodos de imputación.
Para realizar un analisis detallado de las variables se crea un base de datos que contenga solo las variables a analizar.
variables_num <- subset(base1,
select = c( preciom,
areaconst,
banios,
habitaciones)
)
head(variables_num, 3)
| preciom | areaconst | banios | habitaciones |
|---|---|---|---|
| 320 | 150 | 4 | 6 |
| 780 | 380 | 3 | 3 |
| 750 | 445 | 7 | 6 |
knitr::kable(summary (variables_num) )
| preciom | areaconst | banios | habitaciones | |
|---|---|---|---|---|
| Min. : 89.0 | Min. : 30.0 | Min. : 0.000 | Min. : 0.000 | |
| 1st Qu.: 261.2 | 1st Qu.: 140.0 | 1st Qu.: 2.000 | 1st Qu.: 3.000 | |
| Median : 390.0 | Median : 240.0 | Median : 3.000 | Median : 4.000 | |
| Mean : 445.9 | Mean : 264.9 | Mean : 3.555 | Mean : 4.507 | |
| 3rd Qu.: 550.0 | 3rd Qu.: 336.8 | 3rd Qu.: 4.000 | 3rd Qu.: 5.000 | |
| Max. :1940.0 | Max. :1440.0 | Max. :10.000 | Max. :10.000 |
library(GGally)
# coeficiente de correlación de Spearman
matriz_corr <- cor(variables_num,
method = "spearman")
# Visualizar la matriz de correlación
ggpairs(variables_num,
title = " ",
upper = list(continuous = wrap("cor",
method = "spearman",
correlation = matriz_corr)))
preciom - areaconst: Se nota una correlación fuerte positiva de 0.813, es decir, que a medida que el areaconstruida aumenta el precio de la vivienda tiende aumentar.
precio - banios: Tiene una correlación de 0.617 positiva, esto indica que a medida que aumenta el numero de baños tiende aumentar el precio de la vivienda
precio - habitaciones: Tiene una correlación positiva no tan fuerte de 0.398, esto indica que las viviendas con precio alto tienden a tener más habitaciones.
#Creando diagrama de cajas Preciom x Estrato
boxplot(base1$preciom~base1$estrato,
col= c("#FFFFE0", "#E6E6FA", "#87CEEB","#98FB98"),
main= "Precio vivienda en millones por Estrato",
ylab= "Millones",
ylim= c(0,2000),
xlab = "",)
Se evidencia que a medida que el estrato aumenta, la media de los precios de las viviendas en la zona norte de la ciudad de Cali tiende a subir. Este patrón sugiere una relación positiva entre el estrato socioeconómico y los precios de las propiedades, indicando que, en promedio, las viviendas en estratos más altos tienden a tener valores superiores en la zona Norte.
Estime un modelo de regresión lineal múltiple con las variables del punto anterior (precio = f(área construida, estrato, número de cuartos, número de parqueaderos, número de baños ) ) e interprete los coeficientes si son estadísticamente significativos. Las interpretaciones deber están contextualizadas y discutir si los resultados son lógicos. Adicionalmente interprete el coeficiente R2 y discuta el ajuste del modelo e implicaciones (que podrían hacer para mejorarlo).
Dado que solo se utilizara para el modelo las variables de areacont, estrato, habitaciones, parqueaderos, banios y preciom se crea un base de datos con solo estos datos, para evitar errores de datos faltantes.
La variable estrato se encuentra como númerica se convierte a categorica para que el modelo funcione de manera correcta, dado que la variable estrato no se puede manejar como variable continua.
#convirtiendo la variable estrato en factor
base1$estrato <- factor(base1$estrato)
#mostrar variable
str(base1$estrato)
## Factor w/ 4 levels "3","4","5","6": 3 3 4 2 3 2 3 3 1 1 ...
La variable **parqueadero* como se había evidenciado tiene valores faltantes, por consiguiente se decide imputar por la media.
#Calcular la media
media_parqueaderos <- mean(base1$parqueaderos, na.rm = TRUE)
print(media_parqueaderos)
## [1] 2.181609
Dado que la media es 2.1816 un # que no es entero y hablar de 2.18 parqueaderos no es razonable, se imputa por el entero, es decir 2 parqueaderos
#Redondear la media
media_parqueaderos_round <- round(media_parqueaderos)
# Imputar la media a los valores faltantes en la variable parqueadero
base1$parqueaderos[is.na(base1$parqueaderos)] <- media_parqueaderos_round
summary(base1$parqueaderos)
| Min. | 1st Qu. | Median | Mean | 3rd Qu. | Max. |
|---|---|---|---|---|---|
| 1 | 2 | 2 | 2.109418 | 2 | 10 |
base1_dep <- subset(base1,
select = c(
areaconst,
estrato,
banios,
habitaciones,
parqueaderos,
preciom)
)
head(base1_dep, 3)
| areaconst | estrato | banios | habitaciones | parqueaderos | preciom |
|---|---|---|---|---|---|
| 150 | 5 | 4 | 6 | 2 | 320 |
| 380 | 5 | 3 | 3 | 2 | 780 |
| 445 | 6 | 7 | 6 | 2 | 750 |
#Creando modelo regresión lineal múltiple
Modelo_RLM <- lm(preciom ~ areaconst + estrato + habitaciones + parqueaderos + banios, data = base1_dep)
summary(Modelo_RLM)
##
## Call:
## lm(formula = preciom ~ areaconst + estrato + habitaciones + parqueaderos +
## banios, data = base1_dep)
##
## Residuals:
## Min 1Q Median 3Q Max
## -933.85 -71.11 -14.45 48.04 1072.76
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.6492 18.8210 0.460 0.64598
## areaconst 0.8115 0.0429 18.916 < 2e-16 ***
## estrato4 84.8879 17.1771 4.942 9.65e-07 ***
## estrato5 136.2928 16.0700 8.481 < 2e-16 ***
## estrato6 328.4214 26.3102 12.483 < 2e-16 ***
## habitaciones 1.5514 4.0773 0.380 0.70369
## parqueaderos 17.5554 5.6248 3.121 0.00187 **
## banios 23.4042 5.3655 4.362 1.48e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 155.9 on 714 degrees of freedom
## Multiple R-squared: 0.6657, Adjusted R-squared: 0.6624
## F-statistic: 203.1 on 7 and 714 DF, p-value: < 2.2e-16
Ecuación del modelo:
\(\widehat{Precio(Millones COP)}= 8.65 + 0.81\widehat{(areaconst m^{2})} + 84.89\widehat{(estrato4)} + 136.29\widehat{(estrato5)} + 328.42\widehat{(estrato6)} + 1.55\widehat{(Habitaciones)} + 17.55\widehat{(Parqueaderos)} + 23.40\widehat{(Baños)} + error\)
#Coeficientes del modelo
Modelo_RLM$coefficients
## (Intercept) areaconst estrato4 estrato5 estrato6 habitaciones
## 8.6491491 0.8114659 84.8879383 136.2928537 328.4213911 1.5514072
## parqueaderos banios
## 17.5554244 23.4041592
Intercepto: El intercepto se interpreta como el valor de la variable respuesta cuando los demás variables predictoras son iguales a cero, el intercepto es de 8.65, hablar de las demás variables sean 0 es un poco confuso, ya que se estaría suponiendo que todas las variables son estrato3
areaconst: El coeficiente de la variable es de 0.81, esto se podría interpretar como, manteniendo las demás variables predictoras constantes, por cada aumento en la unidad de área construida aumentaría en 0.81 millones el precio de la vivienda.
estrato4: El coeficiente de la variable es de 84.89, esto se podría interpretar como, en comparación con el estrato3 el precio de estas viviendas se encuentran 84.89 millones más altas que las viviendas de estrato3 si se mantienen las demás variables constantes.
estrato5: El coeficiente de la variable es de 136.29, esto se podría interpretar como, en comparación con el estrato3 el precio de estas viviendas se encuentran 136.29 millones más altas que las viviendas de estrato3 si se mantienen las demás variables constantes.
estrato6: El coeficiente de la variable es de 328.42, esto se podría interpretar como, en comparación con el estrato3 el precio de estas viviendas se encuentran 328.42 millones más altas que las viviendas de estrato3 si se mantienen las demás variables constantes.
habitaciones: El coeficiente de la variable es de 1.55, esto se podría interpretar como, manteniendo las demás variables predictoras constantes, por cada aumento en la unidad de número de habitaciones aumentaría en 1.55 millones el precio de la vivienda.
parqueaderos: El coeficiente de la variable es de 17.55, esto se podría interpretar como, manteniendo las demás variables predictoras constantes, por cada aumento en la unidad de número de parqueaderos aumentaría en 17.55 millones el precio de la vivienda.
Baños: El coeficiente de la variable es de 23.4, esto se podría interpretar como, manteniendo las demás variables predictoras constantes, por cada aumento en la unidad de número de baños aumentaría en 23.4 millones el precio de la vivienda.
Todas la variables, excepto Habitaciones se consideran estadisticamente significativos, dado que los valores de P son menores que 0.05. Anteriormente se había notado que la variable menos correlacionada con el precio de las viviendas era la variable habitaciones, lo que nos sugiere esto es que la variable habitaciones puede no ser significativo para el modelo de predecir el precio de las viviendas.
Para mejorar el modelo se podría considerar agregar variables adicionales al modelo, como lo puede ser, la percepción de seguridad de la zona, el estado de la vivienda, lugares cercanos para el esparcimiento (centros comerciales, parques, restaurantes, complejos deportivos).
Realice la validación de supuestos del modelo e interprete los resultados (no es necesario corregir en caso de presentar problemas, solo realizar sugerencias de que se podría hacer).
summary(Modelo_RLM$residuals)
| Min. | 1st Qu. | Median | Mean | 3rd Qu. | Max. |
|---|---|---|---|---|---|
| -933.8461 | -71.10791 | -14.44889 | 0 | 48.04388 | 1072.764 |
#Histograma residuales
hist(Modelo_RLM$residuals,
main = "Histograma Residuales del Modelo ",
xlab = "Residuales Modelo",
ylab = "Frecuencia",
col = "royalblue",
)
Los residuales muestran una media cercana a cero y según el histograma se puede observar que aparentemente la distribución es simétrica. Esto indica que los supuestos de normalidad, homocedasticidad y linealidad de los residuales podrían estar siendo cumplidos. Para una validación más completa se realizan las siguientes pruebas.
#Mostrar gráficas en 2x2
par(mfrow= c(2,2))
plot(Modelo_RLM)
Se observa que:
Residuals vs. Fitted: Los residuales no están dispersos de manera uniforme alrededor de cero. Esto sugiere que existe cierta no linealidad en la relación entre las variables o que la varianza de los residuales no es constante en todos los niveles de los valores ajustados.
Normal Q-Q: El gráfico muestra que en los extremos, los cuantiles de los residuales no siguen una línea diagonal. Esto indica que los residuales no se ajustan perfectamente a una distribución normal. Puede haber desviaciones en las colas de la distribución, lo que podría sugerir problemas con la normalidad de los residuales.
Scale-Location: Se observa una tendencia creciente en los residuales estandarizados, lo que indica que no se distribuyen de manera uniforme alrededor de cero a medida que aumentan los valores ajustados. Esto sugiere una posible violación de la homocedasticidad, es decir, la igualdad de la varianza de los residuales en todos los niveles de los valores ajustados.
Residuals vs Leverage: Se identifican valores atípicos, como el punto 664, que tienen un impacto desproporcionado en el modelo. Estos valores podrían influir significativamente en los resultados del análisis.
A continuación se validan los supuestos con las pruebas de hipótesis.
residualesRML <- residuals(Modelo_RLM)
# Normalidad
shapiro.test(residualesRML)
##
## Shapiro-Wilk normality test
##
## data: residualesRML
## W = 0.82765, p-value < 2.2e-16
#Homoscedasticidad
library(lmtest)
bptest(Modelo_RLM)
##
## studentized Breusch-Pagan test
##
## data: Modelo_RLM
## BP = 134.03, df = 7, p-value < 2.2e-16
#No autocorrelación
dwtest(Modelo_RLM)
##
## Durbin-Watson test
##
## data: Modelo_RLM
## DW = 1.6994, p-value = 1.952e-05
## alternative hypothesis: true autocorrelation is greater than 0
Para abordar este problema de autocorrelación en los residuos se podría agregar nuevas variables al modelo o se podría transformar variables que ya están en el modelo.
Realice una partición en los datos de forma aleatoria donde 70% sea un set para entrenar el modelo y 30% para prueba. Estime el modelo con la muestra del 70%. Muestre los resultados.
Esta partición se realizará con la biblioteca y librería caret
#Instalar y cargar la biblioteca
#install.packages("caret") #solo una vez
#install.packages("recipes") #solo una vez
library(caret)
#semilla
set.seed(38)
#Partir los datos de train
#Variable objetivo
y= base1_dep$preciom
#Indice de partición
Indice_Train <- createDataPartition(
y,
p= 0.7,
list = FALSE
)
Dividir los datos en test y train
#Datos de entrenamiento con el 70% de los datos
X_train <- base1_dep [Indice_Train, ]
#Datos de testeo con el 30% de los datos
X_test <- base1_dep[-Indice_Train, ]
#Datos objetivo de entrenamiento
y_train <- base1_dep$preciom[Indice_Train]
#Datos objetivo de testeo
y_test <- base1_dep$preciom[-Indice_Train]
#Validando partición
dim(X_train)
## [1] 507 6
dim(X_test)
## [1] 215 6
head(X_test, 3)
| areaconst | estrato | banios | habitaciones | parqueaderos | preciom |
|---|---|---|---|---|---|
| 445 | 6 | 7 | 6 | 2 | 750 |
| 355 | 4 | 5 | 5 | 3 | 625 |
| 237 | 5 | 6 | 6 | 2 | 750 |
Ahora, se va a estimar un modelo utilizando los datos de train
Modelo_RLM_Caret <- lm( preciom ~ areaconst + estrato + habitaciones + parqueaderos + banios, data = X_train)
summary(Modelo_RLM_Caret)
##
## Call:
## lm(formula = preciom ~ areaconst + estrato + habitaciones + parqueaderos +
## banios, data = X_train)
##
## Residuals:
## Min 1Q Median 3Q Max
## -754.36 -70.86 -16.71 42.73 1028.82
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.21413 21.51692 0.196 0.84481
## areaconst 0.91164 0.05355 17.025 < 2e-16 ***
## estrato4 87.74462 19.57289 4.483 9.14e-06 ***
## estrato5 126.86541 18.60193 6.820 2.64e-11 ***
## estrato6 283.33587 31.67355 8.946 < 2e-16 ***
## habitaciones 2.97975 4.60208 0.647 0.51762
## parqueaderos 19.34934 6.38951 3.028 0.00259 **
## banios 15.82927 6.15926 2.570 0.01046 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 151.3 on 499 degrees of freedom
## Multiple R-squared: 0.6681, Adjusted R-squared: 0.6635
## F-statistic: 143.5 on 7 and 499 DF, p-value: < 2.2e-16
_** Realice predicciones con el modelo anterior usando los datos de prueba (30%). **_
Predicciones <- predict( Modelo_RLM_Caret, newdata = X_test)
#print(Predicciones)
summary(Predicciones)
| Min. | 1st Qu. | Median | Mean | 3rd Qu. | Max. |
|---|---|---|---|---|---|
| 91.51038 | 267.7025 | 437.0383 | 456.0115 | 552.8144 | 1429.435 |
#Histograma Predicciones
hist(Predicciones,
main = "Histograma Predicciones ",
xlab = "Predicciones",
ylab = "Frecuencia",
col = "mediumseagreen",
)
library(caret)
# Calcular las métricas de rendimiento
metricas <- postResample(pred = Predicciones, obs = X_test$preciom)
# Mostrar las métricas
print(metricas)
## RMSE Rsquared MAE
## 169.4037721 0.6476401 103.1906396
RMSE ( Root Mean Squared Error): Este es el error cuadratico medio de las predicciones, lo que indica este resultado es que en promedio las predicciones del modelo se encuentran 169.40 millones de distancia al valor real.
R cuadrado El R cuadrado indica cuanta varianza en la variable respuesta es explicada por el modelo, en este caso se habla de que el modelo explica un 64.76% de la varianza de los precios de las viviendas.
MAE (Mean Absolute Error): Este indicador representa la media de las diferencias absolutas entre las predicciones y el valor real, lo que nos indica que en promedio el valor predicho se encuentra a 103.19 millones del valor real.