1. Simulasi Data

Data disimulasikan dengan hubungan non-linier (kuadratik) antara x dan y.

set.seed(123)
n <- 100
x <- seq(0, 10, length.out = n)
y <- 5 + 2*x - 0.3*x^2 + rnorm(n, mean = 0, sd = 2)
data <- data.frame(x = x, y = y)

head(data)
##           x        y
## 1 0.0000000 3.879049
## 2 0.1010101 4.738604
## 3 0.2020202 8.509213
## 4 0.3030303 5.719529
## 5 0.4040404 6.017682
## 6 0.5050505 9.363708

Interpretasi:

Berdasarkan output di atas, diperoleh data hasil bangkitan sebanyak 100 data dengan variabel X berada di antara 0 hingga 10 dan variabel Y dibangkitkan dengan persamaan Y=5+2x-0,3x^2 serta nilai error yang berdistribusi normal dengan parameter rata-rata 0 dan standar deviasi

2. Pencocokan Model

Tiga model dicocokkan sebagai pembanding: linier, polinomial derajat 2, dan derajat 3.

model_linier <- lm(y ~ x, data = data)
model_poli2  <- lm(y ~ poly(x, 2, raw = TRUE), data = data)
model_poli3  <- lm(y ~ poly(x, 3, raw = TRUE), data = data)

summary(model_linier)
## 
## Call:
## lm(formula = y ~ x, data = data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -7.4317 -1.4051 -0.0466  2.2166  6.6109 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  9.88171    0.55369  17.847   <2e-16 ***
## x           -0.95028    0.09566  -9.934   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.789 on 98 degrees of freedom
## Multiple R-squared:  0.5017, Adjusted R-squared:  0.4967 
## F-statistic: 98.68 on 1 and 98 DF,  p-value: < 2.2e-16
summary(model_poli2)
## 
## Call:
## lm(formula = y ~ poly(x, 2, raw = TRUE), data = data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.8136 -1.1977 -0.0533  1.3549  4.3891 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)              5.33982    0.53777   9.929   <2e-16 ***
## poly(x, 2, raw = TRUE)1  1.80266    0.24855   7.253    1e-10 ***
## poly(x, 2, raw = TRUE)2 -0.27529    0.02405 -11.447   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.829 on 97 degrees of freedom
## Multiple R-squared:  0.788,  Adjusted R-squared:  0.7837 
## F-statistic: 180.3 on 2 and 97 DF,  p-value: < 2.2e-16
summary(model_poli3)
## 
## Call:
## lm(formula = y ~ poly(x, 3, raw = TRUE), data = data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.7866 -1.1996 -0.0497  1.3386  4.3777 
## 
## Coefficients:
##                          Estimate Std. Error t value Pr(>|t|)    
## (Intercept)              5.282801   0.708433   7.457 3.93e-11 ***
## poly(x, 3, raw = TRUE)1  1.872854   0.616613   3.037  0.00307 ** 
## poly(x, 3, raw = TRUE)2 -0.292931   0.143693  -2.039  0.04424 *  
## poly(x, 3, raw = TRUE)3  0.001176   0.009443   0.125  0.90117    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.838 on 96 degrees of freedom
## Multiple R-squared:  0.7881, Adjusted R-squared:  0.7815 
## F-statistic:   119 on 3 and 96 DF,  p-value: < 2.2e-16

Interpretasi:

Secara keseluruhan, pencocokan ketiga model menunjukkan bahwa model linier belum cukup menggambarkan pola hubungan X dan Y, terlihat dari nilai R^2 yang hanya sebesar 0,5017. Penambahan komponen kuadratik pada model polinomial derajat 2 meningkatkan R^2menjadi 0,7880 dan menurunkan Residual Standard Error menjadi 1,829. Sementara itu, penambahan komponen kubik pada model derajat 3 tidak memberikan peningkatan yang berarti karena koefisien X^3 tidak signifikan dan Residual Standard Error sedikit meningkat. Oleh karena itu, berdasarkan hasil pencocokan model, komponen kuadratik sudah mampu menggambarkan pola hubungan data dengan baik, sedangkan penambahan komponen kubik tidak memberikan kontribusi yang berarti.

3. Perbandingan Model

3.1 Uji ANOVA Bertingkat

anova(model_linier, model_poli2, model_poli3)
## Analysis of Variance Table
## 
## Model 1: y ~ x
## Model 2: y ~ poly(x, 2, raw = TRUE)
## Model 3: y ~ poly(x, 3, raw = TRUE)
##   Res.Df    RSS Df Sum of Sq        F Pr(>F)    
## 1     98 762.42                                 
## 2     97 324.33  1    438.09 129.6931 <2e-16 ***
## 3     96 324.28  1      0.05   0.0155 0.9012    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Interpretasi:

Hasil uji ANOVA bertingkat menunjukkan bahwa penambahan komponen X^2 pada model linier memberikan peningkatan kecocokan model yang signifikan, dengan p-value < 2 × 10⁻¹⁶ < 0,05. Sementara itu, penambahan komponen X^3 pada model polinomial derajat 2 tidak memberikan peningkatan yang signifikan, dengan p-value = 0,9012 > 0,05. Dengan demikian, model polinomial derajat 2 lebih sesuai digunakan karena penambahan komponen kuadratik meningkatkan kecocokan model secara signifikan, sedangkan komponen kubik tidak memberikan tambahan kecocokan yang berarti.

3.2 Perbandingan AIC

AIC(model_linier, model_poli2, model_poli3)
##              df      AIC
## model_linier  3 492.9204
## model_poli2   4 409.4469
## model_poli3   5 411.4307

Interpretasi:

Berdasarkan nilai AIC, model polinomial derajat 2 memiliki nilai AIC paling rendah, yaitu 409,4469, dibandingkan model linier sebesar 492,9204 dan model polinomial derajat 3 sebesar 411,4307. Karena model dengan nilai AIC lebih kecil menunjukkan kecocokan model yang lebih baik dengan mempertimbangkan kompleksitas model, maka model polinomial derajat 2 merupakan model yang paling sesuai berdasarkan kriteria AIC.

4. Visualisasi Perbandingan Model

data$pred_linier <- predict(model_linier)
data$pred_poli2  <- predict(model_poli2)
data$pred_poli3  <- predict(model_poli3)

ggplot(data, aes(x = x, y = y)) +
  geom_point(alpha = 0.5) +
  geom_line(aes(y = pred_linier, color = "Linier"), linewidth = 1) +
  geom_line(aes(y = pred_poli2, color = "Polinomial derajat 2"), linewidth = 1) +
  geom_line(aes(y = pred_poli3, color = "Polinomial derajat 3"), linewidth = 1) +
  labs(title = "Perbandingan Regresi Linier vs Polinomial",
       x = "X", y = "Y", color = "Model") +
  theme_minimal()

#### Interpretasi: Berdasarkan grafik perbandingan regresi linier dan polinomial, terlihat bahwa model linier menghasilkan garis lurus menurun sehingga kurang mampu mengikuti pola data yang melengkung. Sementara itu, model polinomial derajat 2 dan derajat 3 mampu mengikuti pola hubungan X dan Y yang melengkung dengan lebih baik. Kurva polinomial derajat 2 dan derajat 3 juga terlihat sangat berdekatan, menunjukkan bahwa penambahan komponen kubik tidak memberikan perubahan pola yang berarti. Hal ini sejalan dengan hasil ANOVA dan AIC sebelumnya, sehingga model polinomial derajat 2 sudah cukup menggambarkan pola data.

5. Evaluasi dengan RMSE

rmse <- function(actual, predicted) sqrt(mean((actual - predicted)^2))

data.frame(
  Model = c("Linier", "Polinomial derajat 2", "Polinomial derajat 3"),
  RMSE  = c(rmse(data$y, data$pred_linier),
            rmse(data$y, data$pred_poli2),
            rmse(data$y, data$pred_poli3))
)
##                  Model     RMSE
## 1               Linier 2.761195
## 2 Polinomial derajat 2 1.800917
## 3 Polinomial derajat 3 1.800772

Interpretasi:

Berdasarkan hasil evaluasi RMSE, model linier memiliki nilai RMSE sebesar 2,761195, sedangkan model polinomial derajat 2 sebesar 1,800917 dan derajat 3 sebesar 1,800772. Nilai RMSE yang lebih kecil menunjukkan kesalahan prediksi yang lebih rendah. Dengan demikian, model polinomial derajat 2 dan derajat 3 memiliki kinerja prediksi yang jauh lebih baik dibandingkan model linier. Perbedaan RMSE antara model derajat 2 dan derajat 3 juga sangat kecil, sehingga penambahan komponen kubik tidak memberikan peningkatan prediksi yang berarti.

6. Pemilihan Derajat Optimal via Cross-Validation

Skema 10-fold cross-validation diulang 5 kali untuk menguji derajat polinomial 1 sampai 6.

names(data)
## [1] "x"           "y"           "pred_linier" "pred_poli2"  "pred_poli3"
library(caret)

set.seed(42)

kontrol <- trainControl(
  method = "repeatedcv",
  number = 10,
  repeats = 5
)

derajat_max <- 6

hasil_cv <- data.frame(
  derajat = integer(),
  RMSE = numeric()
)

for (d in 1:derajat_max) {
  
  data_cv <- data
  
  # Membuat variabel polynomial
  if (d == 1) {
    data_cv$x_poly <- data_cv$x
  } else {
    poly_x <- poly(data_cv$x, degree = d, raw = TRUE)
    
    data_cv$x_poly <- poly_x[, 1]
    
    for (j in 2:d) {
      data_cv[[paste0("x_poly", j)]] <- poly_x[, j]
    }
  }
  
  # Formula sesuai derajat
  if (d == 1) {
    formula_model <- y ~ x_poly
  } else {
    variabel <- paste0("x_poly", 2:d)
    formula_model <- as.formula(
      paste("y ~ x_poly +", paste(variabel, collapse = " + "))
    )
  }
  
  model_cv <- train(
    formula_model,
    data = data_cv,
    method = "lm",
    trControl = kontrol
  )
  
  hasil_cv <- rbind(
    hasil_cv,
    data.frame(
      derajat = d,
      RMSE = min(model_cv$results$RMSE)
    )
  )
}

hasil_cv
##   derajat     RMSE
## 1       1 2.757805
## 2       2 1.825068
## 3       3 1.859968
## 4       4 1.877628
## 5       5 1.871701
## 6       6 1.880719
derajat_optimal <- hasil_cv$derajat[which.min(hasil_cv$RMSE)]
cat("Derajat optimal berdasarkan CV:", derajat_optimal, "\n")
## Derajat optimal berdasarkan CV: 2
ggplot(hasil_cv, aes(x = derajat, y = RMSE)) +
  geom_line(
    color = "steelblue",
    linewidth = 1
  ) +
  geom_point(
    size = 3,
    color = "steelblue"
  ) +
  geom_vline(
    xintercept = derajat_optimal,
    linetype = "dashed",
    color = "firebrick"
  ) +
  scale_x_continuous(
    breaks = 1:derajat_max
  ) +
  labs(
    title = "Pemilihan Derajat Optimal via 10-Fold Cross-Validation",
    x = "Derajat Polinomial",
    y = "RMSE Rata-rata (Validasi Silang)"
  ) +
  theme_minimal()

#### Interpretasi: Berdasarkan grafik 10-Fold Cross-Validation, nilai RMSE validasi paling tinggi terdapat pada polinomial derajat 1, yaitu sekitar 2,76. Nilai RMSE kemudian turun drastis pada derajat 2 menjadi sekitar 1,80 dan setelah itu cenderung stabil hingga derajat 6. Dengan demikian, derajat 2 dipilih sebagai derajat polinomial optimal, karena memberikan RMSE validasi paling rendah. Penambahan derajat di atas 2 tidak memberikan peningkatan kinerja yang berarti.

7. Model Final

model_final <- lm(y ~ poly(x, derajat_optimal, raw = TRUE), data = data)
summary(model_final)
## 
## Call:
## lm(formula = y ~ poly(x, derajat_optimal, raw = TRUE), data = data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.8136 -1.1977 -0.0533  1.3549  4.3891 
## 
## Coefficients:
##                                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                            5.33982    0.53777   9.929   <2e-16 ***
## poly(x, derajat_optimal, raw = TRUE)1  1.80266    0.24855   7.253    1e-10 ***
## poly(x, derajat_optimal, raw = TRUE)2 -0.27529    0.02405 -11.447   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.829 on 97 degrees of freedom
## Multiple R-squared:  0.788,  Adjusted R-squared:  0.7837 
## F-statistic: 180.3 on 2 and 97 DF,  p-value: < 2.2e-16

Interpretasi:

Berdasarkan hasil model final, diketahui bahwa model yang terbentuk signifikan dengan p-value < 2 × 10⁻¹⁶. Koefisien X dan X^2 juga signifikan. Nilai R^2 = 0,788 menunjukkan bahwa 78,8% variasi Y dapat dijelaskan oleh model, sedangkan 21,2% sisanya dijelaskan oleh faktor lain atau galat. Jadi, model final yang digunakan adalah regresi polinomial derajat 2.

Catatan

  • poly(x, derajat, raw = TRUE) digunakan agar koefisien dapat diinterpretasi langsung sebagai \(\beta_1 x + \beta_2 x^2 + \dots\)
  • Untuk derajat tinggi (\(\geq 4\)), sebaiknya x distandardisasi terlebih dahulu (scale(x)) agar model lebih stabil secara numerik.
  • Pertimbangkan aturan one-standard-error rule: pilih derajat terkecil dengan RMSE \(\leq\) (RMSE minimum + 1 SD) untuk model yang lebih sederhana namun performanya setara.