Creamos los datos necesarios para los ejercicios

set.seed(1)
n <- 200

# Generamos predictores X1 a X5 de una distribución Uniforme en [0, 1]
X1 <- runif(n, 0, 1)
X2 <- runif(n, 0, 1)
X3 <- runif(n, 0, 1)
X4 <- runif(n, 0, 1)
X5 <- runif(n, 0, 1)

# generamos el error
epsilon <- rnorm(n, mean = 0, sd = 2)


Y <- 10 * sin(pi * X1 * X2) + 20 * (X3 - 0.5)^2 + 10 * X4 + 5 * X5 + epsilon

#dataframe
data <- data.frame(Y, X1, X2, X3, X4, X5)
head(data)
##          Y        X1        X2        X3        X4         X5
## 1 15.30827 0.2655087 0.2675082 0.6588776 0.8142518 0.85868745
## 2 13.37820 0.3721239 0.2186453 0.1850700 0.9287772 0.03443876
## 3 16.10909 0.5728534 0.5167968 0.9543781 0.1474810 0.97099715
## 4 21.35450 0.9082078 0.2689506 0.8978485 0.7498217 0.74511014
## 5 18.18876 0.2016819 0.1811683 0.9436971 0.9756573 0.27325524
## 6 27.26476 0.8983897 0.5185761 0.7236908 0.9747925 0.67710610

Crear muestras de entrenamiento y test

indices <- sample(1:n, size = 150)

# Crear conjuntos de entrenamiento y prueba
df_entrenamiento <- data[indices, ]
df_test <- data[indices, ]


head(df_entrenamiento)
##             Y        X1         X2          X3        X4        X5
## 87  16.092105 0.7111212 0.32149212 0.189779918 0.3860996 0.8031190
## 198 15.066600 0.8405070 0.31177237 0.495064106 0.5161239 0.2821363
## 197  2.943658 0.1103606 0.65982103 0.479835794 0.2372649 0.6359498
## 146  9.404216 0.4525708 0.26929378 0.115342325 0.3816321 0.4158602
## 113 16.837547 0.3567269 0.84135172 0.008128455 0.1202823 0.6947533
## 85  10.707153 0.7570871 0.04795913 0.588419190 0.5972405 0.4927158
head(df_test)
##             Y        X1         X2          X3        X4        X5
## 87  16.092105 0.7111212 0.32149212 0.189779918 0.3860996 0.8031190
## 198 15.066600 0.8405070 0.31177237 0.495064106 0.5161239 0.2821363
## 197  2.943658 0.1103606 0.65982103 0.479835794 0.2372649 0.6359498
## 146  9.404216 0.4525708 0.26929378 0.115342325 0.3816321 0.4158602
## 113 16.837547 0.3567269 0.84135172 0.008128455 0.1202823 0.6947533
## 85  10.707153 0.7570871 0.04795913 0.588419190 0.5972405 0.4927158
#librerias relevantes
library(rpart)    
library(randomForest)  
## randomForest 4.7-1.2
## Type rfNews() to see new features/changes/bug fixes.
library(Metrics) 
# Arbol inicial
tree_model <- rpart(Y ~ ., data = df_entrenamiento, control = rpart.control(mincut = 2, minsize = 4))
print(tree_model)
## n= 150 
## 
## node), split, n, deviance, yval
##       * denotes terminal node
## 
##  1) root 150 4466.30600 14.054460  
##    2) X4< 0.4338369 70 1808.77700 11.228340  
##      4) X1< 0.3364993 20  449.76180  6.584919  
##        8) X3>=0.3476537 13  195.94100  4.346209 *
##        9) X3< 0.3476537 7   67.66724 10.742520 *
##      5) X1>=0.3364993 50  755.29710 13.085710  
##       10) X5< 0.5615711 28  413.40920 11.321690  
##         20) X2< 0.2246898 7   58.50565  7.656643 *
##         21) X2>=0.2246898 21  229.53300 12.543370 *
##       11) X5>=0.5615711 22  143.86660 15.330820 *
##    3) X4>=0.4338369 80 1609.24100 16.527320  
##      6) X2< 0.3623184 33  429.81700 13.931780  
##       12) X2< 0.08600476 9  100.04970 11.270720 *
##       13) X2>=0.08600476 24  242.13700 14.929680  
##         26) X5< 0.3346459 9   79.27230 12.912090 *
##         27) X5>=0.3346459 15  104.24710 16.140230 *
##      7) X2>=0.3623184 47  801.01740 18.349710  
##       14) X5< 0.391192 24  291.51940 15.761700  
##         28) X2< 0.8030876 16  151.86260 14.240070 *
##         29) X2>=0.8030876 8   28.51836 18.804980 *
##       15) X5>=0.391192 23  181.01570 21.050240 *
plot(tree_model)
text(tree_model, pretty = 0)

# Podamos
cv_tree <- prune(tree_model, cp = tree_model$cptable[which.min(tree_model$cptable[, "xerror"]), "CP"])

plot(cv_tree)
text(cv_tree, pretty = 0)

#prediccion
tree_pred <- predict(cv_tree, df_test)
#error cuadratico medio
ecm_tree <- mse(df_test$Y, tree_pred)
print(ecm_tree)
## [1] 8.936528

Primero creamos un árbol con 11 hojas. Procedemos a reducir el número a través de la poda. El número final de hojas es 5.

bagging_model <- randomForest(Y ~ ., data = df_entrenamiento, num_predictores = 5, importance = TRUE)
bagging_pred <- predict(bagging_model, df_test)


importance(bagging_model)
##     %IncMSE IncNodePurity
## X1 21.50778      903.8092
## X2 20.40600      795.4912
## X3 12.12068      606.8724
## X4 30.66373     1032.6205
## X5 18.13428      783.8423
varImpPlot(bagging_model)

# ECM
ecm_bagging <- mse(df_test$Y, bagging_pred)
print(ecm_bagging)
## [1] 2.531317

Los predictores más importantes son X4, X1, X2

rf_model <- randomForest(Y ~ ., data = df_entrenamiento, raiz2predic = 2, importance = TRUE)
rf_pred <- predict(rf_model, df_test)

# Importancia de las variables
importance(rf_model)
##     %IncMSE IncNodePurity
## X1 22.55500      926.9187
## X2 23.08535      842.5105
## X3 11.89933      606.6161
## X4 30.30531     1017.9638
## X5 16.86161      757.2018
varImpPlot(rf_model)

# ECM
ecm_rf <- mse(df_test$Y, rf_pred)
print(ecm_rf)
## [1] 2.441371

Los predictores más importantes son X4, X2, X1

Comparanción de los ECMs

ecm_results <- data.frame(
  Método = c("Árbol de Regresión", "Bagging", "Bosques Aleatorios"),
  ECM = c(ecm_tree, ecm_bagging, ecm_rf)
)
print(ecm_results)
##               Método      ECM
## 1 Árbol de Regresión 8.936528
## 2            Bagging 2.531317
## 3 Bosques Aleatorios 2.441371

Como era de esperarse, el modelo con el error más alto es el arbol de regresión. Bagging y Bosques aleatorios quedan igual entre ellos.