Regresión-optimización de hiperpárametros
Daniel Felipe Pérez Grajales
Universidad Nacional de Colombia - Sede Medellín
Efraín Galvis Amaya
Universidad Nacional de Colombia - Sede Medellín
Profesor: Juan David Ospina Arango
Universidad Nacional de Colombia - Sede Medellín
Departamento de Ciencias de la Computación y de la Decisión
Decisiones bajo incertidumbre (Optimización para aprendizaje de máquina)
11 de junio de 2021
Predicción del precio de la vivienda en Taiwan
Utilizando el conjunto de datos Real estate valuation data set Data Set lleve a un ejercicio del predicción del precio de la vivienda en Taiwan utilizando los siguientes métodos de predicción:
Regresión lineal (clásica y elastic net)
Ensambles de árboles: bosques aleatorios y XGBoost
Máquinas de soporte vectorial
Redes neuronales
Es posible dividir el conjunto de datos en entrenamiento y validación o uzar un esquema de validación cruzada.
El conjunto de datos se referncia en la siguiente cita: Citación: Yeh, I. C., & Hsu, T. K. (2018). Building real estate valuation models with comparative approach through case-based reasoning. Applied Soft Computing, 65, 260-271. La definición de las variables es la siguiente:
\(X1\): fecha de la transacción (por ejemplo 2013.250=2013 Marzo, * 2013.500=2013 Junio, etc.)
\(X2\): edad de la casa en años
\(X3\): distancial al MRT (transporte masivo) más cercano en metros
\(X4\): número de tiendas de conveniencia en el vecindario (entero)
\(X5\): latitud (unidad: grados)
\(X6\): longitude (unidad: grados)
\(Y\): precio por unidad de área (10000 Nuevos dólares taiwaneses/ 3.3 \(m^2\))
Data
Descriptivo general
## [1] "Identificación de NA'S por variable"## x1 x2 x3 x4 x5 x6 y
## 0 0 0 0 0 0 0## [1] "Descriptivo1 de variables"## x1 x2 x3 x4
## Min. :2013 Min. : 0.000 Min. : 23.38 Min. : 0.000
## 1st Qu.:2013 1st Qu.: 9.025 1st Qu.: 289.32 1st Qu.: 1.000
## Median :2013 Median :16.100 Median : 492.23 Median : 4.000
## Mean :2013 Mean :17.713 Mean :1083.89 Mean : 4.094
## 3rd Qu.:2013 3rd Qu.:28.150 3rd Qu.:1454.28 3rd Qu.: 6.000
## Max. :2014 Max. :43.800 Max. :6488.02 Max. :10.000
## x5 x6 y
## Min. :24.93 Min. :121.5 Min. : 7.60
## 1st Qu.:24.96 1st Qu.:121.5 1st Qu.: 27.70
## Median :24.97 Median :121.5 Median : 38.45
## Mean :24.97 Mean :121.5 Mean : 37.98
## 3rd Qu.:24.98 3rd Qu.:121.5 3rd Qu.: 46.60
## Max. :25.01 Max. :121.6 Max. :117.50## [1] "Descriptivo2 de variables"## df
##
## 7 Variables 414 Observations
## --------------------------------------------------------------------------------
## x1
## n missing distinct Info Mean Gmd .05 .10
## 414 0 12 0.991 2013 0.3239 2013 2013
## .25 .50 .75 .90 .95
## 2013 2013 2013 2014 2014
##
## lowest : 2012.667 2012.750 2012.833 2012.917 2013.000
## highest: 2013.250 2013.333 2013.417 2013.500 2013.583
##
## Value 2012.667 2012.750 2012.833 2012.917 2013.000 2013.083 2013.167
## Frequency 30 27 31 38 28 46 25
## Proportion 0.072 0.065 0.075 0.092 0.068 0.111 0.060
##
## Value 2013.250 2013.333 2013.417 2013.500 2013.583
## Frequency 32 29 58 47 23
## Proportion 0.077 0.070 0.140 0.114 0.056
## --------------------------------------------------------------------------------
## x2
## n missing distinct Info Mean Gmd .05 .10
## 414 0 236 1 17.71 12.93 1.100 3.500
## .25 .50 .75 .90 .95
## 9.025 16.100 28.150 34.670 37.735
##
## lowest : 0.0 1.0 1.1 1.5 1.7, highest: 40.9 41.3 41.4 42.7 43.8
## --------------------------------------------------------------------------------
## x3
## n missing distinct Info Mean Gmd .05 .10
## 414 0 259 1 1084 1205 90.46 157.61
## .25 .50 .75 .90 .95
## 289.32 492.23 1454.28 2697.66 4082.01
##
## lowest : 23.38284 49.66105 56.47425 57.58945 82.88643
## highest: 4605.74900 5512.03800 6306.15300 6396.28300 6488.02100
## --------------------------------------------------------------------------------
## x4
## n missing distinct Info Mean Gmd .05 .10
## 414 0 11 0.986 4.094 3.371 0 0
## .25 .50 .75 .90 .95
## 1 4 6 8 9
##
## lowest : 0 1 2 3 4, highest: 6 7 8 9 10
##
## Value 0 1 2 3 4 5 6 7 8 9 10
## Frequency 67 46 24 46 31 67 37 31 30 25 10
## Proportion 0.162 0.111 0.058 0.111 0.075 0.162 0.089 0.075 0.072 0.060 0.024
## --------------------------------------------------------------------------------
## x5
## n missing distinct Info Mean Gmd .05 .10
## 414 0 234 1 24.97 0.01382 24.95 24.95
## .25 .50 .75 .90 .95
## 24.96 24.97 24.98 24.98 24.99
##
## lowest : 24.93207 24.93293 24.93363 24.93885 24.94155
## highest: 24.99156 24.99176 24.99800 25.00115 25.01459
## --------------------------------------------------------------------------------
## x6
## n missing distinct Info Mean Gmd .05 .10
## 414 0 232 1 121.5 0.01601 121.5 121.5
## .25 .50 .75 .90 .95
## 121.5 121.5 121.5 121.5 121.5
##
## lowest : 121.4735 121.4752 121.4788 121.4846 121.4951
## highest: 121.5539 121.5548 121.5596 121.5617 121.5663
## --------------------------------------------------------------------------------
## y
## n missing distinct Info Mean Gmd .05 .10
## 414 0 270 1 37.98 15.13 16.49 21.02
## .25 .50 .75 .90 .95
## 27.70 38.45 46.60 54.94 59.17
##
## lowest : 7.6 11.2 11.6 12.2 12.8, highest: 71.0 73.6 78.0 78.3 117.5
## --------------------------------------------------------------------------------Análisis de correlación de variables
Del análisis preliminar se pueden extraer las siguientes conclusiones:
La variable que tienen una mayor relación lineal con el precio por unidad de área: \(x5\): Lat y \(x6\): Long es decir la ubicación del inmueble (\(r= 0.823\) y \(r= 0.760\)).
La variable que tienen medianamente una relación lineal con el precio por unidad de área: \(X3\):distancial al MRT (transporte masivo) más cercano en metros (\(r= -0.507\)).
\(X3\):distancial al MRT (transporte masivo) más cercano en metros y \(X4\): número de tiendas de conveniencia en el vecindario (entero) están medianamente correlacionados (\(r = -0.729\)) y existe una dependencia lineal por lo que posiblemente no sea útil introducir ambos predictores en los modelos.
\(X2\): edad de la casa en años y \(X4\):número de tiendas de conveniencia en el vecindario (entero) están correlacionados (\(r = 0.894\)) y existe una dependencia lineal por lo que posiblemente no sea útil introducir ambos predictores en los modelos.
La variable \(X3\): distancial al MRT (transporte masivo) más cercano en metros se debería hacer una transformación que posiblemente haría más normal su distribución.
Preparación de la información para ser modelada:
selección de la muestra
Escalado de datos:
Advertencia : si se escala la variable respuesta restándole la media entonces en los modelos lineales no se incluye un intercepto
se eliminan las variables lat y long porque no son faciles de interpretar.
creamos nuevas variables
X_scaled$x1_2<-X_scaled$x1^2
X_scaled$x2_2<-X_scaled$x2^2
X_scaled$x3_2<-X_scaled$x3^2
X_scaled$x4_2<-X_scaled$x4^2
datatable(head(X_scaled))Base de entrenamiento y base de validación:
X_tr<-X_scaled[-ix_val,]
X_vl<-X_scaled[ix_val,]
Y_tr<-Y[-ix_val,]
Y_vl<-Y[ix_val,]Modelo lineal
- Calculemos un modelo lineal:
##
## Call:
## lm(formula = y ~ ., data = df_tr)
##
## Residuals:
## Min 1Q Median 3Q Max
## -17.578 -4.408 -0.309 3.477 72.548
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 33.4094 1.0067 33.186 < 2e-16 ***
## x1 2.1350 0.4687 4.555 7.43e-06 ***
## x2 -4.1302 0.5024 -8.221 5.01e-15 ***
## x3 -12.3288 1.0989 -11.219 < 2e-16 ***
## x4 2.2466 0.6555 3.427 0.000689 ***
## x1_2 0.1819 0.5392 0.337 0.736111
## x2_2 2.2387 0.5156 4.342 1.89e-05 ***
## x3_2 2.5533 0.3911 6.529 2.58e-10 ***
## x4_2 -0.1174 0.5026 -0.234 0.815495
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.415 on 322 degrees of freedom
## Multiple R-squared: 0.6553, Adjusted R-squared: 0.6467
## F-statistic: 76.51 on 8 and 322 DF, p-value: < 2.2e-16- Calculo del error cuadrático medio (rmse)
mide cuanto se degrada la predicción del modelo en este caso fue de:
## [1] "RSMS de entrenamiento 8.29998259853576"## [1] "RSMS de validación 8.26858483125769"Se observa que el modelo tuvo un mejor desempeño en la base de validación, se puede concluir que no hay presencia de overfitting
precio por unidad de área
- Evaluemos la importancia de una variable por su impacto en la predicción:
## [1] 8.548125Aporte de la variable x1 a la predicción:
## [1] -2.989677Este valor me determina la importancia de distorcionar una variable predictora y se identifica cual aporta más al crecimiento del error cuadrático medio,distorción de la predicción
- Calculamos la importancia de las ocho variables:
## [1] -3.0724300 -9.8746673 -17.6612211 -1.7915234 0.1372202 -2.6262241
## [7] -6.3976581 0.5863169- Graficar importancia de la variables datos perturbados
Estabilizar importancia media
## [1] -3.0144705 -9.8296065 -17.7590190 -1.6213060 0.1443075 -2.7041869
## [7] -6.2308683 0.1618022- Importancia por coeficientes
- importancia por estadístico t
La variable más importante para la predicción del precio por unidad de área es \(x3\): distancial al MRT (transporte masivo) más cercano en metros
Regresión lineal (elastic net)
- Definición de los parámetros de optimización:
hiperparametros <- expand.grid(
alpha = (1:100)/100,
lambda = 10^(-5:2)
)
control_optimizacion <- caret::trainControl(method = "LGOCV", number = 10,p=0.2)- optimización Modelo:
m02 <- train(
y~-1+.,data = df_tr,
method = "glmnet",
trControl = control_optimizacion,
tuneGrid = hiperparametros
)- Mejor modelo elastic net con el RMSE más bajo
Mejor rmse óptimizado del modelo elastic net:
## [1] "RSMS de entrenamiento 9.15709186946298"## [1] "RSMS de validación 8.02864462129587"Se observa que el modelo tuvo un mejor desempeño en la base de validación, se puede concluir que no hay presencia de overfitting
Optimización de árboles CART
Optimizacion 1 modelo 1
# M1 #####
set.seed(123)
model1 <- train(
y ~., data = train, method = "rpart",
trControl = trainControl("cv", number = 20),
tuneLength = 10
)
#plot(model)
print("mejores parametros")## [1] "mejores parametros"paste(print(model1$bestTune))## cp
## 2 0.006874159## [1] "0.00687415868359284"rmse modelo 1 arból optimizado
tree <- rpart::rpart(y ~ .,
data = train,
method="anova",
cp = as.numeric(model1$bestTune))
#train
pred_tr <- predict(tree, X_tr)
M1_arb_Op_rmse_train <- RMSE(pred_tr, Y_tr$y)
paste("RSMS de entrenamiento", M1_arb_Op_rmse_train)## [1] "RSMS de entrenamiento 7.00210342571251"#val
pred <- predict(tree, X_vl)
M1_arb_Op_rmse_val <- RMSE(pred, Y_vl$y)
paste("RSMS de validación", M1_arb_Op_rmse_val)## [1] "RSMS de validación 9.00577220513169"Se observa que el modelo tuvo un peor desempeño en la base de validación, se puede concluir que hay presencia de overfitting
Planteamos una segunda alternativa de optimizacion para encontrar el mejor modelo, que contiene un grid con de 4235 combinaciones de parametros maxdepth y mincriterion.
# M2####
set.seed(123)
model2 <- train(
y ~., data = train, method = "ctree2",
trControl = trainControl("cv", number = 10),
tuneGrid = expand.grid(maxdepth = c(1:35),
mincriterion = seq(0.85, 0.97, 0.001) )
)
print("mejores parametros para maxdepth y mincriterion")## [1] "mejores parametros para maxdepth y mincriterion"paste(print(model2$bestTune))## maxdepth mincriterion
## 303 3 0.91## [1] "3" "0.91"rmse modelo 2 arból optimizado
## [1] "RSMS de entrenamiento 7.88264081197781"## [1] "RSMS de validación 8.83101760719523"Se observa que el modelo tuvo un peor desempeño en la base de validación, se puede concluir que hay presencia de overfitting
Ahora, planteamos una segunda alternativa de optimizacion para encontrar el mejor modelo, que contiene un grid con de 2001 con el parámetro cp
# M3 ####
numFolds = trainControl( method = "cv", number = 10 )
cpGrid = expand.grid( .cp = seq(0.01,0.03,0.00001))
mod3 <- train(y ~ ., data = train,
method = "rpart",
trControl = numFolds,
tuneGrid = cpGrid )
print("mejores parametros")## [1] "mejores parametros"paste(print(mod3$bestTune))## cp
## 7 0.01006## [1] "0.01006"rmse modelo 3 arból optimizado
## [1] "RSMS de entrenamiento 7.00210342571251"## [1] "RSMS de validación 9.00577220513169"Se observa que el modelo tuvo un peor desempeño en la base de validación, se puede concluir que hay presencia de overfitting
Finalmente se recopilan los distintas salidas de RMSE de los 3 modelos optimizados
## Modelo1 Modelo2 Modelo3
## RSME_Train 7.002103 7.882641 7.002103
## RSME_Val 9.005772 8.831018 9.005772Apesar que los tres metodologías implementadas contienen overfitting se puede concluir que el mejor modelo es optimización del modelo 2, ya que da un RMSE en el set de validación de:
## [1] 8.831018Ensambles de árboles
bosques aleatorios
##
## Call:
## randomForest(formula = y ~ ., data = df_tr, importance = TRUE)
## Type of random forest: regression
## Number of trees: 500
## No. of variables tried at each split: 2
##
## Mean of squared residuals: 59.97074
## % Var explained: 69.99
se concluye que la variable mas importante es \(X3\):distancia al MRT (transporte masivo) más cercano en metros, seguido de la misma variable elevada al cuadrado,en tercer lugar la variable \(X4\): número de tiendas de conveniencia en el vecindario (entero)
El rmse del Random Forest mide cuanto se degrada la predicción del modelo en este caso fue de:
## [1] "RSMS de entrenamiento 7.74407771799045"## [1] "RSMS de validación 7.26003121425599"Se observa que el modelo tuvo un mejor desempeño en la base de validación, se puede concluir que no hay presencia de overfitting
XGBoost
## [1] train-rmse:28.889256
## [2] train-rmse:21.122471
## [3] train-rmse:15.914741
## [4] train-rmse:12.490872
## [5] train-rmse:10.300000
## [6] train-rmse:9.003424
## [7] train-rmse:8.219088
## [8] train-rmse:7.715497
## [9] train-rmse:7.393971
## [10] train-rmse:7.200536El rmse del XGboost mide cuanto se degrada la predicción del modelo en este caso fue de:
## [1] "RSMS de entrenamiento 7.200536"## [1] "RSMS de validación 7.98919336169158"Se observa que el modelo tuvo un peor desempeño en la base de validación, se puede concluir que hay presencia de overfitting
optimización XGboost
- Definición de los parámetros de optimización:
gbmGrid2 <- expand.grid(interaction.depth = c(1, 5, 9),
n.trees = (1:30)*50,
shrinkage = c(0.0001, 0.001, 0.01, 0.1, 0.2, 0.3),
n.minobsinnode = 20)
fitControl <- trainControl(## 10-fold CV
method = "repeatedcv",
number = 10,
repeats = 10)- optimización Modelo XGboost:
set.seed(825)
gbmFit2 <- train(y ~ ., data = df_tr,
method = "gbm",
trControl = fitControl,
verbose = FALSE,
tuneGrid = gbmGrid2)- hiperparámetros mejor modelo XGboost con el RMSE más bajo
Mejor rmse óptimizado del modelo XGboost es:
## [1] "RSMS de entrenamiento 7.43678346920549"## [1] "RSMS de validación 8.00634164537246"Se observa que el modelo tuvo un peor desempeño en la base de validación, se puede concluir que hay presencia de overfitting
Máquinas de soporte vectorial
m_svm <- svm(y ~ ., data = df_tr,type = 'eps-regression')
y_svm_pred <- predict(m_svm)
error_svm <- Y_tr$y - y_svm_pred
rmse_svm_tr <- sqrt(mean(error_svm^2))
#(RMSEsvm=rmse(pred_valid_svm,test$y))el modelo de máquinas de soporte vectorial tiene error cuadrático medio de:
## [1] "RSMS de entrenamiento 7.77926557823923"## [1] "RSMS de validación 7.93779143809829"Se observa que el modelo tuvo un peor desempeño en la base de validación, se puede concluir que hay presencia de overfitting
Redes neuronales
El rmse del Red neuronal que mide cuanto se degrada la predicción del modelo en este caso fue de:
## [1] "RSMS de entrenamiento 5.60482633894204"## [1] "RSMS de validación 19.6932238351513"Se observa que el modelo tuvo un peor desempeño en la base de validación, se puede concluir que hay presencia de overfitting
Preguntas adicionales
1. ¿Qué variables tienen el mayor impacto en el precio de la vivienda?
Las variables \(X3\):distancia al MRT (transporte masivo) más cercano en metros en el análisis de correlación se evidencia una dependecia lineal medianamente al precio por unidad de área de la vivienda en Taiwan.en los modelos regresión lineal se identificó como la variable más importante, al igual que en el modelo de bosque aleatorios.
¿Cómo aporta cada modelo al conocimiento de este impacto?
En la identificación de la importancia de las variables al modelo, el impacto se nota en el error cuadratico medio(RMSE) al hacer distorsión en las variables más importantes como lo fue \(x3\) era significativo el aporte al error de predicción del modelo.
2. ¿Cuál es el mejor modelo entre los usados para resolver este problema?
El mejor modelo con un RMSE de:
## [1] 5.604826precio por unidad de área de la vivienda en Taiwan y fue la red neuronal con el menor error cuadrático medio.
¿Qué criterios se pueden utilizar para responder a esta pregunta?
para este caso por los modelos de: Regresión lineal (clásica y elastic net),Ensambles de árboles: bosques aleatorios y XGBoost,Máquinas de soporte vectorial, Redes neuronales, se busca una medida del error de ajuste de la predicción que se pueda calcular en todos los modelos como lo es el error cuadrático medio(RMSE) y se pueda comparar entre todos.El de menor RMSE será el que tiene menor error al predecir el precio por unidad de área de la vivienda en Taiwan.
De acuerdo a los siguientes resultados:
## Modelo1_rmse_regresion Modelo2_rmse_elastic_net_opt
## RSME Train 8.299983 9.157092
## RSME Validación 8.268585 8.028645
## Modelo3_rmse_arbol_opt Modelo4_rmse_randomForest
## RSME Train 7.882641 7.744078
## RSME Validación 8.831018 7.260031
## Modelo5_rmse_XGBoost_opt Modelo6_rmse_svm Modelo7_rmse_nnet
## RSME Train 7.436783 7.779266 5.604826
## RSME Validación 8.006342 7.937791 19.693224Despues de validar el mejor RMSE para el set de validación es el modelo de Random Forest, ya que nos da un RSME de 7.198296, y se observa que no hay presencia de overfitting.