Regresi Linear Berganda adalah model regresi linear dengan melibatkan lebih dari satu variable bebas atau predictor. Dalam bahasa inggris, istilah ini disebut dengan multiple linear regression.
Disini data diambil dari Excel. Excel adalah salah satu jenis file eksternal yang sering digunakan untuk menyimpan data. Kita dapat menggunakan package {readxl} dengan fungsi read_excel() untuk import data dari file Excel. Argumen path = adalah lokasi dan nama file Excel yang akan kita gunakan.
library(readxl)
## Warning: package 'readxl' was built under R version 4.1.2
datacovidmobility <- read_excel(path = "dataPositifMobilityJuni2020.xlsx")
datacovidmobility
summary() atau ringkasan data digunakan untuk mencari nilai statistik lima serangkai (minimum, Q1, Q2 atau median, Q3, maksimum).
summary(datacovidmobility)
## Tanggal POSITIF Dirawat Sembuh
## Min. :2020-06-01 00:00:00 Min. : 7383 Min. : 951 Min. :2246
## 1st Qu.:2020-06-08 06:00:00 1st Qu.: 8097 1st Qu.:1338 1st Qu.:3246
## Median :2020-06-15 12:00:00 Median : 9176 Median :1396 Median :4264
## Mean :2020-06-15 12:00:00 Mean : 9393 Mean :1416 Mean :4230
## 3rd Qu.:2020-06-22 18:00:00 3rd Qu.:10423 3rd Qu.:1446 3rd Qu.:5203
## Max. :2020-06-30 00:00:00 Max. :13040 Max. :1794 Max. :6512
## Meninggal Self Isolation
## Min. :521.0 Min. :2764
## 1st Qu.:540.2 1st Qu.:2796
## Median :581.5 Median :2868
## Mean :581.0 Mean :2916
## 3rd Qu.:618.8 3rd Qu.:2952
## Max. :641.0 Max. :3299
## retail_and_recreation_percent_change_from_baseline
## Min. :-52.00
## 1st Qu.:-44.00
## Median :-38.00
## Mean :-38.63
## 3rd Qu.:-33.00
## Max. :-27.00
## grocery_and_pharmacy_percent_change_from_baseline
## Min. :-18.00
## 1st Qu.:-15.00
## Median :-13.00
## Mean :-12.97
## 3rd Qu.:-12.25
## Max. : -3.00
## parks_percent_change_from_baseline
## Min. :-91.00
## 1st Qu.:-88.75
## Median :-85.50
## Mean :-86.03
## 3rd Qu.:-84.00
## Max. :-77.00
## transit_stations_percent_change_from_baseline
## Min. :-69.00
## 1st Qu.:-54.00
## Median :-49.50
## Mean :-50.97
## 3rd Qu.:-47.00
## Max. :-40.00
## workplaces_percent_change_from_baseline
## Min. :-68.00
## 1st Qu.:-35.75
## Median :-34.00
## Mean :-32.47
## 3rd Qu.:-28.25
## Max. :-15.00
## residential_percent_change_from_baseline
## Min. :10.00
## 1st Qu.:14.00
## Median :16.00
## Mean :16.17
## 3rd Qu.:18.00
## Max. :27.00
Untuk membuat matriks scatterplot kita hanya perlu memasukkan objek datacovidmobility kedalam fungsi pairs(). Berikut adalah sintaks yang digunakan dan output yang dihasilkan
pairs(datacovidmobility)
Kita juga dapat melakukan drop terhadap panel bawah grafik tersebut, yaitu dengan memasukkan argumen lower.panel=NULL.
pairs(datacovidmobility, lower.panel=NULL)
Fungsi plot() merupakan fungsi umum yang digunakan untuk membuat plot pada R.
plot(datacovidmobility$POSITIF ~ datacovidmobility$retail_and_recreation_percent_change_from_baseline, data=datacovidmobility)
Memvisualisasikan Data dengan Positif sebagai variabel Y dan Google Mobility Index sebagai variabel X
plot(datacovidmobility$POSITIF ~ datacovidmobility$
retail_and_recreation_percent_change_from_baseline+datacovidmobility$
grocery_and_pharmacy_percent_change_from_baseline+datacovidmobility$
parks_percent_change_from_baseline+datacovidmobility$
transit_stations_percent_change_from_baseline+datacovidmobility$
workplaces_percent_change_from_baseline+datacovidmobility$
residential_percent_change_from_baseline, data = datacovidmobility)
Korelasi merupakan keterhubungan antar variabel. Untuk mengukur seberapa jauh hubungan antara satu variabel dengan variabel yang lain kita dapat menggunakan fungsi cor().
a. Korelasi variabel y dengan x1
cor(datacovidmobility$POSITIF,datacovidmobility$
retail_and_recreation_percent_change_from_baseline)
## [1] 0.7612531
Korelasi variabel y dengan x2
cor(datacovidmobility$POSITIF,datacovidmobility$
grocery_and_pharmacy_percent_change_from_baseline)
## [1] 0.557789
Korelasi variabel y dengan x3
cor(datacovidmobility$POSITIF,datacovidmobility$
parks_percent_change_from_baseline)
## [1] 0.7128829
Korelasi variabel y dengan x4
cor(datacovidmobility$POSITIF,datacovidmobility$
transit_stations_percent_change_from_baseline)
## [1] 0.7593921
Korelasi variabel y dengan x5
cor(datacovidmobility$POSITIF,datacovidmobility$
workplaces_percent_change_from_baseline)
## [1] 0.3162402
Korelasi variabel y dengan x6
cor(datacovidmobility$POSITIF,datacovidmobility$
residential_percent_change_from_baseline)
## [1] -0.542098
Dari hasil seluruh output di atas dapat disimpulkan bahwa untuk tingkat keterhubungan antara variabel y dengan x1, x2, x3, x4, dan x5 memiliki hubungan karena nilai yang dihasilkan berjumlah lebih dari 0.3. Sedangkan untuk tingkat keterhubungan antara variabel y dan x6 tidak sangat terhubung karena nilai yang dihasilkan berjumlah kurang dari 0.
Berikut cara melakukan permodelan regresi linier berganda.
model <- lm(datacovidmobility$POSITIF ~ datacovidmobility$Tanggal, data = datacovidmobility)
Selanjutnya dengan model yang telah dibuat di atas, kita akan menggunakan fungsi summary() untuk menjelaskan atau mereview hasil dari model tersebut. Dan dengan ringkasan summary(model) kita dapat melihat informasi terperinci tentang kinerja dan koefisian model.
summary(model)
##
## Call:
## lm(formula = datacovidmobility$POSITIF ~ datacovidmobility$Tanggal,
## data = datacovidmobility)
##
## Residuals:
## Min 1Q Median 3Q Max
## -358.0 -316.5 -245.1 -164.1 3437.0
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.570e+06 3.637e+05 -7.065 1.10e-07 ***
## datacovidmobility$Tanggal 1.620e-03 2.285e-04 7.091 1.03e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 935.7 on 28 degrees of freedom
## Multiple R-squared: 0.6423, Adjusted R-squared: 0.6295
## F-statistic: 50.28 on 1 and 28 DF, p-value: 1.028e-07
Di atas merupakan rincian dari model yang telah dibuat.
Di posisi paling atas terdapat lm formula adalah datacovidmobility\(POSITIF ~ datacovidmobility\)Tanggal, data = datacovidmobility.
Lalu di bawahnya terdapat 5 nilai residual, sebelumnya kita perlu tahu bahwa Residual adalah perbedaan antara nilai nyata dan nilai prediksi. Yang mana semakin kecil nilai residual maka semakin baik atau benar model yang kita buat. Berikut nilai-nilai residual yang dihasilkan:
Nilai minimum = -358.0
Nilai maximum = 3437.0
Nilai median = -245.1
Nilai quartil 1 = -316.5
Nilai quartil 3 = -164.1
Dari nilai-nilai tersebut dapat kita lihat bahwa dalam konteks ini berupa nilai minimum, maximum, median, quartil 1 dan quartil 3. Dapat kita simpulkan bahwa model yang telah kita buat belum bisa dikatakan baik atau benar karena nilai-nilai yang dihasilkan tidak mendekati nol.
Di bawah nilai residual terdapat koefisien, yang mana dalam koefisien tersebut terdapat nilai intersep, dan tanggal. Selain itu juga terdapat nilai-p dari koefisien
Selanjutnya terdapat dua R2 yaitu:
Multiple R-squared (pelipatan daripada variabel): 0.6423. Hal ini menunjukkan bahwa 0.006423% variasi variabel respon, y, dapat dijelaskan oleh variabel prediktor, x.
Multiple R-squared tidak dapat berkurang saat kita menambahkan lebih banyak variabel independen ke model yang kita buat.
Adjusted R-squared (tidak pelipatan pada variabel): 0.6295.
Adjusted R-squared lebih baik ada penambahan variabel. Jadi jika kita menambahkan lebih dari satu variabel ke model, itu hanya meningkat jika itu mengurangi kesalahan prediksi secara keseluruhan
ANOVA (analysis of variance) adalah pengujian yang dilakukan dengan membandingkan varians. Dengan membandingkan varians tersebut, dapat diketahui ada tidaknya perbedaan rata-rata dari tiga atau lebih kelompok. Asumsi normalitas pada ANOVA adalah pada residual yaitu selisih antara Y Prediksi dengan Y Aktual. Tepatnya residual dapat dihitung sebagai berikut: Y Aktual – Y Prediksi. Dimana Y Aktual adalah Y sesungguhnya atau kenyataan. Sedangkan Y prediksi adalah Y hasil persamaan ANOVA.
anova(model)
plot(datacovidmobility$POSITIF ~ datacovidmobility$
Tanggal,
data = datacovidmobility, col = "dodgerblue", pch = 20, cex = 1.5,
main = "Data Covid POSITIF di DKI Jakarta dan Google Mobility Index")
abline(model)
#Add a regression line
Dari Plot di atas perlu kita ketahui bahwa titik-titik hijau yang ada pada grafik tersebut adalah data real dan garis hitam di dalam kotak adalah data prediksi.
plot(cooks.distance(model), pch = 16, col = "dodgerblue")
#Plot the Cooks Distances.
plot(model)
AIC dan BIC banyak digunakan dalam kriteria pemilihan model. AIC berarti Kriteria Informasi Akaike dan BIC berarti Kriteria Informasi Bayesian. Meskipun kedua istilah ini membahas pemilihan model, keduanya tidak sama. Seseorang dapat menemukan perbedaan antara dua pendekatan pemilihan model.
AIC(model)
## [1] 499.547
BIC(model)
## [1] 503.7506
head(predict(model), n = 15)
## 1 2 3 4 5 6 7 8
## 7363.645 7503.608 7643.570 7783.532 7923.495 8063.457 8203.420 8343.382
## 9 10 11 12 13 14 15
## 8483.344 8623.307 8763.269 8903.232 9043.194 9183.156 9323.119
plot(head(predict(model), n = 20))
head(resid(model), n = 11)
## 1 2 3 4 5 6 7
## 19.35484 -44.60756 -104.56997 -183.53237 -239.49477 -277.45717 -257.41958
## 8 9 10 11
## -306.38198 -207.34438 -200.30679 -211.26919
coef(model)
## (Intercept) datacovidmobility$Tanggal
## -2.569904e+06 1.619935e-03
Tabel di bawah merupakan semua proses yang telah dilakukan, sehingga terbuat tabel yang ada nilai residuals dan nilai protected.
datacovidmobility$residuals <- model$residuals
datacovidmobility$predicted <- model$fitted.values
datacovidmobility
scatter.smooth(x=datacovidmobility$Tanggal, y=datacovidmobility$POSITIF,
main="Tanggal ~ POSITIF")
boxplot(datacovidmobility$POSITIF, main="POSITIF",
boxplot.stats(datacovidmobility$POSITIF)$out)
plot(density(datacovidmobility$POSITIF), main="Google Mobility Index: POSITIF",
ylab="Frequency")
coefs <- coef(model)
plot(POSITIF ~ Tanggal, data = datacovidmobility)
abline(coefs)
text(x = 12, y = 10, paste('expression = ', round(coefs[1], 2), '+',
round(coefs[2], 2), '*POSITIF'))
Adanya korelasi antar variabel dapat dilakukan melalui visualisasi menggunakan scatterplot dan perhitungan matematis menggunakan metode Pearson untuk metode parametrik dan metode rangking Spearman dan Kendall untuk metode non-parametrik. Pada R uji korelasi dapat dilakukan dengan menggunakan fungsi cor.test(). Format fungsi tersebut adalah sebagai berikut:
Uji Korelasi Variabel y dengan x1
cor.test(datacovidmobility$
retail_and_recreation_percent_change_from_baseline,
datacovidmobility$POSITIF)
##
## Pearson's product-moment correlation
##
## data: datacovidmobility$retail_and_recreation_percent_change_from_baseline and datacovidmobility$POSITIF
## t = 6.212, df = 28, p-value = 1.037e-06
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.5525143 0.8801389
## sample estimates:
## cor
## 0.7612531
Uji Korelasi Variabel y dengan x2
cor.test(datacovidmobility$
grocery_and_pharmacy_percent_change_from_baseline,
datacovidmobility$POSITIF)
##
## Pearson's product-moment correlation
##
## data: datacovidmobility$grocery_and_pharmacy_percent_change_from_baseline and datacovidmobility$POSITIF
## t = 3.5561, df = 28, p-value = 0.001362
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.2471945 0.7644406
## sample estimates:
## cor
## 0.557789
Uji Korelasi Variabel y dengan x3
cor.test(datacovidmobility$
parks_percent_change_from_baseline,
datacovidmobility$POSITIF)
##
## Pearson's product-moment correlation
##
## data: datacovidmobility$parks_percent_change_from_baseline and datacovidmobility$POSITIF
## t = 5.379, df = 28, p-value = 9.849e-06
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.4744723 0.8538564
## sample estimates:
## cor
## 0.7128829
Uji Korelasi Variabel y dengan x4
cor.test(datacovidmobility$
transit_stations_percent_change_from_baseline,
datacovidmobility$POSITIF)
##
## Pearson's product-moment correlation
##
## data: datacovidmobility$transit_stations_percent_change_from_baseline and datacovidmobility$POSITIF
## t = 6.176, df = 28, p-value = 1.141e-06
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.5494424 0.8791410
## sample estimates:
## cor
## 0.7593921
Uji Korelasi Variabel y dengan x5
cor.test(datacovidmobility$
workplaces_percent_change_from_baseline,
datacovidmobility$POSITIF)
##
## Pearson's product-moment correlation
##
## data: datacovidmobility$workplaces_percent_change_from_baseline and datacovidmobility$POSITIF
## t = 1.7639, df = 28, p-value = 0.08866
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.04969027 0.60731687
## sample estimates:
## cor
## 0.3162402
Uji Korelasi Variabel y dengan x6
cor.test(datacovidmobility$
residential_percent_change_from_baseline,
datacovidmobility$POSITIF)
##
## Pearson's product-moment correlation
##
## data: datacovidmobility$residential_percent_change_from_baseline and datacovidmobility$POSITIF
## t = -3.4136, df = 28, p-value = 0.001972
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.7549287 -0.2259588
## sample estimates:
## cor
## -0.542098
Berdasarkan seluruh output yang dihasilkan, metode Pearson menghasilkan output berupa nilai t uji, derajat kebebasan, nilai p-value, rentang estimasi nilai korelasi berdasarkan tingkat kepercayaan, dan estimasi nilai korelasi.