Perbandingan Metode Variable Selection dengan Classical Linear Regression, Ridge Regression, dan Lasso Regression dalam Memilih Model dengan Peubah Terbaik.

PACKAGE

library(glmnet)
## Loading required package: Matrix
## Loaded glmnet 4.1-8

INPUT DATA

Data yang digunakan adalah data kandungan persentase lemak tubuh dari 252 amatan (Y). Peubah X yang digunakan merupakan 14 pengukuran kondisi tubuh.

data<- read.csv("C:\\Users\\ACER\\Downloads\\dataModel.csv", sep = ";")
dim(data)
## [1] 252  16
str(data)
## 'data.frame':    252 obs. of  16 variables:
##  $ No : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ Y  : num  12.6 6.9 24.6 10.9 27.8 20.6 19 12.8 5.1 12 ...
##  $ X1 : int  23 22 22 26 24 24 26 25 25 23 ...
##  $ X2 : num  70 78.6 69.9 83.8 83.6 95.4 82.1 79.8 86.6 89.9 ...
##  $ X3 : num  172 184 168 184 181 ...
##  $ X4 : num  23.7 23.4 24.7 24.9 25.6 26.5 26.2 23.6 24.6 25.8 ...
##  $ X5 : num  36.2 38.5 34 37.4 34.4 39 36.4 37.8 38.1 42.1 ...
##  $ X6 : num  93.1 93.6 95.8 101.8 97.3 ...
##  $ X7 : num  85.2 83 87.9 86.4 100 94.4 90.7 88.5 82.5 88.6 ...
##  $ X8 : num  94.5 98.7 99.2 101.2 101.9 ...
##  $ X9 : num  59 58.7 59.6 60.1 63.2 66 58.4 60 62.9 63.1 ...
##  $ X10: num  37.3 37.3 38.9 37.3 42.2 42 38.3 39.4 38.3 41.7 ...
##  $ X11: num  21.9 23.4 24 22.8 24 25.6 22.9 23.2 23.8 25 ...
##  $ X12: num  32 30.5 28.8 32.4 32.2 35.7 31.9 30.5 35.9 35.6 ...
##  $ X13: num  27.4 28.9 25.2 29.4 27.7 30.6 27.8 29 31.1 30 ...
##  $ X14: num  17.1 18.2 16.6 18.2 17.7 18.8 17.7 18.8 18.2 19.2 ...

Menghapus Kolom Nomor

Penghapusan terhadap kolom yang tidak digunakan dalam analisis, yaitu kolom Nomor, sehingga tersisa peubah Y,X1,X2,X3,X4,X5,X6,X7,X8,X9,X10,X11,X12,X13,X14.

data <- data[,-1]
View(data)

Ekplorasi Data

# Sebaran peubah Y 
hist(data$Y, col = "blue")

Statistika Deskriptif

summary(data)
##        Y               X1              X2               X3       
##  Min.   : 0.00   Min.   :22.00   Min.   : 53.80   Min.   : 74.9  
##  1st Qu.:12.80   1st Qu.:35.75   1st Qu.: 72.10   1st Qu.:173.4  
##  Median :19.00   Median :43.00   Median : 80.05   Median :177.8  
##  Mean   :18.94   Mean   :44.88   Mean   : 81.16   Mean   :178.2  
##  3rd Qu.:24.60   3rd Qu.:54.00   3rd Qu.: 89.40   3rd Qu.:183.5  
##  Max.   :45.10   Max.   :81.00   Max.   :164.70   Max.   :197.5  
##        X4              X5              X6               X7        
##  Min.   :18.10   Min.   :31.10   Min.   : 79.30   Min.   : 69.40  
##  1st Qu.:23.10   1st Qu.:36.40   1st Qu.: 94.35   1st Qu.: 84.58  
##  Median :25.05   Median :38.00   Median : 99.65   Median : 90.95  
##  Mean   :25.44   Mean   :37.99   Mean   :100.82   Mean   : 92.56  
##  3rd Qu.:27.32   3rd Qu.:39.42   3rd Qu.:105.38   3rd Qu.: 99.33  
##  Max.   :48.90   Max.   :51.20   Max.   :136.20   Max.   :148.10  
##        X8              X9             X10             X11            X12       
##  Min.   : 85.0   Min.   :47.20   Min.   :33.00   Min.   :19.1   Min.   :24.80  
##  1st Qu.: 95.5   1st Qu.:56.00   1st Qu.:36.98   1st Qu.:22.0   1st Qu.:30.20  
##  Median : 99.3   Median :59.00   Median :38.50   Median :22.8   Median :32.05  
##  Mean   : 99.9   Mean   :59.41   Mean   :38.59   Mean   :23.1   Mean   :32.27  
##  3rd Qu.:103.5   3rd Qu.:62.35   3rd Qu.:39.92   3rd Qu.:24.0   3rd Qu.:34.33  
##  Max.   :147.7   Max.   :87.30   Max.   :49.10   Max.   :33.9   Max.   :45.00  
##       X13             X14       
##  Min.   :21.00   Min.   :15.80  
##  1st Qu.:27.30   1st Qu.:17.60  
##  Median :28.70   Median :18.30  
##  Mean   :28.66   Mean   :18.23  
##  3rd Qu.:30.00   3rd Qu.:18.80  
##  Max.   :34.90   Max.   :21.40

PEMODELAN REGRESI

Pembentukan Model Regresi Awal

model <- lm(Y ~., data = data)
coefficients(model)
##  (Intercept)           X1           X2           X3           X4           X5 
## -15.08665241   0.05690320  -0.17844335  -0.02089591   0.06094166  -0.44552271 
##           X6           X7           X8           X9          X10          X11 
##  -0.03113825   0.87895193  -0.20380813   0.22750425  -0.00151003   0.15667072 
##          X12          X13          X14 
##   0.14848174   0.42908539  -1.47922503
summary(model)
## 
## Call:
## lm(formula = Y ~ ., data = data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -10.262  -2.590  -0.106   2.906   9.281 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -15.08665   16.09640  -0.937  0.34957    
## X1            0.05690    0.03004   1.894  0.05938 .  
## X2           -0.17844    0.10989  -1.624  0.10573    
## X3           -0.02090    0.04070  -0.513  0.60816    
## X4            0.06094    0.27788   0.219  0.82660    
## X5           -0.44552    0.21842  -2.040  0.04248 *  
## X6           -0.03114    0.09780  -0.318  0.75048    
## X7            0.87895    0.08547  10.283  < 2e-16 ***
## X8           -0.20381    0.13702  -1.487  0.13822    
## X9            0.22750    0.13558   1.678  0.09466 .  
## X10          -0.00151    0.22980  -0.007  0.99476    
## X11           0.15667    0.20756   0.755  0.45110    
## X12           0.14848    0.16006   0.928  0.35453    
## X13           0.42908    0.18487   2.321  0.02114 *  
## X14          -1.47923    0.49678  -2.978  0.00321 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.996 on 237 degrees of freedom
## Multiple R-squared:  0.749,  Adjusted R-squared:  0.7342 
## F-statistic: 50.52 on 14 and 237 DF,  p-value: < 2.2e-16
AIC(model)
## [1] 1429.899

Pada metode regresi linier berganda model yang didapatkan yaitu \[y = -15.08665+ 0.05690X_{1}-0.17844X_{2}+ -0.02090X_{3}+0.06094X_{4}-0.44552X_{5}-0.03114 X_{6}+ 0.87895X_{7}-0.20381X_{8}+ 0.22750X_{9}-0.00151X_{10}+ 0.15667X_{11}+0.14848X_{12}+0.42908X_{13}-1.47923 X_{14}\] Dapat dikethui juga nilai \[R-squared=0.7342\] dan nilai AIC sebesar \[AIC=1429.899\]

PEMERIKSAAN ASUMSI

Kemudian dilakukan pemeriksaan asumsi : # Multikolinearitas

car::vif(model)
##        X1        X2        X3        X4        X5        X6        X7        X8 
##  2.251927 33.722299  2.253757 16.151978  4.430814 10.685381 13.350737 15.143865 
##        X9       X10       X11       X12       X13       X14 
##  7.962503  4.827959  1.945075  3.675477  2.193314  3.380745

Terlihat bahwa nilai VIF> 10 memiliki potensi untuk terjadi multikolinearitas yaitu peubah X2,X4,X6,X7,dan X8.

Pengujian Asumsi

Uji Asumsi Normalitas

##Kolmogorov-Smirnov Test
ks.test(model$residuals, "pnorm", mean=mean(model$residuals), sd=sd(model$residuals))
## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  model$residuals
## D = 0.042642, p-value = 0.7492
## alternative hypothesis: two-sided

Berdasarkan Kolmogorov-Smirnov test, residual data menyebar normal dengan p-value > 5%.

##Shapiro-Wilk Test
shapiro.test(model$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  model$residuals
## W = 0.9928, p-value = 0.2618

Berdasarkan Shapiro-Wilk test, residual data menyebar normal dengan p-value > 5%.

Uji Asumsi Homoskedastisitas (Gauss Markov)

lmtest::bptest(model)
## 
##  studentized Breusch-Pagan test
## 
## data:  model
## BP = 20.134, df = 14, p-value = 0.126

Karena p-value > 0.05 maka ragam sisaan homogen atau tidak terdapat masalah heteroskedastisitas

Uji Kebebasan Sisaan (Gauss Markov)

library(randtests)
runs.test(model$residuals)
## 
##  Runs Test
## 
## data:  model$residuals
## statistic = -0.37872, runs = 124, n1 = 126, n2 = 126, n = 252, p-value
## = 0.7049
## alternative hypothesis: nonrandomness

Karena p-value > 0.05 maka sisaan saling bebas.

Nilai harapan sisaan sama dengan nol (Gauss Markov)

t.test(model$residuals,
       mu = 0,
       conf.level = 0.95)
## 
##  One Sample t-test
## 
## data:  model$residuals
## t = 3.4818e-16, df = 251, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.481766  0.481766
## sample estimates:
##   mean of x 
## 8.51701e-17

Karena p-value > 0.05 maka nilai harapan sisaan sama dengan nol

SELEKSI PEUBAH PENJELAS

Forward Selection

olsrr::ols_step_forward_p(model)
## 
##                             Selection Summary                              
## --------------------------------------------------------------------------
##         Variable                  Adj.                                        
## Step    Entered     R-Square    R-Square     C(p)         AIC        RMSE     
## --------------------------------------------------------------------------
##    1    X7            0.6621      0.6608    71.0322    1478.8012    4.5144    
##    2    X2            0.7187      0.7164    19.6026    1434.6119    4.1273    
##    3    X14           0.7275      0.7242    13.2630    1428.5726    4.0702    
##    4    X13           0.7351      0.7308     8.1586    1423.5156    4.0217    
##    5    X5            0.7380      0.7326     7.4209    1422.7425    4.0078    
##    6    X1            0.7404      0.7341     7.0826    1422.3496    3.9969    
##    7    X9            0.7444      0.7371     5.3236    1420.4546    3.9743    
##    8    X8            0.7466      0.7383     5.2258    1420.2544    3.9651    
## --------------------------------------------------------------------------

Berdasarkan metode forward, model terbaik adalah model dengan peubah X1,X2,X5,X7,X8,X9,X13, dan X14 dengan nilai \[R-squared=0.7383 \] dan nilai AIC sebesar \[AIC= 1420.2544\]

Backward Elimination

olsrr::ols_step_backward_p(model)
## 
## 
##                            Elimination Summary                             
## --------------------------------------------------------------------------
##         Variable                  Adj.                                        
## Step    Removed     R-Square    R-Square     C(p)         AIC        RMSE     
## --------------------------------------------------------------------------
##    1    X10            0.749      0.7353    13.0000    1427.8988    3.9878    
##    2    X4            0.7489      0.7363    11.0510    1425.9530    3.9799    
##    3    X6            0.7489      0.7374     9.1164    1424.0225    3.9722    
##    4    X3            0.7484      0.7379     7.5960    1422.5317    3.9679    
##    5    X11           0.7476      0.7382     6.2921    1421.2689    3.9655    
##    6    X12           0.7466      0.7383     5.2258    1420.2544    3.9651    
## --------------------------------------------------------------------------

Berdasarkan metode backward, model terbaik adalah model dengan peubah X1, X2, X5, X7,X8,X9, X13,dan X14 dengan nilai \[R-squared=0.7383\] dan nilai AIC sebesar \[AIC= 1420.2544\]. Hasil ini sama dengan hasil dengan metode forward.

Stepwise

olsrr::ols_step_both_p(model)
## 
##                               Stepwise Selection Summary                               
## --------------------------------------------------------------------------------------
##                      Added/                   Adj.                                        
## Step    Variable    Removed     R-Square    R-Square     C(p)         AIC        RMSE     
## --------------------------------------------------------------------------------------
##    1       X7       addition       0.662       0.661    71.0320    1478.8012    4.5144    
##    2       X2       addition       0.719       0.716    19.6030    1434.6119    4.1273    
##    3      X14       addition       0.728       0.724    13.2630    1428.5726    4.0702    
##    4      X13       addition       0.735       0.731     8.1590    1423.5156    4.0217    
## --------------------------------------------------------------------------------------

Berdasarkan metode stepwise, model terbaik adalah model dengan peubah X2,X7,X13, dan X14 dengan nilai \[R-squared= 0.731\] dan nilai AIC sebesar \[AIC= 1423.5156\].

REGRESI RIDGE (glmnet)

Regresi Ridge dapat digunakan untuk mencocokkan modelregresi ketika data mengandung multikolinieritas. Regresi Ridge meminimumkan Jumlah Kuadrat Residual (JKR) prediktor dalam model. Regresi Ridge cenderung menyusutkan estimasi koefisien menuju nol.

Matriks untuk regresi ridge

library(lmridge)
matrix_X <- data.matrix(data[, -1])
matrix_Y <- matrix(data$Y)  # Assuming 'y' is a single-column vector
alpha_ridge = 0
model_Ridge <- glmnet::cv.glmnet(matrix_X,matrix_Y,alpha=alpha_ridge)
summary(model_Ridge)
##            Length Class  Mode     
## lambda     100    -none- numeric  
## cvm        100    -none- numeric  
## cvsd       100    -none- numeric  
## cvup       100    -none- numeric  
## cvlo       100    -none- numeric  
## nzero      100    -none- numeric  
## call         4    -none- call     
## name         1    -none- character
## glmnet.fit  12    elnet  list     
## lambda.min   1    -none- numeric  
## lambda.1se   1    -none- numeric  
## index        2    -none- numeric
# Hasil Regresi Ridge
print(model_Ridge)
## 
## Call:  glmnet::cv.glmnet(x = matrix_X, y = matrix_Y, alpha = alpha_ridge) 
## 
## Measure: Mean-Squared Error 
## 
##     Lambda Index Measure    SE Nonzero
## min 0.6294   100   19.01 1.948      14
## 1se 1.7515    89   20.89 2.191      14
# Memilih nilai lambda terbaik
best_lambda <- model_Ridge$lambda.min
cat("Lambda terbaik:", best_lambda, "\n")
## Lambda terbaik: 0.6294393

Model regresi Ridge menggunakan nilai alpha sebesar 0 dan nilai lambda sebesar 0.6294393 .

# Melakukan prediksi dengan model Ridge terbaik
predictions <- predict(model_Ridge, s = best_lambda, newx = matrix_X)

# koefisien Ridge
coefficients(model_Ridge, s = best_lambda)
## 15 x 1 sparse Matrix of class "dgCMatrix"
##                       s1
## (Intercept) -9.163583273
## X1           0.099932050
## X2          -0.026664933
## X3          -0.042415969
## X4           0.301842944
## X5          -0.354963579
## X6           0.089771994
## X7           0.462324392
## X8          -0.028707737
## X9           0.156789637
## X10         -0.008099125
## X11         -0.035499258
## X12          0.048119989
## X13          0.288150828
## X14         -1.526172688

Nilai koefisien dari setiap peubah bebas terhadap peubah respons (Y). Dapat terlihat bahwa Peubah X7 merupakan peubah dengan koefisien positif tertinggi sebesar \[0.462324392\] yang berpengaruh terhadap model. Sehingga dapat disimpulkan bahwa peubah ini merupakan peubah yang paling berpengaruh terhadap peubah respons jika dibandingkan peubah bebas lainnya.

Pada metode regresi Ridge mempertahankan semua peubah penjelas sehingga peubah regresinya yaitu

\[y = -9.163583273+ 0.099932050X_{1} -0.026664933X_{2}-0.0424159690X_{3}+0.301842944X_{4}-0.354963579X_{5}+0.089771994X_{6}+0.462324392X_{7}-0.028707737X_{8}+ 0.156789637X_{9}-0.008099125X_{10}-0.035499258X_{11}+0.048119989X_{12}+0.288150828X_{13}-1.526172688X_{14}\]

R-SQUARE RIDGE

# R-squared untuk model Ridge
r_squared_ridge <- 1 - sum((matrix_Y - predictions)^2) / sum((matrix_Y - mean(matrix_Y))^2)
r_squared_ridge
## [1] 0.7199482

Pada regresi Ridge, tidak terdapat peubah yang dihilangkan atau semua peubah dimasukkan dalam model. R-Square yang diperoleh dengan regresi ridge adalah

\[R-squared=0.7199482\]

REGRESI LASSO (glmnet)

alpha_Lasso = 1
model_Lasso <- glmnet::cv.glmnet(matrix_X,matrix_Y,alpha=alpha_Lasso)
model_Lasso
## 
## Call:  glmnet::cv.glmnet(x = matrix_X, y = matrix_Y, alpha = alpha_Lasso) 
## 
## Measure: Mean-Squared Error 
## 
##     Lambda Index Measure    SE Nonzero
## min 0.0499    53   17.66 1.556      11
## 1se 0.5106    28   19.20 1.295       4

Hasil Regresi Lasso

print(model_Lasso)
## 
## Call:  glmnet::cv.glmnet(x = matrix_X, y = matrix_Y, alpha = alpha_Lasso) 
## 
## Measure: Mean-Squared Error 
## 
##     Lambda Index Measure    SE Nonzero
## min 0.0499    53   17.66 1.556      11
## 1se 0.5106    28   19.20 1.295       4
# Memilih nilai lambda terbaik
best_lambda1 <- model_Lasso$lambda.min
cat("Lambda terbaik:", best_lambda1, "\n")
## Lambda terbaik: 0.04988199
# Melakukan prediksi dengan model Ridge terbaik
predictions_lasso <- predict(model_Lasso, s = best_lambda1, newx = matrix_X)

# Koefisien Ridge
coefficients(model_Lasso, s = best_lambda1)
## 15 x 1 sparse Matrix of class "dgCMatrix"
##                       s1
## (Intercept) -11.61302871
## X1            0.05392134
## X2           -0.11678382
## X3           -0.03671724
## X4            .         
## X5           -0.38244630
## X6            .         
## X7            0.80387185
## X8           -0.12185029
## X9            0.13279766
## X10           .         
## X11           0.03911244
## X12           0.09130579
## X13           0.36506536
## X14          -1.41260621

Nilai koefisien dari setiap peubah bebas terhadap peubah respons (Y). Dapat terlihat bahwa Peubah X7 merupakan peubah dengan koefisien positif tertinggi yang berpengaruh terhadap model. Sehingga dapat disimpulkan bahwa peubah ini merupakan peubah yang paling berpengaruh terhadap peubah respons jika dibandingkan peubah bebas lainnya.

Pada metode regresi Lasso hanya membuang peubah penjelas X4,X6,dan X10 sehingga peubah regresinya yaitu

\[y = -14.26747068+0.05545516X_{1}-0.15240537X_{2}-0.03017653X_{3} -0.41189328X_{5}+ 0.84016734X_{7}-0.15933006X_{8}+0.18622077X_{9}+ 0.10306501X_{11}+0.12103655X_{12}+ 0.39607291X_{13} -1.45014940X_{14}\]

# Menghitung R-squared untuk model Lasso
r_squared_lasso <- 1 - sum((matrix_Y - predictions_lasso)^2) / sum((matrix_Y - mean(matrix_Y))^2)
r_squared_lasso
## [1] 0.7460467

Melalui regresi Lasso terlihat bahwa peubah penjelas X4,X6,dan X10 yang dihapus.

Perbandingan Metode Klasik, Ridge, dan Lasso

# Membuat data untuk tabel
metode <- c("Klasik","Seleksi Peubah", "Ridge", "Lasso")
Adj.R2 <- c(0.7342 , 0.7383  ,r_squared_ridge, r_squared_lasso)

# Membuat data frame
data_tabel <- data.frame(Metode = metode, `Adj R^2` = Adj.R2)

# Mencetak tabel menggunakan knitr
library(knitr)
kable(data_tabel, format = "markdown")
Metode Adj.R.2
Klasik 0.7342000
Seleksi Peubah 0.7383000
Ridge 0.7199482
Lasso 0.7460467

Interpretasi Model Terbaik

Jika dibandingkan dari ketiga model, berdasarkan nilai R-Square, model terbaik adalah menggunakan regresi Lasso dengan \[R−squared:0.7482014\] Model regresi Lasso yang diperoleh adalah

\[y = -14.26747068+0.05545516X_{1}-0.15240537X_{2}-0.03017653X_{3} -0.41189328X_{5}+ 0.84016734X_{7}-0.15933006X_{8}+0.18622077X_{9}+ 0.10306501X_{11}+0.12103655X_{12}+ 0.39607291X_{13} -1.45014940X_{14}\]