Tengku-Hanis01

Aufgabe

Bearbeiten Sie diese Fallstudie von Tengku Hanis!

Lösung

Die folgende Lösung basiert auf der oben angegebenen Fallstudie.

Pakete laden:

library(tidyverse)

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.1.3     ✔ readr     2.1.4
✔ forcats   1.0.0     ✔ stringr   1.5.0
✔ ggplot2   3.4.4     ✔ tibble    3.2.1
✔ lubridate 1.9.3     ✔ tidyr     1.3.0
✔ 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.0
✔ dials        1.2.0     ✔ tune         1.1.2
✔ infer        1.0.5     ✔ workflows    1.1.3
✔ modeldata    1.2.0     ✔ workflowsets 1.0.1
✔ parsnip      1.1.1     ✔ yardstick    1.2.0
✔ recipes      1.0.8     
── 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()
• Dig deeper into tidy modeling with R at https://www.tmwr.org

library(finetune)

Daten importieren:

data(income, package = "kernlab")

Datensatz vereinfachen:

set.seed(2021)
income2 <- 
  income %>% 
  filter(INCOME == "[75.000-" | INCOME == "[50.000-75.000)") %>% 
  slice_sample(n = 600) %>% 
  mutate(INCOME = fct_drop(INCOME), 
         INCOME = fct_recode(INCOME, 
                             rich = "[75.000-",
                             less_rich = "[50.000-75.000)"), 
         INCOME = factor(INCOME, ordered = F)) %>% 
  mutate(across(-INCOME, fct_drop))

Check:

DataExplorer::plot_missing(income)

{DataExplorer} sieht nach einem nützlichen Paket aus. Check it out hier!

Daten aufteilen (“Spending our data budget”):

set.seed(2021)
dat_index <- initial_split(income2, strata = INCOME)
dat_train <- training(dat_index)
dat_test <- testing(dat_index)

Kreuzvalidierung:

set.seed(2021)
dat_cv <- vfold_cv(dat_train, v = 10, repeats = 1, strata = INCOME)

Rezept:

dat_rec <- 
  recipe(INCOME ~ ., data = dat_train) %>% 
  step_impute_mode(all_predictors()) %>% 
  step_ordinalscore(AGE, EDUCATION, AREA, HOUSEHOLD.SIZE, UNDER18)

Als Modell (im engeren Sinne) nutzen wir ein Random-Forest-Modell:

rf_mod <- 
  rand_forest(mtry = tune(),
              trees = tune(),
              min_n = tune()) %>% 
  set_mode("classification") %>% 
  set_engine("ranger")

Wie man sieht, geben wir 3 Tuningparameter an.

Modell und Rezept zum Workflow zusammenfassen:

rf_wf <- 
  workflow() %>% 
  add_recipe(dat_rec) %>% 
  add_model(rf_mod)

Tuning Grids definieren:

Wichtig ist, dass wir genau die Parameter angeben im Grid, die wir auch zum Tunen getaggt haben. Das kann man händisch erledigen:

# Regular grid:
reg_grid <- grid_regular(mtry(c(1, 13)), 
                         trees(), 
                         min_n(), 
                         levels = 3)

# Random grid mit 100 Kandidaten:
rand_grid <- grid_random(mtry(c(1, 13)), 
                         trees(), 
                         min_n(), 
                         size = 100)

Wir speichern die Vorhersagen aller Folds im Train-Sample, um die Modellgüte im Train- bzw. Validierungssample anschauen zu können:

ctrl <- control_grid(save_pred = T,
                     extract = extract_model)
measure <- metric_set(roc_auc)

Außerdem haben wir als Gütemaß roc_auc definiert.

In der Fallstudie wurde noch extract = extract_model bei control_grid() ergänzt. Das lassen wir der Einfachheit halber mal weg.

Parallelisieren auf mehreren Kernen, um Rechenzeit zu sparen:

library(doParallel)

Loading required package: foreach

Attaching package: 'foreach'

The following objects are masked from 'package:purrr':

    accumulate, when

Loading required package: iterators

Loading required package: parallel

# Create a cluster object and then register: 
cl <- makePSOCKcluster(4)
registerDoParallel(cl)

Wie viele CPUs hat mein Computer?

detectCores(logical = FALSE)

[1] 4

Jetzt geht’s ab: Tuning und Fitting!

Hier das “reguläre Gitter” an Tuningkandidaten:

set.seed(2021)
tune_regular <- 
  rf_wf %>% 
  tune_grid(
    resamples = dat_cv, 
    grid = reg_grid,         
    control = ctrl, 
    metrics = measure)

stopCluster(cl)

Die Modellgüte im Vergleich zwischen den Tuning-Kandidaten kann man sich schön ausgeben lassen:

autoplot(tune_regular)

Geht aber nur, wenn man oben gesagt hat, dass man die Predictions speichern möchte.

Welche Kandidatin war am besten:

show_best(tune_regular)

# A tibble: 5 × 9
   mtry trees min_n .metric .estimator  mean     n std_err .config              
  <int> <int> <int> <chr>   <chr>      <dbl> <int>   <dbl> <chr>                
1     7  2000     2 roc_auc binary     0.690    10  0.0178 Preprocessor1_Model08
2     7  2000    40 roc_auc binary     0.689    10  0.0184 Preprocessor1_Model26
3     7  1000     2 roc_auc binary     0.688    10  0.0162 Preprocessor1_Model05
4    13  1000    21 roc_auc binary     0.687    10  0.0155 Preprocessor1_Model15
5     7  2000    21 roc_auc binary     0.687    10  0.0161 Preprocessor1_Model17

So kann man sich die beste Kandidatin anschauen:

show_best(tune_regular) %>% 
  arrange(-mean) %>% 
  slice(1)

# A tibble: 1 × 9
   mtry trees min_n .metric .estimator  mean     n std_err .config              
  <int> <int> <int> <chr>   <chr>      <dbl> <int>   <dbl> <chr>                
1     7  2000     2 roc_auc binary     0.690    10  0.0178 Preprocessor1_Model08

Aber man kann sich auch von Tidymodels einfach die beste Kandidatin sagen lassen:

best_rf <-
  select_best(tune_regular, "roc_auc")

Auf dieser Basis können wir jetzt den Workflow finalisieren, also die Tuningparameter einfüllen:

final_wf <- 
  rf_wf %>% 
  finalize_workflow(best_rf)
final_wf

══ Workflow ════════════════════════════════════════════════════════════════════
Preprocessor: Recipe
Model: rand_forest()

── Preprocessor ────────────────────────────────────────────────────────────────
2 Recipe Steps

• step_impute_mode()
• step_ordinalscore()

── Model ───────────────────────────────────────────────────────────────────────
Random Forest Model Specification (classification)

Main Arguments:
  mtry = 7
  trees = 2000
  min_n = 2

Computational engine: ranger 

Und mit diesen Werten den ganzen Train-Datensatz fitten:

test_fit <- 
  final_wf %>%
  last_fit(dat_index) 

Wie gut ist das jetzt?

test_fit %>%
  collect_metrics()

# A tibble: 2 × 4
  .metric  .estimator .estimate .config             
  <chr>    <chr>          <dbl> <chr>               
1 accuracy binary         0.576 Preprocessor1_Model1
2 roc_auc  binary         0.599 Preprocessor1_Model1

---

Categories:

• tidymodels
• prediction
• yacsda
• statlearning
• trees
• speed
• string