Materi 2 : Implementasi GLM (Y Normal) dalam R
Akbar Rizki, S.Stat. M.Si
5 Februari 2021 14.00-15.00 WIB
Data
Berikut ini adalah Data Contoh dengan n = 30
## Rows: 30
## Columns: 4
## $ y <dbl> 51, 57, 53, 64, 53, 45, 44, 52, 53, 50, 56, 51, 45, 38, 50, 55, ...
## $ x1 <dbl> 35, 32, 31, 36, 26, 26, 27, 32, 36, 33, 33, 30, 30, 26, 29, 29, ...
## $ x2 <dbl> 15, 10, 11, 10, 9, 10, 6, 13, 10, 14, 11, 9, 11, 8, 7, 11, 9, 10...
## $ x3 <dbl> 14, 18, 17, 19, 18, 15, 12, 15, 14, 14, 16, 15, 12, 12, 15, 19, ...
Data yang terdiri dari 30 baris dan 4 kolom tersebut bisa Bapak/Ibu unduh dari http://ipb.link/datalatihan, dengan nama data_02.csv, data_02.sav, data_02.txt, data_02.xlsx lalu menyimpannya pada folder kerja, misal E:\Pelatihan
Data tersebut juga terdapat pada github https://github.com/Nr5D/data
- y = penghasilan per bulan (dalam ratus ribu rupiah)
- x1 = usia (dalam tahun)
- x2 = lama pendidikan (dalam tahun)
- x3 = lama bekerja (dalam bulan)
Eksplorasi Data (1) : Histogram
hist(normal$y, main="Histogram y", ylab="Frekuensi", xlab="y", col="lightblue")
Keterangan Syntax: 1. membuat histogram melalui fungsi hist dari variabel y di objek normal dengan perintah normal$y
2. main="Histogram y" membuat judul dari histogram
3. xlab="y" memberikan nama pada sumbu x
4. ylab="Frekuensi", memberikan nama pada sumbu y
5. col="lighblue" ini memberikan warna biru muda pada Histogram. Referensi lebih lanjut terkait warna: http://www.sthda.com/english/wiki/colors-in-r
Histogram menunjukkan bahwa Y cenderung menyebar normal karena pola histogram menyerupai lonceng.
Eksplorasi Data (2) : Korelasi
cor(normal, method="pearson")
## y x1 x2 x3
## y 1.0000000 0.8921240 0.2982323 0.9012543
## x1 0.8921240 1.0000000 0.4311743 0.6397021
## x2 0.2982323 0.4311743 1.0000000 0.1977660
## x3 0.9012543 0.6397021 0.1977660 1.0000000
Implementasi GLM (1)
Bentuk umum dari Syntax glm di R adalah :
dimana dalam kasus ini
formula = y ~ x1 + x2 + x3 atau \(y_{duga} = b_0 + b_1 x_1 + b_2 x_2 + b_3 x_3\)
family = gaussian (karena y menyebar normal)
data berasal dari data frame normal
Sehingga syntax pada R adalah
gfit<-glm(y~x1+x2+x3, family = gaussian(link = "identity"), data=normal)
syntax ini, juga menyimpan hasil dari fungsi glm ke objek gfit.
Implementasi GLM (2) : Summary
##
## Call:
## glm(formula = y ~ x1 + x2 + x3, family = gaussian(link = "identity"),
## data = normal)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.13177 -0.78508 -0.03781 0.95502 1.74386
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.26774 1.63634 -1.997 0.0564 .
## x1 0.95734 0.06179 15.494 1.20e-14 ***
## x2 -0.23966 0.12806 -1.871 0.0726 .
## x3 1.78549 0.10782 16.560 2.49e-15 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for gaussian family taken to be 1.434158)
##
## Null deviance: 2208.800 on 29 degrees of freedom
## Residual deviance: 37.288 on 26 degrees of freedom
## AIC: 101.66
##
## Number of Fisher Scoring iterations: 2
Keterangan Output
- Pada \(\alpha\) = 5%, diketahui bahwa
– Intercept dan \(x_2\) tidak nyata
– \(x_1\) dan \(x_3\) nyata
- degree of freedom = derajat bebas, angka ini beserta nilai deviance akan digunakan dalam Uji Deviance
- AIC merupakan ukuran kebaikan model
Implementasi GLM (3)
par(mfrow=c(1,2))
car::qqPlot(gfit$residuals)
## [1] 14 11
urutan<-seq(1:30)
plot(urutan,gfit$residuals, type="b")
Implementasi GLM (4) : Uji Asumsi
Kenormalan Sisaan : Shapiro Test
shapiro.test(gfit$residuals)
##
## Shapiro-Wilk normality test
##
## data: gfit$residuals
## W = 0.95399, p-value = 0.216
Hipotesis pada Uji Shapiro ini adalah
\(H_0\) : sisaan menyebar normal
\(H_1\) : sisaan tidak menyebar normal
Harapan : sisaan menyebar normal
Harapan : tidak tolak \(H_0\) yaitu nilai p-value > \(\alpha\)
misal \(\alpha =\) 5%
dari output p-value = 0.216 ( 21.6% ) lebih besar dari \(\alpha =\) 5%, maka tidak tolak \(H_0\), dapat disimpulkan bahwa sisaan menyebar normal
Implementasi GLM (5)
Prediksi
Misal diketahui data baru berupa :
\(x_1 = 0\)
\(x_2 = 1\)
\(x_3 = 0\)
maka untuk prediksi nilai y, menggunakan fungsi predict, dengan syntax sebagai berikut :
# Memasukkan data baru tersebut menjadi suatu data frame
data_baru<-data.frame(
x1=0,
x2=1,
x3=0
)
data_baru
## x1 x2 x3
## 1 0 1 0
predict(gfit, newdata = data_baru)
## 1
## -3.5074
Implementasi GLM (6)
Prediksi
Misal diketahui data yang akan diprediksi adalah 20 baris data yang terdiri dari x1,x2 dan x3 yang terletak pada folder pribadi atau penyimpanan online, maka
data_baru_banyak<-read.csv("https://raw.githubusercontent.com/Nr5D/data/main/data_03.csv", sep=",", header = TRUE)
dplyr::glimpse(data_baru_banyak)
## Rows: 20
## Columns: 3
## $ x1 <int> 30, 30, 37, 20, 29, 30, 36, 27, 35, 32, 26, 36, 36, 36, 34, 29, ...
## $ x2 <int> 10, 9, 11, 7, 11, 7, 10, 9, 9, 10, 6, 12, 7, 9, 8, 8, 14, 10, 7, 7
## $ x3 <int> 14, 17, 17, 10, 15, 14, 19, 13, 16, 14, 14, 17, 15, 18, 12, 12, ...
hasil_prediksi<-predict(gfit, newdata = data_baru_banyak)
hasil_prediksi
## 1 2 3 4 5 6 7 8
## 48.05279 53.64893 59.87100 32.05639 48.64128 48.77178 62.72431 43.63494
## 9 10 11 12 13 14 15 16
## 56.65015 49.96748 45.18207 58.67400 56.30132 61.17848 48.79050 44.00379
## 17 18 19 20
## 51.36410 44.35262 44.81322 49.59993
Implementasi GLM (7)
Model tanpa x2
gfit2<-glm(y~x1+x3, family = gaussian(link = "identity"), data=normal)
Summary
##
## Call:
## glm(formula = y ~ x1 + x3, family = gaussian(link = "identity"),
## data = normal)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.7666 -0.9149 -0.1130 1.0770 2.3228
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.60703 1.53824 -2.995 0.00582 **
## x1 0.91058 0.05907 15.415 6.64e-15 ***
## x3 1.80820 0.11199 16.146 2.14e-15 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for gaussian family taken to be 1.567076)
##
## Null deviance: 2208.800 on 29 degrees of freedom
## Residual deviance: 42.311 on 27 degrees of freedom
## AIC: 103.45
##
## Number of Fisher Scoring iterations: 2
Implementasi GLM (8)
par(mfrow=c(1,2))
car::qqPlot(gfit2$residuals)
## [1] 14 7
urutan<-seq(1:30)
plot(urutan,gfit2$residuals, type="b")
Implementasi GLM (9) : Uji Asumsi
Kenormalan Sisaan
shapiro.test(gfit2$residuals)
##
## Shapiro-Wilk normality test
##
## data: gfit2$residuals
## W = 0.97113, p-value = 0.5706
Implementasi GLM (10)
Prediksi
Misal diketahui data baru berupa :
\(x_1 = 1\)
\(x_2 = 0\)
\(x_3 = 0\)
maka untuk prediksi nilai y, menggunakan fungsi predict, dengan syntax sebagai berikut :
# Memasukkan data baru tersebut menjadi suatu data frame
data_baru<-data.frame(
x1=1,
x3=0
)
data_baru
## x1 x3
## 1 1 0
predict(gfit2, newdata = data_baru)
## 1
## -3.69645
Implementasi GLM (12)
Prediksi
Misal diketahui data yang akan diprediksi adalah 20 baris data yang terdiri dari x1,x2 dan x3 yang terletak pada folder pribadi atau penyimpanan online, maka
data_baru_banyak<-read.csv("https://raw.githubusercontent.com/Nr5D/data/main/data_03.csv", sep=",", header = TRUE)
dplyr::glimpse(data_baru_banyak)
## Rows: 20
## Columns: 3
## $ x1 <int> 30, 30, 37, 20, 29, 30, 36, 27, 35, 32, 26, 36, 36, 36, 34, 29, ...
## $ x2 <int> 10, 9, 11, 7, 11, 7, 10, 9, 9, 10, 6, 12, 7, 9, 8, 8, 14, 10, 7, 7
## $ x3 <int> 14, 17, 17, 10, 15, 14, 19, 13, 16, 14, 14, 17, 15, 18, 12, 12, ...
hasil_prediksi<-predict(gfit2, newdata = data_baru_banyak)
hasil_prediksi
## 1 2 3 4 5 6 7 8
## 48.02535 53.44996 59.82405 31.68669 48.92297 48.02535 62.52987 43.48539
## 9 10 11 12 13 14 15 16
## 56.19467 49.84651 44.38301 58.91346 55.29705 60.72167 48.05127 43.49835
## 17 18 19 20
## 52.52641 44.39597 44.37005 48.92297
Implementasi GLM (13)
Analysis of Deviance
anova(gfit2,gfit, test = "Chisq")
## Analysis of Deviance Table
##
## Model 1: y ~ x1 + x3
## Model 2: y ~ x1 + x2 + x3
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 27 42.311
## 2 26 37.288 1 5.0229 0.06128 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Hipotesis pada Analysis of Deviance
\(H_0\) : \(x_2\) tidak berpengaruh terhadap y
\(H_1\) : \(x_2\) berpengaruh terhadap y
misal \(\alpha =\) 5%
dari output p-value = 0.06128 ( 6.128% ) lebih besar dari \(\alpha =\) 5%, maka tidak tolak \(H_0\), dapat disimpulkan bahwa \(x_2\) tidak berpengaruh terhadap y