library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.0 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.1
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(tidymodels)
## ── Attaching packages ────────────────────────────────────── tidymodels 1.1.1 ──
## ✔ broom 1.0.5 ✔ rsample 1.2.1
## ✔ dials 1.2.1 ✔ tune 1.1.2
## ✔ infer 1.0.6 ✔ workflows 1.1.4
## ✔ modeldata 1.3.0 ✔ workflowsets 1.0.1
## ✔ parsnip 1.2.0 ✔ yardstick 1.3.0
## ✔ recipes 1.0.10
## ── Conflicts ───────────────────────────────────────── tidymodels_conflicts() ──
## ✖ scales::discard() masks purrr::discard()
## ✖ dplyr::filter() masks stats::filter()
## ✖ recipes::fixed() masks stringr::fixed()
## ✖ dplyr::lag() masks stats::lag()
## ✖ yardstick::spec() masks readr::spec()
## ✖ recipes::step() masks stats::step()
## • Search for functions across packages at https://www.tidymodels.org/find/
El conjunto de datos Air Quality del repositorio de UCI (https://archive.ics.uci.edu/dataset/360/air+quality) contiene mediciones de la calidad del aire tomadas en una estación de monitoreo de calidad del aire en una zona suburbana de Milán, Italia, en el año 2004. Este dataset es ampliamente utilizado para experimentos de análisis de series temporales, modelos de predicción y problemas de regresión.
Variables del dataset:
Date (Fecha): La fecha de la medición en el formato DD/MM/YYYY.
Time (Hora): La hora de la medición en el formato de 24 horas HH.MM.SS.
CO(GT): Concentración de monóxido de carbono (CO) en mg/m³. Este valor fue registrado por un sensor de gas (GT).
PT08.S1(CO): Señal del sensor 1, que es sensible al CO. Este sensor responde principalmente al CO y está representado en unidades arbitrarias de la señal del sensor.
NMHC(GT): Hidrocarburos no metánicos (NMHC) en microgramos por metro cúbico (µg/m³).
C6H6(GT): Concentración de benceno en microgramos por metro cúbico (µg/m³).
PT08.S2(NMHC): Señal del sensor 2, que responde principalmente a hidrocarburos no metánicos.
NOx(GT): Concentración de óxidos de nitrógeno (NOx) en partes por mil millones (ppb).
PT08.S3(NOx): Señal del sensor 3, que responde principalmente a los óxidos de nitrógeno.
NO2(GT): Concentración de dióxido de nitrógeno (NO₂) en partes por mil millones (ppb).
PT08.S4(NO2): Señal del sensor 4, sensible al dióxido de nitrógeno.
PT08.S5(O3): Señal del sensor 5, sensible al ozono (O₃).
T (Temperatura): Temperatura ambiental en grados Celsius (°C).
RH (Humedad Relativa): Humedad relativa en porcentaje (%).
AH (Humedad Absoluta): Humedad absoluta en gramos por metro cúbico (g/m³).
data <- read.csv("https://raw.githubusercontent.com/data-datum/datasets_ml/refs/heads/main/AirQualityUCI.csv", sep = ";", dec = ",", na.strings = c("NA", "-200"))
# Verificar nombres de columnas
names(data)
## [1] "Date" "Time" "CO.GT." "PT08.S1.CO."
## [5] "NMHC.GT." "C6H6.GT." "PT08.S2.NMHC." "NOx.GT."
## [9] "PT08.S3.NOx." "NO2.GT." "PT08.S4.NO2." "PT08.S5.O3."
## [13] "T" "RH" "AH" "X"
## [17] "X.1"
# Seleccionar las variables de interés y eliminar valores faltantes
data <- data %>%
select(PT08.S5.O3., T, RH, NOx.GT., NO2.GT., NMHC.GT.) %>%
drop_na()
set.seed(123) # Para reproducibilidad
p_split <- initial_split(data, prop = 0.8, strata = PT08.S5.O3.)
train_data <- training(p_split)
test_data <- testing(p_split)
Los datos de TRAIN voy a dividirlos en 3-folds
p_folds <- vfold_cv(train_data, v=3)
# Crear la receta de preprocesamiento
rf_recipe <- recipe(PT08.S5.O3. ~ T + RH + NOx.GT. + NO2.GT. + NMHC.GT., data = train_data) %>%
step_normalize(all_predictors()) # Estandarizar las variables predictoras
rf_tune <- rand_forest(
mtry = tune(),
trees = 1000,
min_n = tune()
) %>%
set_mode("regression") %>%
set_engine("ranger", importance = "permutation")
tune_wf <- workflow() %>%
add_recipe(rf_recipe) %>%
add_model(rf_tune)
# Define the metrics, ensuring "mae" is included
metrics <- metric_set(mae, rmse, mape)
set.seed(8577)
doParallel::registerDoParallel()
ranger_tune <-tune_grid(tune_wf,
resamples = p_folds,
metrics=metrics,
grid = 11
)
## i Creating pre-processing data to finalize unknown parameter: mtry
show_best(ranger_tune, metric = "mae")
## # A tibble: 5 × 8
## mtry min_n .metric .estimator mean n std_err .config
## <int> <int> <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 5 6 mae standard 104. 3 1.37 Preprocessor1_Model06
## 2 4 10 mae standard 106. 3 1.36 Preprocessor1_Model05
## 3 2 4 mae standard 106. 3 0.945 Preprocessor1_Model10
## 4 3 18 mae standard 110. 3 1.50 Preprocessor1_Model01
## 5 4 24 mae standard 111. 3 1.67 Preprocessor1_Model03
ranger_tune %>%
collect_metrics() %>%
filter(.metric == "mae") %>%
dplyr::select(mean, min_n, mtry) %>%
pivot_longer(min_n:mtry,
values_to = "value",
names_to = "parameter"
) %>%
ggplot(aes(value, mean, color = parameter)) +
geom_point(show.legend = FALSE) +
facet_wrap(~parameter, scales = "free_x") +
labs(x = NULL, y = "mae")
show_best(ranger_tune, "mae")
## # A tibble: 5 × 8
## mtry min_n .metric .estimator mean n std_err .config
## <int> <int> <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 5 6 mae standard 104. 3 1.37 Preprocessor1_Model06
## 2 4 10 mae standard 106. 3 1.36 Preprocessor1_Model05
## 3 2 4 mae standard 106. 3 0.945 Preprocessor1_Model10
## 4 3 18 mae standard 110. 3 1.50 Preprocessor1_Model01
## 5 4 24 mae standard 111. 3 1.67 Preprocessor1_Model03
best_mae <- select_best(ranger_tune, "mae")
best_mae
## # A tibble: 1 × 3
## mtry min_n .config
## <int> <int> <chr>
## 1 5 6 Preprocessor1_Model06
final_rf <- finalize_workflow(
tune_wf,
best_mae
)
final_rf
## ══ Workflow ════════════════════════════════════════════════════════════════════
## Preprocessor: Recipe
## Model: rand_forest()
##
## ── Preprocessor ────────────────────────────────────────────────────────────────
## 1 Recipe Step
##
## • step_normalize()
##
## ── Model ───────────────────────────────────────────────────────────────────────
## Random Forest Model Specification (regression)
##
## Main Arguments:
## mtry = 5
## trees = 1000
## min_n = 6
##
## Engine-Specific Arguments:
## importance = permutation
##
## Computational engine: ranger
# Define the metrics, ensuring "mae" is included
metrics <- metric_set(mae, rmse, mape)
final_res <- last_fit(final_rf, p_split, metrics=metrics)
## Warning: package 'ranger' was built under R version 4.3.3
collect_metrics(final_res)
## # A tibble: 3 × 4
## .metric .estimator .estimate .config
## <chr> <chr> <dbl> <chr>
## 1 mae standard 102. Preprocessor1_Model1
## 2 rmse standard 127. Preprocessor1_Model1
## 3 mape standard 11.6 Preprocessor1_Model1
# Obtener los valores verdaderos y predichos
predictions <- final_res %>%
collect_predictions()
# Crear el gráfico de valores reales vs valores predichos
ggplot(predictions, aes(x = PT08.S5.O3., y = .pred)) +
geom_point(color = "blue", alpha = 0.6) + # Puntos azules para las observaciones
geom_abline(slope = 1, intercept = 0, color = "red", linetype = "dashed") + # Línea de igualdad
labs(
title = "Predicted vs Real Values of PT08.S5(O3)",
x = "Real Values (PT08.S5(O3))",
y = "Predicted Values"
) +
theme_minimal()
library(vip)
# Fit the model to the entire training data to get the final model
final_fit <- final_rf %>%
fit(data = train_data)
# Plot variable importance
vip(pull_workflow_fit(final_fit)$fit, geom = "point")