load("Valencia_Sale2.RData")
Valencia_Sale2 <- as.data.frame(Valencia_Sale)

#MODELOS PREDICTIVOS

El primer paso en la construcción de los modelos predictivos para el mercado inmobiliario de Valencia en 2018 consistió en la carga y adecuación del conjunto de datos Valencia_Sale. Este fue convertido a un data.frame denominado Valencia_Sale2 para facilitar su manipulación en R.

##Selección de variables para el modelo Variables en el dataframe de 2025 y las análogas en el de 2018:

floor –> FLOORCLEAN price –> PRICE size –> CONSTRUCTEDAREA exterior –> Eliminar propertyType –> Eliminar rooms –> ROOMNUMBER bathrooms –> BATHNUMBER district –> DISTRITO distance –> DISTANCE_TO_CITY_CENTER status –> BUILTTYPEID_1, BUILTTYPEID_2, BUILTTYPEID_3 hasLift –> HASLIFT priceByArea –> UNITPRICE distancia_min_estacion_m –> DISTANCE_TO_METRO RHabitacion_Banyo –> Hay que crearla

Finalmente, variables como RHabitacion_Banyo y priceByArea las descartamos a la hora de hacer el modelo ya que al ser variables derivadas directamente de price/size o rooms/bathrooms generaría fuga de información. No añadirían información nueva y pueden generar multicolinealidad. La multicolinealidad es un problema que puede surgir en los modelos de regresión lineal cuando dos o más variables independientes (predictoras) están altamente correlacionadas entre sí. Es decir, una variable independiente puede predecirse en gran medida a partir de otra.

variables_usar <- c("ASSETID", "FLOORCLEAN", "PRICE", "CONSTRUCTEDAREA", "ROOMNUMBER", "BATHNUMBER", "DISTRITO", "BUILTTYPEID_1", "BUILTTYPEID_2", "BUILTTYPEID_3","DISTANCE_TO_CITY_CENTER", "HASLIFT", "DISTANCE_TO_METRO")
datosmodelo <- Valencia_Sale[, variables_usar]

Creamos una nueva columna llamada status en datosmodelo. Esta nueva columna sintetiza la información que antes estaba repartida en tres variables booleanas: BUILTTYPEID_1, BUILTTYPEID_2, BUILTTYPEID_3

El valor de status se determina así: Si BUILTTYPEID_1==1, entonces status = “Nueva”. Si BUILTTYPEID_1 no es 1 pero BUILTTYPEID_2==1, entonces status = “Restaurar”. Si ninguna de las anteriores (es decir, solo queda la posibilidad de que BUILTTYPEID_3==1), entonces status = “BuenEstado”.

Esto simplifica el modelo, evita multicolinealidad y hace más fácil la interpretación y el tratamiento posterior en modelos estadísticos o de machine learning.

Conviertimos la columna status en un factor de R, es decir, una variable categórica y eliminamos las columnas que ya no nos interesan.

datosmodelo$status <- ifelse(datosmodelo$BUILTTYPEID_1 == 1, "Nueva",
                        ifelse(datosmodelo$BUILTTYPEID_2 == 1, "Restaurar", "BuenEstado"))
datosmodelo$status <- as.factor(datosmodelo$status)
# Ahora puedes eliminar las 3 columnas dummy
datosmodelo <- datosmodelo[, !(names(datosmodelo) %in% c("BUILTTYPEID_1", "BUILTTYPEID_2", "BUILTTYPEID_3"))]

Calculamos la matriz de correlaciones entre las varibles numéricas para evitar problemas de multicolinealidad.

Dado que nuestro dataset tenía información espacial (columna geometry típica de objetos sf en R), procedemos a eliminarla usando st_drop_geometry(), ya que no aporta información relevante a nuestro análisis predictivo y podía generar errores con funciones estándar de R. Posteriormente, convertimos el resultado en un data.frame convencional para asegurar la compatibilidad con procedimientos estadísticos clásicos.

datosmodelo <- as.data.frame(st_drop_geometry(datosmodelo))

Procedemos a un análisis de correlación entre las variables numéricas. Este paso es fundamental para comprender las interrelaciones lineales entre los predictores y, crucialmente, para detectar posibles problemas de multicolinealidad que podrían afectar la estabilidad y la interpretabilidad de los modelos estadísticos que queremos desarrollar.

library(corrplot)
## corrplot 0.95 loaded
num_vars <- sapply(datosmodelo, is.numeric)
datos_num <- datosmodelo[, num_vars, drop=FALSE]
if (ncol(datos_num) > 1) {
  cor_matrix <- cor(datos_num, use = "complete.obs")
  corrplot(cor_matrix, method = "color")
} else {
  print("No hay suficientes variables numéricas para calcular correlación.")
}

En resumen, el análisis de correlación confirma relaciones esperadas en el mercado inmobiliario, como la influencia positiva del tamaño y el número de dependencias en el precio, y la influencia negativa de la distancia a puntos de interés. Si bien se observan correlaciones altas entre algunas variables predictoras (especialmente CONSTRUCTEDAREA, ROOMNUMBER y BATHNUMBER), la decisión previa de no incluir variables directamente derivadas como priceByArea o RHabitacion_Banyo ya ha contribuido a mitigar problemas severos de multicolinealidad.

Comenzamos la etapa de construcción de modelos predictivos dividiendo nuestro conjunto de datos datosmodelo. Para ello, utilizamos la función createDataPartition del paquete caret, estableciendo una semilla (set.seed(129)) para garantizar la reproducibilidad de los resultados. Esta función separó los datos en un conjunto de entrenamiento (trainData), que comprende el 70% de las observaciones, y un conjunto de prueba (testData) con el 30% restante.

set.seed(129)
trainIndex <- createDataPartition(datosmodelo$PRICE, p = 0.7, list = FALSE)
trainData <- datosmodelo[trainIndex, ]
testData  <- datosmodelo[-trainIndex, ]

# Elimina ASSETID de ambos sets
trainData <- dplyr::select(trainData, -ASSETID)
testData  <- dplyr::select(testData, -ASSETID)

Una vez creados los conjuntos de entrenamiento y prueba, eliminamos la variable ASSETID de ambos. Consideramos que ASSETID, al ser un identificador único de cada inmueble, no posee capacidad predictiva y su inclusión podría generar un sobreajuste del modelo a los datos de entrenamiento o interpretaciones erróneas de su importancia.

Antes de ajustar el modelo revisamos si tenemos datos faltantes:

colSums(is.na(trainData))
##              FLOORCLEAN                   PRICE         CONSTRUCTEDAREA 
##                     937                       0                       0 
##              ROOMNUMBER              BATHNUMBER                DISTRITO 
##                       0                       0                       0 
## DISTANCE_TO_CITY_CENTER                 HASLIFT       DISTANCE_TO_METRO 
##                       0                       0                       0 
##                  status 
##                       0

Como la variable número de piso (FLOORCLEAN) es importante para el modelo y sólo el 3.5% aproximadamente de los valores estaban ausentes, decidimos imputar los valores faltantes con la mediana. Esto permite conservar la mayor cantidad de información posible sin introducir sesgos extremos, ya que la mediana es robusta ante valores atípicos.

# Calculamos la mediana (ignorando los NA)
mediana_floor <- median(datosmodelo$FLOORCLEAN, na.rm = TRUE)

floorclean_original <- Valencia_Sale$FLOORCLEAN

# Imputamos los NA con la mediana
datosmodelo$FLOORCLEAN[is.na(datosmodelo$FLOORCLEAN)] <- mediana_floor

par(mfrow = c(1, 2))

# Histograma de FLOORCLEAN antes de la imputación
hist(floorclean_original, 
     main = "FLOORCLEAN (Sin Imputación)",
     xlab = "Número de Piso (FLOORCLEAN)",
     ylab = "Frecuencia",
     col = "lightblue",
     border = "black",
     breaks = 20) # Puedes ajustar el número de 'breaks'

# Histograma de FLOORCLEAN después de la imputación con la mediana
hist(datosmodelo$FLOORCLEAN, 
     main = "FLOORCLEAN (Imputación Mediana)",
     xlab = "Número de Piso (FLOORCLEAN)",
     ylab = "Frecuencia",
     col = "lightgreen",
     border = "black",
     breaks = 20) # Puedes ajustar el número de 'breaks'

# Restauramos la configuración gráfica por defecto
par(mfrow = c(1, 1))

La comparación visual mediante estos histogramas nos permite confirmar que la imputación con la mediana, si bien, como es esperado, incrementa la frecuencia del valor de la mediana (en este caso, el tercer piso), no altera drásticamente la forma general de la distribución de FLOORCLEAN. Las viviendas siguen concentrándose mayoritariamente en pisos bajos. Esta técnica nos permite retener la variable y las observaciones asociadas para el análisis predictivo sin introducir distorsiones extremas en las características subyacentes de la variable.

##REGRESIÓN LINEAL MÚLTIPLE

Usaremos 10-fold cross-validation. El modelo incluirá todas las variables numéricas y factores relevantes. Utilizamos PRICE como variable dependiente y todas las demás como predictoras.

'
# Configuración de la validación cruzada
ctrl <- trainControl(method = "cv", number = 10)

# Ajuste del modelo
modelo_lm <- train(
  PRICE ~ ., 
  data = trainData, 
  method = "lm",
  trControl = ctrl
)

# Evaluación en test y métricas
predicciones <- predict(modelo_lm, newdata = testData)
rmse <- RMSE(predicciones, testData$PRICE)
mae  <- MAE(predicciones, testData$PRICE)
r2   <- R2(predicciones, testData$PRICE)
cat("RMSE:", rmse, "\nMAE:", mae, "\nR2:", r2, "\n")
'
## [1] "\n# Configuración de la validación cruzada\nctrl <- trainControl(method = \"cv\", number = 10)\n\n# Ajuste del modelo\nmodelo_lm <- train(\n  PRICE ~ ., \n  data = trainData, \n  method = \"lm\",\n  trControl = ctrl\n)\n\n# Evaluación en test y métricas\npredicciones <- predict(modelo_lm, newdata = testData)\nrmse <- RMSE(predicciones, testData$PRICE)\nmae  <- MAE(predicciones, testData$PRICE)\nr2   <- R2(predicciones, testData$PRICE)\ncat(\"RMSE:\", rmse, \"\nMAE:\", mae, \"\nR2:\", r2, \"\n\")\n"

Nos da error, nos aseguramos por si acaso de si hay columnas constantes y no hay, nos aseguramos de que todos los nombres de las columnas son válidos para fórmulas en R

names(trainData) <- make.names(names(trainData))
names(testData)  <- make.names(names(testData))

Tenemos variables que son en realidad factores (por ejemplo, HASLIFT como 0/1), probaremos en convertirlas explícitamente a factor:

trainData$HASLIFT <- as.factor(trainData$HASLIFT)
testData$HASLIFT  <- as.factor(testData$HASLIFT)

Las variables categóricas con múltiples niveles, como DISTRITO y status, necesitan una transformación para ser incluidas adecuadamente en un modelo de regresión lineal. Optamos por una codificación one-hot (también conocida como creación de variables dummy). Este proceso se realizó utilizando dummyVars sobre el conjunto datosmodelo (excluyendo ASSETID y PRICE). La función predict de dummyVars generó un nuevo conjunto de columnas binarias, donde cada columna representa una categoría de las variables originales (por ejemplo, DISTRITO.AlgunDistrito, status.Nueva, etc.). Estos datos transformados, junto con la variable PRICE, conformaron un nuevo dataframe datos_listo. A partir de datos_listo, volvimos a generar los conjuntos trainData y testData con la misma semilla (set.seed(129)) y proporción (70/30).

variables_categoricas <- setdiff(names(datosmodelo), c("ASSETID", "PRICE"))
dummies <- dummyVars(~ ., data = datosmodelo[, variables_categoricas, drop=FALSE])
datos_dummies <- predict(dummies, newdata = datosmodelo[, variables_categoricas, drop=FALSE])

datos_listo <- data.frame(
  datos_dummies,
  PRICE = datosmodelo$PRICE
)

# División train/test 
set.seed(129)
trainIndex <- createDataPartition(datos_listo$PRICE, p = 0.7, list = FALSE)
trainData <- datos_listo[trainIndex, ]
testData  <- datos_listo[-trainIndex, ]

# Entrenamiento del modelo (especifica que ASSETID NO es predictor)
ctrl <- trainControl(method = "cv", number = 10)
modelo_lm <- train(
  PRICE ~ .,     
  data = trainData,
  method = "lm",
  trControl = ctrl
)

# Evaluación en test y métricas
predicciones <- predict(modelo_lm, newdata = testData)
rmse <- RMSE(predicciones, testData$PRICE)
mae  <- MAE(predicciones, testData$PRICE)
r2   <- R2(predicciones, testData$PRICE)
cat("RMSE:", rmse, "\nMAE:", mae, "\nR2:", r2, "\n")
## RMSE: 91516.08 
## MAE: 52919.66 
## R2: 0.7270535

El RMSE de 91516.08 € indica que, en promedio, las predicciones del modelo se desvían unos 91,516 € del precio real de las viviendas en el conjunto de prueba. El MAE de 52919.66 € nos dice que el error absoluto promedio de las predicciones es de aproximadamente 52,920 €. Ambas métricas nos dan una medida de la magnitud del error de predicción en las unidades originales de la variable PRICE.

El R² de 0.727 (o 72.7%) es una medida de la bondad de ajuste del modelo. Nos indica que aproximadamente el 72.7% de la variabilidad en los precios de las viviendas en el conjunto de prueba puede ser explicada por las variables predictoras incluidas en nuestro modelo lineal. Un R² de esta magnitud sugiere que el modelo tiene una capacidad explicativa considerable, aunque todavía hay un 27.3% de la varianza que no se explica por el modelo.

# Para ver cómo se comporta el modelo

library(scales)
## Warning: package 'scales' was built under R version 4.4.3
## 
## Adjuntando el paquete: 'scales'
## The following object is masked from 'package:readr':
## 
##     col_factor
comparacion <- data.frame(Real = testData$PRICE, Predicho = predicciones)

ggplot(comparacion, aes(x = Real, y = Predicho)) +
  geom_point(aes(color = "Predicción"), alpha = 0.6) +
  geom_abline(slope = 1, intercept = 0, color = "red", linetype = "dashed", size = 1) +
  scale_x_continuous(labels = comma, name = "Precio real (€)") +
  scale_y_continuous(labels = comma, name = "Precio predicho (€)") +
  scale_color_manual(name = "Leyenda", values = c("Predicción" = "black")) +
  labs(
    title = "Comparación entre Precio Real y Precio Predicho",
    subtitle = "Modelo de Regresión Lineal",
    caption = "Línea roja = ajuste perfecto (y = x)"
  ) +
  theme_minimal(base_size = 14) +
  theme(
    legend.position = "bottom",
    plot.title = element_text(face = "bold"),
    plot.subtitle = element_text(size = 12)
  )
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

# Para ver qué variables son más relevantes para el modelo

importancia <- varImp(modelo_lm, scale = TRUE)
importancia_df <- importancia$importance
importancia_df$variable <- rownames(importancia_df)
library(dplyr)
importancia_df <- arrange(importancia_df, desc(Overall))

ggplot(importancia_df, aes(x = reorder(variable, Overall), y = Overall)) +
  geom_bar(stat = "identity", fill = "skyblue") +
  coord_flip() +
  labs(title = "Importancia de las variables (Regresión Lineal)",
       x = "Variable",
       y = "Importancia (relativa)") +
  theme_minimal()

Observamos que la mayoría de los puntos se agrupan alrededor de la línea de ajuste perfecto, especialmente para los precios más bajos y medios (hasta aproximadamente 750,000 € - 1,000,000 €). Esto indica que el modelo tiene un desempeño razonable en este rango. A medida que los precios reales aumentan (por encima de 1,000,000 €), notamos una mayor dispersión de los puntos y una tendencia del modelo a subestimar los precios más altos. Es decir, para las viviendas más caras, las predicciones tienden a ser inferiores a los valores reales. También hay algunos casos donde el modelo sobreestima. La heterocedasticidad (varianza no constante de los errores) podría ser un factor aquí, donde la variabilidad de las predicciones aumenta con el nivel del precio.

Para identificar qué variables tienen mayor influencia en las predicciones del modelo lineal, utilizamos la función varImp del paquete caret. La Figura 3 muestra la importancia relativa de cada variable predictora en el modelo. El gráfico presenta las variables predictoras ordenadas de mayor a menor importancia (de arriba hacia abajo). La longitud de la barra azul claro indica la importancia relativa de cada variable, con la más influyente (en este caso, CONSTRUCTEDAREA) escalada a un valor de 100. CONSTRUCTEDAREA (Superficie Construida): Es, con diferencia, la variable más importante en nuestro modelo lineal. Esto es coherente con la intuición y los análisis previos de correlación, indicando que el tamaño de la vivienda es el principal factor que el modelo utiliza para predecir el precio. Variables de Distrito: La segunda variable más importante es DISTRITO.L.Eixample. Esto sugiere que pertenecer al distrito de L’Eixample (comparado con el distrito de referencia que se haya omitido en la codificación dummy) tiene un impacto significativo en el precio según el modelo. Otras variables de distrito como DISTRITO.Ciutat.Vella, DISTRITO.El.Pla.del.Real y DISTRITO.Camins.al.Grau también aparecen con una importancia considerable, aunque menor que L’Eixample. Esto resalta la relevancia de la ubicación geográfica a nivel de distrito. Características de la Vivienda: ROOMNUMBER (Número de Habitaciones) y BATHNUMBER (Número de Baños) le siguen en importancia, lo cual es lógico, ya que, junto con la superficie, definen las características básicas de la vivienda. FLOORCLEAN (Número de Piso) también muestra una contribución relevante.

# Transforma la variable objetivo en entrenamiento
trainData$LOG_PRICE <- log(trainData$PRICE)

# Entrena el modelo con validación cruzada
ctrl <- trainControl(method = "cv", number = 5)
modelo_lm_log <- train(
  LOG_PRICE ~ . -PRICE, 
  data = trainData,
  method = "lm",
  trControl = ctrl
)

# Predice (en escala logarítmica)
pred_log <- predict(modelo_lm_log, newdata = testData)

# Vuelve a la escala original
predicciones <- exp(pred_log)

# Métricas de evaluación
rmse <- RMSE(predicciones, testData$PRICE)
mae  <- MAE(predicciones, testData$PRICE)
r2   <- R2(predicciones, testData$PRICE)
cat("RMSE:", rmse, "\nMAE:", mae, "\nR2:", r2, "\n")
## RMSE: 301750.4 
## MAE: 52096.88 
## R2: 0.2511464

Hemos probado a transformar la variable objetivo PRICE aplicándole el logaritmo natural, con la expectativa de mejorar el comportamiento del modelo lineal, especialmente en lo referente a la heterocedasticidad y la asimetría de la variable precio.

Se exploró la transformación logarítmica de la variable objetivo PRICE con el objetivo de mejorar la linealidad de las relaciones y estabilizar la varianza de los errores. Contrariamente a las expectativas iniciales para este tipo de transformación, el modelo resultante, al ser evaluado en la escala original de precios, mostró un deterioro en métricas clave como el RMSE (que aumentó de 91,516 € a 301,750 €) y el R² (que disminuyó de 0.727 a 0.251). El MAE, sin embargo, se mantuvo en un nivel similar (aproximadamente 52,000 €).

Este resultado sugiere que la retrotransformación exponencial de las predicciones logarítmicas pudo haber amplificado ciertos errores, o que la relación subyacente entre los predictores y el precio se ajustaba mejor mediante un modelo lineal directo en la escala original para este conjunto de datos.

##RANDOM FOREST

trainData <- dplyr::select(trainData, -LOG_PRICE)
set.seed(129)

ctrl <- trainControl(method = "cv", number = 3)

modelo_rf <- train(
  PRICE ~ .,
  data = trainData,
  method = "rf",
  trControl = ctrl,
  importance = TRUE,
  ntree = 50
)
# 
# Evaluación
pred_rf <- predict(modelo_rf, newdata = testData)
rmse_rf <- RMSE(pred_rf, testData$PRICE)
mae_rf  <- MAE(pred_rf, testData$PRICE)
r2_rf   <- R2(pred_rf, testData$PRICE)
cat("Random Forest - RMSE:", rmse_rf, "\nMAE:", mae_rf, "\nR2:", r2_rf, "\n")
## Random Forest - RMSE: 67660.69 
## MAE: 32814.38 
## R2: 0.8508601

El modelo Random Forest ha logrado una mejora sustancial en todas las métricas de evaluación:

El RMSE ha disminuido de aproximadamente 91,516 € a 67,661 €. Esto indica que, en promedio, los errores cuadráticos de predicción son considerablemente menores. El MAE también ha disminuido notablemente, pasando de unos 52,920 € a 32,814 €. El error absoluto promedio de las predicciones se ha reducido en casi 20,000 €. El R² ha aumentado de 0.727 a 0.851. Esto significa que el modelo Random Forest es capaz de explicar aproximadamente el 85.1% de la variabilidad en los precios de las viviendas del conjunto de prueba, una mejora de más del 12% en la varianza explicada en comparación con el modelo lineal.

df_pred <- data.frame(
  PrecioReal = testData$PRICE,
  PrecioPredicho = pred_rf
)

library(scales) # Para las etiquetas de los ejes

ggplot(df_pred, aes(x = PrecioReal, y = PrecioPredicho)) +
  geom_point(aes(color = "Predicción"), alpha = 0.5) +
  geom_abline(slope = 1, intercept = 0, color = "red", linetype = "dashed", linewidth = 1) +
  scale_x_continuous(name = "Precio Real (€)", labels = comma) +
  scale_y_continuous(name = "Precio Predicho (€)", labels = comma) +
  scale_color_manual(name = "Leyenda", values = c("Predicción" = "black")) +
  labs(
    title = "Comparación entre Precio Real y Precio Predicho",
    subtitle = "Modelo Random Forest (ntree = 50, cv = 3)",
    caption = "Línea roja discontinua = Ajuste perfecto (Predicho = Real)"
  ) +
  theme_minimal(base_size = 14) +
  theme(
    legend.position = "bottom",
    plot.title = element_text(face = "bold", hjust = 0.5),
    plot.subtitle = element_text(hjust = 0.5),
    axis.title.x = element_text(margin = ggplot2::margin(t = 10, r = 0, b = 0, l = 0, unit = "pt")), # Usamos ggplot2::margin
    axis.title.y = element_text(margin = ggplot2::margin(t = 0, r = 10, b = 0, l = 0, unit = "pt"))  # Usamos ggplot2::margin
  )

var_imp <- varImp(modelo_rf)
plot(var_imp, top = 20, main = "Importancia de las variables (Random Forest)")

Los puntos (predicciones individuales) se muestran más agrupados alrededor de la línea roja discontinua (ajuste perfecto) en comparación con el gráfico del modelo de regresión lineal. Esto es especialmente notable en el rango de precios medios y altos. Aunque todavía existe cierta dispersión, y el modelo puede subestimar o sobrestimar algunas viviendas, la tendencia general es una mejora en la precisión a lo largo de todo el espectro de precios.

El problema de subestimación sistemática de los precios más altos, que era evidente en el modelo lineal, parece haberse mitigado considerablemente, aunque las viviendas de precios extremadamente altos (superiores a 1,500,000 € - 2,000,000 €) siguen presentando mayores errores, lo cual es común dada su escasez y variabilidad inherente.

De manera similar al modelo lineal, CONSTRUCTEDAREA emerge como la variable más importante para el modelo Random Forest, con una puntuación de importancia de 100. Esto reafirma que el tamaño de la vivienda es el factor individual más determinante del precio. Variables de Distancia: A diferencia del modelo lineal donde los distritos específicos tenían una alta importancia justo después de la superficie, aquí las variables DISTANCE_TO_CITY_CENTER (Distancia al Centro de la Ciudad) y DISTANCE_TO_METRO (Distancia al Metro) ocupan el segundo y tercer lugar, respectivamente, con puntuaciones de importancia alrededor de 83 y 78. Esto sugiere que el modelo Random Forest encuentra un fuerte poder predictivo en la proximidad a puntos clave de la ciudad, posiblemente capturando gradientes de precios más continuos que las categorías de distrito por sí solas. Distritos Específicos: DISTRITO.Benimaclet aparece como el distrito individual más importante (cuarto en la general con una importancia cercana a 45), seguido de DISTRITO.Camins.al.Grau y DISTRITO.Poblats.Marítims. Es notable que DISTRITO.L.Eixample, que era la segunda variable más importante en el modelo lineal, ahora tiene una importancia mucho menor en el modelo Random Forest

Vamos a probar a mejorar el modelo cambiando el número de árboles y el número de cross validation anque esto va a suponer un cambio en el tiempo de cómputo. Además, aplicaremos grid search de mtry. Con doParallel para aprovechar todos los núcleos de nuestro CPU:

'
set.seed(129)

tunegrid <- expand.grid(.mtry = c(2, 6, 10))
ctrl <- trainControl(method = "repeatedcv", number = 3, repeats = 1)

library(doParallel)
cl <- makePSOCKcluster(parallel::detectCores() - 1)
registerDoParallel(cl)

modelo_rf <- train(
  PRICE ~ .,
  data = trainData,
  method = "rf",
  trControl = ctrl,
  tuneGrid = tunegrid,
  ntree = 150,         # con 500 no funciona
  importance = TRUE
)

# Evaluación
pred_rf <- predict(modelo_rf, newdata = testData)
rmse_rf <- RMSE(pred_rf, testData$PRICE)
mae_rf  <- MAE(pred_rf, testData$PRICE)
r2_rf   <- R2(pred_rf, testData$PRICE)
cat("Random Forest - RMSE:", rmse_rf, "\nMAE:", mae_rf, "\nR2:", r2_rf, "\n")
'
## [1] "\nset.seed(129)\n\ntunegrid <- expand.grid(.mtry = c(2, 6, 10))\nctrl <- trainControl(method = \"repeatedcv\", number = 3, repeats = 1)\n\nlibrary(doParallel)\ncl <- makePSOCKcluster(parallel::detectCores() - 1)\nregisterDoParallel(cl)\n\nmodelo_rf <- train(\n  PRICE ~ .,\n  data = trainData,\n  method = \"rf\",\n  trControl = ctrl,\n  tuneGrid = tunegrid,\n  ntree = 150,         # con 500 no funciona\n  importance = TRUE\n)\n\n# Evaluación\npred_rf <- predict(modelo_rf, newdata = testData)\nrmse_rf <- RMSE(pred_rf, testData$PRICE)\nmae_rf  <- MAE(pred_rf, testData$PRICE)\nr2_rf   <- R2(pred_rf, testData$PRICE)\ncat(\"Random Forest - RMSE:\", rmse_rf, \"\nMAE:\", mae_rf, \"\nR2:\", r2_rf, \"\n\")\n"

Los resultados mostraron mejoras marginales en el rendimiento: el RMSE se redujo ligeramente a 66,644.44 € y el R² aumentó a 0.8556, mientras que el MAE se mantuvo estable en torno a los 32,860 €.

Cuando un modelo está casi optimizado, pequeños cambios en hiperparámetros (como mtry, ntree o el número de folds) no suelen mejorar mucho el error, y a veces la variación se debe solo al azar del split de datos o al propio Random Forest. El rango de error se mantiene estable, lo que indica que el modelo ya está capturando casi toda la información relevante de los datos.

(Lo comentamos ya que es similar al anterior, para ahorrar tiempo al cargar el fichero)

##XGBoost

Suele ofrecer mejor precisión que Random Forest. Es más rápidos para entrenar. Permite controlar el overfitting fácilmente Tiene muchas opciones de tuning

set.seed(129)

# Grid de hiperparámetros básico
tunegrid <- expand.grid(
  nrounds = 100,           # número de iteraciones
  max_depth = c(3, 6, 9),  # profundidad de los árboles
  eta = c(0.05, 0.1, 0.3), # learning rate
  gamma = 0,
  colsample_bytree = 1,
  min_child_weight = 1,
  subsample = 1
)

ctrl <- trainControl(method = "repeatedcv", number = 5, repeats = 1)

xgb_model <- train(
  PRICE ~ .,
  data = trainData,
  method = "xgbTree",
  trControl = ctrl,
  tuneGrid = tunegrid,
  verbose = FALSE
)

# Evaluación
pred_xgb <- predict(xgb_model, newdata = testData)
rmse_xgb <- RMSE(pred_xgb, testData$PRICE)
mae_xgb  <- MAE(pred_xgb, testData$PRICE)
r2_xgb   <- R2(pred_xgb, testData$PRICE)
cat("XGBoost - RMSE:", rmse_xgb, "\nMAE:", mae_xgb, "\nR2:", r2_xgb, "\n")
## XGBoost - RMSE: 70004.63 
## MAE: 34858.91 
## R2: 0.8402484

En esta primera iteración con un grid de hiperparámetros básico, el modelo XGBoost no ha superado el rendimiento del modelo Random Forest optimizado. Sin embargo, su rendimiento sigue siendo muy bueno y significativamente mejor que el del modelo de regresión lineal inicial. El RMSE del XGBoost (70,004.63 €) es ligeramente superior al del Random Forest optimizado (66,644.44 €). El MAE del XGBoost (34,858.91 €) también es un poco más alto que el del Random Forest (32,860.48 €). El R² del XGBoost (0.840) es ligeramente inferior al del Random Forest (0.856).

df_pred_xgb <- data.frame(
  PrecioReal = testData$PRICE,
  PrecioPredicho = pred_xgb
)

ggplot(df_pred_xgb, aes(x = PrecioReal, y = PrecioPredicho)) +
  geom_point(aes(color = "Predicción"), alpha = 0.5) + # Añadido aes() para la leyenda
  geom_abline(slope = 1, intercept = 0, color = "red", linetype = "dashed", linewidth = 1) + # Usar linewidth en lugar de size
  scale_x_continuous(name = "Precio Real (€)", labels = comma) + # Formato con comas
  scale_y_continuous(name = "Precio Predicho (€)", labels = comma) + # Formato con comas
  scale_color_manual(name = "Leyenda", values = c("Predicción" = "black")) + # Define color y etiqueta para leyenda
  labs(
    title = "Comparación entre Precio Real y Precio Predicho",
    subtitle = "Modelo XGBoost", # Actualiza el subtítulo
    caption = "Línea roja discontinua = Ajuste perfecto (Predicho = Real)"
  ) +
  theme_minimal(base_size = 14) +
  theme(
    legend.position = "bottom",
    plot.title = element_text(face = "bold", hjust = 0.5),
    plot.subtitle = element_text(hjust = 0.5),
    # Intentamos de nuevo con ggplot2::margin, si da error, puedes comentar estas dos líneas
    axis.title.x = element_text(margin = ggplot2::margin(t = 10, r = 0, b = 0, l = 0, unit = "pt")),
    axis.title.y = element_text(margin = ggplot2::margin(t = 0, r = 10, b = 0, l = 0, unit = "pt"))
  )

var_imp_xgb <- varImp(xgb_model)
plot(var_imp_xgb, top = 20, main = "Importancia de las variables (XGBoost)")

Al igual que con Random Forest, se espera que los puntos se agrupen razonablemente bien alrededor de la línea de ajuste perfecto (roja discontinua). La dispersión de los puntos, especialmente en los rangos de precios más altos, nos dará una indicación visual de su precisión y de si presenta patrones de error similares a los modelos anteriores.

Una vez más, y de forma consistente con todos los modelos anteriores (Regresión Lineal y Random Forest), CONSTRUCTEDAREA es, con diferencia, la variable más importante para el modelo XGBoost, alcanzando el máximo de 100 en la escala de importancia. Esto subraya de manera inequívoca que el tamaño de la vivienda es el principal factor predictivo del precio. DISTANCE_TO_CITY_CENTER (Distancia al Centro de la Ciudad): Ocupa el segundo lugar con una importancia considerable, alrededor de 18-20. Esto es similar a lo que vimos con el modelo Random Forest, donde las variables de distancia continua también eran muy relevantes, aunque en Random Forest su importancia relativa era mayor en comparación con XGBoost. BATHNUMBER (Número de Baños): Esta variable aparece en tercer lugar, con una importancia alrededor de 10-12. En el modelo Random Forest, BATHNUMBER no estaba tan arriba en el ranking (estaba por debajo de FLOORCLEAN y varios distritos). Esto sugiere que XGBoost le da un peso relativamente mayor al número de baños que el Random Forest en esta configuración particular.

##LIGHTGBM

No fue posible entrenar LightGBM utilizando caret, ya que este modelo no está incluido en la librería nativa de caret. Se optó por entrenar LightGBM directamente con su propio paquete, lo que permite igualmente comparar sus resultados con los de otros modelos y analizar su desempeño en la predicción del precio de las viviendas.

library(lightgbm)
## Warning: package 'lightgbm' was built under R version 4.4.3
# Prepara los datos en formato matrix y vector
train_x <- as.matrix(subset(trainData, select = -PRICE))
train_y <- trainData$PRICE

test_x <- as.matrix(subset(testData, select = -PRICE))
test_y <- testData$PRICE

# Crea datasets de lightgbm
dtrain <- lgb.Dataset(data = train_x, label = train_y)

# Ajusta el modelo
params <- list(
  objective = "regression",
  metric = "rmse"
)

modelo_lgb <- lgb.train(
  params = params,
  data = dtrain,
  nrounds = 100
)
## [LightGBM] [Info] Auto-choosing row-wise multi-threading, the overhead of testing was 0.001636 seconds.
## You can set `force_row_wise=true` to remove the overhead.
## And if memory is not enough, you can set `force_col_wise=true`.
## [LightGBM] [Info] Total Bins 838
## [LightGBM] [Info] Number of data points in the train set: 19125, number of used features: 26
## [LightGBM] [Info] Start training from score 194440.627451
# Predicciones
pred_lgb <- predict(modelo_lgb, test_x)

# Métricas de evaluación
rmse_lgb <- sqrt(mean((test_y - pred_lgb)^2))
mae_lgb <- mean(abs(test_y - pred_lgb))
r2_lgb <- 1 - sum((test_y - pred_lgb)^2) / sum((test_y - mean(test_y))^2)

cat("LightGBM - RMSE:", rmse_lgb, "\nMAE:", mae_lgb, "\nR2:", r2_lgb, "\n")
## LightGBM - RMSE: 72355.27 
## MAE: 39014.97 
## R2: 0.8293307

El RMSE del LightGBM (72,355.27 €) es superior al del Random Forest optimizado y también al del XGBoost. El MAE del LightGBM (39,014.97 €) es el más alto entre estos tres modelos (RF, XGBoost, LightGBM). El R² del LightGBM (0.829) es el más bajo de los tres. Con esta configuración inicial y sin un ajuste de hiperparámetros más profundo, el modelo LightGBM no ha logrado superar el rendimiento de Random Forest ni de XGBoost. Sin embargo, su rendimiento sigue siendo considerablemente mejor que el modelo de regresión lineal original.

df_disp <- data.frame(
  Precio_Real = test_y,
  Precio_Predicho = pred_lgb
)

ggplot(df_disp, aes(x = Precio_Real, y = Precio_Predicho)) +
  geom_point(alpha = 0.5) +
  geom_abline(slope = 1, intercept = 0, col = "red", linetype = "dashed", size = 1.2) +
  labs(
    title = "Comparación entre Precio Real y Precio Predicho\nLightGBM",
    x = "Precio real (€)",
    y = "Precio predicho (€)"
  ) +
  theme_minimal() +
  scale_x_continuous(labels = scales::comma) +
  scale_y_continuous(labels = scales::comma)

# Obtener importancia de variables del modelo LightGBM
importancia <- lgb.importance(modelo_lgb)
importancia$Feature <- factor(importancia$Feature, levels = importancia$Feature[order(importancia$Gain, decreasing = TRUE)])

ggplot(importancia, aes(x = Gain, y = Feature)) +
  geom_point(size = 3, colour = "blue") +
  labs(
    title = "Importancia de las variables (LightGBM)",
    x = "Importancia",
    y = NULL
  ) +
  theme_minimal()

El gráfico muestra una agrupación de puntos similar a la que vimos con XGBoost. Los puntos se distribuyen alrededor de la línea de ajuste perfecto (roja discontinua). La dispersión parece ser comparable a la de XGBoost. Al igual que con los otros modelos basados en árboles, parece manejar mejor la predicción en rangos de precios bajos y medios, con una mayor dispersión (errores más grandes) en los precios más altos. Visualmente, no parece haber una mejora drástica o un empeoramiento significativo en comparación con el gráfico de XGBoost.

Vamos a intentar mejorar un poco más el modelo del Random Forest. La eliminación de outliers es una de las formas más efectivas de mejorar la precisión de un modelo de regresión en el sector inmobiliario, ya que permite reducir la influencia de anuncios irreales o sesgados en el cálculo del MAE. Al fin y al cabo, hay que tener en cuenta que cualquiera puede subir un anuncio a Idealista y que la gente quiere ganar el máximo dinero posible así que muchas veces tienden a sobreestimar su casa. Vamos a probar con quedarnos con el percentil 95% y ver si mejora el valor del MAE:

library(randomForest)
## randomForest 4.7-1.2
## Type rfNews() to see new features/changes/bug fixes.
## 
## Adjuntando el paquete: 'randomForest'
## The following object is masked from 'package:ggplot2':
## 
##     margin
## The following object is masked from 'package:dplyr':
## 
##     combine
# 1. Calcular el percentil 95 de precio en train y test
p95_train <- quantile(trainData$PRICE, probs = 0.95)
p95_test  <- quantile(testData$PRICE,  probs = 0.95)

# 2. Filtrar los datos (solo eliminamos los precios por encima del percentil 95)
trainData_filtrado <- trainData[trainData$PRICE <= p95_train, ]
testData_filtrado  <- testData[testData$PRICE <= p95_test, ]

# 3. Transformación logarítmica del precio
trainData_filtrado$LOG_PRICE <- log(trainData_filtrado$PRICE)
testData_filtrado$LOG_PRICE  <- log(testData_filtrado$PRICE)

set.seed(129)
modelo_rf_log <- randomForest(
  LOG_PRICE ~ . -PRICE, 
  data = trainData_filtrado
)

# 5. Predicciones en test (en escala logarítmica)
pred_rf_log <- predict(modelo_rf_log, newdata = testData_filtrado)

# 6. Llevar las predicciones de vuelta a escala original
pred_rf_log_exp <- exp(pred_rf_log)

# 7. Calcular métricas en escala original
rmse_log <- sqrt(mean((testData_filtrado$PRICE - pred_rf_log_exp)^2))
mae_log  <- mean(abs(testData_filtrado$PRICE - pred_rf_log_exp))
r2_log   <- 1 - sum((testData_filtrado$PRICE - pred_rf_log_exp)^2) / sum((testData_filtrado$PRICE - mean(testData_filtrado$PRICE))^2)

cat("Random Forest (con filtrado y log) - RMSE:", rmse_log, "\nMAE:", mae_log, "\nR2:", r2_log, "\n")
## Random Forest (con filtrado y log) - RMSE: 38747.33 
## MAE: 25343.29 
## R2: 0.8360346

Tras aplicar un proceso completo de filtrado de outliers (eliminando el 5% superior de precios) y una transformación logarítmica sobre la variable objetivo, el modelo Random Forest ha experimentado una mejora significativa en la predicción de precios inmobiliarios:

El RMSE ha disminuido espectacularmente de ~66,644 € a ~38,747 €. Esto representa una reducción del error cuadrático medio de aproximadamente el 41.9%. El MAE también ha mejorado significativamente, pasando de ~32,860 € a ~25,343 €. Esto es una reducción del error absoluto medio de aproximadamente el 22.9%. Este era tu objetivo principal y el resultado es excelente. El modelo, en promedio, se equivoca en unos 25,343 € en la predicción de los precios de las viviendas dentro del 95% del rango de precios.

El R² ha disminuido ligeramente de ~0.856 a ~0.836. Esto puede parecer contraintuitivo al principio, dado que los errores absolutos y cuadráticos han bajado. Sin embargo, es importante recordar que el R² se calcula sobre los datos de prueba filtrados. Al eliminar los precios más altos (el 5% superior), la varianza total en testData_filtradoPRICE es probablemente menor que la varianza en el testDataPRICE original. El modelo ahora explica el 83.6% de esta nueva varianza (más pequeña) de los precios en el rango del 0 al 95 percentil.

df_disp <- data.frame(
  Precio_Real = testData_filtrado$PRICE,
  Precio_Predicho = pred_rf_log_exp
)

ggplot(df_disp, aes(x = Precio_Real, y = Precio_Predicho)) +
  geom_point(alpha = 0.4) +
  geom_abline(slope = 1, intercept = 0, col = "red", linetype = "dashed", size = 1.2) +
  labs(
    title = "Random Forest: Precio Real vs Precio Predicho (filtrado + log)",
    x = "Precio real (€)",
    y = "Precio predicho (€)"
  ) +
  theme_minimal() +
  scale_x_continuous(labels = comma) +
  scale_y_continuous(labels = comma)

# Gráfico de Residuos

df_disp$residuo <- df_disp$Precio_Real - df_disp$Precio_Predicho

ggplot(df_disp, aes(x = Precio_Real, y = residuo)) +
  geom_point(alpha = 0.4) +
  geom_hline(yintercept = 0, col = "red", linetype = "dashed") +
  labs(
    title = "Residuos del modelo Random Forest (filtrado + log)",
    x = "Precio real (€)",
    y = "Error de predicción (€)"
  ) +
  theme_minimal() +
  scale_x_continuous(labels = comma) +
  scale_y_continuous(labels = comma)

# Importancia de las variables

importancia <- importance(modelo_rf_log)
importancia_df <- data.frame(
  Variable = rownames(importancia),
  Importancia = importancia[, 1]
)
importancia_df <- importancia_df[order(importancia_df$Importancia, decreasing = TRUE), ]

ggplot(importancia_df[1:15, ], aes(x = reorder(Variable, Importancia), y = Importancia)) +
  geom_col(fill = "blue") +
  coord_flip() +
  labs(
    title = "Importancia de las variables en Random Forest (filtrado + log)",
    x = NULL,
    y = "Importancia"
  ) +
  theme_minimal()

Tras aplicar un proceso completo de filtrado de outliers (eliminando el 5% superior de precios) y una transformación logarítmica sobre la variable objetivo, el modelo Random Forest ha experimentado una mejora significativa en la predicción de precios inmobiliarios:

Los puntos parecen agruparse más densamente alrededor de la línea de ajuste perfecto (roja discontinua), especialmente en los rangos de precios más bajos y medios. Esto es consistente con la mejora en RMSE y MAE. La dispersión sigue aumentando para los precios más altos dentro de este rango filtrado, lo cual es un comportamiento típico. Visualmente, el modelo parece tener un buen rendimiento en este subconjunto de datos.

Se observa que los residuos se distribuyen alrededor de la línea de error cero. Sin embargo, se evidencia un patrón de heterocedasticidad: la dispersión de los residuos aumenta a medida que el precio real de la vivienda se incrementa. Esto indica que, si bien el modelo es bastante preciso para las viviendas de menor precio, la magnitud de los errores tiende a ser mayor para las viviendas más caras dentro del rango analizado. No se observan otros sesgos sistemáticos pronunciados. Este comportamiento es común en la modelización de precios inmobiliarios y sugiere que, aunque el modelo ha mejorado significativamente en términos de MAE y RMSE, su precisión absoluta varía con el nivel de precio.

El análisis de importancia de variables para este modelo destaca CONSTRUCTEDAREA como el factor más influyente, seguido por BATHNUMBER, DISTANCE_TO_CITY_CENTER y HASLIFT. Esto indica que, dentro del rango de precios más común, el tamaño, el número de baños, la proximidad al centro y la disponibilidad de ascensor son determinantes clave.

Dadas sus métricas de error significativamente más bajas para el grueso del mercado, este modelo Random Forest (con filtrado y transformación logarítmica) se considera el más adecuado para la predicción de precios de viviendas en este estudio.

##MODELO CON CLUSTERS

A continuación, haremos uso de los clusters previamente calculados para separar el municipio de Valencia y aplicar un modelo(sólo aplicaremos Random Forest por falta de tiempo y por ser el que mejor resultados nos ha dado anteriormente) a cada clúster para ver si de esta forma, es capaz de crear mejores predicciones con un error mayor. Los clústers que tenemos son:

Clúster 1 -> Benimaclet, Jesús Clúster 2 -> Camins al Grau, Campanar Clúster 3 -> Ciutat Vella, El Pla del Real, Extramurs, L’Eixample Clúster 4 -> Algirós, L’Olivereta, Patraix Clúster 5 -> Benicalap, La Saïdia, Poblats Marítims, Quatre Carreres, Rascanya

# 1. Definir los clústers
cluster_list <- list(
  cluster1 = c("Benimaclet", "Jesús"),
  cluster2 = c("Camins al Grau", "Campanar"),
  cluster3 = c("Ciutat Vella", "El Pla del Real", "Extramurs", "L'Eixample"),
  cluster4 = c("Algirós", "L'Olivereta", "Patraix"),
  cluster5 = c("Benicalap", "La Saïdia", "Poblats Marítims", "Quatre Carreres", "Rascanya")
)

# 2. Añadir columna CLUSTER en base a DISTRITO
datosmodelo$CLUSTER <- NA
for (i in seq_along(cluster_list)) {
  datosmodelo$CLUSTER[datosmodelo$DISTRITO %in% cluster_list[[i]]] <- paste0("cluster", i)
}
datosmodelo$CLUSTER <- as.factor(datosmodelo$CLUSTER)

# 3. Bucle para entrenar y evaluar un modelo por clúster
resultados <- list()

for (cl in levels(datosmodelo$CLUSTER)) {
  cat("\nProcesando", cl, "...\n")
  
  # Filtrar datos del clúster
  datos_cluster <- subset(datosmodelo, CLUSTER == cl)
  
  # Dividir en train y test (70/30)
  set.seed(129)
  trainIndex <- createDataPartition(datos_cluster$PRICE, p = 0.7, list = FALSE)
  trainData <- datos_cluster[trainIndex, ]
  testData  <- datos_cluster[-trainIndex, ]
  
  # Entrenar Random Forest (ajusta la fórmula según tus variables, aquí un ejemplo)
  modelo_rf <- randomForest(PRICE ~ . -ASSETID -CLUSTER -DISTRITO, data = trainData)
  
  # Predecir en test
  pred_rf <- predict(modelo_rf, newdata = testData)
  
  # Calcular métricas
  rmse <- sqrt(mean((testData$PRICE - pred_rf)^2))
  mae  <- mean(abs(testData$PRICE - pred_rf))
  r2   <- 1 - sum((testData$PRICE - pred_rf)^2) / sum((testData$PRICE - mean(testData$PRICE))^2)
  
  resultados[[cl]] <- list(
    modelo = modelo_rf,
    rmse = rmse,
    mae = mae,
    r2 = r2
  )
  
  cat(sprintf("Cluster: %s | RMSE: %.2f | MAE: %.2f | R²: %.3f\n", cl, rmse, mae, r2))
}
## 
## Procesando cluster1 ...
## Cluster: cluster1 | RMSE: 33382.22 | MAE: 21486.05 | R²: 0.845
## 
## Procesando cluster2 ...
## Cluster: cluster2 | RMSE: 82146.57 | MAE: 44865.74 | R²: 0.767
## 
## Procesando cluster3 ...
## Cluster: cluster3 | RMSE: 115554.00 | MAE: 64770.64 | R²: 0.768
## 
## Procesando cluster4 ...
## Cluster: cluster4 | RMSE: 40152.91 | MAE: 23122.47 | R²: 0.594
## 
## Procesando cluster5 ...
## Cluster: cluster5 | RMSE: 41114.23 | MAE: 25644.19 | R²: 0.775

Clusters 1, 4 y 5: Tienen métricas de error (RMSE, MAE) más bajas y un R² aceptable o bueno frente al modelo global. Esto indica que la segmentación ha ayudado a captar particularidades de esos barrios y mejorar la predicción. Clusters 2 y 3: Presentan errores mucho más altos y un R² más bajo. Es probable que: Sean zonas más heterogéneas internamente (mucha variabilidad de precios, tipologías, etc.). Los datos de estos clústers sean más escasos o tengan más outliers. El modelo no tenga suficientes datos/potencia predictiva para captar bien sus particularidades.

Vamos a probar a englobar en dos grandes clusters que van a ver centro y periféria para así tener más datos en cada clúster y que el modelo pueda entrenarse mejor. Los clusters son:

Clúster 1 -> Ciutat Vella, L’Eixample, El Pla del Real, Extramurs Clúster 2 -> Poblats Marítims, Benicalap, Rascanya, La Saïdia, Quatre Carreres, Benimaclet, Jesús, Camins al Grau, Patraix, Campanar, Algirós, L’Olivereta

# 1. Definir los dos grandes clústers
cluster_centro <- c("Ciutat Vella", "L'Eixample", "El Pla del Real", "Extramurs")
cluster_periferia <- c(
  "Poblats Marítims", "Benicalap", "Rascanya", "La Saïdia", "Quatre Carreres",
  "Benimaclet", "Jesús", "Camins al Grau", "Patraix", "Campanar", "Algirós", "L'Olivereta"
)

# 2. Crear columna CLUSTER en base a DISTRITO
datosmodelo$CLUSTER <- NA
datosmodelo$CLUSTER[datosmodelo$DISTRITO %in% cluster_centro] <- "Centro"
datosmodelo$CLUSTER[datosmodelo$DISTRITO %in% cluster_periferia] <- "Periferia"
datosmodelo$CLUSTER <- as.factor(datosmodelo$CLUSTER)

# 3. Bucle para entrenar y evaluar un modelo por clúster
resultados <- list()

for (cl in levels(datosmodelo$CLUSTER)) {
  cat("\nProcesando", cl, "...\n")
  
  # Filtrar datos del clúster
  datos_cluster <- subset(datosmodelo, CLUSTER == cl)
  
  # Dividir en train y test (75/25)
  set.seed(129)
  trainIndex <- createDataPartition(datos_cluster$PRICE, p = 0.75, list = FALSE)
  trainData <- datos_cluster[trainIndex, ]
  testData  <- datos_cluster[-trainIndex, ]
  
  # Entrenar Random Forest 
  modelo_rf <- randomForest(PRICE ~ . -ASSETID -CLUSTER -DISTRITO, data = trainData)
  
  # Predecir en test
  pred_rf <- predict(modelo_rf, newdata = testData)
  
  # Calcular métricas
  rmse <- sqrt(mean((testData$PRICE - pred_rf)^2))
  mae  <- mean(abs(testData$PRICE - pred_rf))
  r2   <- 1 - sum((testData$PRICE - pred_rf)^2) / sum((testData$PRICE - mean(testData$PRICE))^2)
  
  resultados[[cl]] <- list(
    modelo = modelo_rf,
    rmse = rmse,
    mae = mae,
    r2 = r2
  )
  
  cat(sprintf("Cluster: %s | RMSE: %.2f | MAE: %.2f | R²: %.3f\n", cl, rmse, mae, r2))
  # ----- GRAFICO: Precio real vs predicho -----
  p <- ggplot(data.frame(Real=testData$PRICE, Predicho=pred_rf), aes(x=Real, y=Predicho)) +
    geom_point(alpha=0.4) +
    geom_abline(slope=1, intercept=0, linetype="dashed", color="red", size=1.2) +
    labs(
      title = paste("Random Forest:", cl, "- Precio Real vs Precio Predicho"),
      x = "Precio real (€)",
      y = "Precio predicho (€)"
    ) +
    theme_minimal()
  
  print(p)
}
## 
## Procesando Centro ...
## Cluster: Centro | RMSE: 105834.78 | MAE: 63127.57 | R²: 0.796

## 
## Procesando Periferia ...
## Cluster: Periferia | RMSE: 47913.14 | MAE: 29731.10 | R²: 0.752

Los resultados nos muestran que el modelo por clústeres (Centro y Periferia) no mejora el error respecto al modelo global, especialmente en el clúster Centro. ¿Por qué los resultados no mejoran?

Tamaño de muestra reducido Al dividir en dos clústeres, cada modelo se entrena con menos datos, especialmente en el centro (barrios como Ciutat Vella, L’Eixample, etc. tienen menos viviendas que toda la ciudad). Esto afecta la capacidad del Random Forest de generalizar y aumenta el error, sobre todo en presencia de muchos outliers o alta dispersión de precios.

Mayor heterogeneidad interna en el “Centro” El centro suele tener una mayor variedad de tipos de vivienda y precios (pisos de lujo, viviendas antiguas, áticos, etc.), lo que complica la predicción. Si no tienes suficientes variables que expliquen bien esa variedad, el modelo lo sufre.

Presencia de outliers/extremos En el gráfico se ve que hay bastantes puntos alejados de la diagonal, sobre todo para precios altos. El centro suele tener propiedades muy caras (“outliers”) que el modelo no predice bien. Si no aplicaste filtrado de outliers ni transformación logarítmica (como sí hiciste en el modelo global), el error se dispara.

Menos capacidad para aprender patrones generales Cuando entrenas el modelo global, puede aprender relaciones que son válidas para toda la ciudad (por ejemplo, el tamaño, el estado, etc.). Al segmentar demasiado, puedes perder esa generalidad.

Vamos a aplicar la transformación logarítmica y usar sólo el percentil 95 a ver si mejoran los resultados.

# 1. Definir los dos grandes clústers
cluster_centro <- c("Ciutat Vella", "L'Eixample", "El Pla del Real", "Extramurs")
cluster_periferia <- c(
  "Poblats Marítims", "Benicalap", "Rascanya", "La Saïdia", "Quatre Carreres",
  "Benimaclet", "Jesús", "Camins al Grau", "Patraix", "Campanar", "Algirós", "L'Olivereta"
)

# 2. Crear columna CLUSTER en base a DISTRITO
datosmodelo$CLUSTER <- NA
datosmodelo$CLUSTER[datosmodelo$DISTRITO %in% cluster_centro] <- "Centro"
datosmodelo$CLUSTER[datosmodelo$DISTRITO %in% cluster_periferia] <- "Periferia"
datosmodelo$CLUSTER <- as.factor(datosmodelo$CLUSTER)

# 3. Bucle para entrenar y evaluar un modelo por clúster
resultados <- list()
resumen <- data.frame(
  Cluster = character(),
  Precio_Medio = numeric(),
  MAE = numeric(),
  Error_Relativo = numeric(),
  stringsAsFactors = FALSE
)

for (cl in levels(datosmodelo$CLUSTER)) {
  cat("\nProcesando", cl, "...\n")
  
  # Filtrar datos del clúster
  datos_cluster <- subset(datosmodelo, CLUSTER == cl)
  
  # Dividir en train y test (80/20)
  set.seed(129)
  trainIndex <- createDataPartition(datos_cluster$PRICE, p = 0.8, list = FALSE)
  trainData <- datos_cluster[trainIndex, ]
  testData  <- datos_cluster[-trainIndex, ]
  
  # Entrenar Random Forest (ajusta la fórmula según tus variables)
  modelo_rf <- randomForest(PRICE ~ . -ASSETID -CLUSTER -DISTRITO, data = trainData)
  
  # Predecir en test
  pred_rf <- predict(modelo_rf, newdata = testData)
  
  # Calcular métricas
  rmse <- sqrt(mean((testData$PRICE - pred_rf)^2))
  mae  <- mean(abs(testData$PRICE - pred_rf))
  r2   <- 1 - sum((testData$PRICE - pred_rf)^2) / sum((testData$PRICE - mean(testData$PRICE))^2)
  
  # Calcular precio medio y error relativo
  precio_medio <- mean(testData$PRICE)
  error_relativo <- mae / precio_medio
  
  # Guardar resultados
  resultados[[cl]] <- list(
    modelo = modelo_rf,
    rmse = rmse,
    mae = mae,
    r2 = r2,
    precio_medio = precio_medio,
    error_relativo = error_relativo
  )
  
  # Añadir a la tabla resumen
  resumen <- rbind(
    resumen,
    data.frame(
      Cluster = cl,
      Precio_Medio = round(precio_medio, 2),
      MAE = round(mae, 2),
      Error_Relativo = round(error_relativo, 4)
    )
  )
  
  cat(sprintf("Cluster: %s | RMSE: %.2f | MAE: %.2f | R²: %.3f | Precio medio: %.2f | Error relativo: %.4f\n",
              cl, rmse, mae, r2, precio_medio, error_relativo))
  
  # ----- GRAFICO: Precio real vs predicho -----
  p <- ggplot(data.frame(Real=testData$PRICE, Predicho=pred_rf), aes(x=Real, y=Predicho)) +
    geom_point(alpha=0.4) +
    geom_abline(slope=1, intercept=0, linetype="dashed", color="red", size=1.2) +
    labs(
      title = paste("Random Forest:", cl, "- Precio Real vs Precio Predicho"),
      x = "Precio real (€)",
      y = "Precio predicho (€)"
    ) +
    theme_minimal()
  
  print(p)
}
## 
## Procesando Centro ...
## Cluster: Centro | RMSE: 114154.94 | MAE: 65161.87 | R²: 0.764 | Precio medio: 325286.01 | Error relativo: 0.2003

## 
## Procesando Periferia ...
## Cluster: Periferia | RMSE: 47733.20 | MAE: 28921.02 | R²: 0.745 | Precio medio: 144519.49 | Error relativo: 0.2001

cat("\n----- RESUMEN DE CLÚSTERES -----\n")
## 
## ----- RESUMEN DE CLÚSTERES -----
print(resumen)
##     Cluster Precio_Medio      MAE Error_Relativo
## 1    Centro     325286.0 65161.87         0.2003
## 2 Periferia     144519.5 28921.02         0.2001

Vamos a probar con XGBoost para ver si obtenemos mejores resultados

library(xgboost)
## Warning: package 'xgboost' was built under R version 4.4.3
## 
## Adjuntando el paquete: 'xgboost'
## The following object is masked from 'package:dplyr':
## 
##     slice
library(Matrix)
## 
## Adjuntando el paquete: 'Matrix'
## The following objects are masked from 'package:tidyr':
## 
##     expand, pack, unpack
for (cl in levels(datosmodelo$CLUSTER)) {
  cat("\nProcesando", cl, "con XGBoost...\n")
  
  # 1. Filtrar datos del clúster
  datos_cluster <- subset(datosmodelo, CLUSTER == cl)
  
  # 2. Eliminar outliers por percentil 95
  p95 <- quantile(datos_cluster$PRICE, 0.95)
  datos_cluster <- datos_cluster[datos_cluster$PRICE <= p95, ]
  
  # 3. Transformación logarítmica
  datos_cluster$LOG_PRICE <- log(datos_cluster$PRICE)
  
  # 4. Train/test split
  set.seed(129)
  trainIndex <- createDataPartition(datos_cluster$LOG_PRICE, p = 0.8, list = FALSE)
  trainData <- datos_cluster[trainIndex, ]
  testData  <- datos_cluster[-trainIndex, ]
  
  # 5. Preparar matrices para XGBoost
  # Quita columnas que no quieras usar
  X_train <- model.matrix(LOG_PRICE ~ . -PRICE -ASSETID -CLUSTER -DISTRITO, data = trainData)[, -1]
  X_test <- model.matrix(LOG_PRICE ~ . -PRICE -ASSETID -CLUSTER -DISTRITO, data = testData)[, -1]
  
  y_train <- trainData$LOG_PRICE
  y_test <- testData$LOG_PRICE
  
  # 6. Entrenar XGBoost
  dtrain <- xgb.DMatrix(data = X_train, label = y_train)
  dtest <- xgb.DMatrix(data = X_test)
  
  params <- list(
    objective = "reg:squarederror",
    eval_metric = "rmse"
  )
  
  modelo_xgb <- xgboost(
    data = dtrain,
    params = params,
    nrounds = 150,
    verbose = 0
  )
  
  # 7. Predicciones y métricas en escala original
  pred_xgb_log <- predict(modelo_xgb, dtest)
  pred_xgb_exp <- exp(pred_xgb_log)
  
  rmse <- sqrt(mean((testData$PRICE - pred_xgb_exp)^2))
  mae  <- mean(abs(testData$PRICE - pred_xgb_exp))
  r2   <- 1 - sum((testData$PRICE - pred_xgb_exp)^2) / sum((testData$PRICE - mean(testData$PRICE))^2)
  
  cat(sprintf("XGBoost - Cluster: %s | RMSE: %.2f | MAE: %.2f | R²: %.3f\n", cl, rmse, mae, r2))
  
  # 8. Gráfico real vs predicho
  library(ggplot2)
  p <- ggplot(data.frame(Real=testData$PRICE, Predicho=pred_xgb_exp), aes(x=Real, y=Predicho)) +
    geom_point(alpha=0.4) +
    geom_abline(slope=1, intercept=0, linetype="dashed", color="red", size=1.2) +
    labs(
      title = paste("XGBoost:", cl, "- Precio Real vs Precio Predicho"),
      x = "Precio real (€)",
      y = "Precio predicho (€)"
    ) +
    theme_minimal()
  print(p)
}
## 
## Procesando Centro con XGBoost...
## XGBoost - Cluster: Centro | RMSE: 77897.29 | MAE: 53804.80 | R²: 0.742

## 
## Procesando Periferia con XGBoost...
## XGBoost - Cluster: Periferia | RMSE: 33992.29 | MAE: 24416.76 | R²: 0.680

Al igual que pasaba en los modelos sin clusters, el resultado del modelo con XGBoost no mejora al del Random Forest.

testData$Predicho <- pred_xgb_exp  # ajusta el nombre de tu vector de predicción
testData$ErrorAbs <- abs(testData$PRICE - testData$Predicho)
testData$ErrorRel <- abs(testData$PRICE - testData$Predicho) / testData$PRICE * 100

# Gráfico de residuos absolutos vs precio real
ggplot(testData, aes(x=PRICE, y=ErrorAbs)) +
  geom_point(alpha=0.2, size=1, color="black") +
  geom_smooth(method="loess", color="blue", size=1.2, se=FALSE) +
  labs(title="Error absoluto vs Precio real",
       x="Precio real (€)", y="Error absoluto (€)") +
  theme_minimal() +
  theme(plot.title = element_text(size=18, face="bold"),
        axis.title = element_text(size=14))
## `geom_smooth()` using formula = 'y ~ x'

# Boxplot por distrito
ggplot(testData, aes(x=DISTRITO, y=ErrorAbs)) +
  geom_boxplot() +
  labs(title="Error absoluto por distrito",
       x="Distrito", y="Error absoluto (€)") +
  theme_minimal() +
  coord_flip()

# Top Errores
head(testData[order(-testData$ErrorAbs), c("PRICE", "Predicho", "ErrorAbs", "ErrorRel", "DISTRITO")], 10)

Gráfico Error absoluto vs Precio real

Para viviendas de bajo precio, el error absoluto tiende a ser mucho menor y se concentra cerca de cero, lo que indica que el modelo predice con mayor precisión en este rango. Para viviendas de precio elevado, los errores absolutos aumentan considerablemente, mostrando una mayor dispersión y algunos errores extremos. Esto sugiere que el modelo tiene más dificultad para predecir correctamente los precios altos, probablemente debido a una mayor heterogeneidad de este tipo de propiedades y/o a la presencia de outliers.

El modelo es más preciso para viviendas de precio medio-bajo y tiende a cometer errores más grandes en las viviendas más caras. Esto es habitual en problemas de precios inmobiliarios, ya que las viviendas de alto valor suelen presentar características singulares que los modelos no siempre capturan bien.

Gráfico Error Absoluto por Distrito

Hay distritos como Campanar, Camins al Grau y Benimaclet donde el error absoluto es más alto y más disperso (cajas y bigotes más largos, más outliers). Otros distritos como Algirós y Jesús presentan menor error absoluto y menos dispersión, lo que indica que el modelo predice mejor en estas zonas. Presencia de outliers

Top Errores

En algunos casos el modelo sobreestima (predicho mucho mayor que real, como en Patraix y L’Olivereta) y en otros subestima (predicho mucho menor que real, como en Benimaclet y Camins al Grau). Esto sugiere que el modelo no logra capturar bien las características particulares de ciertas viviendas excepcionales (por ejemplo, viviendas muy reformadas, áticos, bajos, o con características poco frecuentes en el distrito).

Winsorizar es una técnica estadística para limitar el impacto de los valores extremos (outliers) en una variable (por ejemplo, el precio). En vez de eliminar los outliers, lo que hace es “recortar” sus valores y sustituirlos por un percentil determinado.

Reduce la influencia de los valores extremos (outliers) en el entrenamiento de modelos. Hace los modelos más robustos y menos sensibles a anomalías. Ayuda a que métricas como el RMSE y el MAE no se disparen innecesariamente por unos pocos valores atípicos.

Eliminar 5% atípico de cada distrito

A continuación, desarrollaremos otro modelo pero esta vez teniendo en cuenta lso valores atípicos de cada distrito y eliminandolos por arriba para ver si conseguimos obtener un mejor resultado

library(dplyr)
library(randomForest)
library(ggplot2)
library(caret)

# 1. Elimina el 5% superior de precios por distrito
datos_filtrados <- datosmodelo %>%
  group_by(DISTRITO) %>%
  mutate(pct_95 = quantile(PRICE, 0.95, na.rm = TRUE)) %>%
  filter(PRICE <= pct_95) %>%
  ungroup() %>%
  dplyr::select(-pct_95)

# 2. Aplica logaritmo natural al precio
datos_filtrados <- datos_filtrados %>%
  mutate(LOG_PRICE = log(PRICE))

# 3. Define los 5 clústeres
cluster_list <- list(
  cluster1 = c("Benimaclet", "Jesús"),
  cluster2 = c("Camins al Grau", "Campanar"),
  cluster3 = c("Ciutat Vella", "El Pla del Real", "Extramurs", "L'Eixample"),
  cluster4 = c("Algirós", "L'Olivereta", "Patraix"),
  cluster5 = c("Benicalap", "La Saïdia", "Poblats Marítims", "Quatre Carreres", "Rascanya")
)
datos_filtrados$CLUSTER <- NA
for (i in seq_along(cluster_list)) {
  datos_filtrados$CLUSTER[datos_filtrados$DISTRITO %in% cluster_list[[i]]] <- paste0("cluster", i)
}
datos_filtrados$CLUSTER <- as.factor(datos_filtrados$CLUSTER)

# 4. Bucle para entrenar y evaluar un modelo por clúster
resultados <- list()
errores_top10 <- data.frame()
resumen <- data.frame(
  Cluster = character(),
  Precio_Medio = numeric(),
  MAE = numeric(),
  RMSE = numeric(),
  R2 = numeric(),
  Error_Relativo = numeric(),
  stringsAsFactors = FALSE
)
errores_distritos <- data.frame()

for (cl in levels(datos_filtrados$CLUSTER)) {
  cat("\nProcesando", cl, "...\n")
  datos_cluster <- subset(datos_filtrados, CLUSTER == cl)

  set.seed(129)
  trainIndex <- createDataPartition(datos_cluster$LOG_PRICE, p = 0.7, list = FALSE)
  trainData <- datos_cluster[trainIndex, ]
  testData  <- datos_cluster[-trainIndex, ]

  # Entrena Random Forest sobre LOG_PRICE
  modelo_rf <- randomForest(LOG_PRICE ~ . -ASSETID -CLUSTER -DISTRITO -PRICE, data = trainData)

  # Predice en test
  pred_log <- predict(modelo_rf, newdata = testData)
  pred_orig <- exp(pred_log)

  # Errores
  abs_error <- abs(testData$PRICE - pred_orig)
  resultados[[cl]] <- data.frame(
    CLUSTER = cl,
    DISTRITO = testData$DISTRITO,
    REAL = testData$PRICE,
    PREDICHO = pred_orig,
    ABS_ERROR = abs_error
  )
  errores_distritos <- bind_rows(errores_distritos, resultados[[cl]])

  # Métricas del cluster
  precio_medio <- mean(testData$PRICE)
  mae <- mean(abs_error)
  rmse <- sqrt(mean((testData$PRICE - pred_orig)^2))
  r2 <- 1 - sum((testData$PRICE - pred_orig)^2) / sum((testData$PRICE - mean(testData$PRICE))^2)
  error_relativo <- mae / precio_medio

  resumen <- rbind(
    resumen,
    data.frame(
      Cluster = cl,
      Precio_Medio = round(precio_medio, 2),
      MAE = round(mae, 2),
      RMSE = round(rmse, 2),
      R2 = round(r2, 3),
      Error_Relativo = round(error_relativo, 4)
    )
  )

  cat(sprintf("Cluster: %s | Precio medio: %.2f | MAE: %.2f | RMSE: %.2f | R²: %.3f | Error relativo: %.4f\n",
              cl, precio_medio, mae, rmse, r2, error_relativo))

  # Acumula los top 10 errores de este cluster
  top10 <- resultados[[cl]] %>% arrange(desc(ABS_ERROR)) %>% head(10)
  errores_top10 <- bind_rows(errores_top10, top10)
}
## 
## Procesando cluster1 ...
## Cluster: cluster1 | Precio medio: 131801.00 | MAE: 19183.74 | RMSE: 30146.38 | R²: 0.792 | Error relativo: 0.1456
## 
## Procesando cluster2 ...
## Cluster: cluster2 | Precio medio: 241767.34 | MAE: 35535.67 | RMSE: 51832.93 | R²: 0.787 | Error relativo: 0.1470
## 
## Procesando cluster3 ...
## Cluster: cluster3 | Precio medio: 298298.70 | MAE: 55597.11 | RMSE: 85156.73 | R²: 0.746 | Error relativo: 0.1864
## 
## Procesando cluster4 ...
## Cluster: cluster4 | Precio medio: 125259.84 | MAE: 19260.01 | RMSE: 25655.97 | R²: 0.690 | Error relativo: 0.1538
## 
## Procesando cluster5 ...
## Cluster: cluster5 | Precio medio: 121546.59 | MAE: 21259.27 | RMSE: 30350.35 | R²: 0.711 | Error relativo: 0.1749
# 5. Boxplot del error absoluto por distrito
ggplot(errores_distritos, aes(x = reorder(DISTRITO, ABS_ERROR, FUN = median), y = ABS_ERROR)) +
  geom_boxplot(fill = "orange") +
  labs(title = "Boxplot del error absoluto por Distrito",
       x = "Distrito", y = "Error absoluto (€)") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

# 6. Top 10 errores absolutos más grandes (en escala original)
errores_top10 <- errores_top10 %>% arrange(desc(ABS_ERROR)) %>% head(10)
cat("\nTOP 10 de errores absolutos más grandes (precio real, predicho, distrito y cluster):\n")
## 
## TOP 10 de errores absolutos más grandes (precio real, predicho, distrito y cluster):
print(errores_top10[, c("CLUSTER", "DISTRITO", "REAL", "PREDICHO", "ABS_ERROR")])
##       CLUSTER        DISTRITO   REAL PREDICHO ABS_ERROR
## 753  cluster3      L'Eixample 953000 433753.0  519247.0
## 225  cluster3      L'Eixample 949000 443670.0  505330.0
## 802  cluster3      L'Eixample 857000 409244.7  447755.3
## 20   cluster3 El Pla del Real 814000 378769.9  435230.1
## 1562 cluster3      L'Eixample 695000 274605.7  420394.3
## 758  cluster3      L'Eixample 686000 269457.4  416542.6
## 1930 cluster3      L'Eixample 907000 492615.5  414384.5
## 222  cluster3      L'Eixample 869000 457862.0  411138.0
## 762  cluster3      L'Eixample 808000 415916.1  392083.9
## 1861 cluster3      L'Eixample 961000 586230.9  374769.1
# 7. Tabla resumen por cluster
cat("\n----- RESUMEN DE CLÚSTERES -----\n")
## 
## ----- RESUMEN DE CLÚSTERES -----
print(resumen)
##    Cluster Precio_Medio      MAE     RMSE    R2 Error_Relativo
## 1 cluster1     131801.0 19183.74 30146.38 0.792         0.1456
## 2 cluster2     241767.3 35535.67 51832.93 0.787         0.1470
## 3 cluster3     298298.7 55597.11 85156.73 0.746         0.1864
## 4 cluster4     125259.8 19260.01 25655.97 0.690         0.1538
## 5 cluster5     121546.6 21259.27 30350.35 0.711         0.1749

Los resultados mejoran un poco pero el error en el clúster 3 sigue siendo muy alto. Como observamos en el boxplot del error absoluto por distrito que barrios como Extramurs, Ciutat Vella, El Pla del Real y L´Eixample, todos dercanos al centro, poseen valores atípicos mucho más elevados.

Por último, quitaremos el 5% de los valores por clúster, de tal forma que los dos modelos que se entrenen no tendrán en cuenta los valores atípicos de cada clúster.

library(dplyr)
library(randomForest)
library(ggplot2)
library(caret)

# 1. Aplica logaritmo natural al precio (lo adelanto para evitar recalcular después de filtrar)
# Puedes dejar este paso después del filtrado si prefieres

# 2. Define los dos clústeres
cluster_centro <- c("Ciutat Vella", "L'Eixample", "El Pla del Real", "Extramurs")
cluster_periferia <- setdiff(unique(datosmodelo$DISTRITO), cluster_centro)

datosmodelo$CLUSTER <- NA
datosmodelo$CLUSTER[datosmodelo$DISTRITO %in% cluster_centro] <- "Centro"
datosmodelo$CLUSTER[datosmodelo$DISTRITO %in% cluster_periferia] <- "Periferia"
datosmodelo$CLUSTER <- as.factor(datosmodelo$CLUSTER)

# 3. Elimina el 5% superior de precios por clúster
datos_filtrados <- datosmodelo %>%
  group_by(CLUSTER) %>%
  mutate(pct_95 = quantile(PRICE, 0.75, na.rm = TRUE)) %>%
  filter(PRICE <= pct_95) %>%
  ungroup() %>%
  dplyr::select(-pct_95)

# 4. Aplica logaritmo natural al precio
datos_filtrados <- datos_filtrados %>%
  mutate(LOG_PRICE = log(PRICE))

# 5. Bucle para entrenar y evaluar modelo por clúster
resultados <- list()
errores_top10 <- data.frame()
errores_distritos <- data.frame()
resumen <- data.frame(
  Cluster = character(),
  Precio_Medio = numeric(),
  MAE = numeric(),
  RMSE = numeric(),
  R2 = numeric(),
  Error_Relativo = numeric(),
  stringsAsFactors = FALSE
)

for (cl in levels(datos_filtrados$CLUSTER)) {
  cat("\nProcesando", cl, "...\n")
  datos_cluster <- subset(datos_filtrados, CLUSTER == cl)

  set.seed(129)
  trainIndex <- createDataPartition(datos_cluster$LOG_PRICE, p = 0.8, list = FALSE)
  trainData <- datos_cluster[trainIndex, ]
  testData  <- datos_cluster[-trainIndex, ]

  # Entrena Random Forest sobre LOG_PRICE
  modelo_rf <- randomForest(LOG_PRICE ~ . -ASSETID -CLUSTER -DISTRITO -PRICE, data = trainData)

  # Predice en test
  pred_log <- predict(modelo_rf, newdata = testData)
  pred_orig <- exp(pred_log)

  # Errores
  abs_error <- abs(testData$PRICE - pred_orig)
  resultados[[cl]] <- data.frame(
    CLUSTER = cl,
    DISTRITO = testData$DISTRITO,
    REAL = testData$PRICE,
    PREDICHO = pred_orig,
    ABS_ERROR = abs_error
  )
  errores_distritos <- bind_rows(errores_distritos, resultados[[cl]])
  
  # Métricas evaluación cluster
  precio_medio <- mean(testData$PRICE)
  mae <- mean(abs_error)
  rmse <- sqrt(mean((testData$PRICE - pred_orig)^2))
  r2 <- 1 - sum((testData$PRICE - pred_orig)^2) / sum((testData$PRICE - mean(testData$PRICE))^2)
  error_relativo <- mae / precio_medio

  cat(sprintf(
    "Cluster: %s | Precio medio: %.2f | MAE: %.2f | RMSE: %.2f | R²: %.3f | Error relativo: %.4f\n",
    cl, precio_medio, mae, rmse, r2, error_relativo
  ))
  
  resumen <- rbind(
    resumen,
    data.frame(
      Cluster = cl,
      Precio_Medio = round(precio_medio, 2),
      MAE = round(mae, 2),
      RMSE = round(rmse, 2),
      R2 = round(r2, 3),
      Error_Relativo = round(error_relativo, 4)
    )
  )

  # Acumula los top 10 errores de este cluster
  top10 <- resultados[[cl]] %>% arrange(desc(ABS_ERROR)) %>% head(10)
  errores_top10 <- bind_rows(errores_top10, top10)
}
## 
## Procesando Centro ...
## Cluster: Centro | Precio medio: 223477.15 | MAE: 35612.88 | RMSE: 47232.59 | R²: 0.603 | Error relativo: 0.1594
## 
## Procesando Periferia ...
## Cluster: Periferia | Precio medio: 105677.54 | MAE: 16572.46 | RMSE: 21548.59 | R²: 0.644 | Error relativo: 0.1568
# 6. Boxplot del error absoluto por distrito
ggplot(errores_distritos, aes(x = reorder(DISTRITO, ABS_ERROR, FUN = median), y = ABS_ERROR)) +
  geom_boxplot(fill = "orange") +
  labs(title = "Boxplot del error absoluto por Distrito",
       x = "Distrito", y = "Error absoluto (€)") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

# 7. Top 10 errores absolutos más grandes
errores_top10 <- errores_top10 %>% arrange(desc(ABS_ERROR)) %>% head(10)
cat("\nTOP 10 de errores absolutos más grandes (precio real, predicho, distrito y cluster):\n")
## 
## TOP 10 de errores absolutos más grandes (precio real, predicho, distrito y cluster):
print(errores_top10[, c("CLUSTER", "DISTRITO", "REAL", "PREDICHO", "ABS_ERROR")])
##     CLUSTER        DISTRITO   REAL PREDICHO ABS_ERROR
## 941  Centro    Ciutat Vella  66000 249954.4  183954.4
## 221  Centro    Ciutat Vella 346000 180625.4  165374.6
## 838  Centro       Extramurs 117000 281879.0  164879.0
## 687  Centro El Pla del Real 364000 208325.4  155674.6
## 866  Centro       Extramurs 134000 289377.5  155377.5
## 402  Centro      L'Eixample 380000 227383.1  152616.9
## 511  Centro El Pla del Real 356000 204983.5  151016.5
## 771  Centro      L'Eixample  99000 249373.6  150373.6
## 496  Centro       Extramurs 119000 268962.2  149962.2
## 133  Centro      L'Eixample 386000 236750.6  149249.4
# 8. Tabla resumen por cluster
cat("\n----- RESUMEN DE CLÚSTERES -----\n")
## 
## ----- RESUMEN DE CLÚSTERES -----
print(resumen)
##     Cluster Precio_Medio      MAE     RMSE    R2 Error_Relativo
## 1    Centro     223477.1 35612.88 47232.59 0.603         0.1594
## 2 Periferia     105677.5 16572.46 21548.59 0.644         0.1568

Si bien los resultados pueden parecer aceptables a primera vista, al observar el valor de R² ,que ronda 0.6, queda claro que el modelo explica una proporción relativamente baja de la variabilidad en los datos, lo que limita su capacidad predictiva.