Današnje predavanje se odnosi na pregled najvažnijih principa i tehnika strojnog učenja (ML - Machine Learning). Cilj je dati sveobuhvatan pregled osnovnih koncepata i smjestiti ML metodologiju u širi kontekst kvantitativnih medodoloških pristupa koje već poznajemo. U predavanju ćemo objasniti:

  1. što je strojno učenje i koja mu je svrha,
  2. koje su najčešće korištene metode,
  3. kako odrediti kvalitetu procijenjenog modela.

Što je strojno učenje?

ML je dio umjetne inteligencije (Artificial Inteligence - AI) i spada u područje računalnih znanosti koje koriste statističke tehnike kako bi se računalu omogućilo “učenje” na osnovi podataka bez eksplicitnog programiranja. “Učenje” se pri tome odnosi na progresivno poboljšanje performansi u provedbi nekog zadatka. ML je usko povezano sa područjima rudarenja podataka (Data Mining), prepoznavanja obrazaca (Pattern Recognition), inferencijalnom statistikom i statističkim učenjem, a često je riječ o istim “stvarima” sa različitim nazivom. Neki tipični ML primjeri uključuju diskriminaciju spam elektroničke pošte, automatsko tagiranje Facebook slika, preporuke kupnje, identifikaciju bolesti na osnovi simptoma, predviđanje tečaja…

Strojno učenje uključuje izradu i korištenje algoritama s ciljem učenja na osnovi podataka. Sposobnost stroja (računala) da nešto nauči se odražava u poboljšanju peformansi - outputa pri povećanju informacijskog inputa. Informacijski input se zapravo odnosi na podatke slično data.frame objektima koje smo do sada stalno koristili. Jednostavnije rečeno, ukoliko računalu damo više podataka, ono će više (bolje) naučiti (model). Strojno učenje i statistika (kakvu smo do sada upoznali) se u bitnom ne razlikuju mnogo, dapače među njima ima više sličnosti nego razlika. Razlog za nerazumijevanje i nepristupačnost strojnog učenja (statističarima) često predstavlja specifičan žargon strojnog učenja (ali i naprednije informatičke vještine koje zahtijeva ML). ML ima specifičnu terminologiju koju koriste analitičari i inženjeri strojnog učenja, a često se pritom misli na pojmove koje već poznajemo iz statistike…varijable se u ML tako nazivaju karakteristike (feature, label), procjena je trening modela i sl.

Osnovni pojmovi

Strojno učenje se često dijeli na nadgledano i ne-nadgledano:

  • U nadgledanom (supervised) ML cilj je pronaći funkciju koja pripisuje vrijednost ili klasu novim (neviđenim) opservacijama uz dani skup varijabli (i.e. labeled opservacija).

  • Ne-nadgledano (unsupervised) ML ne zahtijeva definirane varijable (i.e. labeld opservacije) već algoritam samostalno pronalazi strukturu u podatcima.

Najvažniji zadatci koje obavlja ML su:

  • Smanjenje dimenzionalnosti označava pojednostavljenje skupa varijabli u niže dimenzionalni prostor. Klasični primjer je redukcija velikog teksta (mnoštva dokumenata) na najbinije teme (grupe u koje dokumenti pripadaju).

  • Klastering je procedura grupiranja skupa inputa pri čemu grupe nisu unaprijed poznate.

  • Klasifikacija je procedura selekcije inputa u dvije ili više klasa pri čemu je cilj ispravno selektirati nove, neviđene opservacije. Klasični primjer je selekcija pacijenta u kategorije zdrav ili bolestan.

  • Regresija se odnosi na sličan problem uz razliku što je output kontinuirana varijabla (ne diskretna/kategorijalna).

Neki od novijih i popularnih pristupa u ML-u su neuralne mreže (Neural Network), pojačano učenje (Reinforcement Learning) i duboko učenje (Deep Learning). Detaljniji opis svakog pojedinačnog pristupa je dovoljno opsežan da bi se mogao smjestiti u zaseban kolegij. Ovdje kratko ističemo duboko učenje kao dio ML-a koji se koristi za modeliranje visoko nelinearnih podataka i doživljava veliki uspjeh (popularnost) zbog koristi u mnoštvu praktičnih primjena poput prepoznavanja govora (speech recognition), autonomnih vozila, procesuiranja prirodnog jezika (natural language processing) i dr. Većina koncepata dubokog učenja je (u matematici) poznata već duže vrijeme no zbog razvoja u računalnim kapacitetima postaje sve šire zastupljena u prakski.

Uvodni primjeri ML pristupa

squares <- data.frame(
  size = c("small", "big","medium"),
  edge = c("dotte", "stripped", "normal"),
  color = c("green", "yellow","green")
)

Dimenzije ovog podatkovnog skupa su:

dim(squares)
## [1] 3 3

Struktura podatakovnog skupa:

str(squares)
## 'data.frame':    3 obs. of  3 variables:
##  $ size : chr  "small" "big" "medium"
##  $ edge : chr  "dotte" "stripped" "normal"
##  $ color: chr  "green" "yellow" "green"

Pregled varijabli:

summary(squares)
##      size               edge              color          
##  Length:3           Length:3           Length:3          
##  Class :character   Class :character   Class :character  
##  Mode  :character   Mode  :character   Mode  :character

Formulirajmo problem labeling-a nepoznatog kvadrata (square) na osnovi prethodnog znanja koje računalo (stroj) ima ukoliko mu input-amo squares podatke. Prikaz procedure bi izgledao otprilike ovako:

Procedura strojnog učenja (ML - machine learning) rješava ovakav (labeling) problem na osnovi naučenog znanja koje je formirano rješavanjem sličnih problema (i.e. inputu square podataka). Pojednostavljeno rečeno, računalu damo podatke, ono nauči strukturu i na osnovi toga ima sposobnost dati odgovor na novim, neviđenim podatcima. Glavna svrha ML je konstrukcija prediktivnih modela na osnovi podataka koji su sposobni dati odgovore na slična, ali ipak nova i nepostavljena pitanja.

Jednostavan i intuitivan primjer ML je također i korištenje regresijskog modela (koji već poznajemo) na podatcima o visini i težini ljudi kako bismo predvidjeli visinu nove osobe na osnovi njene težine:

Realniji primjer predikcijskog ML modela

Za izgradnju jednostavnog predikcijskog modela ćemo koristiti wage podatkovni skup iz ISLR paketa. Modelom ćemo pokušati opisati odnos između starosti radnika i visine plaće. Očekujemo da će stariji radnici imati više iskustva pa će prema tome biti i bolje plaćeni. Prvo ćemo učitati i pregledati podatke:

library(ISLR)
data(Wage)
summary(Wage)
##       year           age                     maritl           race     
##  Min.   :2003   Min.   :18.00   1. Never Married: 648   1. White:2480  
##  1st Qu.:2004   1st Qu.:33.75   2. Married      :2074   2. Black: 293  
##  Median :2006   Median :42.00   3. Widowed      :  19   3. Asian: 190  
##  Mean   :2006   Mean   :42.41   4. Divorced     : 204   4. Other:  37  
##  3rd Qu.:2008   3rd Qu.:51.00   5. Separated    :  55                  
##  Max.   :2009   Max.   :80.00                                          
##                                                                        
##               education                     region               jobclass   
##  1. < HS Grad      :268   2. Middle Atlantic   :3000   1. Industrial :1544  
##  2. HS Grad        :971   1. New England       :   0   2. Information:1456  
##  3. Some College   :650   3. East North Central:   0                        
##  4. College Grad   :685   4. West North Central:   0                        
##  5. Advanced Degree:426   5. South Atlantic    :   0                        
##                           6. East South Central:   0                        
##                           (Other)              :   0                        
##             health      health_ins      logwage           wage       
##  1. <=Good     : 858   1. Yes:2083   Min.   :3.000   Min.   : 20.09  
##  2. >=Very Good:2142   2. No : 917   1st Qu.:4.447   1st Qu.: 85.38  
##                                      Median :4.653   Median :104.92  
##                                      Mean   :4.654   Mean   :111.70  
##                                      3rd Qu.:4.857   3rd Qu.:128.68  
##                                      Max.   :5.763   Max.   :318.34  
##

Potom ćemo procijeniti odnos između plaće i iskustva (starosti radnika) u linearnom regresijkom modelu:

lm_wage <- lm(wage ~ age, data = Wage)

Za primjer ćemo pokušati predvidjeti zaradu 60 godišnjeg radnika. Za to je potrebno stvoriti novi podatkovni skup:

unseen <- data.frame(age = 60)

Na kraju ćemo korištenjem fukcije predict uz kombinaciju regresijskog objekta lm_wage i novog podatkovnosg skupa unseen provesti predikcijsku proceduru:

predict(lm_wage, unseen)
##        1 
## 124.1413

Na osnovi linearnog modela i Wage podataka smo predvidjeli da će 60 godišnji radnik zarađivati 124 $ na dan.

ML procedure

U najširem je smislu moguće izdvojiti tri vrste ML procedura:

  1. Klasifikacija

  2. Regresija

  3. Klastering

Klasifikacija

Svrha klasifikacije je predviđanje kategorije novih opservacija na osnovi prethodnog znanja (podataka, informacija). ML definira klasifikator (funkciju, model) koji će odlučiti kojoj kategoriji pripada novi objekt. Klasifikacijske tehnike imaju kvalitativni output, a kategorije moraju uvijek biti unaprijed definirane. Primjeri klasifikacije su medicinska dijagoza, prepoznavanje životinja, lica…

Pogledajmo jedan praktični primjer klasifikacijske procedure elektronske pošte koja je prethodno klasificirana kao spam/no spam. Feature koji želimo predvidjeti je avg_capital_seq, a u klasifikatoru se referiramo na njega sa x. Kriterij za odluku da li je pošta spam ili nije ćemo definirati kao prosječni broj velikih tiskanih slova u svakom pojedinačnom mail-u.

# Učitaj podatke
> emails <- read.table("../Podatci/emails.dat", header = T, sep = "")

# Pogledaj dimenzije podataka
> dim(emails)
[1] 13  2

# Definiraj klasifikator spam_classifier()
> spam_classifier <- function(x){
    prediction <- rep(NA, length(x)) # predikcijiski vektor
    prediction[x > 4] <- 1
    prediction[x >= 3 & x <= 4] <- 0
    prediction[x >= 2.2 & x < 3] <- 1
    prediction[x >= 1.4 & x < 2.2] <- 0
    prediction[x > 1.25 & x < 1.4] <- 1
    prediction[x <= 1.25] <- 0
    return(prediction) # predikcije je 0 ili 1
  }

# primijeni klasifikator na  avgCapitalSeq kolonu
> spamPred <- sapply(emails$avgCapitalSeq, spamClassifier)

# Usporedi predikcije sa stvarnim podatcima
> spam_pred == emails$spam
 [1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE

U ovoj proceduri je spam_classifier točno filtrirao spam poštu u svih 13 slučajeva. Valja istaknuti da je ovo bio samo demonstrativni primjer te da klasifikacija na novim podatcima gotovo sigurno ne bi bila tako precizna!

Regresija

Regresijska analiza se koristi za predviđanje kontinuirane ili kvantitativne vrijednosti na osnovi prethodnih podataka (informacija).

Regresija je slična klasifikaciji uz bitnu razliku da želimo predvidjeti vrijednost, a ne klasu kao kod klasifikacije. Tipični regresijski problemi se odnose na predviđanje kreditnog rejtinga na osnovi redovitosti plaćanja u prošlosti, analize pretplata na časopis kroz vrijeme, vjerojatnost zaposlenja na osnovi ocjena… Regresijske tehnike imaju kvantitativni output i zahtijevaju input-output podatkovni skup.

Za praktični primjer (još jedan!) analiziramo broj pogleda na LinkedIn profilu u regresijskom kontekstu. Cilj je predvidjeti broj pregleda profila u sljedeća 3 dana na osnovi podataka o pregledima u prethodna 3 tjedna. Za to je potrebno definirati vektor LinkedIn pregleda:

linkedin <- c(5,7,4,9,11,10,14,17,13,11,18,17,21,21,24,23,28,35,21,27,23)

Potom valja napraviti vektor dana:

days <- rep(NA, length(linkedin))
days <- seq(1, 21, by = 1)

Sada možemo procijeniti regresijski odnos:

linkedin_lm <- lm(linkedin ~ days)

Za predviđanje je potrebno specificirati dane za koje želimo rezultate te potom provesti predikciju:

future_days <- data.frame(days = 22:24)
linkedin_pred <- predict(linkedin_lm, future_days)

Na kraju pogledajmo jednu vizualizaciju predviđanja:

plot(linkedin ~ days, xlim = c(1, 24))
points(22:24, linkedin_pred, col = "green")

Klastering

Procedura klasteringa ima za cilj grupirati objekte u klastere na način da su klasteri međusobno različiti, a objekti unutar svakog klastera slični. Klastering je sličan klasifikaciji uz razliku da ne postoje prethodno definirane kategorije već klaster procedura sortira objekte prema sličnosti (ne prema kategorijama). Za klastering nije potrebno definirati kategorije (labels) i ne postoji točan ili netočan output…postoji beskonačno mnogo mogućih klastera.

Za primjer klastering procedure ćemo koristiti otprije poznati podatkovni skup iris. Pri tome je cilj grupirati iris cvijetove u 3 zasebna klastera, prema karakteristikama cvjetova koje su dane u podatcima. Pogledajmo prvo podatke:

library(datasets)
data(iris)

Pregled podataka:

str(iris)
## 'data.frame':    150 obs. of  5 variables:
##  $ Sepal.Length: num  5.1 4.9 4.7 4.6 5 5.4 4.6 5 4.4 4.9 ...
##  $ Sepal.Width : num  3.5 3 3.2 3.1 3.6 3.9 3.4 3.4 2.9 3.1 ...
##  $ Petal.Length: num  1.4 1.4 1.3 1.5 1.4 1.7 1.4 1.5 1.4 1.5 ...
##  $ Petal.Width : num  0.2 0.2 0.2 0.2 0.2 0.4 0.3 0.2 0.2 0.1 ...
##  $ Species     : Factor w/ 3 levels "setosa","versicolor",..: 1 1 1 1 1 1 1 1 1 1 ...

Pregled dimenzija podatkovnog skupa:

dim(iris)
## [1] 150   5

Pregled opservacija:

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

Deskriptivna statistika podataka:

summary(iris)
##   Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
##  Min.   :4.300   Min.   :2.000   Min.   :1.000   Min.   :0.100  
##  1st Qu.:5.100   1st Qu.:2.800   1st Qu.:1.600   1st Qu.:0.300  
##  Median :5.800   Median :3.000   Median :4.350   Median :1.300  
##  Mean   :5.843   Mean   :3.057   Mean   :3.758   Mean   :1.199  
##  3rd Qu.:6.400   3rd Qu.:3.300   3rd Qu.:5.100   3rd Qu.:1.800  
##  Max.   :7.900   Max.   :4.400   Max.   :6.900   Max.   :2.500  
##        Species  
##  setosa    :50  
##  versicolor:50  
##  virginica :50  
##                 
##                 
##

Za procjenu je potebno maknuti zadnju kolonu jer sadrži string varijable i pospremiti ju u zasebni objekt:

my_iris <- iris[-5]
species <- iris$Species

Potom ćemo provesti klastering pomoću kmeans() funkcije:

kmeans_iris <- kmeans(my_iris, 3)

Usporedimo sada rezultate klastering procedure sa orginalnim podatcima pomoću table funkcije. Procijenjene grupe se nalaze u cluster dijelu objekta kmeans_iris:

table(kmeans_iris$cluster, species)
##    species
##     setosa versicolor virginica
##   1     50          0         0
##   2      0          2        36
##   3      0         48        14

Rezultati imaju visok stupanj preciznosti. Za kraj je uvijek korisno i vizualno prikazati rezultate:

plot(Petal.Length ~ Petal.Width,
     data = my_iris,
     col = kmeans_iris$cluster)

Nakon ovih primjera ćemo pogledati kako izgleda ML procedura “od početka”…

Eksplorativna analiza i nenadgledano strojno učenje

Prvo ćemo učitati podatke Breast Cancer Wisconsin (Diagnostic) Data Set iz UCI Machine Learning repozitorija.

library(tidyverse) 

url1 <-"https://archive.ics.uci.edu/ml/machine-learning-databases/breast-cancer-wisconsin/wdbc.data"
url2 <- "https://raw.githubusercontent.com/fpsom/2020-07-machine-learning-sib/master/data/wdbc.colnames.csv"

breastCancerData <- read_csv(url1, col_names = FALSE)

breastCancerDataColNames <- read_csv(url2, col_names = FALSE)

colnames(breastCancerData) <- breastCancerDataColNames$X1

Nakon učitavanja podataka ćemo pogledati kako izgledaju podatci:

breastCancerData %>% head(10)
## # A tibble: 10 x 32
##          ID Diagnosis Radius.Mean Texture.Mean Perimeter.Mean Area.Mean
##       <dbl> <chr>           <dbl>        <dbl>          <dbl>     <dbl>
##  1   842302 M                18.0         10.4          123.      1001 
##  2   842517 M                20.6         17.8          133.      1326 
##  3 84300903 M                19.7         21.2          130       1203 
##  4 84348301 M                11.4         20.4           77.6      386.
##  5 84358402 M                20.3         14.3          135.      1297 
##  6   843786 M                12.4         15.7           82.6      477.
##  7   844359 M                18.2         20.0          120.      1040 
##  8 84458202 M                13.7         20.8           90.2      578.
##  9   844981 M                13           21.8           87.5      520.
## 10 84501001 M                12.5         24.0           84.0      476.
## # ... with 26 more variables: Smoothness.Mean <dbl>, Compactness.Mean <dbl>,
## #   Concavity.Mean <dbl>, Concave.Points.Mean <dbl>, Symmetry.Mean <dbl>,
## #   Fractal.Dimension.Mean <dbl>, Radius.SE <dbl>, Texture.SE <dbl>,
## #   Perimeter.SE <dbl>, Area.SE <dbl>, Smoothness.SE <dbl>,
## #   Compactness.SE <dbl>, Concavity.SE <dbl>, Concave.Points.SE <dbl>,
## #   Symmetry.SE <dbl>, Fractal.Dimension.SE <dbl>, Radius.Worst <dbl>,
## #   Texture.Worst <dbl>, Perimeter.Worst <dbl>, Area.Worst <dbl>, ...

U okviru eksplorativne analize ćemo pogledati kako izgledaju podatci. U tu svrhu je potrebno napraviti nekoliko prilagodbi podataka:

# napravi faktorsku varijablu
breastCancerData$Diagnosis <- as.factor(breastCancerData$Diagnosis)
# ukloni prvu kolonu 
breastCancerDataNoID <- breastCancerData[2:ncol(breastCancerData)]
# vizualizacija
library(GGally)  
ggpairs(breastCancerDataNoID[1:5], aes(color=Diagnosis, alpha=0.4))

Sljedeća stvar koju moramo napraviti je prilagodba podataka u smislu centriranja i skaliranja. To je potrebno napraviti jer prediktori imaju vrlo varijablne prosjeke i standardne devijacije pa nisu usporedivi.

# koristimo funkcije iz ML paketa caret
library(caret)
ppv <- preProcess(breastCancerDataNoID, method = c("center", "scale"))
breastCancerDataNoID_tr <- predict(ppv, breastCancerDataNoID)
# pregledaj podatke prije centriranja i skaliranja
breastCancerDataNoID[1:5] %>% summary()
##  Diagnosis  Radius.Mean      Texture.Mean   Perimeter.Mean     Area.Mean     
##  B:357     Min.   : 6.981   Min.   : 9.71   Min.   : 43.79   Min.   : 143.5  
##  M:212     1st Qu.:11.700   1st Qu.:16.17   1st Qu.: 75.17   1st Qu.: 420.3  
##            Median :13.370   Median :18.84   Median : 86.24   Median : 551.1  
##            Mean   :14.127   Mean   :19.29   Mean   : 91.97   Mean   : 654.9  
##            3rd Qu.:15.780   3rd Qu.:21.80   3rd Qu.:104.10   3rd Qu.: 782.7  
##            Max.   :28.110   Max.   :39.28   Max.   :188.50   Max.   :2501.0

Podatci nakon predprocesnih prilagodbi centriranja i skaliranja izgledaju ovako:

breastCancerDataNoID_tr[1:5] %>% summary()
##  Diagnosis  Radius.Mean       Texture.Mean     Perimeter.Mean   
##  B:357     Min.   :-2.0279   Min.   :-2.2273   Min.   :-1.9828  
##  M:212     1st Qu.:-0.6888   1st Qu.:-0.7253   1st Qu.:-0.6913  
##            Median :-0.2149   Median :-0.1045   Median :-0.2358  
##            Mean   : 0.0000   Mean   : 0.0000   Mean   : 0.0000  
##            3rd Qu.: 0.4690   3rd Qu.: 0.5837   3rd Qu.: 0.4992  
##            Max.   : 3.9678   Max.   : 4.6478   Max.   : 3.9726  
##    Area.Mean      
##  Min.   :-1.4532  
##  1st Qu.:-0.6666  
##  Median :-0.2949  
##  Mean   : 0.0000  
##  3rd Qu.: 0.3632  
##  Max.   : 5.2459

Vidljivo je da varijable imaju prosjek 0 i da su zadržali distribucijska svojstva kao originalne varijable. Pogledajmo nove varijable i grafički:

ggpairs(breastCancerDataNoID_tr[1:5], aes(color=Diagnosis, alpha=0.4))

Smanjenje dimenzionalnosti i PCA

Smanjenje dimenzionalnosti jedan je od važnih ciljeva strojnog učenja pri čemu je analiza glavnih komponenti (PCA) jedna od najzastupljenijih metoda za tu svrhu. Bez da ulazimo u metodološke detalje PCA, pogledajmo kako ova procedura izgleda u praksi:

# provedi PCA
ppv_pca <- prcomp(breastCancerData[3:ncol(breastCancerData)],
                  center = TRUE,
                  scale. = TRUE)
# pregledaj rezultate
summary(ppv_pca)
## Importance of components:
##                           PC1    PC2     PC3     PC4     PC5     PC6     PC7
## Standard deviation     3.6444 2.3857 1.67867 1.40735 1.28403 1.09880 0.82172
## Proportion of Variance 0.4427 0.1897 0.09393 0.06602 0.05496 0.04025 0.02251
## Cumulative Proportion  0.4427 0.6324 0.72636 0.79239 0.84734 0.88759 0.91010
##                            PC8    PC9    PC10   PC11    PC12    PC13    PC14
## Standard deviation     0.69037 0.6457 0.59219 0.5421 0.51104 0.49128 0.39624
## Proportion of Variance 0.01589 0.0139 0.01169 0.0098 0.00871 0.00805 0.00523
## Cumulative Proportion  0.92598 0.9399 0.95157 0.9614 0.97007 0.97812 0.98335
##                           PC15    PC16    PC17    PC18    PC19    PC20   PC21
## Standard deviation     0.30681 0.28260 0.24372 0.22939 0.22244 0.17652 0.1731
## Proportion of Variance 0.00314 0.00266 0.00198 0.00175 0.00165 0.00104 0.0010
## Cumulative Proportion  0.98649 0.98915 0.99113 0.99288 0.99453 0.99557 0.9966
##                           PC22    PC23   PC24    PC25    PC26    PC27    PC28
## Standard deviation     0.16565 0.15602 0.1344 0.12442 0.09043 0.08307 0.03987
## Proportion of Variance 0.00091 0.00081 0.0006 0.00052 0.00027 0.00023 0.00005
## Cumulative Proportion  0.99749 0.99830 0.9989 0.99942 0.99969 0.99992 0.99997
##                           PC29    PC30
## Standard deviation     0.02736 0.01153
## Proportion of Variance 0.00002 0.00000
## Cumulative Proportion  1.00000 1.00000

Detaljniji pregled PCA objekta izgleda ovako:

str(ppv_pca)
## List of 5
##  $ sdev    : num [1:30] 3.64 2.39 1.68 1.41 1.28 ...
##  $ rotation: num [1:30, 1:30] -0.219 -0.104 -0.228 -0.221 -0.143 ...
##   ..- attr(*, "dimnames")=List of 2
##   .. ..$ : chr [1:30] "Radius.Mean" "Texture.Mean" "Perimeter.Mean" "Area.Mean" ...
##   .. ..$ : chr [1:30] "PC1" "PC2" "PC3" "PC4" ...
##  $ center  : Named num [1:30] 14.1273 19.2896 91.969 654.8891 0.0964 ...
##   ..- attr(*, "names")= chr [1:30] "Radius.Mean" "Texture.Mean" "Perimeter.Mean" "Area.Mean" ...
##  $ scale   : Named num [1:30] 3.524 4.301 24.299 351.9141 0.0141 ...
##   ..- attr(*, "names")= chr [1:30] "Radius.Mean" "Texture.Mean" "Perimeter.Mean" "Area.Mean" ...
##  $ x       : num [1:569, 1:30] -9.18 -2.39 -5.73 -7.12 -3.93 ...
##   ..- attr(*, "dimnames")=List of 2
##   .. ..$ : NULL
##   .. ..$ : chr [1:30] "PC1" "PC2" "PC3" "PC4" ...
##  - attr(*, "class")= chr "prcomp"

Rezultati analize izgledaju ovako:

# vizualizacija
library(ggbiplot)

ggbiplot(ppv_pca, choices=c(1, 2),
         labels=rownames(breastCancerData),
         ellipse=TRUE,
         groups = breastCancerData$Diagnosis,
         obs.scale = 1,
         var.axes=TRUE, var.scale = 1) +
  ggtitle("PCA of Breast Cancer Dataset")+
  theme_minimal()+
  theme(legend.position = "bottom")

Klastering

Ovo je također jedna od vrlo popularnih metoda nenadgledanog strojnog učenja koja prepoznaje strukturu u podatcima na način da organizira podatke (opservacije) u interno koherentne i eksterno heterogene klastere. Drugačije rečeno, opservacije unutar klastera trebaju biti jako slične, a opservacije između klastera jako različite. Klasterski algoritmi se najčešće dijele na aglomerativne i divizivne. Aglomerativni započinju sa svakom opservacijom u pojedinačnom klasteru i sukcesivno ih spajaju do ispunjenja kriterija zaustavljanja. Divizivni počinju od jednog jedinstvenog klastera i dijele ga do ispunjenja kriterija zaustavljanje. I jedna i druga vrsta su zapravo procedure grupiranja opservacija prema nekom kriteriju, najčešće blizine u prostoru.

Sada ćemo provesti klastering na podatcima o tumoru koristeći algoritam k-means. Ovaj algoritam provodi grupiranje po principu minimizacije varijance unutar klastera. Postoji više varijanti k-means klastera, a u ovom primjeru koristimo standardni algoritam koji funkcionira na principu ukupne unutar-klasterske varijacije definirane kao sume kvadrata ekuklidovskih udaljenosti između opservacije i pripadajućeg centroida.

Provedimo k-means algoritam:

set.seed(1)
km.out <- kmeans(breastCancerData[3:ncol(breastCancerData)],
                 centers=2,  # broj klastera
                 nstart=20)  #  broj višestrukih inicijalnih konfiguracija

Rezultati provedene procedure izgledaju ovako:

str(km.out)
## List of 9
##  $ cluster     : int [1:569] 2 2 2 1 2 1 2 1 1 1 ...
##  $ centers     : num [1:2, 1:30] 12.6 19.4 18.6 21.7 81.1 ...
##   ..- attr(*, "dimnames")=List of 2
##   .. ..$ : chr [1:2] "1" "2"
##   .. ..$ : chr [1:30] "Radius.Mean" "Texture.Mean" "Perimeter.Mean" "Area.Mean" ...
##  $ totss       : num 2.57e+08
##  $ withinss    : num [1:2] 28559677 49383423
##  $ tot.withinss: num 77943100
##  $ betweenss   : num 1.79e+08
##  $ size        : int [1:2] 438 131
##  $ iter        : int 1
##  $ ifault      : int 0
##  - attr(*, "class")= chr "kmeans"

Sada možemo vizualizirati dobivene rezultate:

ggplot(as.data.frame(ppv_pca$x),
       aes(x=PC1, y=PC2,
           color=as.factor(km.out$cluster),
           shape = breastCancerData$Diagnosis)) +
  geom_point( alpha = 0.6, size = 3) +
  theme_minimal()+
  theme(legend.position = "bottom") +
  labs(title = "K-Means clusteri i PCA", x = "PC1", y = "PC2", color = "Cluster", shape = "Diagnosis")

Kod klasteringa je potrebno prethodno kalibrirati algoritam na način da se unaprijed specificira broj klastera. Taj broj se određuje tzv. elbow metodom koja funkcionira u nekoliko koraka. Prvo je potrebno izračunati ukupnu sumu kvadrata odstupanja (varijabilnost) putem sljedeće funkcije:

kmean_withinss <- function(k) {
  cluster <- kmeans(breastCancerData[3:ncol(breastCancerData)], k)
  return (cluster$tot.withinss)
}

Nakon toga možemo probati sa 2 klastera:

kmean_withinss(2)
## [1] 77943100

Potom ćemo probati za sekvencu klastera od 1 do 20:

max_k <-20
# provedi funkciju na klasterskom k rasponu  
wss <- sapply(2:max_k, kmean_withinss)

Na kraju ćemo prikazati rezultate grafički:

elbow <- data.frame(2:max_k, wss)

ggplot(elbow, aes(x = X2.max_k, y = wss)) +
  geom_point() +
  geom_line() +
  scale_x_continuous(breaks = seq(1, 20, by = 1))

Dodatni primjer

U ovom primjeru želimo grupirati automobile na osnovi snage motora i težine. Pri tome koristimo mtcars podatkovni skup. Prvo ćemo učitati, izabrati i pregledati podatke:

data(mtcars)
cars <- subset(mtcars, select = c(wt,hp))
head(cars)
##                      wt  hp
## Mazda RX4         2.620 110
## Mazda RX4 Wag     2.875 110
## Datsun 710        2.320  93
## Hornet 4 Drive    3.215 110
## Hornet Sportabout 3.440 175
## Valiant           3.460 105

Upoznajmo se dodatno s podatcima:

str(cars)
## 'data.frame':    32 obs. of  2 variables:
##  $ wt: num  2.62 2.88 2.32 3.21 3.44 ...
##  $ hp: num  110 110 93 110 175 105 245 62 95 123 ...

Pogledajmo i deskriptivnu statistiku:

summary(cars)
##        wt              hp       
##  Min.   :1.513   Min.   : 52.0  
##  1st Qu.:2.581   1st Qu.: 96.5  
##  Median :3.325   Median :123.0  
##  Mean   :3.217   Mean   :146.7  
##  3rd Qu.:3.610   3rd Qu.:180.0  
##  Max.   :5.424   Max.   :335.0

Za provedbu klasteringa ćemo koristiti kmeans() funkciju koju već poznajemo:

km_cars <- kmeans(cars,2)

Pogledajmo rezultate:

km_cars$cluster
##           Mazda RX4       Mazda RX4 Wag          Datsun 710      Hornet 4 Drive 
##                   1                   1                   1                   1 
##   Hornet Sportabout             Valiant          Duster 360           Merc 240D 
##                   1                   1                   2                   1 
##            Merc 230            Merc 280           Merc 280C          Merc 450SE 
##                   1                   1                   1                   1 
##          Merc 450SL         Merc 450SLC  Cadillac Fleetwood Lincoln Continental 
##                   1                   1                   2                   2 
##   Chrysler Imperial            Fiat 128         Honda Civic      Toyota Corolla 
##                   2                   1                   1                   1 
##       Toyota Corona    Dodge Challenger         AMC Javelin          Camaro Z28 
##                   1                   1                   1                   2 
##    Pontiac Firebird           Fiat X1-9       Porsche 914-2        Lotus Europa 
##                   1                   1                   1                   1 
##      Ford Pantera L        Ferrari Dino       Maserati Bora          Volvo 142E 
##                   2                   1                   2                   1

Broj ispod modela predstavlja pripadnost klasteru!

Pogledajmo sada klasterske centroide. Centroidi predstavljaju “centar” svakog klastera i mogu se dohvatiti kao centers element procijenjenog klasterskog objekta na sljedeći način:

km_cars$centers
##         wt       hp
## 1 2.911320 118.2000
## 2 4.309857 248.4286

Za kraj ćemo još vizualno prikazati rezultate:

plot(cars, col = km_cars$cluster)
points(km_cars$centers, pch = 22, bg = c(1, 2), cex = 2)

Nadgledano strojno učenje

Nadgledano ML se bavi predikcijom ishoda (labels) poput preživio/nije preživio na osnovi učenja sa označenih (labelled) podataka. Cilj je izvršiti “trening” modela na postojećim podatcima kako bi model mogao predvidjeti nove, nepoznate podatke. Ti novi podatci se najčešće odnose na dio (80/20,70/30,60/40) postojećeg uzorka koji se izdvaja iz “trening” dijela uzorka kako bi se na njima provela procjena kvalitete predikcije. Procjena kvalitete se provodi na osnovi različitih kriterija kojima se procjenjuje točnost modela. Standardni koraci u provedbi nadgledanog strojnog učenja su: provedba eksplorativne podatkovne analize, provedba nekog početnog modela, poboljšanje i preuređenje varijabli, procijena boljeg modela. Neki od često korištenih modela su stabla odlučivanja (decision trees), bayesov klasifikator, knn klasifikator, support vector machines i dr.

Stabla odlučivanja

Bez da ulazimo u detalje oko ovog algoritma, pogledajmo kako on funkcionira u praksi. Prvo ćemo podijeliti podake na trening i test uzorak u odnosu 70/30.

# podijeli uzorak na training i test 
set.seed(1000)
ind <- sample(2, nrow(breastCancerData), replace=TRUE, prob=c(0.7, 0.3))
breastCancerData.train <- breastCancerDataNoID[ind==1,]
breastCancerData.test <- breastCancerDataNoID[ind==2,]

Nakon toga ćemo specificirati model i definirati parametre:

# učitaj knjužnice
library(rpart)
library(rpart.plot)
# specificiraj formulu
myFormula <- Diagnosis ~ Radius.Mean + Area.Mean + Texture.SE
# specificiraj algoritam
breastCancerData.model <- rpart(myFormula,
                                method = "class",
                                data = breastCancerData.train,
                                minsplit = 10, # minimalni broj instanci u čvoru
                                minbucket = 1, # minimalni broj instanci u listu
                                maxdepth = 3,  # maksimalna dubina drveta
                                cp = -1)       # parametar kompleksnosti;intuitivno
# pregledaj rezultate
print(breastCancerData.model$cptable)
##             CP nsplit rel error    xerror       xstd
## 1  0.709459459      0 1.0000000 1.0000000 0.06469979
## 2  0.016891892      1 0.2905405 0.3040541 0.04262349
## 3  0.006756757      3 0.2567568 0.3310811 0.04421819
## 4  0.000000000      4 0.2500000 0.3310811 0.04421819
## 5 -1.000000000      6 0.2500000 0.3310811 0.04421819
rpart.plot(breastCancerData.model)

Izaberimo model sa najmanjom predikcijskom pogreškom:

opt <- which.min(breastCancerData.model$cptable[, "xerror"])
cp <- breastCancerData.model$cptable[opt, "CP"]
# očisti drvo
breastCancerData.pruned.model <- prune(breastCancerData.model, cp = cp)
# vizualiziraj
rpart.plot(breastCancerData.pruned.model)

# pogledaj rezultate
table(predict(breastCancerData.pruned.model, type="class"), breastCancerData.train$Diagnosis)
##    
##       B   M
##   B 233  35
##   M   8 113

Nakon provedenog modela valja pogledati kvalitetu predikcije na testnom dijelu uzorka:

# provedi predikciju
BreastCancer_pred <- predict(breastCancerData.pruned.model, newdata = breastCancerData.test, type="class")
# vizualizacija
plot(BreastCancer_pred ~ Diagnosis, data = breastCancerData.test,
     xlab = "Observed",
     ylab = "Prediction")

# pregled rezultata
table(BreastCancer_pred, breastCancerData.test$Diagnosis)
##                  
## BreastCancer_pred   B   M
##                 B 114  15
##                 M   2  49

Dodatni primjer

U prethodnom primjeru smo koristili kmeans() funkciju za klastering na iris podatcima. U ovom primjeru ćemo koristiti rekurzivno particioniranje na istim podatcima (no nećemo maknuti zadnju kolonu).

Pogledajmo za početak strukturu podataka:

summary(iris)
##   Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
##  Min.   :4.300   Min.   :2.000   Min.   :1.000   Min.   :0.100  
##  1st Qu.:5.100   1st Qu.:2.800   1st Qu.:1.600   1st Qu.:0.300  
##  Median :5.800   Median :3.000   Median :4.350   Median :1.300  
##  Mean   :5.843   Mean   :3.057   Mean   :3.758   Mean   :1.199  
##  3rd Qu.:6.400   3rd Qu.:3.300   3rd Qu.:5.100   3rd Qu.:1.800  
##  Max.   :7.900   Max.   :4.400   Max.   :6.900   Max.   :2.500  
##        Species  
##  setosa    :50  
##  versicolor:50  
##  virginica :50  
##                 
##                 
##

U sljedećem koraku ćemo napraviti model nadgledanog strojnog učenja koristeći funkciju rpart() iz istoimenog paketa. Specifikacija decision tree modela izgleda ovako:

library(rpart)
tree <- rpart(Species ~ Sepal.Length + Sepal.Width + Petal.Length + Petal.Width,
              data = iris,
              method = "class")

Potom je potrebno specificirati podatke za koje želimo predviđanje:

unseen <- data.frame(Sepal.Length = c(5.3, 7.2),
                       Sepal.Width = c(2.9, 3.9),
                       Petal.Length = c(1.7, 5.4),
                       Petal.Width = c(0.8, 2.3))

U zadnjem koraku izvršimo predviđanje pomoću predict funkcije:

predict(tree, unseen, type = "class")
##         1         2 
##    setosa virginica 
## Levels: setosa versicolor virginica

Random Forest algoritam

Riječ je o ansambl algoritmu koji konstruira višestruka stabla odlučivanja tako da je konačna predikcija izvedena kao kombinacija predikcija svakog pojedinačnog stabla odlučivanja. Ovaj algoritam često daje bolje rezultate ali pod cijenu lošije interpretabilnosti. Pogledajmo kako algoritam funkcionira u praksi. Prvo ćemo specificirati model:

library(randomForest)
set.seed(1000)
# specifikacija model
rf <- randomForest(Diagnosis ~ ., data = breastCancerData.train,
                   ntree=100,
                   proximity=T)
# pogledaj rezultate
table(predict(rf), breastCancerData.train$Diagnosis)
##    
##       B   M
##   B 233  11
##   M   8 137

Pogledajmo modelski objekt:

print(rf)
## 
## Call:
##  randomForest(formula = Diagnosis ~ ., data = breastCancerData.train,      ntree = 100, proximity = T) 
##                Type of random forest: classification
##                      Number of trees: 100
## No. of variables tried at each split: 5
## 
##         OOB estimate of  error rate: 4.88%
## Confusion matrix:
##     B   M class.error
## B 233   8  0.03319502
## M  11 137  0.07432432

Također je moguće pregledati važnost varijabli:

importance(rf)
##                         MeanDecreaseGini
## Radius.Mean                    7.5456172
## Texture.Mean                   2.8935249
## Perimeter.Mean                 6.1184541
## Area.Mean                     10.7099799
## Smoothness.Mean                1.0152042
## Compactness.Mean               2.1720691
## Concavity.Mean                10.7662569
## Concave.Points.Mean           17.1142616
## Symmetry.Mean                  0.7182040
## Fractal.Dimension.Mean         0.8412011
## Radius.SE                      2.8843242
## Texture.SE                     0.8093371
## Perimeter.SE                   2.7452789
## Area.SE                        5.9382151
## Smoothness.SE                  0.7157598
## Compactness.SE                 0.5104771
## Concavity.SE                   0.6647103
## Concave.Points.SE              0.9539497
## Symmetry.SE                    0.7490675
## Fractal.Dimension.SE           0.6424414
## Radius.Worst                  22.1621219
## Texture.Worst                  4.2069603
## Perimeter.Worst               23.4641800
## Area.Worst                    22.1238001
## Smoothness.Worst               2.3129352
## Compactness.Worst              2.1932002
## Concavity.Worst                5.2577815
## Concave.Points.Worst          20.8737773
## Symmetry.Worst                 1.8072152
## Fractal.Dimension.Worst        1.4227531
varImpPlot(rf)

Sada ćemo provesti predikciju na testnom dijelu uzorka:

BreastCancer_pred_RD <- predict(rf, newdata = breastCancerData.test)
table(BreastCancer_pred_RD, breastCancerData.test$Diagnosis)
##                     
## BreastCancer_pred_RD   B   M
##                    B 116   4
##                    M   0  60
plot(margin(rf, breastCancerData.test$Diagnosis))

Provedi predikciju sa smanjenim brojem varijabli rangiranima prema važnosti:

result <- rfcv(breastCancerData.train, breastCancerData.train$Diagnosis, cv.fold=3)
with(result, plot(n.var, error.cv, log="x", type="o", lwd=2))

Regresija

Prije provedbe regresijskog modela je korisno pogledati korelaciju između varijabli:

# korelacija između Radius.Mean i Concave.Points.Mean / Area.Mean
cor(breastCancerData$Radius.Mean, breastCancerData$Concave.Points.Mean)
## [1] 0.8225285
cor(breastCancerData$Concave.Points.Mean, breastCancerData$Area.Mean)
## [1] 0.8232689
# selektiraj varijable koje se koriste u modelu
bc <- select(breastCancerData,Radius.Mean,Concave.Points.Mean,Area.Mean)

Napravi regresijski model:

bc_model_full <- lm(Radius.Mean ~ Concave.Points.Mean + Area.Mean, data=bc)
bc_model_full
## 
## Call:
## lm(formula = Radius.Mean ~ Concave.Points.Mean + Area.Mean, data = bc)
## 
## Coefficients:
##         (Intercept)  Concave.Points.Mean            Area.Mean  
##             7.68087              2.72493              0.00964

Napravimo sada predikciju na osnovi provedenog modela:

# predikcija
preds <- predict(bc_model_full)
# vizualni prikaz
plot(preds, bc$Radius.Mean, xlab = "Predikcija", ylab = "Opservacija")
abline(a = 0, b = 1)

# pregled punog modela
summary(bc_model_full)
## 
## Call:
## lm(formula = Radius.Mean ~ Concave.Points.Mean + Area.Mean, data = bc)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.8307 -0.1827  0.1497  0.3608  0.7411 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)         7.6808702  0.0505533 151.936   <2e-16 ***
## Concave.Points.Mean 2.7249328  1.0598070   2.571   0.0104 *  
## Area.Mean           0.0096400  0.0001169  82.494   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5563 on 566 degrees of freedom
## Multiple R-squared:  0.9752, Adjusted R-squared:  0.9751 
## F-statistic: 1.111e+04 on 2 and 566 DF,  p-value: < 2.2e-16

Ovo je procjena na cijelom podatkovnom skupu pa ne znamo kaka će biti rezultat na novim, neviđenim podatcima. Zbog toga ćemo napraviti podjelu na trening/test i provjeriti kvalitetu modela, prvo na trening uzorku:

set.seed(123)

ind <- sample(2, nrow(bc), replace=TRUE, prob=c(0.75, 0.25))
bc_train <- bc[ind==1,]
bc_test <- bc[ind==2,]

# linearna regresija na trening podatcima
(bc_model <- lm(Radius.Mean ~ Concave.Points.Mean + Area.Mean, data=bc_train))
## 
## Call:
## lm(formula = Radius.Mean ~ Concave.Points.Mean + Area.Mean, data = bc_train)
## 
## Coefficients:
##         (Intercept)  Concave.Points.Mean            Area.Mean  
##             7.69215              3.18929              0.00957
# pregled rezultata
summary(bc_model)
## 
## Call:
## lm(formula = Radius.Mean ~ Concave.Points.Mean + Area.Mean, data = bc_train)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.7465 -0.1879  0.1736  0.3828  0.6830 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)         7.6921484  0.0606712 126.784   <2e-16 ***
## Concave.Points.Mean 3.1892857  1.3200486   2.416   0.0161 *  
## Area.Mean           0.0095705  0.0001455  65.787   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5802 on 428 degrees of freedom
## Multiple R-squared:  0.9729, Adjusted R-squared:  0.9727 
## F-statistic:  7670 on 2 and 428 DF,  p-value: < 2.2e-16
# predikcija modela
bc_train$pred <- predict(bc_model)

# vizalizacija stvarnih i predikcijskih vrijednosti
ggplot(bc_train, aes(x = pred, y = Radius.Mean)) +
  geom_point() +
  geom_abline(color = "blue")

Predikcija na testnim podatcima:

# provedi predikciju
bc_test$pred <- predict(bc_model , newdata=bc_test)
# vizualiziraj rezultate
# plot the ground truths vs predictions for test set and examine the plot. Does it look as good with the predictions on the training set?
ggplot(bc_test, aes(x = pred, y = Radius.Mean)) +
  geom_point() +
  geom_abline(color = "blue")

Procjena kvalitete modela

Procjena kvalitete modela kod nadgledanog ML se vrši usporedbom stvarnih label-a sa onima koje je generirao model. U dobroj procjeni je razlika između istih vrlo mala. Kod ne-nadgledanog ML je procjena kvalitete modela nešto složenija pošto ne postoje label-i za usporedbu. Vrlo često se susrećemo i sa kvazi nadgledanim ML gdje dio opservacija ima pripisane label-e , a dio ne. Procedura procjene kvalitete modela u ML ima određene specifičnosti i ponešto se razlikuje od procjene kvalitete modela u statistici (model fit). Zbog toga ćemo kasnije pogledati detaljnije o kojim procedurama se radi i kako se provode.

Procjena prediktivne sposobnosti ML modela uvelike ovisi o empirisjkom problemu i korištenoj metodi. Tri su najčešće korištena kriterija za procjenu kvalitete modela:

  • Preciznost (accuracy)
  • Brzina (computation time)
  • Mogućnost interpretacije rezultata (interpretability)

Svaki od tri ML pristupa (klasifikacija, regresija, klaster) ima specifične kriterije za procjenu kvalitete modela.

Klasifikacija

U klasifikacijskim procedurama su glavni kriteriji kvalitete modela preciznost (accuracy) i greška (error). Ove mjere odažavaju broj pogrešaka i točnih klasifikacija. Definiramo ih na sljedeći način:

\[\mbox{preciznost} = \frac{točne Klasifikacije}{ukupne Klasifikacije}\] i

\[\mbox{greška} = 1 - preciznost\] U klasifikaciji se često koristi i confusion matrica koja u svakoj od ćelija sadrži broj instanci koje su klasificirane na određeni način. Za primjer jednostavne binarne klasifikacije bi confusion matrica izgledala ovako:

Kratice znače sljedeće: p-pozitivno, n-negativno, tp-točno pozitivno, fp-netočno pozitivno, tn-točno negativno, fn-netočno negativno. Iz elemenata confusion matrice se mogu izračunati mjere preciznosti poput: točnost (accuracy), preciznost (precision) i recall:

\[\mbox{točnost} = \frac{tp+tn}{tp+fp+tn+fn}\]

\[\mbox{preciznost} = \frac{tp}{tp+fp}\]

\[\mbox{recall} = \frac{tp}{tp+fn}\]

Praktični primjer

U ovom primjeru ćemo pokušati predvidjeti vjerojatnost preživljavanja nesreće na Titaniku s obzirom na karakteristike putnika (osobe) koje uključuju dob, spol i putničku klasu. Prvo učitavamo i izabiremo potrebne podatke:

titanic = read.csv("../Dta/Titanic_train.csv", na.strings = "") # dostupno na Kaggle!
titanic_dta <- subset(titanic, !is.na(Age), select = c(Survived, Pclass, Sex, Age))
str(titanic_dta)
## 'data.frame':    714 obs. of  4 variables:
##  $ Survived: int  0 1 1 1 0 0 0 1 1 1 ...
##  $ Pclass  : int  3 1 3 1 3 1 3 3 2 3 ...
##  $ Sex     : chr  "male" "female" "female" "female" ...
##  $ Age     : num  22 38 26 35 35 54 2 27 14 4 ...
summary(titanic_dta)
##     Survived          Pclass          Sex                 Age       
##  Min.   :0.0000   Min.   :1.000   Length:714         Min.   : 0.42  
##  1st Qu.:0.0000   1st Qu.:1.000   Class :character   1st Qu.:20.12  
##  Median :0.0000   Median :2.000   Mode  :character   Median :28.00  
##  Mean   :0.4062   Mean   :2.237                      Mean   :29.70  
##  3rd Qu.:1.0000   3rd Qu.:3.000                      3rd Qu.:38.00  
##  Max.   :1.0000   Max.   :3.000                      Max.   :80.00

Potom ćemo procijeniti decision tree model:

titanic_model <- rpart(Survived ~ ., data = titanic_dta, method = "class")

Sada možemo procijeniti preciznost modela:

pred <- predict(titanic_model, titanic_dta, type = "class")
table(titanic_dta$Survived, pred)
##    pred
##       0   1
##   0 371  53
##   1  78 212

Confusion matrica pokazuje da je model točno klasificirao 212 od 265 preživjelih i 371 od 449 preminulih.

Izračunajmo pokazatelje preciznosti procjene:

cfmtx <- table(titanic_dta$Survived, pred)

TP <- cfmtx[1, 1] #
FN <- cfmtx[1, 2] # 
FP <- cfmtx[2, 1] # 
TN <- cfmtx[2, 2] # 

# Točnost: acc
acc <- (TP + TN) / (TP + FN + FP + TN)
acc
## [1] 0.8165266
# Preciznost: prec
prec <- TP/(TP+FP)
prec
## [1] 0.8262806
# Recall: rec
rec <- TP/(TP+FN)
rec
## [1] 0.875

Točnost i preciznost modela su ~82% dok je model krivo predvidio sudbinu 18% putnika na Titaniku. Nije loše, ali nije ni posebno dobro…

Regresija

Kvaliteta regresijskog odnosa se procijenjuje kriterijem rezidualnih kvadrata odstupanja (RMSE- root mean square error) koje smo već spominjali:

\[RMSE = \sqrt{\frac{1}{n}\Sigma_{i=1}^{n}{({y_i -y\hat{}_i})^2}}\]

RMSE procjenjuje udaljenost između stvarnih i predviđenih opservacija, a cilj je postići što manji broj.

Praktični primjer

U ovom primjeru je cilj predvidjeti pritisak zvuka aviona (dec) na osnovi frekvencije vjetra (freq), kuta krila (angle) i nekih drugih varijabli (ch_length). Prvo ćemo učitati, preimenovati i pregledati podatke:

library(dplyr)
air <- read.table("../Dta/airfoil_self_noise.dat",header = T) %>% # dostupno na UCIMLR
dplyr::rename(freq = "X800",
              angle = "X0" ,
              ch_length = "X0.3048",
              velocity = "X71.3",
              thickness = "X0.00266337",
              dec = "X126.201" )
 
str(air)
## 'data.frame':    1502 obs. of  6 variables:
##  $ freq     : int  1000 1250 1600 2000 2500 3150 4000 5000 6300 8000 ...
##  $ angle    : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ ch_length: num  0.305 0.305 0.305 0.305 0.305 ...
##  $ velocity : num  71.3 71.3 71.3 71.3 71.3 71.3 71.3 71.3 71.3 71.3 ...
##  $ thickness: num  0.00266 0.00266 0.00266 0.00266 0.00266 ...
##  $ dec      : num  125 126 128 127 126 ...

Potom procjenjujemo model i predikciju pa računamo RMSE:

fit <- lm(dec ~ freq + angle + ch_length, data = air)
predd <- predict(fit,air)
rmse <- sqrt((1/nrow(air)) * sum( (air$dec - predd) ^ 2))
rmse
## [1] 5.217465

Rezultat je 5.2 decibela…to nam doduše ne govori mnogo.Za procjenu kvalitete modela je treba usporediti rmse nekog drugog modela. Upravo to ćemo sada učiniti. Prvo procjenjujemo novi model koji sadrži više kontrolnih varijabli:

fit2 <- lm(dec ~ freq + angle + ch_length + velocity + thickness, data = air)
predd2 <- predict(fit2,air)
rmse2 <- sqrt((1/nrow(air)) * sum( (air$dec - predd2) ^ 2))
rmse2
## [1] 4.800694

Novi i bogatiji model ima niži pripadajući rmse od 4.8 pa za možemo zaključiti da se radi o boljem modelu.

Klastering

U slučaju klastringa ne postoje explicitne varijable (labels) pa se kvaliteta procjene radi na osnovi mjera udaljenosti unutar i između klastera. Te se mjere odnose na sličnost unutar klastera (poželjno visoka WSS), sličnost između klastera (poželjno niska BSS) i Dunnov indeks.

Praktični primjer

Sljedeći primjer se koristi podatke o različitim karakteristikama sjemena no pošto ne znamo nazive vrijabli (labels) ali znamo da podatci sadržavaju tri vrste sjemena…koristiti ćemo klastering metodu. Prvo učitamo podatke i pregledamo strukturu:

seeds <- read.table("../Dta/seeds_dataset.txt",header = T)  # dostupno na UCIMLR
str(seeds)
## 'data.frame':    209 obs. of  8 variables:
##  $ X15.26: num  14.9 14.3 13.8 16.1 14.4 ...
##  $ X14.84: num  14.6 14.1 13.9 15 14.2 ...
##  $ X0.871: num  0.881 0.905 0.895 0.903 0.895 ...
##  $ X5.763: num  5.55 5.29 5.32 5.66 5.39 ...
##  $ X3.312: num  3.33 3.34 3.38 3.56 3.31 ...
##  $ X2.221: num  1.02 2.7 2.26 1.35 2.46 ...
##  $ X5.22 : num  4.96 4.83 4.8 5.17 4.96 ...
##  $ X1    : int  1 1 1 1 1 1 1 1 1 1 ...

Potom provedemo kmeans klastering sa tri klastera:

km_seeds <- kmeans(seeds[1:7], 3)

Pogledajmo i vizualizaciju rezultata:

plot(seeds$X5.763 ~ seeds$X0.871, data = seeds, col = km_seeds$cluster) # duljina ~ kompaktnost

Za kraj pogledajmo još i mjeru kvalitete klasterske procjene koja se računa kao omjer unutar- i među- klasterske sličnosti. Podatci su dostupni kao dio procijenjenog klasterskog objekta km_seeds:

km_seeds$tot.withinss/km_seeds$betweenss
## [1] 0.2752817

Unutar- klasterska sličnost je niža od među- klsterske sličnostišto znači da su klasteri dobro odijeljeni i kompaktni. To je jasno vidljivo i na vizualnom prikazu pa možemo zaključiti da klasteri dobro predstavljaju tri vrste sjemena.

Trening vs. test podjela

Bitna razlika između strojnog učenja i statistike je prediktivna vs deskriptivna moć modela. Modeli nadgledanog ML bi trebali imati veliku prediktivnu moć, tj. dobro predvidjeti odnos na novim, neviđenim opservacijama. Klasična statistika se bazira na kvaliteti modela s obzirom na dane (postojeće) podatke tj. model treba dobro opisati podatke (poznate opservacije).

Zbog toga procjena prediktivnog ML modela ne smije uključivati potpuni podatkovni skup nego samo trening dio podataka, a kvaliteta modela se procjenjuje na test dijelu podataka. Za pouzdanu procjenu je bitno napraviti ovakvu podjelu podataka i osigurati da među njima nema preklapanja. Ova podjela je važna samo kod nadgledanog ML. Kod ne-nadgledanog ML nije potrebno dijeliti podatke.

Podjela podataka izgleda ovako:

X varijabla se odnosi na opservacije, F varijable na karakteristike (label), a Y varijabla na klasu opservacija. Na training podatcima radimo procjenu modela, a prediktivnu kvalitetu provjeravamo na test dijelu podataka.

Postupak prati sljedeći hodogram:

Podjela podataka na training i test podatke je proizvoljna no postoji nekoliko standardnih i okvirnih pravila. Najčešće se podatci dijele u omjeru 3/1 imajući na umu da će model dati bolju procjenu ako je “treniran” na više podaka. Također je važno postići distribucijsku sličnost u podjeli podataka za klasifikacijsku analizu te “promiješati” (randomizirati) podatke u regresijskoj analizi. Ovdje valja spomenuti i cross-validation proceduru kojom se dodatno pospješuje efekt podjele podatka.

Praktični primjer

Prethodni primjer na titanic podatkovnom skupu nije uključivao podjelu podatka. To ćemo sada ispraviti s ciljem poboljšanja kvalitete procjene. Podatke ćemo podijeliti u omjeru 70/30. Podjelu radimo na sljedeći način:

# promješaj podatke
n <- nrow(titanic_dta)
shuffled <- titanic_dta[sample(n),]

# podijeli na train i test
train_indices <- 1:round(0.7 * n)
train <- shuffled[train_indices, ]
test_indices <- (round(0.7 * n) + 1):n
test <- shuffled[test_indices, ]

# pregledaj podatke
str(train)
## 'data.frame':    500 obs. of  4 variables:
##  $ Survived: int  0 0 1 0 0 1 1 0 0 0 ...
##  $ Pclass  : int  3 2 1 2 3 1 2 3 2 2 ...
##  $ Sex     : chr  "male" "male" "male" "male" ...
##  $ Age     : num  8 25 26 21 28 52 42 29 28 47 ...
str(test)
## 'data.frame':    214 obs. of  4 variables:
##  $ Survived: int  0 1 0 0 1 1 0 0 1 0 ...
##  $ Pclass  : int  2 1 2 3 1 1 3 3 2 3 ...
##  $ Sex     : chr  "male" "female" "male" "male" ...
##  $ Age     : num  35 62 36 16 28 21 40.5 33 62 9 ...

Razmotrimo ponovno otprije korišteni model na novim, podijeljenim titanic podatcima:

# procjeni model
tree_split <- rpart(Survived ~ ., train, method = "class")
# predviđanje modela
preddd <- predict(tree_split, test, type="class")
# izračunaj confusion matricu
conf <- table(test$Survived, preddd)
conf
##    preddd
##       0   1
##   0 125   5
##   1  42  42

Kvalitetu modela izračunamo na isti način kao u prethodnom slučaju:

TP <- conf[1, 1] #
FN <- conf[1, 2] # 
FP <- conf[2, 1] # 
TN <- conf[2, 2] # 

# Točnost: acc
acc <- (TP + TN) / (TP + FN + FP + TN)
acc
## [1] 0.7803738
# Preciznost: prec
prec <- TP/(TP+FP)
prec
## [1] 0.748503
# Recall: rec
rec <- TP/(TP+FN)
rec
## [1] 0.9615385

Rezultati ukazuju na nešto manju točnost i preciznost modela, ~76% i ~78% u odnosu na ~82% koju smo imali u prvom slučaju. Ovi rezultati predstavljaju realniju procjenu prediktivne moći modela.

Dodatni resursi