Regresi Ridge dengan R |
Halo teman-teman, kali ini kita akan melanjutkan belajar bersamanya mengenai pemodelan karena tadi pagi kita telah belajar bersama bagaimana menerapkan fungsi do.call() di R. Model yang kita bahas kali ini adalah model yang agaknya masih keluarga dekat dengan model regresi linier sederhana atau berganda, namanya adalah regresi ridge (ridge regression model).
Acapkali data yang kita gunakan untuk pemodelan regresi linier tidak berjalan mulus dan lancar-lancar saja. Entah dari hubungannya sebenarnya tidak linier, atau yang lebih sering mengalami gangguan asumsi klasik tertentu. Regresi ridge ini merupakan bentuk lain dari regresi yang mampu mengakomodir adanya bias akibat adanya multikolinearitas di antara variabel independen di dalam model. Adapun sifat dari penduga parameter model ini adalah bias namun konsisten karena memiliki kemampuan untuk menurunkan Mean Square Error (MSE). Kendati demikian, pada praktiknya, model regresi ridge ini masih debatable atau masih menjadi bahan perdebatan, terutama dari aspek efektivitas modelnya.
Selain itu, karakteristik dari model regresi ridge ini adalah uji signifikansi koefisien modelnya tidak dapat dilakukan secara langsung karena di dalam koefisiennya mengandung bias yang dalam formulasinya dinotasikan sebagai c. Demikian halnya dengan uji analisis varians dari model, kita tidak bisa melakukannya secara langsung. Perdebatan terhadap model ini juga terlihat akibat bias yang dimasukkan ke dalam model digunakan untuk meningkatkan prediksi.
Lantas bagaimana pengaplikasian regresi ridge dengan R? Kali ini kita akan mempraktikkan dengan menggunakan data dummy yang kita bangkitkan secara manual. Adapun langkah-langkah pemodelannya dapat mengikuti beberapa code berikut:
Code:
y <- c(1, 2, 3, 5, 1, 3, 7, 4, 9, 8)
x1 <- c(3, 2, 5, 7, 3, 5, 9, 6, 12, 10)
x2 <- c(1, 2, 3, 5, 1, 3, 6, 4, 8, 8)
df = data.frame(x1, x2, y)
#regresi OLS
ols <- lm(y~x1+x2, data = df)
#ringkasan hasil
summary(ols)
Hasil:
Call:
lm(formula = y ~ x1 + x2, data = df)
Residuals:
Min 1Q Median 3Q Max
-0.31623 -0.16148 -0.09776 0.11088 0.49885
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.4905 0.2416 -2.030 0.08191 .
x1 0.2577 0.1290 1.998 0.08587 .
x2 0.7787 0.1633 4.768 0.00204 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.3087 on 7 degrees of freedom
Multiple R-squared: 0.991, Adjusted R-squared: 0.9884
F-statistic: 385.2 on 2 and 7 DF, p-value: 6.925e-08
Code:
#Uji nonmultikol
library(olsrr)
ols_vif_tol(ols)
#uji homos
library(lmtest)
bptest(ols)
#uji nonautokorel
dwtest(ols)
#uji normalitas
ks.test(ols$residuals, ecdf(ols$residuals))
Hasil:
Variables Tolerance VIF
1 x1 0.0586826 17.04083
2 x2 0.0586826 17.04083
studentized Breusch-Pagan test
data: ols
BP = 2.5441, df = 2, p-value = 0.2803
Durbin-Watson test
data: ols
DW = 2.6447, p-value = 0.8619
alternative hypothesis: true autocorrelation is greater than 0
One-sample Kolmogorov-Smirnov test
data: ols$residuals
D = 0.1, p-value = 0.9996
alternative hypothesis: two-sided
Terlihat bahwa seluruh uji asumsi klasik memenuhi kecuali uji non-multikolinearitas, karena nilai dari VIF > 10, maka alternatifnya menggunakan model regresi ridge
Code:
#PEMODELAN RIDGE REGRESSION
library(glmnet)
#Definisikan variabel dependen
y <- df$y
#Definisikan variabel independen
x <- data.matrix(df[1:2])
#Model ridge
ridge <- glmnet(x, y, alpha = 0)
#ringkasan model ridge
summary(ridge)
Hasil:
Length Class Mode
a0 100 -none- numeric
beta 200 dgCMatrix S4
df 100 -none- numeric
dim 2 -none- numeric
lambda 100 -none- numeric
dev.ratio 100 -none- numeric
nulldev 1 -none- numeric
npasses 1 -none- numeric
jerr 1 -none- numeric
offset 1 -none- logical
call 4 -none- call
nobs 1 -none- numeric
Code:
#mencari lamda optimal dengan K-fold
cv_ridge <- cv.glmnet(x, y, nfolds = 3)
lamop <- cv_ridge$lambda.min
#Plot lambda
plot(cv_ridge)
Hasil:
Plot lambda untuk mendapatkan lambda optimal |
Code:
#Model Final
ridgemod <- glmnet(x, y, alpha = 0, lambda = lamop)
coef(ridgemod)
Hasil:
3 x 1 sparse Matrix of class "dgCMatrix"
s0
(Intercept) -0.3570000
x1 0.3762883
x2 0.5668323
#Model regresi ridgenya: -0,3570000 + 0,3762883x1 + 0,5668323x2
Code:
#Menggunakan model fit untuk prediksi
y_predicted <- predict(ridgemod, s = lamop, newx = x)
#menghitung nilai Sum Square of Total dan Sum Square of Error
sst <- sum((y - mean(y))^2)
sse <- sum((y_predicted - y)^2)
#Menghitung nilai R Square model
rsq <- 1 - sse/sst
rsq
#atau bisa juga dilihat dari dev.ratio
ridgemod$dev.ratio
Hasil:
[1] 0.9863052
[1] 0.9863052
Code:
#Nilai MSE OLS
olspred <- predict(ols, newdata = df)
mean((olspred-df$y)^2)
#Nilai MSE Ridge
ridgepred <- predict(ridgemod, newx = x)
mean((ridgepred-df$y)^2)
Hasil:
OLS: [1] 0.06671445
Rigde: [1] 0.1014786
Terlihat bahwa MSE dari model ridge juga tidak lebih rendah daripada regresi OLS. Sedangkan kalau dilihat dari R square-nya relatif sama nilainya. Baik, demikian sekilas belajar bersama kita mengenai regresi ridge dengan R. Jangan lupa untuk terus menyimak unggahan berikutnya dan jangan lupa share. Selamat mempraktikkan!