# Instalacija i učitavanje potrebnih paketa
library(rpart)
library(caTools)
library(caret)
## Loading required package: ggplot2
## Loading required package: lattice
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(Boruta)
library(cvms)
library(tibble)
library(yardstick)
## 
## Attaching package: 'yardstick'
## The following objects are masked from 'package:caret':
## 
##     precision, recall, sensitivity, specificity
library(ggplot2)
library(rpart.plot)
library(randomForest)
## randomForest 4.7-1.2
## Type rfNews() to see new features/changes/bug fixes.
## 
## Attaching package: 'randomForest'
## The following object is masked from 'package:dplyr':
## 
##     combine
## The following object is masked from 'package:ggplot2':
## 
##     margin
library(gbm)
## Loaded gbm 2.2.2
## This version of gbm is no longer under development. Consider transitioning to gbm3, https://github.com/gbm-developers/gbm3
library(DALEX)
## Welcome to DALEX (version: 2.4.3).
## Find examples and detailed introduction at: http://ema.drwhy.ai/
## Additional features will be available after installation of: ggpubr.
## Use 'install_dependencies()' to get all suggested dependencies
## 
## Attaching package: 'DALEX'
## The following object is masked from 'package:dplyr':
## 
##     explain
library(DALEXtra)

1. Učitavanje podataka

Učitat ćemo skup podataka iris pomoću funkcije data(iris), koji sadrži podatke o tri vrste cvijeta irisa (setosa, versicolor i virginica) te četiri značajke: duljina i širina čašičnog lista (sepal length i sepal width) i latice (petal length i petal width).

Funkcija head(iris) prikazuje prvih nekoliko redaka skupa podataka.

Zatim postavljamo slučajnost pomoću set.seed(123) kako bi rezultati bili reproducibilni, odnosno kako bi slučajno odabrani uzorci uvijek bili isti.

Podaci se dijele na trening i testni skup korištenjem funkcije sample(), koja nasumično odabire 70% redaka skupa podataka za trening, dok preostalih 30% čini testni skup.

Trening podaci se spremaju u varijablu train_data, a podaci za testiranje u test_data.

data(iris)
head(iris)
##   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1          5.1         3.5          1.4         0.2  setosa
## 2          4.9         3.0          1.4         0.2  setosa
## 3          4.7         3.2          1.3         0.2  setosa
## 4          4.6         3.1          1.5         0.2  setosa
## 5          5.0         3.6          1.4         0.2  setosa
## 6          5.4         3.9          1.7         0.4  setosa
set.seed(123)  # Reproducibilnost
train_indices <- sample(1:nrow(iris), 0.7 * nrow(iris))  # 70% za trening i 30% za testiranje
train_data <- iris[train_indices, ]
test_data <- iris[-train_indices, ]

2. Modeli

Kreiramo tri različita modela strojnog učenja za predviđanje vrste cvijeta na temelju ostalih značajki iz skupa podataka.

Prvi model je Decision Tree, implementiran pomoću funkcije rpart(), koja gradi hijerarhijski model odlučivanja temeljen na kriteriju particioniranja podataka.

Drugi model je Random Forest, stvoren korištenjem funkcije randomForest(), koja kombinira više stabala odlučivanja (100 stabala u ovom slučaju) kako bi poboljšala točnost i otpornost na prekomjerno prilagođavanje (overfitting).

Treći model je Gradient Boosting, implementiran pomoću funkcije gbm(), koja koristi postupak slijednog treniranja slabih učenika (stabala odlučivanja) kako bi smanjila grešku.

Parametri modela uključuju korištenje distribucije multinomial za klasifikaciju više klasa, 100 stabala i dubinu interakcije postavljenu na 3.

# 1. Decision Tree
dt_model <- rpart(Species ~ ., data = train_data)

# 2. Random Forest
rf_model <- randomForest(Species ~ ., data = train_data, ntree = 100)

# 3. Gradient Boosting
gbm_model <- gbm(Species ~ ., data = train_data, distribution = "multinomial", 
                 n.trees = 100, interaction.depth = 3)
## Warning: Setting `distribution = "multinomial"` is ill-advised as it is
## currently broken. It exists only for backwards compatibility. Use at your own
## risk.

2.1. Decision tree

Ovdje koristimo model stabla odluke kako bi analizirali značajke i predviđanja na skupu podataka iris.

Brafički prikaz stabla odluke pomoću funkcije rpart.plot(dt_model) omogućuje vizualizaciju pravaca donošenja odluka unutar modela. Ovo stablo pokazuje kako se različite značajke, poput duljine i širine čašičnog lista, koriste za klasifikaciju različitih vrsta cvijeta (Species).

Zatim, uz pomoć funkcije Boruta, analiziramo značajke koje najviše utječu na model.

Boruta je algoritam za selekciju značajki koji prepoznaje koje varijable imaju značajan utjecaj na predviđanje ciljne varijable. Možemo vidjeti importance (važnost) značajki koja nam pokazuje koje značajke su najvažnije za donošenje ispravnih odluka u modelu. Rezultati su filtrirani tako da prikazujemo samo one značajke koje nisu označene kao “odbačene” (Rejected), a zatim su sortirani prema njihovoj važnosti. Ovaj postupak omogućuje nam da saznamo koje varijable imaju najveći utjecaj na predviđanje vrste cvijeta.

Nakon analize značajki, predikcija na testnom skupu podataka se koristi za ocjenu performansi modela. Prediktivne vrijednosti generirane iz modela stabla odluke uspoređuju se s pravim vrijednostima u testnom skupu, što rezultira matricom zabune (confusion matrix). Ova matrica daje uvid u to koliko je model bio točan u klasifikaciji svake od vrsta cvijeta, te prikazuje broj točnih i netočnih predikcija za svaku klasu. Na temelju matrice zabune, koristi se funkcija plot_confusion_matrix() za vizualizaciju tih rezultata kako bi lakše napravili interpretaciju klasifikacije različitih vrsta.

# Grafička prezentacija modela
rpart.plot(dt_model)  # Decision Tree

Možemo vidjeti koliko je klasifikacija bila točna (u trening skupu) promatranjem doljnjih čvorova. Setosa je ispravno klasificirana svaki put, versicolor je u 11% slučajeva bila pogrešno klasificirana kao virginica, a virginica jeu 3% slučajeva bila pogrešno klasificirana kao versicolor.

boruta_output <- Boruta(Species ~ ., data = train_data, doTrace = 0)
rough_fix_mod <- TentativeRoughFix(boruta_output)
## Warning in TentativeRoughFix(boruta_output): There are no Tentative attributes!
## Returning original object.
boruta_signif <- getSelectedAttributes(rough_fix_mod)
importances <- attStats(rough_fix_mod)
importances <- importances[importances$decision != "Rejected", c("meanImp", "decision")]
importances[order(-importances$meanImp), ]
##                meanImp  decision
## Petal.Width  29.123083 Confirmed
## Petal.Length 28.371186 Confirmed
## Sepal.Length 15.150757 Confirmed
## Sepal.Width   9.482064 Confirmed
plot(boruta_output, ces.axis = 0.7, las = 2, xlab = "", main = "Feature importance")

Važnost značajki prepoznajemo po boji, zeleno obojane značajke su važne. Crvena bi predstavljala nevažne značajke, a plava predstavlja varijablu koju boruta koristi za određivanje značaja. Što je box plot viši na grafu prema y-osi, to je veća važnost. Možemo vidjeti da su Petal.Length i Petal.Width najznačajnije varijable, te ih prate Sepal.Length i Sepal.Width, ali su otprilike upola važne.

predict_dt <- predict(dt_model, newdata = test_data, type = "class")

confusionMatrix(test_data$Species, predict_dt)
## Confusion Matrix and Statistics
## 
##             Reference
## Prediction   setosa versicolor virginica
##   setosa         14          0         0
##   versicolor      0         18         0
##   virginica       0          1        12
## 
## Overall Statistics
##                                           
##                Accuracy : 0.9778          
##                  95% CI : (0.8823, 0.9994)
##     No Information Rate : 0.4222          
##     P-Value [Acc > NIR] : 8.826e-16       
##                                           
##                   Kappa : 0.9662          
##                                           
##  Mcnemar's Test P-Value : NA              
## 
## Statistics by Class:
## 
##                      Class: setosa Class: versicolor Class: virginica
## Sensitivity                 1.0000            0.9474           1.0000
## Specificity                 1.0000            1.0000           0.9697
## Pos Pred Value              1.0000            1.0000           0.9231
## Neg Pred Value              1.0000            0.9630           1.0000
## Prevalence                  0.3111            0.4222           0.2667
## Detection Rate              0.3111            0.4000           0.2667
## Detection Prevalence        0.3111            0.4000           0.2889
## Balanced Accuracy           1.0000            0.9737           0.9848
cm <- confusionMatrix(test_data$Species, predict_dt)
cfm <- as_tibble(cm$table)
plot_confusion_matrix(cfm, target_col = "Reference", prediction_col = "Prediction", counts_col = "n")
## Warning in plot_confusion_matrix(cfm, target_col = "Reference", prediction_col
## = "Prediction", : 'ggimage' is missing. Will not plot arrows and zero-shading.
## Warning in plot_confusion_matrix(cfm, target_col = "Reference", prediction_col
## = "Prediction", : 'rsvg' is missing. Will not plot arrows and zero-shading.

Pomoću matrice zabune možemo vidjeti “heatmapu” krivih klasifikacija, te ako pogledamo statistike, naš model postiže točnost od 97.78%

2.2. Random forest

Prvo smo napravili grafički prikaz važnosti značajki modela korištenjem funkcije varImpPlot(rf_model). Ovaj grafikon prikazuje koje su značajke najvažnije za donošenje odluka unutar Random Forest modela. Viša važnost označava veći utjecaj značajke na predviđanje ciljne varijable, odnosno vrste cvijeta.

Sljedeće, izračunali smo preciznost modela kao komplement od srednje stope greške na svim stablima korištenjem rf_model$err.rate. Preciznost se izračunava kao 1 - srednja stopa greške, što daje uvid u ukupnu točnost modela.

Također, prikazali smo stopu grešaka za svako od stabala u Random Forest modelu. Ovdje se koristi err_dfcolMeans() računa prosječnu stopu greške za svaki skup podataka (uključujući i Out-Of-Bag (OOB) greške, koje se koriste za procjenu modela bez potrebe za dodatnim testnim skupom). OOB greška je posebice korisna jer daje uvid u učinkovitost modela bez korištenja testnog skupa podataka.

Zatim smo prikazali grafikon OOB greške za svako stablo u modelu. Korištenjem plot(), grafikon prikazuje kako se OOB greška mijenja s povećanjem broja stabala u Random Forest modelu. Ovaj grafikon pomaže u razumijevanju koliko stabala je potrebno za stabilizaciju greške modela.

Nakon procjene modela, koristimo funkciju predikcija na testnim podacima (predict(rf_model, newdata = test_data)), koja generira predikcije za vrste cvijeta u testnom skupu. Predviđene klase se uspoređuju s stvarnim vrijednostima u matrici zabune, koja pokazuje koliko je model bio točan u klasifikaciji svake vrste cvijeta. Na temelju ove matrice, funkcija plot_confusion_matrix() prikazuje grafički prikaz rezultata.

# Plot feature importance for the Random Forest model
varImpPlot(rf_model, main = "Feature Importance - Random Forest")

Vidimo da je važnost značajki gotovo ista kao i kod stabla odluke.

# Preciznost
1 - mean(rf_model$err.rate)
## [1] 0.9386804
# Stopa grešaka

err_df <- as.data.frame(rf_model$err.rate)
err_df %>%
  select(-"OOB") %>%
  colMeans()
##     setosa versicolor  virginica 
##  0.0000000  0.0875510  0.0966383
# Grafikon OOB greške
oob_error <- rf_model$err.rate[, "OOB"]
plot(oob_error, type = "l", col = "red", lwd = 2,
     xlab = "Broj stabala", ylab = "OOB Greška",
     main = "Out-of-Bag Error")

OOB greška je minimalna kada je broj stabala otprilike 20 te nešto manje od 40.

predict_rf <- predict(rf_model, newdata = test_data, type = "class")

cm_rf <- confusionMatrix(test_data$Species, predict_rf)
cfm <- as_tibble(cm_rf$table)
plot_confusion_matrix(cfm, target_col = "Reference", prediction_col = "Prediction", counts_col = "n")
## Warning in plot_confusion_matrix(cfm, target_col = "Reference", prediction_col
## = "Prediction", : 'ggimage' is missing. Will not plot arrows and zero-shading.
## Warning in plot_confusion_matrix(cfm, target_col = "Reference", prediction_col
## = "Prediction", : 'rsvg' is missing. Will not plot arrows and zero-shading.

2.3. Gradient boosting

Na početku smo napravili grafički prikaz procesa učenja modela pomoću funkcije gbm.perf(gbm_model, method = “OOB”), koja prikazuje učinkovitost modela kroz različite metode procjene (OOB, test ili cross-validation). OOB metoda koristi podatke koji nisu uključeni u treniranje svakog stabla za testiranje modela i pomaže u procjeni performansi bez dodatnog testnog skupa.

Nakon toga, važnost značajki modela je prikazali smo pomoću funkcije summary(gbm_model). Ova funkcija daje uvid u relativnu važnost značajki za donošenje odluka u Gradient Boosting modelu. Rezultati su također prikazani u obliku bar grafa pomoću ggplot2, gdje su značajke rangirane prema njihovoj važnosti.

Nakon toga radimo predviđanje vjerojatnosti za svaku klasu koristeći predict() s tipom response. Ovaj korak generira vjerojatnosti za svaku od klasa (vrsta cvijeta), a zatim se koristi funkcija which.max() kako bi se odabrala klasa s najvećom vjerojatnošću za svaku instancu. To omogućuje pretvaranje predviđenih vjerojatnosti u konkretna predviđanja klasa.

Nakon toga radimo matricu zabune, koja prikazuje broj točnih i netočnih predikcija za svaku od vrsta cvijeta.

# Grafikon treniranja
gbm.perf(gbm_model, method = "OOB") 
## OOB generally underestimates the optimal number of iterations although predictive performance is reasonably competitive. Using cv_folds>1 when calling gbm usually results in improved predictive performance.

## [1] 14
## attr(,"smoother")
## Call:
## loess(formula = object$oobag.improve ~ x, enp.target = min(max(4, 
##     length(x)/10), 50))
## 
## Number of Observations: 100 
## Equivalent Number of Parameters: 8.32 
## Residual Standard Error: 0.01602
# Važnost značajki za Gradient Boosting
summary(gbm_model, main = "Važnost značajki")

##                       var   rel.inf
## Petal.Length Petal.Length 69.819432
## Petal.Width   Petal.Width 26.078168
## Sepal.Length Sepal.Length  2.565552
## Sepal.Width   Sepal.Width  1.536847
# Izdvajanje važnosti
importance_gbm <- summary(gbm_model, plotit = FALSE)
importance_gbm_df <- data.frame(Feature = rownames(importance_gbm), Importance = importance_gbm[, "rel.inf"])

# Prikaz važnosti
ggplot(importance_gbm_df, aes(x = reorder(Feature, Importance), y = Importance)) +
  geom_bar(stat = "identity", fill = "darkgreen") +
  coord_flip() +
  labs(title = "Važnost značajki", x = "Značajke", y = "Važnost")

# Predviđanje vjerojatnosti
gbm_probs <- predict(gbm_model, newdata = test_data, n.trees = 100, type = "response")

# Pretvorba u klase
predict_gbm <- colnames(gbm_probs)[apply(gbm_probs, 1, which.max)]

# Confusion matrix
cm_gbm <- confusionMatrix(as.factor(test_data$Species), as.factor(predict_gbm))

# Pretvaranje u tibble, prikaz
cfm <- as_tibble(cm_gbm$table)

# Prikaz cm-a
plot_confusion_matrix(cfm, target_col = "Reference", prediction_col = "Prediction", counts_col = "n")
## Warning in plot_confusion_matrix(cfm, target_col = "Reference", prediction_col
## = "Prediction", : 'ggimage' is missing. Will not plot arrows and zero-shading.
## Warning in plot_confusion_matrix(cfm, target_col = "Reference", prediction_col
## = "Prediction", : 'rsvg' is missing. Will not plot arrows and zero-shading.

3. Predikcije

Prvo, generiramo predikcije za testni skup podataka za svaki model. Predikcije za stablo odluke koriste funkciju predict() s tipom class, te se radi predviđanje klase za svaku instancu u testnom skupu. Za Random Forest model, predikcije se također generiraju pomoću predict(), dok za Gradient Boosting model predikcije predstavljaju vjerojatnosti klasa koje se zatim pretvaraju u konkretne klase pomoću funkcije which.max().

Nakon toga, izračunavamo točnost svakog modela. Točnost se računa kao postotak točnih predikcija u odnosu na sve predikcije. Točnost za svaki model je izračunata kao srednja vrijednost usporedbe predviđenih klasa (dt_pred, rf_pred, gbm_pred) sa stvarnim vrijednostima test_data$Species. Ova točnost daje jasan uvid u to koliko su modeli uspješni u klasifikaciji vrsta cvijeta.

Nakon toga prikazujemo rezultate točnosti za svaki model, gdje se ispisuju specifične vrijednosti točnosti za Decision Tree, Random Forest, i Gradient Boosting modele.

Na kraju, provodimo analizu važnosti prediktora za Random Forest model. Funkcija randomForest::importance(rf_model) prikazuje važnost svake značajke za model.

# Predikcije na testnom skupu
dt_pred <- predict(dt_model, test_data, type = "class")
rf_pred <- predict(rf_model, test_data)
gbm_pred <- predict(gbm_model, test_data, n.trees = 100, type = "response")
gbm_pred <- colnames(gbm_pred)[apply(gbm_pred, 1, which.max)]  # Pretvaranje u klase

# Metričke evaluacije
library(caret)
## 1. Točnost
dt_accuracy <- mean(dt_pred == test_data$Species)
rf_accuracy <- mean(rf_pred == test_data$Species)
gbm_accuracy <- mean(gbm_pred == test_data$Species)

# Rezultati
cat("Točnost:\n")
## Točnost:
cat("Decision Tree: ", dt_accuracy, "\n")
## Decision Tree:  0.9777778
cat("Random Forest: ", rf_accuracy, "\n")
## Random Forest:  0.9777778
cat("Gradient Boosting: ", gbm_accuracy, "\n")
## Gradient Boosting:  0.9777778
# Analiza važnosti prediktora
## Random Forest važnost
rf_importance <- randomForest::importance(rf_model)
print(rf_importance)
##              MeanDecreaseGini
## Sepal.Length        10.180560
## Sepal.Width          2.188665
## Petal.Length        30.365444
## Petal.Width         26.398998

4. DALEX

Ovaj blok koda koristi DALEX paket za interpretaciju modela, te analizu performansi i važnosti značajki.

Prvo koristimo funkciju explain(), kreiramo objekte objašnjenja za svaki od modela (explainer_dt, explainer_rf, i explainer_gbm). Ovi objekti služe nam za objašnjavanje svakog modela na temelju trening podataka. Za svaki model, osim ciljne varijable (Species), svi ostali atributi skupa podataka (train_data[, -5]) se koriste kao ulazni podaci. Ova analiza omogućava bolje razumijevanje što svaki model “nauči” i kako donosi svoje odluke.

Zatim, pomoću funkcije model_performance(), izračunavamo performanse svakog modela. Ovaj korak daje uvid u različite metrikama koje DALEX koristi za procjenu uspješnosti modela, uključujući točnost, preciznost, F1 rezultat i sl. Rezultati su pohranjeni u varijablama mp_dt, mp_rf, i mp_gbm za svaki od modela.

Za daljnju usporedbu performansi modela, koristimo funkciju plot() koja prikazuje metrike za sva tri modela u jednom grafu, kako bi ih lakše usporedili.

Posljednji korak analize uključuje analizu važnosti značajki za Random Forest model koristeći funkciju model_parts(). S ovom funkcijom možemo procijeniti kako svaka značajka doprinosi promjeni u predviđanjima modela. Specifično, koristimo opciju type = “difference”, što znači da se procjenjuje kako promjena u svakoj značajki utječe na razliku u modelovim predviđanjima. Možemo prikazati podatke korištenjem plot(fi_rf).

# DALEX analiza
explainer_dt <- explain(model = dt_model, data = train_data[, -5], y = train_data$Species, label = "Decision Tree")
## Preparation of a new explainer is initiated
##   -> model label       :  Decision Tree 
##   -> data              :  105  rows  4  cols 
##   -> target variable   :  105  values 
##   -> predict function  :  yhat.rpart  will be used (  default  )
##   -> predicted values  :  No value for predict function target column. (  default  )
##   -> model_info        :  package rpart , ver. 4.1.23 , task multiclass (  default  ) 
##   -> predicted values  :  predict function returns multiple columns:  3  (  default  ) 
##   -> residual function :  difference between 1 and probability of true class (  default  )
##   -> residuals         :  numerical, min =  0 , mean =  0.08597039 , max =  0.9705882  
##   A new explainer has been created!
explainer_rf <- explain(model = rf_model, data = train_data[, -5], y = train_data$Species, label = "Random Forest")
## Preparation of a new explainer is initiated
##   -> model label       :  Random Forest 
##   -> data              :  105  rows  4  cols 
##   -> target variable   :  105  values 
##   -> predict function  :  yhat.randomForest  will be used (  default  )
##   -> predicted values  :  No value for predict function target column. (  default  )
##   -> model_info        :  package randomForest , ver. 4.7.1.2 , task multiclass (  default  ) 
##   -> predicted values  :  predict function returns multiple columns:  3  (  default  ) 
##   -> residual function :  difference between 1 and probability of true class (  default  )
##   -> residuals         :  numerical, min =  0 , mean =  0.02866667 , max =  0.39  
##   A new explainer has been created!
explainer_gbm <- explain(model = gbm_model, data = train_data[, -5], y = train_data$Species, label = "Gradient Boosting")
## Preparation of a new explainer is initiated
##   -> model label       :  Gradient Boosting 
##   -> data              :  105  rows  4  cols 
##   -> target variable   :  105  values 
##   -> predict function  :  yhat.gbm  will be used (  default  )
##   -> predicted values  :  No value for predict function target column. (  default  )
##   -> model_info        :  package gbm , ver. 2.2.2 , task multiclass (  default  ) 
##   -> predicted values  :  predict function returns multiple columns:  3  (  default  ) 
##   -> residual function :  difference between 1 and probability of true class (  default  )
##   -> residuals         :  numerical, min =  3.389657e-08 , mean =  0.005540643 , max =  0.2617436  
##   A new explainer has been created!
# Performanse modela
mp_dt <- model_performance(explainer_dt)
mp_rf <- model_performance(explainer_rf)
mp_gbm <- model_performance(explainer_gbm)

# Grafička usporedba
plot(mp_dt, mp_rf, mp_gbm)

Pomoću ovog grafa možemmo vidjeti da je Gradient Boost model najtočniji, dok je Decision tree najmanje precizan.

# Analiza važnosti značajki
fi_rf <- model_parts(explainer_rf, type = "difference")
plot(fi_rf)