Pages

Thursday, 18 November 2021

Naive Bayes Classifier: teori dan contoh R

 Download data:

  1. data link: https://stats.idre.ucla.edu/stat/data/hsbdemo.dta


Baca data
Install library "haven" untuk membaca format data ".dta" yang merupakan format data untuk Stata.

install.packages("haven")
library(haven)

setwd("D:/") # setting directory
data <- read_dta('hsbdemo.dta')

Data

Untuk tujuan demonstrasi, kita akan membuat classifier Niave Bayes di sini. Kita akan menggunakan kumpulan data yang berisi informasi 200 siswa. Skor mereka dalam mata pelajaran yang berbeda dan pilihan pendidikan mereka (umum, akademik atau kejuruan). Ada variabel lain yang menunjukkan status sosial ekonomi dan jenis kelamin mereka. Kami akan membuat classifier Naive Bayer di sini.

Mari kita lihat data. Dalam hal ini pro (general = 1, academic =2, vocational = 3).

data <- as.data.frame(data)
head(data)
## female ses prog read write math science socst honors awards cid ## 1 1 1 vocational 34 35 41 29 26 0 0 1 ## 2 0 2 general 34 33 41 36 36 0 0 1 ## 3 0 3 vocational 39 39 44 26 42 0 0 1 ## 4 0 1 vocational 37 37 42 33 32 0 0 1 ## 5 0 2 vocational 39 31 40 39 51 0 0 1 ## 6 1 3 general 42 36 42 31 39 0 0 1
## Correlation matrix within general
datag <- subset(data, prog == 1)
cor(datag[,6:10])
##              read     write      math   science     socst
## read    1.0000000 0.4739121 0.3945974 0.6586988 0.5418732
## write   0.4739121 1.0000000 0.3586417 0.5629392 0.6505204
## math    0.3945974 0.3586417 1.0000000 0.5752819 0.3787115
## science 0.6586988 0.5629392 0.5752819 1.0000000 0.4222026
## socst   0.5418732 0.6505204 0.3787115 0.4222026 1.0000000
## Correlation matrix within academic
datag <- subset(data, prog == 2)
cor(datag[,6:10])
##              read     write      math   science     socst
## read    1.0000000 0.5608413 0.6917634 0.6250391 0.5851566
## write   0.5608413 1.0000000 0.6130255 0.5128848 0.4538175
## math    0.6917634 0.6130255 1.0000000 0.6410174 0.4591657
## science 0.6250391 0.5128848 0.6410174 1.0000000 0.4383806
## socst   0.5851566 0.4538175 0.4591657 0.4383806 1.0000000
## Correlation matrix within vocational
datag <- subset(data, prog == 3)
cor(datag[,6:10])
##              read     write      math   science     socst
## read    1.0000000 0.4615702 0.4570520 0.5132068 0.4325037
## write   0.4615702 1.0000000 0.5090928 0.5225355 0.4926333
## math    0.4570520 0.5090928 1.0000000 0.5706508 0.3769207
## science 0.5132068 0.5225355 0.5706508 1.0000000 0.3348232
## socst   0.4325037 0.4926333 0.3769207 0.3348232 1.0000000

Kita melihat variabel numerik sangat berkorelasi. Tetapi Naive Bayes mengambil kondisi bebas bersyarat. Jadi untuk tujuan demonstrasi, kita hanya akan memilih variabel sciencedan socst untuk membuat classifier kita.

 Statistika Deskriptif: visualisasi dengan plot

Mari kita lihat apakah prediktor pilihan kita memiliki kekuatan diskriminatif. Di bawah ini kita dapat memplot boxplot dan densityplot dari prediktor yang dipilih untuk setiap kategori target.

# Visualisasi data
# boxplot
# one box per variety

# format data untuk mengubah kategori prog=(1,2,3) 
# menjadi program=(general,academic,covactiona) 
head(data)
data$program <- "0"
data$idx <- seq.int(nrow(data))
data$program[(subset(data, prog ==1))$idx] <- "general"
data$program[(subset(data, prog ==2))$idx] <- "academic"
data$program[(subset(data, prog ==3))$idx] <- "vocational"

#install.packages("ggplot2")
library(ggplot2)
s <- ggplot(data, aes(x=program, y=science, fill=program)) + 
  geom_boxplot(alpha=.3) + ggtitle("Science Score boxplot")

ss <- ggplot(data, aes(x=program, y=socst, fill=program)) + 
  geom_boxplot(alpha=.3) + ggtitle("Socst Score boxplot")

ds <- ggplot(data, aes(x=science, fill=program)) + 
  geom_density(adjust=1.5, alpha=.3) + ggtitle("Science Score density")

dss <- ggplot(data, aes(x=socst, fill=program)) + 
  geom_density(adjust=1.5, alpha=.3) + ggtitle("Socst Score density")

#install.packages("ggpubr")
library(ggpubr)
ggarrange(s, ss, ds, dss + rremove("x.text"), 
          labels = c("A", "B", "C", "D"),
          ncol = 2, nrow = 2)

Kita mengamati beberapa hal dari visualisasi:

  1. Kita memiliki beberapa kekuatan diskriminatif untuk variabel yang dipilih.

  2. Sebagian besar distribusi mendekati normal. Setidaknya, mereka tidak terlalu jauh dari normalitas. Mari kita jalankan tes formal (tes Anderson Darling) untuk normalitas bersyarat.

Tes Anderson-Darling adalah tes kesesuaian yang menentukan seberapa cocok data Anda dengan distribusi tertentu. Tes ini paling sering digunakan untuk melihat apakah data Anda mengikuti distribusi normal atau tidak.

 Tes semacam ini dapat digunakan untuk memeriksa normalitas, yang merupakan asumsi umum dalam banyak uji statistik seperti regresi, ANOVA, dan uji-t.

 Menghitung statistik uji Anderson–Darling untuk sampel yang dipilih dari distribusi tertentu dan menentukan apakah akan menolak atau menerima hipotesis bahwa sampel diambil dari distribusi tersebut.
# Visualisasi data
library(nortest)
GenScience <- ad.test((subset(data, prog ==1))$science)
GenSocst <- ad.test((subset(data, prog ==1))$socst)
AcaScience <- ad.test((subset(data, prog ==2))$science)
AcaSocst <- ad.test((subset(data, prog ==2))$socst)
VocScience <- ad.test((subset(data, prog ==3))$science)
VocSocst <- ad.test((subset(data, prog ==3))$socst)

p <-data.frame(GenScience$p.value,GenSocst$p.value,AcaScience$p.value,AcaSocst$p.value,VocScience$p.value,VocSocst$p.value)
colnames(p) <- c("GenScience","GenSocst","AcaScience","AcaSocst","VocScience","VocSocst")
p
## P values of Anderson Darling Test
##      GenScience GenSoct AcaScience AcaSoct VocScience VocSoct
## [1,]       0.48     0.1       0.01       0       0.38     0.1

Jadi kita melihat skor kelompok akademik secara statistik berbeda dari kelompok lainnya. Tetapi untuk 4 lainnya, kita memiliki nilai p yang tidak signifikan (> 0.05), menunjukkan bahwa kelompok ini tidak terdistribusi normal.

Membuat pengklasifikasi Naive Bayes

Sekarang kita akan membuat classifier Naive Bayes untuk data kita. Kami akan membuat partisi 70/30 untuk set pelatihan dan validasi. Kami menggunakan paket createDataPartitionfrom caretuntuk membuat partisi seimbang. Ini sangat penting dalam konteks ini karena kita akan menghitung probabilitas sebelumnya dari jumlah data pelatihan. Jadi, itu harus mengikuti proporsi dataset induk.

install.packaes(caret)
library(caret)
## Warning: package 'caret' was built under R version 3.4.3
set.seed(7267166)
trainIndex=createDataPartition(data$prog, p=0.7)$Resample1
train=data[trainIndex, ]
test=data[-trainIndex, ]

## check the balance
print(table(data$prog))
## 
##   academic    general vocational 
##        105         45         50
print(table(train$prog))

## 
##   academic    general vocational 
##         74         32         35
install.packages("naivebayes")
library(naivebayes)
## Warning: package 'naivebayes' was built under R version 3.4.3

Sekarang kita akan membuat classifier. Kita akan menggunakan e1071paket.

library(e1071)
## Warning: package 'e1071' was built under R version 3.4.3
NBclassfier=naiveBayes(prog~science+socst, data=train)
print(NBclassfier)
## 
## Naive Bayes Classifier for Discrete Predictors
## 
## Call:
## naiveBayes.default(x = X, y = Y, laplace = laplace)
## 
## A-priori probabilities:
## Y
##   academic    general vocational 
##  0.5248227  0.2269504  0.2482270 
## 
## Conditional probabilities:
##             science
## Y                [,1]     [,2]
##   academic   54.21622 9.360761
##   general    52.18750 8.847954
##   vocational 47.31429 9.969871
## 
##             socst
## Y                [,1]      [,2]
##   academic   56.58108  9.635845
##   general    51.12500  8.377196
##   vocational 44.82857 10.279865

The naiveBayesFungsi diasumsikan distribusi gaussian untuk varibles numerik. Juga, apriori dihitung dari porsi data pelatihan. Prioris ditampilkan saat objek dicetak. Nilai Y adalah mean dan standar deviasi dari prediktor dalam setiap kelas.

Sekarang mari kita membuat prediksi untuk data latih dan untuk data uji.

NBclassifier=naive_bayes(prog~science+socst,usekernel=T, data=train)
printALL=function(model){
  trainPred=predict(model, newdata = train, type = "class")
  trainTable=table(train$prog, trainPred)
  testPred=predict(NBclassfier, newdata=test, type="class")
  testTable=table(test$prog, testPred)
  trainAcc=(trainTable[1,1]+trainTable[2,2]+trainTable[3,3])/sum(trainTable)
  testAcc=(testTable[1,1]+testTable[2,2]+testTable[3,3])/sum(testTable)
  message("Contingency Table for Training Data")
  print(trainTable)
  message("Contingency Table for Test Data")
  print(testTable)
  message("Accuracy")
  print(round(cbind(trainAccuracy=trainAcc, testAccuracy=testAcc),3))
}
printALL(NBclassifier)
##             trainPred
##              academic general vocational
##   academic         61       4          9
##   general          14      12          6
##   vocational       13       4         18
##             testPred
##              academic general vocational
##   academic         27       0          4
##   general           7       0          6
##   vocational        5       0         10
##      trainAccuracy testAccuracy
## [1,]         0.645        0.627


## Contingency Table for Training Data
##             trainPred
##              academic general vocational
##   academic         60       0         14
##   general          24       0          8
##   vocational       15       0         20
## Contingency Table for Test Data
##             testPred
##              academic general vocational
##   academic         27       0          4
##   general           7       0          6
##   vocational        5       0         10
## Accuracy
##      trainAccuracy testAccuracy
## [1,]         0.567        0.627

Jadi dalam hal ini, kami memiliki akurasi pengujian yang lebih tinggi dalam data pengujian! Dan akurasi 62,7% untuk tes ini tidak terlalu buruk, mengingat kami hanya menggunakan dua variabel. Tetapi satu hal yang dilakukan oleh pengklasifikasi ini dengan sangat buruk adalah ia tidak memprediksi kelas "umum" apa pun. Kolom tengah semuanya nol dalam tabel kontingensi.

Mari kita coba lagi dengan menggunakan sesvariabel status sosial ekonomi ( ).

newNBclassifier=naive_bayes(prog~ses+science+socst,usekernel=T, data=train)
printALL=function(model){
  trainPred=predict(model, newdata = train, type = "class")
  trainTable=table(train$prog, trainPred)
  testPred=predict(newNBclassfier, newdata=test, type="class")
  testTable=table(test$prog, testPred)
  trainAcc=(trainTable[1,1]+trainTable[2,2]+trainTable[3,3])/sum(trainTable)
  testAcc=(testTable[1,1]+testTable[2,2]+testTable[3,3])/sum(testTable)
  message("Contingency Table for Training Data")
  print(trainTable)
  message("Contingency Table for Test Data")
  print(testTable)
  message("Accuracy")
  print(round(cbind(trainAccuracy=trainAcc, testAccuracy=testAcc),3))
}
printALL(newNBclassifier)
## Results after including ses
## Contingency Table for Training Data
##             trainPred
##              academic general vocational
##   academic         58       5         11
##   general          16      10          6
##   vocational       16       1         18
## Contingency Table for Test Data
##             testPred
##              academic general vocational
##   academic         27       0          4
##   general           7       0          6
##   vocational        5       0         10
## Accuracy
##      trainAccuracy testAccuracy
## [1,]          0.61        0.627

Jadi dimasukkannya sesvariabel membuat peningkatan yang baik. Baik akurasi pelatihan dan pengujian meningkat. Juga, sekarang ada beberapa pengamatan yang diklasifikasikan sebagai "general" setidaknya untuk data train.

Di atas kita memasukkan sesvariabel sebagai numerik. Tetapi memasukkan sebagai factor memberikan peningkatan yang tidak begitu berarti.

newNBclassifier=naive_bayes(prog~as.factor(ses)+science+socst,usekernel=T, data=train)
printALL=function(model){
  trainPred=predict(model, newdata = train, type = "class")
  trainTable=table(train$prog, trainPred)
  testPred=predict(newNBclassfier, newdata=test, type="class")
  testTable=table(test$prog, testPred)
  trainAcc=(trainTable[1,1]+trainTable[2,2]+trainTable[3,3])/sum(trainTable)
  testAcc=(testTable[1,1]+testTable[2,2]+testTable[3,3])/sum(testTable)
  message("Contingency Table for Training Data")
  print(trainTable)
  message("Contingency Table for Test Data")
  print(testTable)
  message("Accuracy")
  print(round(cbind(trainAccuracy=trainAcc, testAccuracy=testAcc),3))
}
printALL(newNBclassifier)
## Results after converting ses as factor
## Contingency Table for Training Data
##             trainPred
##              academic general vocational
##   academic         60       3         11
##   general          15      10          7
##   vocational       16       1         18
## Contingency Table for Test Data
##             testPred
##              academic general vocational
##   academic         27       0          4
##   general           7       0          6
##   vocational        5       0         10
## Accuracy
##      trainAccuracy testAccuracy
## [1,]         0.624        0.627

Splitting dataset dan k-fold cross validation

Tantangan utama dalam merancang model pembelajaran mesin adalah membuatnya bekerja secara akurat pada data yang tidak terlihat. Untuk menget...