library(mlr)
library(tidyverse)
library(ggplot2)
library(fastDummies)
library(e1071)
library(corrplot)
setwd("D:/Studies/5 II/Duomenu gavyba/mine/BIG HOMEWORK")

Prieš pradėdamas darbą su duomenų rinkiniu, iš originalaus failo padariau du: pirmasis (kurį naudosiu šiame namų darbe) yra be praleistų reikšmių, o antrojo kiekvienoje eilutėje yra bent po vieną praleistą reikšmę.

Adult duomenų rinkinys

1. Duomenų charakterizavimas

Adult duomenų rinkinys yra skirtas klasifikavimo užduočiai spręsti, t.y. nustatyti, ar žmogus uždirba daugiau ar mažiau nei 50 tūkst. dolerių per metus.

https://archive.ics.uci.edu/ml/datasets/Adult

data1 <- read.csv("adult_be.csv", sep = ",")[,-1]
names(data1) <- c("age","workclass","fnlwgt","education","education_num","marital_status","occupation","relationship","race","sex","capital_gain","capital_loss","hours_per_week","native_country","income")
head(data1)
##   age        workclass fnlwgt education education_num        marital_status
## 1  50 Self_emp_not_inc  83311 Bachelors            13    Married_civ_spouse
## 2  38          Private 215646   HS_grad             9              Divorced
## 3  53          Private 234721      11th             7    Married_civ_spouse
## 4  28          Private 338409 Bachelors            13    Married_civ_spouse
## 5  37          Private 284582   Masters            14    Married_civ_spouse
## 6  49          Private 160187       9th             5 Married_spouse_absent
##          occupation  relationship  race    sex capital_gain capital_loss
## 1   Exec_managerial       Husband White   Male            0            0
## 2 Handlers_cleaners Not_in_family White   Male            0            0
## 3 Handlers_cleaners       Husband Black   Male            0            0
## 4    Prof_specialty          Wife Black Female            0            0
## 5   Exec_managerial          Wife White Female            0            0
## 6     Other_service Not_in_family Black Female            0            0
##   hours_per_week native_country income
## 1             13  United_States  <=50K
## 2             40  United_States  <=50K
## 3             40  United_States  <=50K
## 4             40           Cuba  <=50K
## 5             40  United_States  <=50K
## 6             16        Jamaica  <=50K

Šiame duomenų rinkinyje yra 15 atributų: 6-ių iš jų yra tolydžios reikšmės, o likę yra kategoriniai kintamieji. Atributai aprašo šias reikšmes:

  1. age - stebinio amžių.
  2. workclass - kuriai sociologinei padėčiai priklauso (Private, Self-emp-not-inc, Self-emp-inc, Federal-gov, Local-gov, State-gov, Without-pay, Never-worked).
  3. fnlwgt - sampling weight (the number of units in the target population that the responding unit represents).
  4. education - kokį išsilavinimą stebinys įgyjęs (Preschool,1st-4th, 5th-6th, 7th-8th, 9th, 10th, 11th, 12th, HS-grad, Prof-school, Some-college, Bachelors, Assoc-acdm, Assoc-voc, Masters, Doctorate).
  5. education_num - kiek metų praleista mokantis.
  6. marital_status - stebinio vedybinė padėtis (Married-civ-spouse, Divorced, Never-married, Separated, Widowed, Married-spouse-absent, Married-AF-spouse).
  7. occupation - stebinio pareigos darbe (Tech-support, Craft-repair, Other-service, Sales, Exec-managerial, Prof-specialty, Handlers-cleaners, Machine-op-inspct, Adm-clerical, Farming-fishing, Transport-moving, Priv-house-serv, Protective-serv, Armed-Forces).
  8. relationship - stebinio šeimyninė padėtis (Wife, Own-child, Husband, Not-in-family, Other-relative, Unmarried).
  9. race - stebinio rasė (White, Asian-Pac-Islander, Amer-Indian-Eskimo, Other, Black).
  10. sex - stebinio lytis (Female, Male).
  11. capital_gain - pajamos gautos iš kitų šaltinių nei darbo užmokęstis.
  12. capital_loss - pajamos prarastos iš kitų šaltinių nei darbo užmokęstis.
  13. hours_per_week - kiek valandų dirba per savaitę
  14. native_country - stebinio gimtoji šalis (United-States, Cambodia, England, Puerto-Rico, Canada, Germany, Outlying-US(Guam-USVI-etc), India, Japan, Greece, South, China, Cuba, Iran, Honduras, Philippines, Italy, Poland, Jamaica, Vietnam, Mexico, Portugal, Ireland, France, Dominican-Republic, Laos, Ecuador, Taiwan, Haiti, Columbia, Hungary, Guatemala, Nicaragua, Scotland, Thailand, Yugoslavia, El-Salvador, Trinadad&Tobago, Peru, Hong, Holand-Netherlands).
  15. income - ar per metus uždirba daugiau nei 50 tūkst. dolerių (<=50K, >50K)

Kiek metų praleista mokantis (education_num) pakankamai gerai gali atvaizduoti aukščiausią išsilavinimo lygį. Be to, fnlwgt nesusijęs su income kintamuoju. Dėl šių priežasčių pašalinsime education ir fnlwgt kintamuosius prieš kuriant modelį.

data1 <- data1[,-c(3:4)]
data1$income <- as.factor(data1$income)

Duomenų apdorojimas

Primiausia normalizavome visus skaitinius kintamuosius pasinaudodami min-max normalizavimo metodu:

\[\begin{align*} x_{ij_{norm}} = \displaystyle \frac{x_{ij}-x_{j_{min}}}{x_{ij_{max}}-x_{ij_{min}}}. \end{align*}\]

min.max.norm <- function(x, x.max, x.min)
{
  return((x-x.min)/(x.max-x.min))
}
for(i in c(1,9,10,11))
{
  max <- max(data1[,i])
  min <- min(data1[,i])
  for(ii in 1:nrow(data1))
  {
    data1[ii,i] <- min.max.norm(data1[ii,i], max, min)
  }
}

Toliau vizualiai tikriname kintamųjų reikšmes.

ggplot(data1, aes(x=workclass, fill=income)) +
  geom_bar()

temp <- NULL
for(i in 1:nrow(data1))
{
  ifelse(data1[i,2] == "Private", temp <- c(temp,1), temp <- c(temp,0))
}
data1$workclass <- temp
table(data1$workclass)
## 
##     0     1 
##  7875 22286

Kadangi grafike matome labai didelį vienos klasės kiekį kintamojame Private, nusprendžiau šį kintamąjį pakeisti į fiktyvų binarinį kintamąjį, kuris nusakys ar objekto workclass reikšmė yra Private, ar ne.

ggplot(data1, aes(x=marital_status, fill=income)) +
  geom_bar()

temp <- NULL
for(i in 1:nrow(data1))
{
  ifelse(data1[i,4] == "Married_civ_spouse", temp <- c(temp,0), temp <- c(temp,1))
}
data1$marital_status <- temp
table(data1$marital_status)
## 
##     0     1 
## 14065 16096

marital_status kintamąjį keisiu į fiktyvų binarinį kintamąjį, kuris nusakys ar marital_status reikšmė yra Married_civ_spouse. Taip darome, nes likusiuose kintamuosiuose vos ne išskirtinai tik vienos klasės kintamieji.

par(mfrow = c(1,2))
hist(data1$education_num)
hist(data1$hours_per_week)

Iš pateiktų education_num ir hours_per_week histogramų, matome, kad yra labai dažnai pasikartojančių reikšmių, bet su jomis nieko papildomai daryti nereikia.

ggplot(data1, aes(x=occupation, fill=income)) +
  geom_bar()

temp <- NULL
for(i in 1:nrow(data1))
{
  ifelse(data1[i,5] == "Craft_repair" |
         data1[i,5] == "Exec_managerial" |
         data1[i,5] == "Prof_specialty" |
         data1[i,5] == "Sales"
         , temp <- c(temp,1), temp <- c(temp,0))
}
data1$occupation <- temp
table(data1$occupation)
## 
##     0     1 
## 14517 15644

Iš occupation reikšmių dažnių stulpelinės diagramos matome, kad šios kokybinės reikšmės yra ganėtinai tolygiai pasiskirsčiusios. Todėl pamėginsime sugrupuoti šiuos kintamuosius tokia tvarka: keturis kintamuosius, kurių priklausymas mėlynai (>50K) grupei yra didžiausias sujungsime į vieną grupę, o likusius į kitą. Tuomet occupation binarinis kintamasis nusakys ar priklauso pirmąjai grupei ar antrajai.

ggplot(data1, aes(x=relationship, fill=income)) +
  geom_bar()

temp <- NULL
for(i in 1:nrow(data1))
{
  ifelse(data1$relationship[i] == "Husband", temp <- c(temp,0), temp <- c(temp,1))
}
data1$relationship <- temp

relationship kintamojo reikšmių dažnių stulpelinėje diagramoje, taip pat, matome, kad pasiskirstymas yra pakankamai tolygus. Visgi Husband kintamasis vos ne išimtinai turi visas >50K klasei priklausančias reikšmes, todėl relationship kintamąjį pakeisime į binarinį fiktyvų kintamąjį, kuris nusakys ar relationship reikšmė yra Husband.

ggplot(data1, aes(x=race, fill=income)) +
  geom_bar()

for(i in 1:nrow(data1))
{
  if(data1[i,7] != "White" & data1[i,7] != "Black") {data1[i,7] <- "Other"}
}
data1 <- dummy_columns(data1, select_columns = "race")
data1 <- data1[,-ncol(data1)]

Analizuojant race kintamojo reikšmes, nuspręsta palikti tris grupes (White,Black,Other). Tuomet sudarysime \(n-1\) fiktyvių binarinių kintamųjų, kai \(n\) yra skirtingų race kintamojo reikšmių kiekis, t.y. šiuo atveju \(n = 3\).

ggplot(data1, aes(x=sex, fill=sex)) +
  geom_bar()

data1 <- dummy_columns(data1, select_columns = "sex")
data1 <- data1[,-ncol(data1)]

Kintamasis sex gerai pasiskirstęs ir turi tik dvi reikšmes, todėl lengvai pakeičiamas į fiktyvų binarinį kintamąjį.

par(mfrow = c(1,2))
hist(data1$capital_gain)
hist(data1$capital_loss)

Kintamieji captial_gain ir capital_loss yra labai prastai pasiskirstę todėl neša labai mažai informacijos. Dėl to šiuos kintamuosius pašalinsime.

ggplot(data1, aes(x=native_country, fill=income)) +
  geom_bar()

Ne išimtis ir native_country kintamasis, kurio nešama informacija yra praktiškai beprasmė šiame tyrime. Jis taip pat pašalinamas iš duomenų rinkinio.

str(data1)
## 'data.frame':    30161 obs. of  11 variables:
##  $ age           : num  0.452 0.288 0.493 0.151 0.274 ...
##  $ workclass     : num  0 1 1 1 1 1 0 1 1 1 ...
##  $ education_num : Ord.factor w/ 16 levels "1"<"2"<"3"<"4"<..: 13 9 7 13 14 5 9 14 13 10 ...
##  $ marital_status: num  0 1 0 0 0 1 0 1 0 0 ...
##  $ occupation    : num  1 0 0 1 1 0 1 1 1 1 ...
##  $ relationship  : num  0 1 0 1 1 1 0 1 0 0 ...
##  $ hours_per_week: num  0.122 0.398 0.398 0.398 0.398 ...
##  $ income        : Factor w/ 2 levels "<=50K",">50K": 1 1 1 1 1 1 2 2 2 2 ...
##  $ race_Black    : int  0 0 1 1 0 1 0 0 0 1 ...
##  $ race_Other    : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ sex_Female    : int  0 0 0 1 1 1 0 1 0 0 ...
corrplot::corrplot(cor(data1[,-c(3,8)]), type = "upper")

Kaip galima pastebėti iš kintamųjų koreliacijų diagramos, stipresnė koreliacija galima pastebėti tarp marital_status ir relationship kintamųjų, bet kadangi nekilo problemų dėl kolinearumo, kol kas jo nei vieno nešalinsime.

Pasirinktiems duomenims sudarykite LDA klasifikatorių pagal pavyzdį 128-131 psl. iš 5 skyriaus “Machine learning with R, tidyverse and mlr”.

data1_Task <- makeClassifTask(data = data1[,-3], target = "income")
lda <- makeLearner("classif.lda")
ldaModel <- train(lda, data1_Task)

ldaModelData <- getLearnerModel(ldaModel)
ldaPreds <- predict(ldaModelData)$x

head(ldaPreds)

df <- cbind(data1, ldaPreds)

ggplot(df, aes(x=LD1, fill=income, color=income)) +
  geom_histogram()

qda <- makeLearner("classif.qda")
qdaModel <- train(qda, data1_Task)

kFold <- makeResampleDesc(method = "RepCV", folds = 2, reps = 5,
                          stratify = TRUE)
set.seed(123)
ldaCV <- resample(learner = lda, task = data1_Task, resampling = kFold,
                  measures = list(mmce, acc))
qdaCV <- resample(learner = qda, task = data1_Task, resampling = kFold,
                  measures = list(mmce, acc))
ldaCV$aggr
qdaCV$aggr

calculateConfusionMatrix(ldaCV$pred, relative = TRUE)
calculateConfusionMatrix(qdaCV$pred, relative = TRUE)
##          LD1
## 1  0.9223140
## 2 -1.1575162
## 3  0.5133259
## 4  1.1576762
## 5  1.4416336
## 6 -1.7045853

## mmce.test.mean  acc.test.mean 
##      0.2052571      0.7947429
## mmce.test.mean  acc.test.mean 
##      0.2475993      0.7524007
## Relative confusion matrix (normalized by row/column):
##         predicted
## true     <=50K     >50K      -err.-   
##   <=50K  0.87/0.86 0.13/0.41 0.13     
##   >50K   0.42/0.14 0.58/0.59 0.42     
##   -err.-      0.14      0.41 0.21     
## 
## 
## Absolute confusion matrix:
##         predicted
## true      <=50K   >50K -err.-
##   <=50K  982345 150305 150305
##   >50K   159233 216167 159233
##   -err.- 159233 150305 309538
## Relative confusion matrix (normalized by row/column):
##         predicted
## true     <=50K     >50K      -err.-   
##   <=50K  0.75/0.90 0.25/0.50 0.25     
##   >50K   0.25/0.10 0.75/0.50 0.25     
##   -err.-      0.10      0.50 0.25     
## 
## 
## Absolute confusion matrix:
##         predicted
## true      <=50K   >50K -err.-
##   <=50K  852536 280114 280114
##   >50K    93278 282122  93278
##   -err.-  93278 280114 373392

Kaip matome iš rezultatų, po 500 iteracijų LDA metodo tikslumas \(\sim 79\%\), o QDA \(\sim 75\%\). Nors tie keturi procentai atrodo nedaug, iš confusion matrix matome, kad per tuos 500 iteracijų QDA metodas \(\sim 64\) tūkstančiais daugiau kartų neteisingai suklasifikavo nematytus objektus.

Padarykite vieną “hold-out” testą ir išveskite klasifikavimo lenteles “confusion matrix” mokymo ir testavimo duomenims.

data1_Task <- makeClassifTask(data = data1[,-3], target = "income")

lda <- makeLearner("classif.lda")
holdout <- makeResampleDesc(method = "Holdout", split = 4/5, stratify = TRUE)
set.seed(123)
holdoutCV_lda <- resample(learner = lda, task = data1_Task, resampling = holdout, measures = list(mmce, acc))
## Resampling: holdout
## Measures:             mmce      acc
## [Resample] iter 1:    0.2181336 0.7818664
## 
## Aggregated Result: mmce.test.mean=0.2181336,acc.test.mean=0.7818664
## 
calculateConfusionMatrix(holdoutCV_lda$pred, relative = TRUE)
## Relative confusion matrix (normalized by row/column):
##         predicted
## true     <=50K     >50K      -err.-   
##   <=50K  0.86/0.85 0.14/0.44 0.14     
##   >50K   0.45/0.15 0.55/0.56 0.45     
##   -err.-      0.15      0.44 0.22     
## 
## 
## Absolute confusion matrix:
##         predicted
## true     <=50K >50K -err.-
##   <=50K   3897  634    634
##   >50K     682  820    682
##   -err.-   682  634   1316

Holdout metodui padalinome mokymo imtį į penkias dalis ir keturias iš jų naudojome mokymui bei vieną testavimui.

Gauname, kad LDA modelio tikslumas \(\sim 78\%\), kai naudojame \(4/5\) duomenų mokymui. \(1316\) iš \(6033\) objektų klasifikuota neteisingai. Visgi pastebime, kad mažas skirtumas tarp klaidų nustatant klases skaičiaus, t.y. spėja vienodai klaidingai.

Automobiles duomenų rinkinys

1. Duomenų charakterizavimas

Automobiles duomenų rinkinys buyo naudojamas regresijos užduočiai spręsti, t.y. prognozuoti automobilio kainą pagal skaitines ir binarines reikšmes. Kadangi šiame uždavinyje yra nemažai kategorinių reikšmių, galima jų pagalba sudaryti ir klasifikavimo uždavinį.

https://archive.ics.uci.edu/ml/datasets/Automobile

data2 <- read.csv("auto_be.csv", sep = ",")[,-1]
head(data2)
##   symboling normalized_losses make fuel_type aspiration num_of_doors body_style
## 1         2               164 audi       gas        std         four      sedan
## 2         2               164 audi       gas        std         four      sedan
## 3         1               158 audi       gas        std         four      sedan
## 4         1               158 audi       gas      turbo         four      sedan
## 5         2               192  bmw       gas        std          two      sedan
## 6         0               192  bmw       gas        std         four      sedan
##   drive_wheels engine_location wheel_base length width height curb_weight
## 1          fwd           front       99.8  176.6  66.2   54.3        2337
## 2          4wd           front       99.4  176.6  66.4   54.3        2824
## 3          fwd           front      105.8  192.7  71.4   55.7        2844
## 4          fwd           front      105.8  192.7  71.4   55.9        3086
## 5          rwd           front      101.2  176.8  64.8   54.3        2395
## 6          rwd           front      101.2  176.8  64.8   54.3        2395
##   engine_type num_of_cylinders engine_size fuel_system bore stroke
## 1         ohc             four         109        mpfi 3.19    3.4
## 2         ohc             five         136        mpfi 3.19    3.4
## 3         ohc             five         136        mpfi 3.19    3.4
## 4         ohc             five         131        mpfi 3.13    3.4
## 5         ohc             four         108        mpfi 3.50    2.8
## 6         ohc             four         108        mpfi 3.50    2.8
##   compression_ratio horsepower peak_rpm city_mpg highway_mpg price
## 1              10.0        102     5500       24          30 13950
## 2               8.0        115     5500       18          22 17450
## 3               8.5        110     5500       19          25 17710
## 4               8.3        140     5500       17          20 23875
## 5               8.8        101     5800       23          29 16430
## 6               8.8        101     5800       23          29 16925

Šiame duomenų rinkinyje yra 159 objektai su 26 atributai: 15 iš jų yra tolydžios reikšmės, 1 sveikojo skaičiaus ir likę 10 yra kategoriniai kintamieji.

  1. symboling - draudimo rizikos lygis {-3; -2; -1; 0; 1; 2; 3}.
  2. normalized_losses - santykinis vidutinis nuostolių atlyginimas per apdraustos transporto priemonės metus continuous from 65 to 256.
  3. make - automobilio markė (alfa-romero, audi, bmw, chevrolet, dodge, honda, isuzu, jaguar, mazda, mercedes-benz, mercury, mitsubishi, nissan, peugot, plymouth, porsche, renault, saab, subaru, toyota, volkswagen, volvo)
  4. fuel_type - kuro tipas (diesel, gas).
  5. aspiration - oro tiekimo varikliui tipas (std, turbo).
  6. num_of_doors - durų kiekis (four, two).
  7. body_style - kėbulo tipas (hardtop, wagon, sedan, hatchback, convertible).
  8. drive_wheels - kurie ratai varomieji (4wd, fwd, rwd).
  9. engine_location - kur variklis (front, rear).
  10. wheel_base - continuous from 86.6 120.9.
  11. length - automobilio ilgis continuous from 141.1 to 208.1.
  12. width - automobilio plotis continuous from 60.3 to 72.3.
  13. height - automobilio aukštis continuous from 47.8 to 59.8.
  14. curb_weight - continuous from 1488 to 4066.
  15. engine_type - variklio tipas (dohc, dohcv, l, ohc, ohcf, ohcv, rotor).
  16. num_of_cylinders - cilindrų kiekis variklyje (eight, five, four, six, three, twelve, two).
  17. engine_size - variklio dydis continuous from 61 to 326.
  18. fuel_system - kuro padavimo sistema (1bbl, 2bbl, 4bbl, idi, mfi, mpfi, spdi, spfi).
  19. bore - vidinio cilindro skersmuo continuous from 2.54 to 3.94.
  20. stroke - stūmuoklio etapų kiekis reikalingas pasukti alkūninį veleną vienam ciklui continuous from 2.07 to 4.17.
  21. compression_ratio - vidaus degimo variklio degimo kameros ir cilindro tūrio santykis continuous from 7 to 23.
  22. horsepower - automobilio arklio galių kiekis continuous from 48 to 288.
  23. peak_rpm - maksimalus apsukų kiekis per minutę continuous from 4150 to 6600.
  24. city_mpg - automobilio pritaikymo miesto sąlygoms įvertinimas continuous from 13 to 49.
  25. highway_mpg - automobilio pritaikymo greitkelio sąlygoms įvertinimas continuous from 16 to 54.
  26. price - automobilio kaina continuous from 5118 to 45400.

Duomenų apdorojimas

apply(data2[,c(4,5,6,7,8,9,15,16,18)], 2, table)
## $fuel_type
## 
## diesel    gas 
##     15    144 
## 
## $aspiration
## 
##   std turbo 
##   132    27 
## 
## $num_of_doors
## 
## four  two 
##   95   64 
## 
## $body_style
## 
## convertible     hardtop   hatchback       sedan       wagon 
##           2           5          56          79          17 
## 
## $drive_wheels
## 
## 4wd fwd rwd 
##   8 105  46 
## 
## $engine_location
## 
## front 
##   159 
## 
## $engine_type
## 
## dohc    l  ohc ohcf ohcv 
##    8    8  123   12    8 
## 
## $num_of_cylinders
## 
## eight  five  four   six three 
##     1     7   136    14     1 
## 
## $fuel_system
## 
## 1bbl 2bbl  idi  mfi mpfi spdi 
##   11   63   15    1   64    5
data2$num_of_doors <- as.factor(data2$num_of_doors)
data2 <- data2[,-c(1,9)]

Iš kokybinių kintamųjų dažnių lentelės matome, kad geriausiai pasiskirstęs kintamasis yra num_of_doors, todėl siekdami sukurti kokybišką ML modelį, klasifikuosime šį kintamąjį. Be to, iš karto galime pastebėti, kad kintamasis engine_location turi tik vieną reikšmę (neneša informacijos), todėl jau šiame etape jį pašalinsime. Pašalinsime ir symboling kintamąjį, kuris netinka LDA ir QDA metodams, nes yra ranguotas kokybinis kintamasis.

Tuomet atlikome skaitinių kintamųjų normalizavimą, kad išvengti skirtingose skalėse esančių kintamųjų persvaros. Normalizavome pasinaudodami min-max normalizavimo metodu:

\[\begin{align*} x_{ij_{norm}} = \displaystyle \frac{x_{ij}-x_{j_{min}}}{x_{ij_{max}}-x_{ij_{min}}}. \end{align*}\]

par(mfrow = c(2,4))
hist(data2$normalized_losses, main = "normalized_losses")
hist(data2$wheel_base, main = "wheel_base")
hist(data2$length, main = "length")
hist(data2$width, main = "width")

hist(data2$height, main = "height")
hist(data2$curb_weight, main = "curb_weight")
hist(data2$engine_size, main = "engine_size")
hist(data2$bore, main = "bore")

hist(data2$stroke, main = "stroke")
hist(data2$compression_ratio, main = "compression_ratio")
hist(data2$horsepower, main = "horsepower")
hist(data2$peak_rpm, main = "peak_rpm")

hist(data2$city_mpg, main = "city_mpg")
hist(data2$highway_mpg, main = "highway_mpg")
hist(data2$price, main = "price")

Iš historgramų matyti, kad tik compression_ratio kintamasis neneša informacijos, todėl jį pašaliname.

Tiesinei priklausomybei tarp kiekybinių kintamųjų patikrinti naudosime Pirsono koreliacijos koeficientą (angl. Pearson corelation coefficient):

\[\begin{align*} \rho_{X,Y} = \frac{cov(X,Y)}{\sigma_X \sigma_Y}; \end{align*}\]

kur \(cov\) yra kovariacija, \(\sigma_X\) yra \(X\)-o standartinis nuokrypis ir \(\sigma_Y\) - \(Y\) standartinis nuokrypis.

corrplot(cor(data2[,c(1,8:12,15,18,20:24)]), method = "pie", type = "upper")

data2 <- data2[,-c(9,10,12,15,19,20,22)]
corrplot(cor(data2[,c(1,8:9,13:17)]), method = "number", type = "upper")

Siekiant sumažinti tą pačią informaciją nešančių kintamųjų kiekį (tuo pačiu sumažinant modelio šališkumą), iš modelio pašalinome labai aukštą koreliaciją \((\rho_{X,Y} > 0.8)\) turinčius kintamuosius, t.y. pašalinome length, width, curb_weight, engine_size, compression_ratio, horsepower ir city_mpg.

Toliau vizualiai tikriname kokybinius kintamuosius.

Autoriaus pastaba: tik padarius visus kokybinius kintamuosius į fiktyvius binarinius (grupavimo pagalba) pavyko padaryti kad QDA nemestų “rank defficiency” error’o.

ggplot(data = data2,  aes(x=make, fill=num_of_doors)) +
  geom_bar()

data2 <- data2[,-2]

Kadangi automobilių markių yra labai daug, o klasių pasiskirstymas per jas yra labai išsibarstęs - negalime protingai jų sugrupuoti.

ggplot(data = data2,  aes(x=fuel_type, fill=num_of_doors)) +
  geom_bar()

temp <- NULL
for(i in 1:nrow(data2))
{
  ifelse(data2$fuel_type[i] == "gas", 
         temp <- c(temp,1), 
         temp <- c(temp,0))
}
data2$fuel_type <- temp

Kadangi kintamasis fuel_type turi tik dvi reikšmes pakeitėme jas į fiktyvų binarinį kintamąjį, kuris nurodo ar kuro tipas yra gas, ar ne.

ggplot(data = data2,  aes(x=aspiration, fill=num_of_doors)) +
  geom_bar()

temp <- NULL
for(i in 1:nrow(data2))
{
  ifelse(data2$aspiration[i] == "std", 
         temp <- c(temp,1), 
         temp <- c(temp,0))
}
data2$aspiration <- temp

Kaip ir ankstesniam kintamąjam, aspiration kintamąjį pakeitėm į fiktyvų binarinį kintamąjį, kuris nurodo, ar oro tiekimas automobiliui yra std tipo, ar ne.

ggplot(data = data2,  aes(x=num_of_cylinders, fill=num_of_doors)) +
  geom_bar()

for(i in 1:nrow(data2))
{
  if(data2$num_of_cylinders[i] == "three"){data2$num_of_cylinders[i] <- 3}
  if(data2$num_of_cylinders[i] == "four"){data2$num_of_cylinders[i] <- 4}
  if(data2$num_of_cylinders[i] == "five"){data2$num_of_cylinders[i] <- 5}
  if(data2$num_of_cylinders[i] == "six"){data2$num_of_cylinders[i] <- 6}
  if(data2$num_of_cylinders[i] == "eight"){data2$num_of_cylinders[i] <- 8}  
}
data2$num_of_cylinders <- as.numeric(data2$num_of_cylinders)

Pastebėjau, kad kintamojo num_of_cylinders reikšmės yra paprasčiausiai skaičiai užrašyti žodžiais, tad atkeičiau atgal į skaičius.

ggplot(data = data2,  aes(x=body_style, fill=num_of_doors)) +
  geom_bar()

for(i in 1:nrow(data2))
{
  ifelse(data2$body_style[i] == "sedan",
     data2$body_style[i] <- 1,
     data2$body_style[i] <- 0)
}
data2$body_style <- as.numeric(data2$body_style)

body_style kintamąjam nusprendžiau grupuoti, pagal tai ar priklauso ar nepriklauso populiariausiam variantui, t.y. ar reikšmė lygi sedan, ar ne.

ggplot(data = data2,  aes(x=drive_wheels, fill=num_of_doors)) +
  geom_bar()

for(i in 1:nrow(data2))
{
  ifelse(data2$drive_wheels[i] == "fwd",
         data2$drive_wheels[i] <- 1,
         data2$drive_wheels[i] <- 0)
}
data2$drive_wheels <- as.numeric(data2$drive_wheels)

Analogiškai ankstesnio kintamojo atvejui, drive_wheels kintamasis pakeistas į fiktyvų binarinį kintamąjį kuris nurodo ar reikšmė yra fwd, ar ne.

ggplot(data = data2,  aes(x=engine_type, fill=num_of_doors)) +
  geom_bar()

for(i in 1:nrow(data2))
{
  ifelse(data2$engine_type[i] != "ohc",
    data2$engine_type[i] <- 0,
    data2$engine_type[i] <- 1)
}
data2$engine_type <- as.numeric(data2$engine_type)

Kaip ir ankstesniem dviem kintamiesiams, engine_type grupuoju pagal dažniausią reikšmę ohc. Taigi fiktyvus binarinis kintamasis nurodo ar reikšmė yra ohc, ar ne.

ggplot(data = data2,  aes(x=fuel_system, fill=num_of_doors)) +
  geom_bar()

for(i in 1:nrow(data2))
{
  ifelse(data2$fuel_system[i] == "2bbl",
         data2$fuel_system[i] <- 1,
         data2$fuel_system[i] <- 0)
}
data2$fuel_system <- as.numeric(data2$fuel_system)

Dėl problemų su QDA, išbandžiau įvairius grupavimo variantus, ir galiausiai suveikė fuel_system kintamojo reikšmes sugrupavus į dvi grupęs: arba reikšmė lygi 2bbl, arba ne.

Pasirinktiems duomenims sudarykite LDA klasifikatorių pagal pavyzdį 128-131 psl. iš 5 skyriaus “Machine learning with R, tidyverse and mlr”.

data2_Task <- makeClassifTask(data = data2, target = "num_of_doors")
lda <- makeLearner("classif.lda")
ldaModel <- train(lda, data2_Task)


ldaModelData <- getLearnerModel(ldaModel)
ldaPreds <- predict(ldaModelData)$x
head(ldaPreds)

df <- cbind(data2, ldaPreds)

ggplot(df, aes(x=LD1, fill=num_of_doors, color=num_of_doors)) +
  geom_histogram()

qda <- makeLearner("classif.qda")
qdaModel <- train(qda, data2_Task)

kFold <- makeResampleDesc(method = "RepCV", folds = 10, reps = 50,
                          stratify = TRUE)
ldaCV <- resample(learner = lda, task = data2_Task, resampling = kFold,
                  measures = list(mmce, acc))
set.seed(1232131231)
qdaCV <- resample(learner = qda, task = data2_Task, resampling = kFold,
                  measures = list(mmce, acc))
ldaCV$aggr
qdaCV$aggr

calculateConfusionMatrix(ldaCV$pred, relative = TRUE)
calculateConfusionMatrix(qdaCV$pred, relative = TRUE)
##           LD1
## 1 -0.31267527
## 2 -0.03460106
## 3 -1.26614105
## 4 -1.38193057
## 5  0.40131860
## 6  0.40493861

## mmce.test.mean  acc.test.mean 
##      0.1874495      0.8125505
## mmce.test.mean  acc.test.mean 
##      0.2270441      0.7729559
## Relative confusion matrix (normalized by row/column):
##         predicted
## true     four      two       -err.-   
##   four   0.88/0.82 0.12/0.20 0.12     
##   two    0.29/0.18 0.71/0.80 0.29     
##   -err.-      0.18      0.20 0.19     
## 
## 
## Absolute confusion matrix:
##         predicted
## true     four  two -err.-
##   four   4180  570    570
##   two     921 2279    921
##   -err.-  921  570   1491
## Relative confusion matrix (normalized by row/column):
##         predicted
## true     four      two       -err.-   
##   four   0.84/0.79 0.16/0.26 0.16     
##   two    0.32/0.21 0.68/0.74 0.32     
##   -err.-      0.21      0.26 0.23     
## 
## 
## Absolute confusion matrix:
##         predicted
## true     four  two -err.-
##   four   3970  780    780
##   two    1025 2175   1025
##   -err.- 1025  780   1805

Kaip matome iš rezultatų, po 500 iteracijų LDA metodo tikslumas \(\sim 81\%\), o QDA \(\sim 77\%\). Iš confusion matrix matome, kad QDA metodas \(\sim 300\)-ais objektų daugiau suklasifikavo neteisingai, nei LDA metodas.

Padarykite vieną “hold-out” testą ir išveskite klasifikavimo lenteles “confusion matrix” mokymo ir testavimo duomenims.

data2_Task <- makeClassifTask(data = data2, target = "num_of_doors")
lda <- makeLearner("classif.lda")
holdout <- makeResampleDesc(method = "Holdout", split = 4/5, stratify = TRUE)
set.seed(123)
holdoutCV_lda <- resample(learner = lda, 
                          task = data2_Task, 
                          resampling = holdout, 
                          measures = list(mmce, acc))
## Resampling: holdout
## Measures:             mmce      acc
## [Resample] iter 1:    0.1875000 0.8125000
## 
## Aggregated Result: mmce.test.mean=0.1875000,acc.test.mean=0.8125000
## 
calculateConfusionMatrix(holdoutCV_lda$pred, relative = TRUE)
## Relative confusion matrix (normalized by row/column):
##         predicted
## true     four      two       -err.-   
##   four   0.95/0.78 0.05/0.11 0.05     
##   two    0.38/0.22 0.62/0.89 0.38     
##   -err.-      0.22      0.11 0.19     
## 
## 
## Absolute confusion matrix:
##         predicted
## true     four two -err.-
##   four     18   1      1
##   two       5   8      5
##   -err.-    5   1      6

Holdout metodui padalinome mokymo imtį į penkias dalis ir keturias iš jų naudojome mokymui bei vieną testavimui.

Gauname, kad LDA modelio tikslumas \(\sim 81\%\), kai naudojame \(4/5\) duomenų mokymui. \(6\) iš \(32\) objektų klasifikuota neteisingai.