Pertemuan 8 - Kriteria Pemilihan Model Terbaik

A. Metode Seleksi Peubah Penjelas

Metode seleksi peubah penjelas berguna untuk menemukan variabel model terbaik dan signifikan terhadap model. Perlu diingat bahwa model regresi yang disukai adalah model dengan unsur peubah penjelas paling sederhana. backward, forward, dan stepwise selection adalah metode seleksi peubah yang akan dibahas. Sebelum menggunakan ketiga teknik tersebut, kita perlu memahami uji F-parsial dan uji F-sekuensial. Uji F-parsial digunakan untuk backward selection, sedangkan uji F-sekuensial digunakan untuk forward selection.

Uji F-Sekuensial

Uji F-sekuensial digunakan untuk melihat pengaruh r peubah penjelas yang ditambahkan ke dalam model. Contoh penerapan uji F-sekuensial menggunakan data Data Supervisor.csv adalah sebagai berikut:

Data Supervisor terdiri dari 4 peubah penjelas dan 30 observasi

Supervisor <- read.csv("C:/Users/Putera Ramadhan/OneDrive/Documents/IPB University/Semester 4/STA1231 - Analisis Regresi/UAS/Praktikum/Minggu 9/Data Supervisor.csv")
Supervisor <- na.omit(Supervisor)
sum(is.na(Supervisor))
## [1] 0
str(Supervisor)
## 'data.frame':    30 obs. of  5 variables:
##  $ Y : int  43 63 71 61 81 43 58 71 72 67 ...
##  $ X1: int  51 64 70 63 78 55 67 75 82 61 ...
##  $ X2: int  30 51 68 45 56 49 42 50 72 45 ...
##  $ X3: int  39 54 69 47 66 44 56 55 67 47 ...
##  $ X4: int  61 63 76 54 71 54 66 70 71 62 ...
head(Supervisor)
##    Y X1 X2 X3 X4
## 1 43 51 30 39 61
## 2 63 64 51 54 63
## 3 71 70 68 69 76
## 4 61 63 45 47 54
## 5 81 78 56 66 71
## 6 43 55 49 44 54

Penentuan peubah yang akan dimasukkan ke dalam model didasarkan pada nilai korelasi antara peubah penjelas terhadap peubah respon. Peubah penjelas dengan korelasi paling tinggi dimasukkan terlebih dahulu.

## Korelasi all X vs Y
cor(Supervisor[2:5], Supervisor$Y)
##         [,1]
## X1 0.8254176
## X2 0.4261169
## X3 0.6236782
## X4 0.5901390

Peubah X1 dimasukkan ke model 1 (model tanpa peubah penjelas) terlebih dahulu sehingga menjadi model 2

## Model 1: Y = b0
model1 <- lm(Y ~ 1, data = Supervisor)
summary(model1)
## 
## Call:
## lm(formula = Y ~ 1, data = Supervisor)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -24.6333  -5.8833   0.8667   7.1167  20.3667 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   64.633      2.222   29.08   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 12.17 on 29 degrees of freedom
## Model 2: Y = bo + b1X1 (memasukkan X1 ke model)
model2 <- lm(Y ~ X1, data = Supervisor)
summary(model2)
## 
## Call:
## lm(formula = Y ~ X1, data = Supervisor)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -12.8799  -5.9905   0.1783   6.2978   9.6294 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 14.37632    6.61999   2.172   0.0385 *  
## X1           0.75461    0.09753   7.737 1.99e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.993 on 28 degrees of freedom
## Multiple R-squared:  0.6813, Adjusted R-squared:  0.6699 
## F-statistic: 59.86 on 1 and 28 DF,  p-value: 1.988e-08

Peubah X1 berpengaruh nyata terhadap Y pada taraf nyata 5%.

Selanjutnya, memilih peubah penjelas lain yang akan dimasukkan ke dalam model berdasarkan nilai korelasi antara peubah respon dengan peubah penjelas sisanya.

cor(Supervisor[3:5], Supervisor$Y)
##         [,1]
## X2 0.4261169
## X3 0.6236782
## X4 0.5901390

Peubah X3 akan dimasukkan ke dalam model 2 sehingga menjadi model 3

## Model 3: Y = b0 + b1X1 + b3X3
model3 <- lm(Y ~ X1 + X3, data = Supervisor)
summary(model3)
## 
## Call:
## lm(formula = Y ~ X1 + X3, data = Supervisor)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -11.5568  -5.7331   0.6701   6.5341  10.3610 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   9.8709     7.0612   1.398    0.174    
## X1            0.6435     0.1185   5.432 9.57e-06 ***
## X3            0.2112     0.1344   1.571    0.128    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.817 on 27 degrees of freedom
## Multiple R-squared:  0.708,  Adjusted R-squared:  0.6864 
## F-statistic: 32.74 on 2 and 27 DF,  p-value: 6.058e-08

Hasil uji t di atas menunjukkan bahwa peubah X3 tidak berpengaruh nyata terhadap Y pada taraf nyata 5%.

Untuk menentukan apakah X3 memberikan pengaruh ketika dimasukkan ke dalam model, digunakan uji F-sekuensial seperti berikut

## JK(X1) dari model 2
aov(model2)
## Call:
##    aov(formula = model2)
## 
## Terms:
##                       X1 Residuals
## Sum of Squares  2927.584  1369.382
## Deg. of Freedom        1        28
## 
## Residual standard error: 6.993319
## Estimated effects may be unbalanced
JK.b1 <- 2927.584
## JK(X1, X3) dari model 3
aov(model3)
## Call:
##    aov(formula = model3)
## 
## Terms:
##                        X1        X3 Residuals
## Sum of Squares  2927.5843  114.7334 1254.6490
## Deg. of Freedom         1         1        27
## 
## Residual standard error: 6.816779
## Estimated effects may be unbalanced
JK.b1b2 <- 2927.5843 + 114.7334
## Uji F-sekuensial model 2 dan 3
# H0: b2 = 0 vs H1: b2 != 0
# Menghitung F-hitung
JK.b2_b1 <- JK.b1b2 - JK.b1 # JK(b2|b1)
JK.b2_b1
## [1] 114.7337
KT.b2_b1 <- JK.b2_b1 / 1 # KT(b2|b1)
KT.res <- 1254.6490 # KT(b1, b2)

F.hit <- KT.b2_b1 / KT.res
F.tabel <- qf(0.95, 1, (nrow(Supervisor) - 3))
cbind(F.hit, F.tabel)
##           F.hit  F.tabel
## [1,] 0.09144685 4.210008

F hitung dari output diatas lebih kecil dari pada F tabel sehingga tak tolak H0. Artinya, Peubah X3 tidak perlu ditambahkan ke dalam model 2.

All Possible Regression

Pada metode ini, akan diuji semua kemungkinan model subset regresi dari semua kandidat peubah penjelas. Banyaknya model yang mungkin yaitu \(2^k − 1\), di mana \(k\) = banyaknya peubah penjelas. Contoh penerapan metode ini masih akan menggunakan data Data Supervisor.csv dan akan digunakan untuk forward, backward, stepwise, best subsets regression, serta PRESS.

str(Supervisor)
## 'data.frame':    30 obs. of  5 variables:
##  $ Y : int  43 63 71 61 81 43 58 71 72 67 ...
##  $ X1: int  51 64 70 63 78 55 67 75 82 61 ...
##  $ X2: int  30 51 68 45 56 49 42 50 72 45 ...
##  $ X3: int  39 54 69 47 66 44 56 55 67 47 ...
##  $ X4: int  61 63 76 54 71 54 66 70 71 62 ...
head(Supervisor)
##    Y X1 X2 X3 X4
## 1 43 51 30 39 61
## 2 63 64 51 54 63
## 3 71 70 68 69 76
## 4 61 63 45 47 54
## 5 81 78 56 66 71
## 6 43 55 49 44 54
library(olsrr)
## Warning: package 'olsrr' was built under R version 4.4.3
## 
## Attaching package: 'olsrr'
## The following object is masked from 'package:datasets':
## 
##     rivers
# Jumlah kombinasi model = 2^k - 1, k = peubah penjelas
full.model <- lm(Y ~ ., data = Supervisor)
ols_step_all_possible(full.model)
##    Index N  Predictors  R-Square Adj. R-Square Mallow's Cp
## 1      1 1          X1 0.6813142     0.6699325    1.976859
## 3      2 1          X3 0.3889745     0.3671521   27.640838
## 4      3 1          X4 0.3482640     0.3249877   31.214733
## 2      4 1          X2 0.1815756     0.1523461   45.848013
## 6      5 2       X1 X3 0.7080152     0.6863867    1.632823
## 7      6 2       X1 X4 0.6838979     0.6604830    3.750034
## 5      7 2       X1 X2 0.6830639     0.6595872    3.823252
## 10     8 2       X3 X4 0.4506703     0.4099792   24.224669
## 8      9 2       X2 X3 0.4075138     0.3636260   28.013299
## 9     10 2       X2 X4 0.3815018     0.3356872   30.296846
## 11    11 3    X1 X2 X3 0.7150044     0.6821203    3.019249
## 13    12 3    X1 X3 X4 0.7082916     0.6746329    3.608557
## 12    13 3    X1 X2 X4 0.6862100     0.6500034    5.547063
## 14    14 3    X2 X3 X4 0.4587139     0.3962578   25.518541
## 15    15 4 X1 X2 X3 X4 0.7152237     0.6696595    5.000000

Berdasarkan output di atas, model ke 6, 7, 11, dan 12 menghasilkan nilai \(R^2\) dan \(Adj.R^2\) sangat besar serta Mallow’s Cp yang kecil mendekati k + 1. Model 6 terdiri dari dua peubah (X1 dan X4), model 7 terdiri dari dua peubah (X1 dan X2), model 11 terdiri dari tiga peubah(X1, X2, dan X3), dan model 12 terdiri dari tiga peubah (X1, X3, dan X4). Perubahan \(R^2\) dan \(Adj.R^2\) yang paling besar terjadi ketika menambahkan satu X ke model yang terdiri dari 3 peubah bebas, yaitu terjadi pada model 11 dan 12. Model 11 memperbaiki model 6 dan menjadi model terbaik.

Forward Selection

Forward selection diawali dengan model tanpa peubah penjelas, kemudian memasukkan peubah penjelas satu per satu ke dalam model berdasarkan nilai korelasi seperti contoh uji F-sekuensial di atas. Setiap peubah yang masuk ke dalam model tidak akan dikeluarkan. Proses seleksi akan berhenti jika tidak ada lagi peubah penjelas yang signifikan di antara yang tersisa.

full <- lm(Y ~ ., data = Supervisor)
null <- lm(Y ~ 1, data = Supervisor)
## Forward Selection
step(null, scope = list(lower = null, upper = full), direction = "both")
## Start:  AIC=150.93
## Y ~ 1
## 
##        Df Sum of Sq    RSS    AIC
## + X1    1   2927.58 1369.4 118.63
## + X3    1   1671.41 2625.6 138.16
## + X4    1   1496.48 2800.5 140.09
## + X2    1    780.22 3516.7 146.92
## <none>              4297.0 150.93
## 
## Step:  AIC=118.63
## Y ~ X1
## 
##        Df Sum of Sq    RSS    AIC
## + X3    1    114.73 1254.6 118.00
## <none>              1369.4 118.63
## + X4    1     11.10 1358.3 120.38
## + X2    1      7.52 1361.9 120.46
## - X1    1   2927.58 4297.0 150.93
## 
## Step:  AIC=118
## Y ~ X1 + X3
## 
##        Df Sum of Sq    RSS    AIC
## <none>              1254.7 118.00
## - X3    1    114.73 1369.4 118.63
## + X2    1     30.03 1224.6 119.28
## + X4    1      1.19 1253.5 119.97
## - X1    1   1370.91 2625.6 138.16
## 
## Call:
## lm(formula = Y ~ X1 + X3, data = Supervisor)
## 
## Coefficients:
## (Intercept)           X1           X3  
##      9.8709       0.6435       0.2112
fw <- ols_step_forward_p(full.model)
fw
## 
## 
##                              Stepwise Summary                              
## -------------------------------------------------------------------------
## Step    Variable        AIC        SBC       SBIC        R2       Adj. R2 
## -------------------------------------------------------------------------
##  0      Base Model    238.070    240.873    150.751    0.00000    0.00000 
##  1      X1            205.764    209.967    120.906    0.68131    0.66993 
##  2      X3            205.139    210.744    120.967    0.70802    0.68639 
## -------------------------------------------------------------------------
## 
## Final Model Output 
## ------------------
## 
##                          Model Summary                          
## ---------------------------------------------------------------
## R                       0.841       RMSE                 6.467 
## R-Squared               0.708       MSE                 41.822 
## Adj. R-Squared          0.686       Coef. Var           10.547 
## Pred R-Squared          0.640       AIC                205.139 
## MAE                     5.612       SBC                210.744 
## ---------------------------------------------------------------
##  RMSE: Root Mean Square Error 
##  MSE: Mean Square Error 
##  MAE: Mean Absolute Error 
##  AIC: Akaike Information Criteria 
##  SBC: Schwarz Bayesian Criteria 
## 
##                                ANOVA                                 
## --------------------------------------------------------------------
##                 Sum of                                              
##                Squares        DF    Mean Square      F         Sig. 
## --------------------------------------------------------------------
## Regression    3042.318         2       1521.159    32.735    0.0000 
## Residual      1254.649        27         46.468                     
## Total         4296.967        29                                    
## --------------------------------------------------------------------
## 
##                                  Parameter Estimates                                   
## --------------------------------------------------------------------------------------
##       model     Beta    Std. Error    Std. Beta      t       Sig      lower     upper 
## --------------------------------------------------------------------------------------
## (Intercept)    9.871         7.061                 1.398    0.174    -4.618    24.359 
##          X1    0.644         0.118        0.704    5.432    0.000     0.400     0.887 
##          X3    0.211         0.134        0.204    1.571    0.128    -0.065     0.487 
## --------------------------------------------------------------------------------------

Backward Selection

Backward selection diawali dengan model penuh (semua peubah penjelas ada di dalam model), kemudian mengeluarkan peubah penjelas satu per satu berdasarkan signifikansi hasil uji F-parsial. Setiap peubah penjelas yang dikeluarkan dari model tidak akan dimasukkan kembali. Proses seleksi akan berhenti jika tidak ada lagi peubah penjelas yang tidak signifikan di antara yang tersisa.

## Backward selection
step(full, direction = "backward")
## Start:  AIC=121.25
## Y ~ X1 + X2 + X3 + X4
## 
##        Df Sum of Sq    RSS    AIC
## - X4    1      0.94 1224.6 119.28
## - X2    1     29.79 1253.5 119.97
## <none>              1223.7 121.25
## - X3    1    124.67 1348.3 122.16
## - X1    1   1102.21 2325.9 138.52
## 
## Step:  AIC=119.28
## Y ~ X1 + X2 + X3
## 
##        Df Sum of Sq    RSS    AIC
## - X2    1     30.03 1254.7 118.00
## <none>              1224.6 119.28
## - X3    1    137.25 1361.9 120.46
## - X1    1   1321.28 2545.9 139.23
## 
## Step:  AIC=118
## Y ~ X1 + X3
## 
##        Df Sum of Sq    RSS    AIC
## <none>              1254.7 118.00
## - X3    1    114.73 1369.4 118.63
## - X1    1   1370.91 2625.6 138.16
## 
## Call:
## lm(formula = Y ~ X1 + X3, data = Supervisor)
## 
## Coefficients:
## (Intercept)           X1           X3  
##      9.8709       0.6435       0.2112
bw <- ols_step_backward_p(full.model)
bw
## 
## 
##                              Stepwise Summary                              
## -------------------------------------------------------------------------
## Step    Variable        AIC        SBC       SBIC        R2       Adj. R2 
## -------------------------------------------------------------------------
##  0      Full Model    208.389    216.796    125.172    0.71522    0.66966 
##  1      X4            206.412    213.418    122.789    0.71500    0.68212 
##  2      X2            205.139    210.744    120.967    0.70802    0.68639 
## -------------------------------------------------------------------------
## 
## Final Model Output 
## ------------------
## 
##                          Model Summary                          
## ---------------------------------------------------------------
## R                       0.841       RMSE                 6.467 
## R-Squared               0.708       MSE                 41.822 
## Adj. R-Squared          0.686       Coef. Var           10.547 
## Pred R-Squared          0.640       AIC                205.139 
## MAE                     5.612       SBC                210.744 
## ---------------------------------------------------------------
##  RMSE: Root Mean Square Error 
##  MSE: Mean Square Error 
##  MAE: Mean Absolute Error 
##  AIC: Akaike Information Criteria 
##  SBC: Schwarz Bayesian Criteria 
## 
##                                ANOVA                                 
## --------------------------------------------------------------------
##                 Sum of                                              
##                Squares        DF    Mean Square      F         Sig. 
## --------------------------------------------------------------------
## Regression    3042.318         2       1521.159    32.735    0.0000 
## Residual      1254.649        27         46.468                     
## Total         4296.967        29                                    
## --------------------------------------------------------------------
## 
##                                  Parameter Estimates                                   
## --------------------------------------------------------------------------------------
##       model     Beta    Std. Error    Std. Beta      t       Sig      lower     upper 
## --------------------------------------------------------------------------------------
## (Intercept)    9.871         7.061                 1.398    0.174    -4.618    24.359 
##          X1    0.644         0.118        0.704    5.432    0.000     0.400     0.887 
##          X3    0.211         0.134        0.204    1.571    0.128    -0.065     0.487 
## --------------------------------------------------------------------------------------

Stepwise Regression

Stepwise Regression digunakan untuk mencari model terbaik dan dapat mengatasi multikolinearitas. Metode tersebut mengombinasikan forward dan backward selection. Proses seleksi diawali dari model tanpa peubah penjelas, lalu memasukkan peubah X berdasarkan nilai korelasi tertinggi (langkah selanjutnya seperti forward selection). Peubah penjelas yang sudah masuk ke dalam model, dapat dikeluarkan menggunakan backward selection. Pada metode ini, digunakan dua taraf nyata, yaitu taraf nyata masuk dan keluar yang besar nilainya berdasarkan bidang yang diteliti.

## Stepwise
step(null, scope = list(lower = null, upper = full), direction = "both")
## Start:  AIC=150.93
## Y ~ 1
## 
##        Df Sum of Sq    RSS    AIC
## + X1    1   2927.58 1369.4 118.63
## + X3    1   1671.41 2625.6 138.16
## + X4    1   1496.48 2800.5 140.09
## + X2    1    780.22 3516.7 146.92
## <none>              4297.0 150.93
## 
## Step:  AIC=118.63
## Y ~ X1
## 
##        Df Sum of Sq    RSS    AIC
## + X3    1    114.73 1254.6 118.00
## <none>              1369.4 118.63
## + X4    1     11.10 1358.3 120.38
## + X2    1      7.52 1361.9 120.46
## - X1    1   2927.58 4297.0 150.93
## 
## Step:  AIC=118
## Y ~ X1 + X3
## 
##        Df Sum of Sq    RSS    AIC
## <none>              1254.7 118.00
## - X3    1    114.73 1369.4 118.63
## + X2    1     30.03 1224.6 119.28
## + X4    1      1.19 1253.5 119.97
## - X1    1   1370.91 2625.6 138.16
## 
## Call:
## lm(formula = Y ~ X1 + X3, data = Supervisor)
## 
## Coefficients:
## (Intercept)           X1           X3  
##      9.8709       0.6435       0.2112

B. Kriteria Model Terbaik

Best Subset Regression

Memilih model terbaik dari beberapa model subset regresi terbaik (dari hasil all possible regression) berdasarkan beberapa kriteria baku seperti memiliki nilai \(R^2\), \(Adj.R^2\), dan \(Pred.R^2\) terbesar, Mallow’s Cp terkecil, serta AIC terkecil.

bs.hc <- ols_step_best_subset(full.model)
bs.hc
##  Best Subsets Regression  
## --------------------------
## Model Index    Predictors
## --------------------------
##      1         X1          
##      2         X1 X3       
##      3         X1 X2 X3    
##      4         X1 X2 X3 X4 
## --------------------------
## 
##                                                     Subsets Regression Summary                                                     
## -----------------------------------------------------------------------------------------------------------------------------------
##                        Adj.        Pred                                                                                             
## Model    R-Square    R-Square    R-Square     C(p)       AIC         SBIC        SBC         MSEP         FPE       HSP       APC  
## -----------------------------------------------------------------------------------------------------------------------------------
##   1        0.6813      0.6699      0.6379    1.9769    205.7638    120.9063    209.9674    1467.4370    52.1669    1.8114    0.3642 
##   2        0.7080      0.6864      0.6402    1.6328    205.1387    120.9666    210.7435    1396.1991    51.1153    1.7872    0.3569 
##   3        0.7150      0.6821       0.629    3.0192    206.4119    122.7889    213.4179    1417.2894    53.3807    1.8840    0.3727 
##   4        0.7152      0.6697      0.6039    5.0000    208.3888    125.1725    216.7960    1475.2072    57.1048    2.0395    0.3987 
## -----------------------------------------------------------------------------------------------------------------------------------
## AIC: Akaike Information Criteria 
##  SBIC: Sawa's Bayesian Information Criteria 
##  SBC: Schwarz Bayesian Criteria 
##  MSEP: Estimated error of prediction, assuming multivariate normality 
##  FPE: Final Prediction Error 
##  HSP: Hocking's Sp 
##  APC: Amemiya Prediction Criteria

Model terbaik yaitu model 2 dengan peubah X1 dan X3 di dalam model.

Validasi Model

PRESS (Prediction Sum of Square)

PRESS merupakan kombinasi dari all possible regression, analisis sisaan, dan teknik validasi yang digunakan untuk mengukur validitas model. Model yang baik memiliki nilai PRESS yang kecil.

## PRESS
PRESS1 <- ols_press(lm(Y ~ X1, data = Supervisor))
PRESS2 <- ols_press(lm(Y ~ X1 + X3, data = Supervisor))
PRESS3 <- ols_press(lm(Y ~ X1 + X2 + X3, data = Supervisor))
PRESS4 <- ols_press(full.model)
cbind(PRESS1, PRESS2, PRESS3, PRESS4)
##        PRESS1   PRESS2   PRESS3   PRESS4
## [1,] 1556.065 1545.966 1594.106 1701.952

Model ke-2 dengan peubah X1 dan X3 dapat dikatakan sebagai model terbaik karena memiliki PRESS paling kecil.