Machine learning adalah sub bidang kecerdasan buatan (Artificial Intelligence). Secara umum tujuan dari machine learning adalah untuk memahami struktur data dan memasukkan data tersebut ke dalam model yang dapat dipahami dan dimanfaatkan oleh orang-orang. Meskipun machine learning adalah bidang dalam ilmu komputer, teknisnya seringkali menggunakan analisis statistik untuk menghasilkan nilai sesuai tujuan analisisnya. Oleh sebab itu, machine learning memfasilitasi komputer dalam membangun model dari data sampel untuk mengotomatisasi proses pengambilan keputusan berdasarkan input data. Secara umum, perbedaan dari machine learning dan statistics antara lain
Permodelan statistics lebih tentang menemmukan hubungan antara variabel dan signifikasi hubungan tersebut dan dapat juga untuk menghasilkan prediksi, sedangkan machine learning adalah semua tentang hasil.
Statistics menarik kesimpulan dari sampel, sedangkan machine learning menemukan pola prediktif yang dapat digeneralisasikan
Pembelajaran statistik mengacu pada seperangkat alat yang luas untuk memahami data. Alat-alat ini dapat diklasifikasikan sebagai supervised learning dan unsupervised learning. Secara umum, statistical supervised learning melibatkan pembangunan model statistik untuk memprediksi atau memperkirakan output berdasarkan satu atau lebih input. Pengaplikasiannya seringkali dilakukan pada bidang bisnis, kedokteran, astrofisika, dan kebijakan publik. Berbeda dengan supervised learning, statistical unsupervised learning ada input tetapi tidak ada supervising output. Namun, kita dapat mempelajari hubungan dan struktur dari data tersebut.
Supervised learning terdiri dari:
Classification
Regression
Unsupervised learning terdiri dari:
Clustering
Association
Dimesionality Reduction
Beberapa metode dalam regression supervised learning yang dapat digunakan antara lain:
Ridge regression. Metode ini meminimalkan jumlah kuadrat galat yang terikat pada regularisasi \(L_2\) dari koefisiennya. Ridge regression berguna untuk menduga koefisien regresi dengan peubah prediktornya saling berkorelasi (multikolinearitasnya tinggi). Selain itu, metode ini juga digunakan pada permodelan regresi dengan peubah prediktornya yang sangat banyak (\(p>>n\)). Sifat dugaan parameternya adalah berbias dan mempertahankan semua prediktornya.
Lasso regression. Metode ini melakukan penyusutan berkelanjutan (shrinkage) dan pemilihan peubah otomatis secara bersamaan. Lasso regression berguna dalam seleksi peubah prediktor dalam model. Sama halnya dengan Ridge regression, metode ini juga dapat digunakan pada permodelan regresi dengan peubah prediktornya yang sangat banyak (\(p>>n\)). Sifat dugaan parameternya adalah berbias dan memungkinkan tidak semua peubah prediktor dipilih.
Model averaging. Metode ini membangun beberapa kandidat model untuk dikombinasikan menjadi sebuah model final. Umumnya model ini diterapkan dalam rangka untuk memperoleh nilai prediksi peubah respon.
Dalam rangka meningkatkan kemudahan dalam interpretasi model regresi, metode yang digunakan adalah feature selection atau variable selection. Dimana metode ini menyisihkan peubah prediktor yang tidak relevan pada model regresi linear berganda. Beberapa metode dalam variable selection adalah
Subset selection. Pendekatan ini melibatkan pengidentifikasian subset dari \(p\) buah prediktor yang diyakini terkait dengan respon. Kemudian, dugaan modelnya diperoleh menggunakan OLS pada sebagian prediktor yang telah dipilih.
Shrinkage. Pendekatan ini melibatkan pembuatan model yang melibatkan semua \(p\) buah prediktor. Namun, koefisien ini diduga dengan penyusutan menuju nol. Penyusutan ini (juga dikenal sebagai regularisasi) memiliki efek mengurangi keragaman. Oleh karena itu, metode penyusutan juga dapat melakukan pemilihan peubah prediktor.
Dimension reduction. Pendekatan ini memproyeksikan \(p\) buah prediktor ke dalam subruang \(M\)-dimensi, di mana \(M<p\). Hal ini diperoleh dengan menghitung \(M\) buah kombinasi linier yang berbeda, atau proyeksi, dari variabel. Kemudian proyeksi \(M\) ini digunakan sebagai prediktor agar sesuai dengan model regresi linier dengan OLS.
Salah satu metode dalam subset selection adalah metode stepwise. Metode stepwise terbagi lagi menjadi tiga, yaitu
Forward stepwise selection. Forward Stepwise Selection dimulai dengan model yang tidak mengandung prediktor, dan kemudian menambahkan prediktor ke model, satu per satu, sampai semua prediktor berada dalam model. Forward Stepwise Selection dimulai dengan model yang tidak mengandung prediktor, dan kemudian menambahkan prediktor ke model, satu per satu, sampai semua prediktor berada dalam model.
Backward stepwise selection. Dimulai dengan model penuh yang berisi semua \(p\) prediktor, dan kemudian secara iteratifmenghilangkan prediktor yang paling tidak berguna, satu per satu.
Hybrid stepwise selection. Pendekatan seleksi subset terbaik, Forward Stepwise Selection, dan Backward Stepwise Selection umumnya memberikan model yang serupa tetapi tidak identik.
Analisis gerombol merupakan suatu metode peubah ganda untuk mengelompokkan n objek ke dalam m gerombol (\(m\le n\)) berdasarkan karakter-karakternya (Johnsons & Winchern, 2002). Tujuan dari penggerombolan ini adalah menemukan gerombol alamiah dari sekumpulan unit pengamatan. Ada dua metode dalam analisis gerombol satu tahapan, yaitu:
Metode Berhirarki. Metode penggerombolan berhirarki digunakan apabila banyak gerombol yang akan dibentuk belum diketahui sebelumnya. Metode ini dapat dibedakan atas dua metode yaitu metode penggabungan (agglomerative) dan metode pembagian (divisive). Dalam metode berhirarki terdapat beberapa ukuran jarak antar gerombol, antara lain metode pautan tunggal (single linkage), pautan lengkap (complete linkage), pautan rataan (average linkage), metode ward, dan metode centroid. Fungsi jarak yang sering digunakan diantaranya jarak Euclidean dan jarak Mahalanobis.
Metode Tak Berhirarki. Metode ini digunakan apabila banyak gerombol yang akan dibentuk sudah diketahui terlebih dahulu. Salah satu contohnya adalah metode K-means. Pada metode ini harus ditentukan terlebih dahulu besarnya k, yaitu banyaknya gerombol. Pemilihan k dapat ditentukan secara subjektif berdasarkan latar belakang bidang masing-masing. Jarak yang biasanya digunakan adalah jarak Euclidean.
library(broom)
library(glmnet)
library(glmnetUtils)
library(leaps)
library(varbvs)
library(ggplot2)
library(tidyverse)
library(plyr)
library(readr)
library(dplyr)
library(caret)
library(repr)
library(MASS)
library(MuMIn)
library(DT)
library(ggcorrplot)
library(tidyverse)
library(mlr3verse)
library(mlr3fselect)
library(DataExplorer)
library(skimr)
library(leaps)
library(factoextra)
library(FactoMineR)
library(MASS)
set.seed(123)
x <- cbind(1,matrix(rnorm(100*15,1,1),100,15))
e <- matrix(rnorm(100),100,1)
b <- c(1,rep(0:4,each=3))
y <- x%*%b+e
Bangkitan <- data.frame(y,x[, -1])
# Structure of Data Bangkitan
str(Bangkitan)
## 'data.frame': 100 obs. of 16 variables:
## $ y : num 22 13.3 29.5 27.8 31.5 ...
## $ X1 : num 0.44 0.77 2.56 1.07 1.13 ...
## $ X2 : num 0.2896 1.2569 0.7533 0.6525 0.0484 ...
## $ X3 : num 3.199 2.312 0.735 1.543 0.586 ...
## $ X4 : num 0.2848 0.2473 0.0615 -0.0525 0.5628 ...
## $ X5 : num 0.926 -0.169 0.365 0.971 1.671 ...
## $ X6 : num 0.3981 0.0063 2.0268 1.7511 -0.5092 ...
## $ X7 : num 2.074 0.973 0.967 -0.516 1.79 ...
## $ X8 : num 0.272 -0.54 0.307 1.119 -0.365 ...
## $ X9 : num 1.356 0.342 1.855 2.153 1.276 ...
## $ X10: num -0.0141 0.2087 1.2996 2.6391 2.0846 ...
## $ X11: num 0.0042 -0.04 0.982 0.8678 -1.5493 ...
## $ X12: num 1.916 1.8006 0.0634 -0.4008 1.1603 ...
## $ X13: num 1.62 0.242 1.852 0.252 1.63 ...
## $ X14: num 0.25 0.678 -0.148 1.354 1.425 ...
## $ X15: num -0.0861 0.3347 1.7148 0.5683 1.2276 ...
Boston Dataset berasal dari library MASS. Boston dataset
merupakan data atau informasi yang dikumpulkan oleh US Census Service
mengenai perumahan di area Boston MASS.
data("Boston",package = "MASS")
datatable(Boston)
set.seed(100)
index = sample(1:nrow(Boston), 0.7*nrow(Boston))
trainBoston = Boston[index,]
testBoston = Boston[-index,]
dim(trainBoston)
## [1] 354 14
dim(testBoston)
## [1] 152 14
Dataset yang digunakan dalam analisis ini adalah data pengunjung suatu pusat perbelanjaan (mall) yang berformat .csv. Data ini berjumlah 200 observasi (pengunjung) yang pada awalnya memiliki lima peubah, tetapi hanya digunakan tiga peubah yaitu Age, Annual Income, dan Spending Score.
Costumer <- read.csv("D:/Kuliah/Semester 5/STA1381 Pengantar Sains Data/Tugas UAS/Data Mall_Customer.csv",sep=";")
Costumer <- Costumer[,-c(1:2)]
datatable(Costumer)
corBangkitan <- cor(Bangkitan)
ggcorrplot(corBangkitan,lab=TRUE,type="lower")
Plot korelasi di atas menunjukkan hubungan antara peubah respon Y dan peubah penjelas X (X1-X15) serta antar peubah penjelas. Semakin pekat warnanya maka semakin kuat korelasinya. Begitupula sebaliknya, semakin pudar warnanya maka semakin rendah korelasinya. Mayoritas peubah memiliki hubungan yang cenderung rendah satu sama lain. Korelasi tertinggi dimiliki oleh peubah Y dan X14 dengan nilai sebesar 0.57. Terdapat beberapa hubungan yang nilainya 0 antar peubah penjelas.
plot_histogram(data = Bangkitan,nrow=3,ncol = 3,
geom_histogram_args = list(fill="blue"),title="Variables Distribution")
Histogram di atas menujukkan sebaran dari masing-masing peubah yang digunakan. Terlihat bahwa seluruh peubah cenderung memiliki distribusi normal. Hal ini disebabkan oleh pada deklarasi data bangkitan digunakan angka random yang memiliki distribusi normal.
mod1Bangkitan <- lm( y~.,data=Bangkitan)
summary(mod1Bangkitan)
##
## Call:
## lm(formula = y ~ ., data = Bangkitan)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.10684 -0.46276 -0.01661 0.46989 2.28986
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.04427 0.34971 2.986 0.0037 **
## X1 0.11730 0.09988 1.174 0.2435
## X2 -0.06743 0.09221 -0.731 0.4667
## X3 -0.11002 0.09564 -1.150 0.2533
## X4 1.19139 0.09006 13.230 <2e-16 ***
## X5 1.05484 0.09385 11.240 <2e-16 ***
## X6 1.20593 0.09652 12.494 <2e-16 ***
## X7 1.98520 0.08618 23.035 <2e-16 ***
## X8 1.96980 0.09241 21.317 <2e-16 ***
## X9 1.97923 0.08532 23.197 <2e-16 ***
## X10 2.98632 0.08850 33.745 <2e-16 ***
## X11 2.70147 0.09007 29.993 <2e-16 ***
## X12 3.10295 0.09271 33.468 <2e-16 ***
## X13 3.90424 0.09896 39.454 <2e-16 ***
## X14 4.05354 0.09192 44.097 <2e-16 ***
## X15 3.96807 0.09359 42.397 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.8465 on 84 degrees of freedom
## Multiple R-squared: 0.9952, Adjusted R-squared: 0.9943
## F-statistic: 1152 on 15 and 84 DF, p-value: < 2.2e-16
Dari summary model di atas, terlihat peubah-peubah yang berpengaruh secara signifikan dan ukuran kebaikan model. Hanya peubah X1, X2, dan X3 yang tidak berpengaruh secara signifikan terhadap peubah Y (p-value lebih dari 0.05), sedangkan peubah lainnya berpengaruh secara signifikan terhadap peubah Y. Model regresi linear di atas juga memiliki R-Squared 99.52 persen dan Adjusted R-Square 99.43 persen (mendekati 1). Hal ini mengindikasikan bahwa model di atas sudah sangat baik.
mod2Bangkitan <- cv.glmnet(x[,-1], y, alpha=0)
mod2Bangkitan
##
## Call: glmnet::cv.glmnet(x = x[, -1], y = y, alpha = 0)
##
## Measure: Mean-Squared Error
##
## Lambda Index Measure SE Nonzero
## min 0.6390 100 1.214 0.1645 15
## 1se 0.7697 98 1.329 0.1875 15
coef(mod2Bangkitan, s="lambda.min")
## 16 x 1 sparse Matrix of class "dgCMatrix"
## s1
## (Intercept) 2.33592307
## V1 0.13434607
## V2 -0.06508831
## V3 -0.15128275
## V4 1.19243040
## V5 1.02851269
## V6 1.13955222
## V7 1.89575258
## V8 1.87588686
## V9 1.87474594
## V10 2.83582175
## V11 2.67424486
## V12 3.00441230
## V13 3.71498903
## V14 3.89244467
## V15 3.72533796
mod3Bangkitan <- cv.glmnet(x[, -1], y,alpha=1)
mod3Bangkitan
##
## Call: glmnet::cv.glmnet(x = x[, -1], y = y, alpha = 1)
##
## Measure: Mean-Squared Error
##
## Lambda Index Measure SE Nonzero
## min 0.02192 62 0.8601 0.1573 15
## 1se 0.12840 43 0.9976 0.1773 13
coef(mod3Bangkitan, s="lambda.min")
## 16 x 1 sparse Matrix of class "dgCMatrix"
## s1
## (Intercept) 1.24575589
## V1 0.08935937
## V2 -0.04348371
## V3 -0.09603699
## V4 1.16854810
## V5 1.03637141
## V6 1.17953379
## V7 1.95844496
## V8 1.95923371
## V9 1.95533042
## V10 2.97005339
## V11 2.69931151
## V12 3.08595662
## V13 3.88648305
## V14 4.04911077
## V15 3.94717560
mod1Bangkitan <- lm(y~.,data=Bangkitan, na.action = na.fail)
mod4Bangkitan <- dredge(global.model=mod1Bangkitan, m.lim=c(10,20))
mod5Bangkitan <- model.avg(mod4Bangkitan, delta<4)
summary(mod5Bangkitan)
##
## Call:
## model.avg(object = mod4Bangkitan, subset = delta < 4)
##
## Component model call:
## lm(formula = y ~ <7 unique rhs>, data = Bangkitan, na.action = na.fail)
##
## Component models:
## df logLik AICc delta weight
## 2+3+4+5+6+7+10+11+12+13+14+15 14 -118.64 270.21 0.00 0.28
## 1+2+3+4+5+6+7+10+11+12+13+14+15 15 -117.64 270.99 0.77 0.19
## 2+3+4+5+6+7+9+10+11+12+13+14+15 15 -117.69 271.10 0.89 0.18
## 2+3+4+5+6+7+8+10+11+12+13+14+15 15 -118.24 272.19 1.98 0.11
## 1+2+3+4+5+6+7+9+10+11+12+13+14+15 16 -116.83 272.22 2.01 0.10
## 1+2+3+4+5+6+7+8+10+11+12+13+14+15 16 -117.30 273.15 2.93 0.07
## 2+3+4+5+6+7+8+9+10+11+12+13+14+15 16 -117.33 273.21 3.00 0.06
##
## Term codes:
## X1 X10 X11 X12 X13 X14 X15 X2 X3 X4 X5 X6 X7 X8 X9
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
##
## Model-averaged coefficients:
## (full average)
## Estimate Std. Error Adjusted SE z value Pr(>|z|)
## (Intercept) 0.97710 0.32728 0.33130 2.949 0.00319 **
## X10 2.98843 0.08835 0.08961 33.350 < 2e-16 ***
## X11 2.71154 0.08999 0.09127 29.710 < 2e-16 ***
## X12 3.08944 0.09261 0.09392 32.893 < 2e-16 ***
## X13 3.89029 0.09875 0.10015 38.844 < 2e-16 ***
## X14 4.06766 0.09157 0.09286 43.802 < 2e-16 ***
## X15 3.98206 0.09273 0.09403 42.348 < 2e-16 ***
## X4 1.17946 0.08986 0.09114 12.942 < 2e-16 ***
## X5 1.05778 0.09269 0.09399 11.254 < 2e-16 ***
## X6 1.19404 0.09567 0.09703 12.306 < 2e-16 ***
## X7 1.97873 0.08609 0.08732 22.660 < 2e-16 ***
## X8 1.99454 0.09180 0.09309 21.427 < 2e-16 ***
## X9 1.98490 0.08524 0.08646 22.958 < 2e-16 ***
## X1 0.04627 0.08573 0.08634 0.536 0.59205
## X3 -0.04153 0.07992 0.08049 0.516 0.60592
## X2 -0.01732 0.05464 0.05517 0.314 0.75352
##
## (conditional average)
## Estimate Std. Error Adjusted SE z value Pr(>|z|)
## (Intercept) 0.97710 0.32728 0.33130 2.949 0.00319 **
## X10 2.98843 0.08835 0.08961 33.350 < 2e-16 ***
## X11 2.71154 0.08999 0.09127 29.710 < 2e-16 ***
## X12 3.08944 0.09261 0.09392 32.893 < 2e-16 ***
## X13 3.89029 0.09875 0.10015 38.844 < 2e-16 ***
## X14 4.06766 0.09157 0.09286 43.802 < 2e-16 ***
## X15 3.98206 0.09273 0.09403 42.348 < 2e-16 ***
## X4 1.17946 0.08986 0.09114 12.942 < 2e-16 ***
## X5 1.05778 0.09269 0.09399 11.254 < 2e-16 ***
## X6 1.19404 0.09567 0.09703 12.306 < 2e-16 ***
## X7 1.97873 0.08609 0.08732 22.660 < 2e-16 ***
## X8 1.99454 0.09180 0.09309 21.427 < 2e-16 ***
## X9 1.98490 0.08524 0.08646 22.958 < 2e-16 ***
## X1 0.12730 0.09954 0.10097 1.261 0.20738
## X3 -0.11845 0.09543 0.09680 1.224 0.22106
## X2 -0.07376 0.09247 0.09379 0.786 0.43164
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Berdasarkan Model Averaging, diperoleh model terbaik yakni (secara urut) X4, X5, X6, X7, X8, X9, X10, X11, X12, X13, X14, X15. Peubah-peubah tersebut terbukti berpengaruh secara signifikan terhadap peubah Y karena memiliki nilai p-value kurang dari 0.05. Selain itu dapat terlihat pula nilai koefisien tiap-tiap peubah.
stepAIC(mod1Bangkitan, direction = "forward")
## Start: AIC=-18.76
## y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 +
## X12 + X13 + X14 + X15
##
## Call:
## lm(formula = y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + X9 +
## X10 + X11 + X12 + X13 + X14 + X15, data = Bangkitan, na.action = na.fail)
##
## Coefficients:
## (Intercept) X1 X2 X3 X4 X5
## 1.04427 0.11730 -0.06743 -0.11002 1.19139 1.05484
## X6 X7 X8 X9 X10 X11
## 1.20593 1.98520 1.96980 1.97923 2.98632 2.70147
## X12 X13 X14 X15
## 3.10295 3.90424 4.05354 3.96807
Berdasarkan metode Forward Stepwise, didapatkan model terbaik yakni model dengan peubah penjelas lengkap dari X1-X15. Model ini memiliki nilai AIC sebear -18.76. Model ini dapat dituliskan sebagai berikut
\[\begin{aligned} &\hat{y}=1.04427+0.11730*X1-0.06743*X2-0.11002*X3+1.19139*X4+1.05484*X5+1.20593*X6+1.98520*7+1.96980*X8+1.97923*X9+2.98632*X10+2.70147*X11+3.10295*X12+3.90424*X13+4.05354*X14+3.96807*X15\\ \end{aligned}\]
stepAIC(mod1Bangkitan, direction = "backward")
## Start: AIC=-18.76
## y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 +
## X12 + X13 + X14 + X15
##
## Df Sum of Sq RSS AIC
## - X2 1 0.38 60.58 -20.124
## - X3 1 0.95 61.14 -19.196
## - X1 1 0.99 61.18 -19.130
## <none> 60.19 -18.759
## - X5 1 90.53 150.73 71.029
## - X6 1 111.87 172.06 84.270
## - X4 1 125.42 185.62 91.851
## - X8 1 325.63 385.83 165.022
## - X7 1 380.25 440.44 178.261
## - X9 1 385.60 445.80 179.469
## - X11 1 644.66 704.85 225.282
## - X12 1 802.68 862.88 245.510
## - X10 1 816.01 876.21 247.043
## - X13 1 1115.45 1175.65 276.440
## - X15 1 1288.09 1348.29 290.142
## - X14 1 1393.50 1453.69 297.669
##
## Step: AIC=-20.12
## y ~ X1 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 +
## X13 + X14 + X15
##
## Df Sum of Sq RSS AIC
## - X3 1 0.98 61.56 -20.515
## - X1 1 1.05 61.63 -20.401
## <none> 60.58 -20.124
## - X5 1 95.32 155.90 72.403
## - X6 1 113.07 173.65 83.187
## - X4 1 125.09 185.67 89.879
## - X8 1 328.09 388.67 163.756
## - X7 1 380.87 441.45 176.488
## - X9 1 388.08 448.65 178.108
## - X11 1 644.94 705.51 223.376
## - X12 1 802.51 863.09 243.535
## - X10 1 827.30 887.88 246.366
## - X13 1 1115.10 1175.68 274.443
## - X15 1 1300.32 1360.89 289.073
## - X14 1 1403.94 1464.52 296.411
##
## Step: AIC=-20.52
## y ~ X1 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 + X13 +
## X14 + X15
##
## Df Sum of Sq RSS AIC
## - X1 1 1.24 62.80 -20.517
## <none> 61.56 -20.515
## - X5 1 97.86 159.42 72.639
## - X6 1 112.92 174.49 81.667
## - X4 1 125.06 186.62 88.393
## - X8 1 338.82 400.38 164.725
## - X7 1 381.07 442.63 174.757
## - X9 1 394.91 456.47 177.835
## - X11 1 660.86 722.42 223.743
## - X12 1 802.07 863.63 241.598
## - X10 1 826.38 887.94 244.374
## - X13 1 1129.42 1190.98 273.736
## - X15 1 1363.75 1425.31 291.697
## - X14 1 1403.00 1464.56 294.414
##
## Step: AIC=-20.52
## y ~ X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 + X13 + X14 +
## X15
##
## Df Sum of Sq RSS AIC
## <none> 62.80 -20.517
## - X5 1 96.84 159.65 70.779
## - X6 1 112.35 175.15 80.046
## - X4 1 123.84 186.65 86.404
## - X8 1 356.89 419.70 167.436
## - X7 1 379.92 442.72 172.777
## - X9 1 393.81 456.61 175.866
## - X11 1 662.56 725.36 222.150
## - X12 1 806.40 869.21 240.241
## - X10 1 825.69 888.50 242.436
## - X13 1 1128.22 1191.02 271.740
## - X15 1 1383.66 1446.46 291.171
## - X14 1 1455.06 1517.86 295.989
##
## Call:
## lm(formula = y ~ X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 +
## X13 + X14 + X15, data = Bangkitan, na.action = na.fail)
##
## Coefficients:
## (Intercept) X4 X5 X6 X7 X8
## 0.9519 1.1729 1.0569 1.1892 1.9757 2.0091
## X9 X10 X11 X12 X13 X14
## 1.9876 2.9885 2.7176 3.0811 3.8822 4.0757
## X15
## 3.9922
Berdasarkan metode Backward Stepwise, diperoleh model yang berbeda dengan Forward Stepwise. Model ini memiliki nilai AIC yang lebih rendah yakni -20.52. Model tersebut dapat dituliskan sebagai berikut
\[\begin{aligned} &\hat{y}=0.9519+1.1729*X4+1.0569*X5+1.1892*X6+1.9757*7+2.0091*X8+1.9876*X9+2.9885*X10+2.7176*X11+3.0811*X12+ 3.8822*X13+4.075*X14+3.9922*X15\\ \end{aligned}\]
stepAIC(mod1Bangkitan, direction = "both")
## Start: AIC=-18.76
## y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 +
## X12 + X13 + X14 + X15
##
## Df Sum of Sq RSS AIC
## - X2 1 0.38 60.58 -20.124
## - X3 1 0.95 61.14 -19.196
## - X1 1 0.99 61.18 -19.130
## <none> 60.19 -18.759
## - X5 1 90.53 150.73 71.029
## - X6 1 111.87 172.06 84.270
## - X4 1 125.42 185.62 91.851
## - X8 1 325.63 385.83 165.022
## - X7 1 380.25 440.44 178.261
## - X9 1 385.60 445.80 179.469
## - X11 1 644.66 704.85 225.282
## - X12 1 802.68 862.88 245.510
## - X10 1 816.01 876.21 247.043
## - X13 1 1115.45 1175.65 276.440
## - X15 1 1288.09 1348.29 290.142
## - X14 1 1393.50 1453.69 297.669
##
## Step: AIC=-20.12
## y ~ X1 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 +
## X13 + X14 + X15
##
## Df Sum of Sq RSS AIC
## - X3 1 0.98 61.56 -20.515
## - X1 1 1.05 61.63 -20.401
## <none> 60.58 -20.124
## + X2 1 0.38 60.19 -18.759
## - X5 1 95.32 155.90 72.403
## - X6 1 113.07 173.65 83.187
## - X4 1 125.09 185.67 89.879
## - X8 1 328.09 388.67 163.756
## - X7 1 380.87 441.45 176.488
## - X9 1 388.08 448.65 178.108
## - X11 1 644.94 705.51 223.376
## - X12 1 802.51 863.09 243.535
## - X10 1 827.30 887.88 246.366
## - X13 1 1115.10 1175.68 274.443
## - X15 1 1300.32 1360.89 289.073
## - X14 1 1403.94 1464.52 296.411
##
## Step: AIC=-20.52
## y ~ X1 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 + X13 +
## X14 + X15
##
## Df Sum of Sq RSS AIC
## - X1 1 1.24 62.80 -20.517
## <none> 61.56 -20.515
## + X3 1 0.98 60.58 -20.124
## + X2 1 0.42 61.14 -19.196
## - X5 1 97.86 159.42 72.639
## - X6 1 112.92 174.49 81.667
## - X4 1 125.06 186.62 88.393
## - X8 1 338.82 400.38 164.725
## - X7 1 381.07 442.63 174.757
## - X9 1 394.91 456.47 177.835
## - X11 1 660.86 722.42 223.743
## - X12 1 802.07 863.63 241.598
## - X10 1 826.38 887.94 244.374
## - X13 1 1129.42 1190.98 273.736
## - X15 1 1363.75 1425.31 291.697
## - X14 1 1403.00 1464.56 294.414
##
## Step: AIC=-20.52
## y ~ X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 + X13 + X14 +
## X15
##
## Df Sum of Sq RSS AIC
## <none> 62.80 -20.517
## + X1 1 1.24 61.56 -20.515
## + X3 1 1.17 61.63 -20.401
## + X2 1 0.49 62.31 -19.307
## - X5 1 96.84 159.65 70.779
## - X6 1 112.35 175.15 80.046
## - X4 1 123.84 186.65 86.404
## - X8 1 356.89 419.70 167.436
## - X7 1 379.92 442.72 172.777
## - X9 1 393.81 456.61 175.866
## - X11 1 662.56 725.36 222.150
## - X12 1 806.40 869.21 240.241
## - X10 1 825.69 888.50 242.436
## - X13 1 1128.22 1191.02 271.740
## - X15 1 1383.66 1446.46 291.171
## - X14 1 1455.06 1517.86 295.989
##
## Call:
## lm(formula = y ~ X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 +
## X13 + X14 + X15, data = Bangkitan, na.action = na.fail)
##
## Coefficients:
## (Intercept) X4 X5 X6 X7 X8
## 0.9519 1.1729 1.0569 1.1892 1.9757 2.0091
## X9 X10 X11 X12 X13 X14
## 1.9876 2.9885 2.7176 3.0811 3.8822 4.0757
## X15
## 3.9922
Model yang didapatkan dari metode Forward and Backward Stepwise sama dengan yang diperoleh dari metode Backward Stepwise. Model ini memiliki nilai AIC -20.52. Model tersebut dapat dituliskan sebagai berikut
\[\begin{aligned} &\hat{y}=0.9519+1.1729*X4+1.0569*X5+1.1892*X6+1.9757*7+2.0091*X8+1.9876*X9+2.9885*X10+2.7176*X11+3.0811*X12+ 3.8822*X13+4.075*X14+3.9922*X15\\ \end{aligned}\]
Model terbaik yang didapatkan adalah model hasil dari metode Backward Stepwise dan Forward and Backward Stepwise. Model tersebut dapat dituliskan sebagai berikut
\[\begin{aligned} &\hat{y}=0.9519+1.1729*X4+1.0569*X5+1.1892*X6+1.9757*7+2.0091*X8+1.9876*X9+2.9885*X10+2.7176*X11+3.0811*X12+ 3.8822*X13+4.075*X14+3.9922*X15\\ \end{aligned}\]
plot_intro(data = Boston,
geom_label_args = list(size=2.5))
datatable(skim_without_charts(data = Boston))
Terlihat dari plot di atas, tidak terdapat missing values (NA) pada baris, kolom, maupun observasi. Selain itu, terlihat pula dari tabel di atas nilai mean dan standar deviasi dari tiap-tiap peubah.
corBoston <- cor(Boston)
ggcorrplot(corBoston,lab = TRUE,type="lower")
Plot korelasi di atas menunjukkan hubungan antara peubah respon (medv) dan peubah penjelas serta antar peubah penjelas. Semakin pekat warnanya maka semakin kuat korelasinya. Begitupula sebaliknya, semakin pudar warnanya maka semakin rendah korelasinya. Terlihat bahwa mayoritas hubungan antar peubah memiliki warna biru keunguan. Hal ini mengindikasikan bahwa mayoritas peubah memiliki hubungan yang lemah. Korelasi tertinggi dimiliki oleh peubah rad dan tax dengan nilai 0.91 (mendekati 1).
plot_histogram(data = Boston,nrow=3,ncol = 3,
geom_histogram_args = list(fill="blue"),title="Variables Distribution")
Plot histogram digunakan untuk melihat sebaran dari tiap-tiap peubah. Terlihat bahwa tiap-tiap peubah memilli sebaran yang beragam. Terdapat beberapa peubah yang memiliki sebaran normal, menjulur ke kanan, dan menjulur ke kiri.
linearBoston <- lm(medv~.,data=Boston)
summary(linearBoston)
##
## Call:
## lm(formula = medv ~ ., data = Boston)
##
## Residuals:
## Min 1Q Median 3Q Max
## -15.595 -2.730 -0.518 1.777 26.199
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.646e+01 5.103e+00 7.144 3.28e-12 ***
## crim -1.080e-01 3.286e-02 -3.287 0.001087 **
## zn 4.642e-02 1.373e-02 3.382 0.000778 ***
## indus 2.056e-02 6.150e-02 0.334 0.738288
## chas 2.687e+00 8.616e-01 3.118 0.001925 **
## nox -1.777e+01 3.820e+00 -4.651 4.25e-06 ***
## rm 3.810e+00 4.179e-01 9.116 < 2e-16 ***
## age 6.922e-04 1.321e-02 0.052 0.958229
## dis -1.476e+00 1.995e-01 -7.398 6.01e-13 ***
## rad 3.060e-01 6.635e-02 4.613 5.07e-06 ***
## tax -1.233e-02 3.760e-03 -3.280 0.001112 **
## ptratio -9.527e-01 1.308e-01 -7.283 1.31e-12 ***
## black 9.312e-03 2.686e-03 3.467 0.000573 ***
## lstat -5.248e-01 5.072e-02 -10.347 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.745 on 492 degrees of freedom
## Multiple R-squared: 0.7406, Adjusted R-squared: 0.7338
## F-statistic: 108.1 on 13 and 492 DF, p-value: < 2.2e-16
Berdasarkan summary di atas, terlihat bahwa pada model regresi penuh
peubah indus dan age tidak berpengaruh secara
signifikan terhadap peubah respon medv. Model ini memiliki
ukuran kebaikan model R-Squared sebesar 74.06 persen dan Adjusted
R-Squared 73.38 persen.
# Membuat Costum Control Parameter
customBoston <- trainControl(method = "repeatedcv",
number = 10,
repeats = 5,
verboseIter = T)
# Fitting Ridge Regression Model
set.seed(1234)
ridgeBoston <- train(medv~.,trainBoston,
method="glmnet",
tuneGrid=expand.grid(alpha=0,lambda=seq(0.0001,1,length=5)),
trControl=customBoston)
## + Fold01.Rep1: alpha=0, lambda=1
## - Fold01.Rep1: alpha=0, lambda=1
## + Fold02.Rep1: alpha=0, lambda=1
## - Fold02.Rep1: alpha=0, lambda=1
## + Fold03.Rep1: alpha=0, lambda=1
## - Fold03.Rep1: alpha=0, lambda=1
## + Fold04.Rep1: alpha=0, lambda=1
## - Fold04.Rep1: alpha=0, lambda=1
## + Fold05.Rep1: alpha=0, lambda=1
## - Fold05.Rep1: alpha=0, lambda=1
## + Fold06.Rep1: alpha=0, lambda=1
## - Fold06.Rep1: alpha=0, lambda=1
## + Fold07.Rep1: alpha=0, lambda=1
## - Fold07.Rep1: alpha=0, lambda=1
## + Fold08.Rep1: alpha=0, lambda=1
## - Fold08.Rep1: alpha=0, lambda=1
## + Fold09.Rep1: alpha=0, lambda=1
## - Fold09.Rep1: alpha=0, lambda=1
## + Fold10.Rep1: alpha=0, lambda=1
## - Fold10.Rep1: alpha=0, lambda=1
## + Fold01.Rep2: alpha=0, lambda=1
## - Fold01.Rep2: alpha=0, lambda=1
## + Fold02.Rep2: alpha=0, lambda=1
## - Fold02.Rep2: alpha=0, lambda=1
## + Fold03.Rep2: alpha=0, lambda=1
## - Fold03.Rep2: alpha=0, lambda=1
## + Fold04.Rep2: alpha=0, lambda=1
## - Fold04.Rep2: alpha=0, lambda=1
## + Fold05.Rep2: alpha=0, lambda=1
## - Fold05.Rep2: alpha=0, lambda=1
## + Fold06.Rep2: alpha=0, lambda=1
## - Fold06.Rep2: alpha=0, lambda=1
## + Fold07.Rep2: alpha=0, lambda=1
## - Fold07.Rep2: alpha=0, lambda=1
## + Fold08.Rep2: alpha=0, lambda=1
## - Fold08.Rep2: alpha=0, lambda=1
## + Fold09.Rep2: alpha=0, lambda=1
## - Fold09.Rep2: alpha=0, lambda=1
## + Fold10.Rep2: alpha=0, lambda=1
## - Fold10.Rep2: alpha=0, lambda=1
## + Fold01.Rep3: alpha=0, lambda=1
## - Fold01.Rep3: alpha=0, lambda=1
## + Fold02.Rep3: alpha=0, lambda=1
## - Fold02.Rep3: alpha=0, lambda=1
## + Fold03.Rep3: alpha=0, lambda=1
## - Fold03.Rep3: alpha=0, lambda=1
## + Fold04.Rep3: alpha=0, lambda=1
## - Fold04.Rep3: alpha=0, lambda=1
## + Fold05.Rep3: alpha=0, lambda=1
## - Fold05.Rep3: alpha=0, lambda=1
## + Fold06.Rep3: alpha=0, lambda=1
## - Fold06.Rep3: alpha=0, lambda=1
## + Fold07.Rep3: alpha=0, lambda=1
## - Fold07.Rep3: alpha=0, lambda=1
## + Fold08.Rep3: alpha=0, lambda=1
## - Fold08.Rep3: alpha=0, lambda=1
## + Fold09.Rep3: alpha=0, lambda=1
## - Fold09.Rep3: alpha=0, lambda=1
## + Fold10.Rep3: alpha=0, lambda=1
## - Fold10.Rep3: alpha=0, lambda=1
## + Fold01.Rep4: alpha=0, lambda=1
## - Fold01.Rep4: alpha=0, lambda=1
## + Fold02.Rep4: alpha=0, lambda=1
## - Fold02.Rep4: alpha=0, lambda=1
## + Fold03.Rep4: alpha=0, lambda=1
## - Fold03.Rep4: alpha=0, lambda=1
## + Fold04.Rep4: alpha=0, lambda=1
## - Fold04.Rep4: alpha=0, lambda=1
## + Fold05.Rep4: alpha=0, lambda=1
## - Fold05.Rep4: alpha=0, lambda=1
## + Fold06.Rep4: alpha=0, lambda=1
## - Fold06.Rep4: alpha=0, lambda=1
## + Fold07.Rep4: alpha=0, lambda=1
## - Fold07.Rep4: alpha=0, lambda=1
## + Fold08.Rep4: alpha=0, lambda=1
## - Fold08.Rep4: alpha=0, lambda=1
## + Fold09.Rep4: alpha=0, lambda=1
## - Fold09.Rep4: alpha=0, lambda=1
## + Fold10.Rep4: alpha=0, lambda=1
## - Fold10.Rep4: alpha=0, lambda=1
## + Fold01.Rep5: alpha=0, lambda=1
## - Fold01.Rep5: alpha=0, lambda=1
## + Fold02.Rep5: alpha=0, lambda=1
## - Fold02.Rep5: alpha=0, lambda=1
## + Fold03.Rep5: alpha=0, lambda=1
## - Fold03.Rep5: alpha=0, lambda=1
## + Fold04.Rep5: alpha=0, lambda=1
## - Fold04.Rep5: alpha=0, lambda=1
## + Fold05.Rep5: alpha=0, lambda=1
## - Fold05.Rep5: alpha=0, lambda=1
## + Fold06.Rep5: alpha=0, lambda=1
## - Fold06.Rep5: alpha=0, lambda=1
## + Fold07.Rep5: alpha=0, lambda=1
## - Fold07.Rep5: alpha=0, lambda=1
## + Fold08.Rep5: alpha=0, lambda=1
## - Fold08.Rep5: alpha=0, lambda=1
## + Fold09.Rep5: alpha=0, lambda=1
## - Fold09.Rep5: alpha=0, lambda=1
## + Fold10.Rep5: alpha=0, lambda=1
## - Fold10.Rep5: alpha=0, lambda=1
## Aggregating results
## Selecting tuning parameters
## Fitting alpha = 0, lambda = 0.5 on full training set
ridgeBoston
## glmnet
##
## 354 samples
## 13 predictor
##
## No pre-processing
## Resampling: Cross-Validated (10 fold, repeated 5 times)
## Summary of sample sizes: 318, 319, 318, 319, 319, 318, ...
## Resampling results across tuning parameters:
##
## lambda RMSE Rsquared MAE
## 0.000100 4.340939 0.7712224 3.102308
## 0.250075 4.340939 0.7712224 3.102308
## 0.500050 4.340939 0.7712224 3.102308
## 0.750025 4.346518 0.7708322 3.100617
## 1.000000 4.365780 0.7693706 3.097781
##
## Tuning parameter 'alpha' was held constant at a value of 0
## RMSE was used to select the optimal model using the smallest value.
## The final values used for the model were alpha = 0 and lambda = 0.50005.
# Mean Validation Score
mean(ridgeBoston$resample$RMSE)
## [1] 4.340939
# Model Plot
plot(ridgeBoston, main = "Ridge Regression")
# Variables Importance Plot
plot(varImp(ridgeBoston,scale=T), main="Variables Importance")
Berdasarkan plot Variables Importance di atas, terlihat bahwa peubah
nox merupakan peubah paling penting diantara peubah-peubah
yang lain terhadap model.
lambdasBoston <- 10^seq(2, -3, by = -.1)
ridge_reg = glmnet(as.matrix(trainBoston[,-14]), trainBoston[,14], nlambda = 25, alpha = 0,
family = 'gaussian', lambda = lambdasBoston)
summary(ridge_reg)
## Length Class Mode
## a0 51 -none- numeric
## beta 663 dgCMatrix S4
## df 51 -none- numeric
## dim 2 -none- numeric
## lambda 51 -none- numeric
## dev.ratio 51 -none- numeric
## nulldev 1 -none- numeric
## npasses 1 -none- numeric
## jerr 1 -none- numeric
## offset 1 -none- logical
## call 7 -none- call
## nobs 1 -none- numeric
ridge_reg
##
## Call: glmnet::glmnet(x = as.matrix(trainBoston[, -14]), y = trainBoston[, 14], nlambda = 25, alpha = 0, family = "gaussian", lambda = lambdasBoston)
##
## Df %Dev Lambda
## 1 13 30.01 100.000
## 2 13 34.00 79.430
## 3 13 38.04 63.100
## 4 13 42.05 50.120
## 5 13 45.95 39.810
## 6 13 49.70 31.620
## 7 13 53.25 25.120
## 8 13 56.59 19.950
## 9 13 59.69 15.850
## 10 13 62.52 12.590
## 11 13 65.09 10.000
## 12 13 67.36 7.943
## 13 13 69.34 6.310
## 14 13 71.03 5.012
## 15 13 72.44 3.981
## 16 13 73.62 3.162
## 17 13 74.57 2.512
## 18 13 75.34 1.995
## 19 13 75.94 1.585
## 20 13 76.42 1.259
## 21 13 76.79 1.000
## 22 13 77.08 0.794
## 23 13 77.30 0.631
## 24 13 77.46 0.501
## 25 13 77.58 0.398
## 26 13 77.67 0.316
## 27 13 77.73 0.251
## 28 13 77.77 0.200
## 29 13 77.81 0.158
## 30 13 77.83 0.126
## 31 13 77.84 0.100
## 32 13 77.85 0.079
## 33 13 77.86 0.063
## 34 13 77.86 0.050
## 35 13 77.87 0.040
## 36 13 77.87 0.032
## 37 13 77.87 0.025
## 38 13 77.87 0.020
## 39 13 77.87 0.016
## 40 13 77.87 0.013
## 41 13 77.87 0.010
## 42 13 77.87 0.008
## 43 13 77.87 0.006
## 44 13 77.87 0.005
## 45 13 77.87 0.004
## 46 13 77.87 0.003
## 47 13 77.87 0.003
## 48 13 77.87 0.002
## 49 13 77.87 0.002
## 50 13 77.87 0.001
## 51 13 77.87 0.001
# Tuning Parameter
set.seed(100)
cv_ridgeBoston <- cv.glmnet(as.matrix(trainBoston[,-14]), trainBoston[,14], alpha = 0, lambda = lambdasBoston)
optimal_lambdaBoston <- cv_ridgeBoston$lambda.min
optimal_lambdaBoston
## [1] 0.1258925
# Coefficients
coef(cv_ridgeBoston, s="lambda.min")
## 14 x 1 sparse Matrix of class "dgCMatrix"
## s1
## (Intercept) 33.512762193
## crim -0.116541536
## zn 0.045303128
## indus 0.015047943
## chas 2.864206507
## nox -10.422379184
## rm 3.659253006
## age -0.012198505
## dis -1.384865675
## rad 0.227556069
## tax -0.008280337
## ptratio -0.902912189
## black 0.008342041
## lstat -0.623529573
# Compute R^2 from true and predicted values
eval_results <- function(true, predicted, df) {
SSE <- sum((predicted - true)^2)
SST <- sum((true - mean(true))^2)
R_square <- 1 - SSE / SST
RMSE = sqrt(SSE/nrow(df))
# Model performance metrics
data.frame(RMSE = RMSE,
Rsquare = R_square)
}
# Prediction and evaluation on train data
predictions_trainBoston <- predict(ridge_reg, s = optimal_lambdaBoston, newx = as.matrix(trainBoston[,-14]))
eval_results(trainBoston[,14], predictions_trainBoston, trainBoston)
## RMSE Rsquare
## 1 4.186191 0.7782747
# Prediction and evaluation on test data
predictions_testBoston <- predict(ridge_reg, s = optimal_lambdaBoston, newx = as.matrix(testBoston[,-14]))
eval_results(testBoston[,14], predictions_testBoston, testBoston)
## RMSE Rsquare
## 1 5.830073 0.6410751
Terlihat dari tabel prediksi dan evaluasi di atas, hasil yang ditunjukkan sudah cukup baik. Nilai RMSE pada data training dan data testing sudah cukup rendah yakni sekitar 4.19 dan 5.83. Sedangkan untuk nilai R-Square berturut-turut sebesar 77.83 persen dan 64.11 persen.
# Setting alpha = 1 implements lasso regression
lasso_regBoston <- cv.glmnet(as.matrix(trainBoston[,-14]), trainBoston[,14], alpha = 1,
lambda = lambdasBoston, standardize = TRUE, nfolds = 5)
# Best
lambda_bestBoston <- lasso_regBoston$lambda.min
lambda_bestBoston
## [1] 0.01995262
lasso_modelBoston <- glmnet(as.matrix(trainBoston[,-14]), trainBoston[,14], alpha = 1,
lambda = lambda_bestBoston, standardize = TRUE)
# Prediction and evaluation on train data
predictions_trainBoston <- predict(lasso_modelBoston, s = lambda_bestBoston, newx = as.matrix(trainBoston[,-14]))
eval_results(trainBoston[,14], predictions_trainBoston, trainBoston)
## RMSE Rsquare
## 1 4.185119 0.7783883
# Prediction and evaluation on test data
predictions_testBoston <- predict(lasso_modelBoston, s = lambda_bestBoston, newx = as.matrix(testBoston[,-14]))
eval_results(testBoston[,14], predictions_testBoston, testBoston)
## RMSE Rsquare
## 1 5.844554 0.6392898
Terlihat pada tabel prediksi dan evaluasi dari model regresi Lasso, Hasil yang ditunjukkan sudah cukup baik dan tidak jauh berbeda dengan hasil pada model regresi Ridge.
# Model Declare
mod1Boston <- lm(medv~.,trainBoston, na.action = na.fail)
mod2Boston <- dredge(global.model=mod1Boston, m.lim=c(10,14))
mod3Boston <- model.avg(mod2Boston, delta<4)
summary(mod3Boston)
##
## Call:
## model.avg(object = mod2Boston, subset = delta < 4)
##
## Component model call:
## lm(formula = medv ~ <4 unique rhs>, data = trainBoston, na.action =
## na.fail)
##
## Component models:
## df logLik AICc delta weight
## 2+3+4+5+7+8+9+10+11+12+13 13 -1009.26 2045.59 0.00 0.49
## 1+2+3+4+5+7+8+9+10+11+12+13 14 -1008.94 2047.12 1.53 0.23
## 2+3+4+5+6+7+8+9+10+11+12+13 14 -1009.12 2047.48 1.89 0.19
## 1+2+3+4+5+6+7+8+9+10+11+12+13 15 -1008.80 2049.02 3.43 0.09
##
## Term codes:
## age black chas crim dis indus lstat nox ptratio rad
## 1 2 3 4 5 6 7 8 9 10
## rm tax zn
## 11 12 13
##
## Model-averaged coefficients:
## (full average)
## Estimate Std. Error Adjusted SE z value Pr(>|z|)
## (Intercept) 35.650212 5.767077 5.787586 6.160 < 2e-16 ***
## black 0.008197 0.002781 0.002791 2.937 0.003318 **
## chas 2.769576 0.974071 0.977528 2.833 0.004608 **
## crim -0.122049 0.030904 0.031014 3.935 8.31e-05 ***
## dis -1.442410 0.210692 0.211425 6.822 < 2e-16 ***
## lstat -0.647875 0.055948 0.056144 11.540 < 2e-16 ***
## nox -11.525187 3.940808 3.954705 2.914 0.003565 **
## ptratio -0.921584 0.143636 0.144146 6.393 < 2e-16 ***
## rad 0.265188 0.067213 0.067451 3.932 8.44e-05 ***
## rm 3.517981 0.473696 0.475375 7.400 < 2e-16 ***
## tax -0.009769 0.003694 0.003707 2.636 0.008401 **
## zn 0.049457 0.014677 0.014729 3.358 0.000785 ***
## age -0.003714 0.010012 0.010037 0.370 0.711326
## indus 0.009484 0.037995 0.038109 0.249 0.803472
##
## (conditional average)
## Estimate Std. Error Adjusted SE z value Pr(>|z|)
## (Intercept) 35.650212 5.767077 5.787586 6.160 < 2e-16 ***
## black 0.008197 0.002781 0.002791 2.937 0.003318 **
## chas 2.769576 0.974071 0.977528 2.833 0.004608 **
## crim -0.122049 0.030904 0.031014 3.935 8.31e-05 ***
## dis -1.442410 0.210692 0.211425 6.822 < 2e-16 ***
## lstat -0.647875 0.055948 0.056144 11.540 < 2e-16 ***
## nox -11.525187 3.940808 3.954705 2.914 0.003565 **
## ptratio -0.921584 0.143636 0.144146 6.393 < 2e-16 ***
## rad 0.265188 0.067213 0.067451 3.932 8.44e-05 ***
## rm 3.517981 0.473696 0.475375 7.400 < 2e-16 ***
## tax -0.009769 0.003694 0.003707 2.636 0.008401 **
## zn 0.049457 0.014677 0.014729 3.358 0.000785 ***
## age -0.011705 0.014911 0.014964 0.782 0.434076
## indus 0.033960 0.065866 0.066101 0.514 0.607419
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
stepAIC(linearBoston, direction = "forward")
## Start: AIC=1589.64
## medv ~ crim + zn + indus + chas + nox + rm + age + dis + rad +
## tax + ptratio + black + lstat
##
## Call:
## lm(formula = medv ~ crim + zn + indus + chas + nox + rm + age +
## dis + rad + tax + ptratio + black + lstat, data = Boston)
##
## Coefficients:
## (Intercept) crim zn indus chas nox
## 3.646e+01 -1.080e-01 4.642e-02 2.056e-02 2.687e+00 -1.777e+01
## rm age dis rad tax ptratio
## 3.810e+00 6.922e-04 -1.476e+00 3.060e-01 -1.233e-02 -9.527e-01
## black lstat
## 9.312e-03 -5.248e-01
Berdasarkan metode Forward Stepwise, terlihat bahwa model terbaik adalah model penuh. Model ini memiliki nilai AIC 1589.64 dan dapat dituliskan sebagai berikut:
\[\begin{aligned} &\hat{medv}=(3.646e+01)-(1.080e-01)*crim+(4.642e-02)*zn+(2.056e-02)*indus+(2.687e+00)*chas-(1.777e+01)*nox+(3.810e+00)*rm+(6.922e-04)*age-(1.476e+00)*dis+(3.060e-01)*rad -(1.233e-02)*tax-(9.527e-01)*ptratio+(9.312e-03)*black-(5.248e-01)*lstat\\ \end{aligned}\]
stepAIC(linearBoston, direction = "backward")
## Start: AIC=1589.64
## medv ~ crim + zn + indus + chas + nox + rm + age + dis + rad +
## tax + ptratio + black + lstat
##
## Df Sum of Sq RSS AIC
## - age 1 0.06 11079 1587.7
## - indus 1 2.52 11081 1587.8
## <none> 11079 1589.6
## - chas 1 218.97 11298 1597.5
## - tax 1 242.26 11321 1598.6
## - crim 1 243.22 11322 1598.6
## - zn 1 257.49 11336 1599.3
## - black 1 270.63 11349 1599.8
## - rad 1 479.15 11558 1609.1
## - nox 1 487.16 11566 1609.4
## - ptratio 1 1194.23 12273 1639.4
## - dis 1 1232.41 12311 1641.0
## - rm 1 1871.32 12950 1666.6
## - lstat 1 2410.84 13490 1687.3
##
## Step: AIC=1587.65
## medv ~ crim + zn + indus + chas + nox + rm + dis + rad + tax +
## ptratio + black + lstat
##
## Df Sum of Sq RSS AIC
## - indus 1 2.52 11081 1585.8
## <none> 11079 1587.7
## - chas 1 219.91 11299 1595.6
## - tax 1 242.24 11321 1596.6
## - crim 1 243.20 11322 1596.6
## - zn 1 260.32 11339 1597.4
## - black 1 272.26 11351 1597.9
## - rad 1 481.09 11560 1607.2
## - nox 1 520.87 11600 1608.9
## - ptratio 1 1200.23 12279 1637.7
## - dis 1 1352.26 12431 1643.9
## - rm 1 1959.55 13038 1668.0
## - lstat 1 2718.88 13798 1696.7
##
## Step: AIC=1585.76
## medv ~ crim + zn + chas + nox + rm + dis + rad + tax + ptratio +
## black + lstat
##
## Df Sum of Sq RSS AIC
## <none> 11081 1585.8
## - chas 1 227.21 11309 1594.0
## - crim 1 245.37 11327 1594.8
## - zn 1 257.82 11339 1595.4
## - black 1 270.82 11352 1596.0
## - tax 1 273.62 11355 1596.1
## - rad 1 500.92 11582 1606.1
## - nox 1 541.91 11623 1607.9
## - ptratio 1 1206.45 12288 1636.0
## - dis 1 1448.94 12530 1645.9
## - rm 1 1963.66 13045 1666.3
## - lstat 1 2723.48 13805 1695.0
##
## Call:
## lm(formula = medv ~ crim + zn + chas + nox + rm + dis + rad +
## tax + ptratio + black + lstat, data = Boston)
##
## Coefficients:
## (Intercept) crim zn chas nox rm
## 36.341145 -0.108413 0.045845 2.718716 -17.376023 3.801579
## dis rad tax ptratio black lstat
## -1.492711 0.299608 -0.011778 -0.946525 0.009291 -0.522553
Pada metode Backward Stepwise, terlihat hasil yang berbeda daripada metode Forward Stepwise. Model ini memiliki nilai AIC yang lebih rendah yakni 1585.76 dan dapat dituliskan sebagai berikut
\[\begin{aligned} &\hat{medv}=36.341145-0.108413*crim+0.045845*zn+2.718716*indus-17.376023*nox+3.801579*rm-1.492711*dis+0.299608*rad-0.011778*tax-0.946525*ptratio+0.009291*black-0.522553*lstat\\ \end{aligned}\]
stepAIC(linearBoston, direction = "both")
## Start: AIC=1589.64
## medv ~ crim + zn + indus + chas + nox + rm + age + dis + rad +
## tax + ptratio + black + lstat
##
## Df Sum of Sq RSS AIC
## - age 1 0.06 11079 1587.7
## - indus 1 2.52 11081 1587.8
## <none> 11079 1589.6
## - chas 1 218.97 11298 1597.5
## - tax 1 242.26 11321 1598.6
## - crim 1 243.22 11322 1598.6
## - zn 1 257.49 11336 1599.3
## - black 1 270.63 11349 1599.8
## - rad 1 479.15 11558 1609.1
## - nox 1 487.16 11566 1609.4
## - ptratio 1 1194.23 12273 1639.4
## - dis 1 1232.41 12311 1641.0
## - rm 1 1871.32 12950 1666.6
## - lstat 1 2410.84 13490 1687.3
##
## Step: AIC=1587.65
## medv ~ crim + zn + indus + chas + nox + rm + dis + rad + tax +
## ptratio + black + lstat
##
## Df Sum of Sq RSS AIC
## - indus 1 2.52 11081 1585.8
## <none> 11079 1587.7
## + age 1 0.06 11079 1589.6
## - chas 1 219.91 11299 1595.6
## - tax 1 242.24 11321 1596.6
## - crim 1 243.20 11322 1596.6
## - zn 1 260.32 11339 1597.4
## - black 1 272.26 11351 1597.9
## - rad 1 481.09 11560 1607.2
## - nox 1 520.87 11600 1608.9
## - ptratio 1 1200.23 12279 1637.7
## - dis 1 1352.26 12431 1643.9
## - rm 1 1959.55 13038 1668.0
## - lstat 1 2718.88 13798 1696.7
##
## Step: AIC=1585.76
## medv ~ crim + zn + chas + nox + rm + dis + rad + tax + ptratio +
## black + lstat
##
## Df Sum of Sq RSS AIC
## <none> 11081 1585.8
## + indus 1 2.52 11079 1587.7
## + age 1 0.06 11081 1587.8
## - chas 1 227.21 11309 1594.0
## - crim 1 245.37 11327 1594.8
## - zn 1 257.82 11339 1595.4
## - black 1 270.82 11352 1596.0
## - tax 1 273.62 11355 1596.1
## - rad 1 500.92 11582 1606.1
## - nox 1 541.91 11623 1607.9
## - ptratio 1 1206.45 12288 1636.0
## - dis 1 1448.94 12530 1645.9
## - rm 1 1963.66 13045 1666.3
## - lstat 1 2723.48 13805 1695.0
##
## Call:
## lm(formula = medv ~ crim + zn + chas + nox + rm + dis + rad +
## tax + ptratio + black + lstat, data = Boston)
##
## Coefficients:
## (Intercept) crim zn chas nox rm
## 36.341145 -0.108413 0.045845 2.718716 -17.376023 3.801579
## dis rad tax ptratio black lstat
## -1.492711 0.299608 -0.011778 -0.946525 0.009291 -0.522553
Pada metode Forward and Backward Stepwise, diperoleh hasil yang sama dengan metode Backward Stepwise. Model ini memiliki nilai AIC 1585.76 dan dapat dituliskan sebagai berikut
\[\begin{aligned} &\hat{medv}=36.341145-0.108413*crim+0.045845*zn+2.718716*indus-17.376023*nox+3.801579*rm-1.492711*dis+0.299608*rad-0.011778*tax-0.946525*ptratio+0.009291*black-0.522553*lstat\\ \end{aligned}\]
Model terbaik yang didapatkan adalah model hasil dari metode Backward Stepwise dan Forward and Backward Stepwise. Model tersebut dapat dituliskan sebagai berikut
\[\begin{aligned} &\hat{medv}=36.341145-0.108413*crim+0.045845*zn+2.718716*indus-17.376023*nox+3.801579*rm-1.492711*dis+0.299608*rad-0.011778*tax-0.946525*ptratio+0.009291*black-0.522553*lstat\\ \end{aligned}\]
corCostumer <- cor(Costumer)
ggcorrplot(corCostumer,lab = TRUE,type="lower")
Berdasarkan output di atas, ketiga peubah cenderung memiliki korelasi yang rendah satu sama lain. Korelasi terbesar dimiliki oleh peubah Age dengan Spending Score, yakni sebesar -0.327. Sedangkan korelasi paling rendah dimiliki oleh peubah Annual Income dengan Spending Score, yakni sebesar 0.010.
data.mean <- data.frame("Variable"=c("Age","Annual Income","Spending Score"),"Mean"=c(mean(Costumer$Age),mean(Costumer$Annual.Income),mean(Costumer$Spending.Score)))
ggplot(data.mean,aes(x=Variable,y=Mean))+geom_bar(stat="identity",fill=c("red","green","blue"))+
theme_classic()+ggtitle(label = "Barplot of Each Variable")
Berdasarkan output di atas, peubah Annual Income memiliki mean paling tinggi diantara kedua peubah yang lain, yakni sebesar 60.56. Sedangkan mean paling rendah dimiliki oleh peubah Age, yakni sebesar 38.85. Mean yang ditampilkan pada output di atas merupakan nilai mean dari data yang belum terstandarisasi.
plot_histogram(data = Costumer,ncol = 3,
geom_histogram_args = list(fill="blue"),title="Variables Distribution")
Standarisasi peubah merupakan proses transformasi peubah menjadi peubah yang memiliki rata-rata nol dan simpangan baku satu. Proses standarisasi ini dilakukan jika kita melihat perbedaan satuan pengukuran peubah-peubah yang digunakan. Standarisasi dilakukan karena metode k-means menggunakan konsep jarak antara objek/amatan, yang mana sensitif terhadap satuan pengukuran.
Costumer.stdz = scale(Costumer)
apply(Costumer.stdz, 2, mean)
## Age Annual.Income Spending.Score
## -1.016906e-16 -8.144310e-17 -1.096708e-16
apply(Costumer.stdz, 2, sd)
## Age Annual.Income Spending.Score
## 1 1 1
Terlihat pada output di atas, nilai rata-rata (mean) dari tiap-tiap peubah mendekati nilai nol dan simpangan baku-nya bernilai satu.
# Optimum Cluster
fviz_nbclust(Costumer.stdz, FUNcluster = hcut, method = "silhouette", hc_method = "complete", hc_metric="euclidean")
Terlihat berdasarkan output di atas, dengan menggunakan metode complete linkage diperoleh jumlah klaster optimal yakni sebanyak lima klaster.
# Dendogram
fviz_dend(hclust(dist(Costumer.stdz, method = "euclidean"), method = "complete"),k=5)
# Cluster's Group
hc.data<- eclust(Costumer, stand = TRUE, FUNcluster = "hclust", k=5, hc_method = "complete", hc_metric = "euclidean", graph = F)
datatable(aggregate(Costumer,by=list(cluster=hc.data$cluster), mean))
cek1 <- data.frame("Observasi"=c(seq(1:200)),"Klaster"=hc.data$cluster)
datatable(aggregate(Observasi~Klaster,cek1,paste,collapse=", "))
# Cluster's Plot
fviz_cluster(hc.data)
Terlihat berdasarkan output di atas, plot klaster terbagi menjadi lima dan terdapat beberapa amatan yang masih tumpang tindih dalam klasternya.
# Optimum Cluster
fviz_nbclust(Costumer.stdz, FUNcluster = hcut, method = "silhouette", hc_method = "average", hc_metric="euclidean")
Terlihat berdasarkan output di atas, dengan menggunakan metode average linkage diperoleh jumlah klaster optimal yakni sebanyak lima klaster. Hasil ini sama dengan metode complete linkage.
# Dendogram
fviz_dend(hclust(dist(Costumer.stdz, method = "euclidean"), method = "average"),k=5)
# Cluster's Group
hc.data<- eclust(Costumer, stand = TRUE, FUNcluster = "hclust", k=5, hc_method = "average", hc_metric = "euclidean", graph = F)
datatable(aggregate(Costumer,by=list(cluster=hc.data$cluster), mean))
cek2 <- data.frame("Observasi"=c(seq(1:200)),"Klaster"=hc.data$cluster)
datatable(aggregate(Observasi~Klaster,cek2,paste,collapse=", "))
# Cluster's Plot
fviz_cluster(hc.data)
Terlihat berdasarkan output di atas, plot klaster terbagi menjadi lima dan terdapat beberapa amatan yang masih tumpang tindih dalam klasternya. Hal ini juga tidak jauh berbeda dengan penggunaan metode complete linkage.
# Optimum Cluster
fviz_nbclust(Costumer.stdz, FUNcluster = hcut, method = "silhouette", hc_method = "centroid", hc_metric="euclidean")
Terlihat dari output, jumlah klaster optimal dari metode centroid linkage adalah sebanyak dua klaster saja.
# Dendogram
fviz_dend(hclust(dist(Costumer.stdz, method = "euclidean"), method = "centroid"),k=2)
# Cluster's Group
hc.data<- eclust(Costumer, stand = TRUE, FUNcluster = "hclust", k=2, hc_method = "centroid", hc_metric = "euclidean", graph = F)
datatable(aggregate(Costumer,by=list(cluster=hc.data$cluster), mean))
cek3 <- data.frame("Observasi"=c(seq(1:200)),"Klaster"=hc.data$cluster)
datatable(aggregate(Observasi~Klaster,cek3,paste,collapse=", "))
# Cluster's Plot
fviz_cluster(hc.data)
Hasil yang berbeda terlihat pada output di atas apabila dibandingkan dengan metode average dan complete linkage. Pembagian amatan/observasi yang masuk ke dalam tiap-tiap klaster yang tidak seimbang. Hanya terdapat tiga amatan yang masuk ke dalam klaster 2.
# Optimum Cluster
fviz_nbclust(Costumer.stdz, FUNcluster = hcut, method = "silhouette", hc_method = "single", hc_metric="euclidean")
Terlihat dari output di atas, jumlah klaster optimal dari metode single linkage adalah sebanyak dua klaster saja. Hal ini sama dengan hasil menggunakan metode centroid linkage.
# Dendogram
fviz_dend(hclust(dist(Costumer.stdz, method = "euclidean"), method = "single"),k=2)
# Cluster's Group
hc.data<- eclust(Costumer, stand = TRUE, FUNcluster = "hclust", k=2, hc_method = "single", hc_metric = "euclidean", graph = F)
datatable(aggregate(Costumer,by=list(cluster=hc.data$cluster), mean))
cek4 <- data.frame("Observasi"=c(seq(1:200)),"Klaster"=hc.data$cluster)
datatable(aggregate(Observasi~Klaster,cek4,paste,collapse=", "))
# Cluster's Plot
fviz_cluster(hc.data)
Hasil yang cukup unik terlihat pada output di atas. Mayoritas amatan/observasi masuk ke dalam klaster 1. Hanya terdapat satu amatan yang masuk ke dalam klaster 2.
Dalam menentukan metode terbaik dalam analisis klaster berhirarki dapat dilakukan dengan beberapa pendekatan, antara lain dengan secara visual/eksploratif atau perhitungan menggunakan koefisien korelasi cophenetic. Dalam analisis kali ini, hanya digunakan pedekatan secara visual atau eksploratif dalam menentukan metode terbaik bagi klaster berhirarki. Metode terbaik tersebut adalah complete linkage. Dapat terlihat pada output metode complete linkage, hanya terdapat sedikit amatan yang terdeteksi tumpang tindih apabila dibandingkan dengan metode average linkage. Selain itu terlihat pula pada cluster’s plot-nya, klaster yang terbentuk pembagiannya lebih merata apabila dibandingkan dengan metode-metode yang lain. Hal ini dapat dilihat pada tabel anggota klaster di atas, perbandingan anggota antara klaster 1 hingga 5 cukup seimbang.
#Penentuan k dengan within sum square
fviz_nbclust(Costumer.stdz, FUNcluster = kmeans, method = "wss")
Penentuan jumlah klaster berdasarkan plot optimal diambil pada garis siku (siku turun untuk koefisien WSS). Berdasarkan output di atas, terlihat bahwa titik jumlah klaster = 6 paling membentuk siku dibandingkan pada jumlah yang lain. Dengan demikian, kandidat jumlah klaster untuk metode K-Means yakni k=6.
#Ditetapkan k=6
kmeans.data <- eclust(Costumer, stand = TRUE, FUNcluster = "kmeans", k=6, graph = F)
datatable(aggregate(Costumer, by=list(cluster=kmeans.data$cluster), FUN = min))
datakmeans <- data.frame("Observasi"=c(seq(1:200)),"Klaster"=kmeans.data$cluster)
datatable(aggregate(Observasi~Klaster,datakmeans,paste,collapse=", "))
fviz_cluster(kmeans.data)
Terlihat pada output metode complete linkage, gerombol 1 merupakan gerombol dari pengunjung/pelanggan remaja-dewasa dengan annual income dan spending score yang dinilai cukup. Gerombol 2 merupakan gerombol dari pengunjung/pelanggan remaja-dewasa yang mempunyai annual income rendah tetapi spending score yang tinggi. Gerombol 3 merupakan gerombol dari pengunjung/pelanggan berusia dewasa-lanjut yang mempunyai annual income dan spending score yang dinilai cukup. Gerombol 4 merupakan gerombol dari pengunjung/pelanggan berusia matang yang memiliki annual income dan spending score yang tinggi. Gerombol 5 merupakan gerombol dari pengunjung/pelanggan berusia matang yang memiliki annual income tinggi tetapi spending score yang rendah.
Berdasarkan output metode K-means di atas, 6 gerombol yang terbentuk merepresentasikan 6 kelompok dari pelanggan/pengunjung dari suatu pusat perbelajaan. Gerombol 1 merupakan pengunjung/pelanggan dengan usia dewasa-matang yang memiliki annual income tinggi tetapi spending score yang rendah. Gerombol 2 merupakan pengunjung/pelanggan dengan usia matang-lanjut yang mempunyai annual income dan spending score yang cenderung seimbang atau tetapi di bawah rata-rata. Gerombol 3 merupakan pengunjung/pelanggan berusia remaja-dewasa yang memiliki annual income rendah tetapi spending score yang tinggi. Gerombol 4 merupakan pengunjung/pelanggan berusia dewasa-lanjut yang memiliki annual income tinggi tetapi spending score rendah. Gerombol 5 merupakan pengunjung/pelanggan berusia dewasa-matang yang memiliki annual income dan spending score yang tinggi. Gerombol 6 merupakan pengunjung/pelanggan berusia dewasa-matang yang memiliki annual income dan spending score yang cukup seimbang di bawah rata-rata.
Dari dua kandidat metode analisis gerombol di atas, hanya dipilih satu metode terbaik yakni analisis gerombol berhirarki dengan menggunakan complete linkage. Metode tersebut dipilih karena secara plot lebih sedikit terjadi “tumpang tindih” apabila dibandingkan dengan metode gerombol tak-berhirarki K-Means k=6. Selain itu, 5 gerombol yang terbentuk berdasarkan metode complete linkage sudah dapat merepresentasikan setiap kelompok dari pengunjung/pelanggan suatu pusat perbelanjaan.
Manajemen suatu pusat perbelanjaan dapat menargetkan pengunjung/pelanggan yang berusia antara 20-29 tahun atau 30-40 tahun. Pada usia 20-29 tahun, pengunjung/pelanggan umumnya merupakan remaja-dewasa yang pergi ke suatu pusat perbelanjaan untuk sekedar refreshing, “cuci mata”, atau membeli keperluan yang diinginkan. Sedangkan, pada usia 30-40 tahun umumnya mereka merupakan pengunjung/pelanggan yang pergi ke suatu pusat perbelanjaan untuk memenuhi kebutuhan hidup/rumah tangga mereka. Hal ini juga didukung dengan annual income yang cenderung tinggi pada rentang usia tersebut.
Dalam analisis suatu kasus dapat digunakan dua analisis gerombol, yakni analisis gerombol berhirarki dan tak berhirarki. Pada studi kasus costumer mall, data lebih cocok untuk dianalisis menggunakan metode analisis gerombol berhirarki complete linkage. Hal ini disebabkan oleh jumlah klaster yang optimal sehingga gerombolnya dapat merepresentasikan tiap-tiap kelompok dengan baik dan minimnya terjadi “tumpang tindih” pada klasternya. Pada studi kasus tersebut juga manajemen suatu pusat perbelanjaan dapat lebih menargetkan penjualan kepada pengunjung/pelanggan remaja-dewasa yang berusia 20-29 tahun. Selain itu, manajemen juga perlu untuk menjaga kualitas dari produk/barang supaya dapat menjaga kepercayaan pengunjung/pelanggan dari kalangan dewasa-matang yang berusia 30-40 tahun.