1 Pendahuluan

Uji Ramsey RESET digunakan untuk mendeteksi kesalahan spesifikasi model regresi linier. Dalam simulasi ini, kita akan mengamati power uji Ramsey RESET berdasarkan kombinasi dari:

Kita juga menguji seberapa sering uji Ramsey mendeteksi ketidaksesuaian model jika hanya satu variabel digunakan dari dua prediktor.

2 Persiapan Packages

library(lmtest)
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(sandwich)
library(tibble)
library(purrr)
library(ggplot2)
library(plotly)
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
library(tidyr)

3 Parameter Simulasi

n_values <- c(20, 40, 60, 80, 100)
sigma_values <- c(1, 5, 10, 15, 20)
r2_values <- c(0.1, 0.3, 0.5, 0.7, 0.9)

4 Fungsi Simulasi

simulasi_once <- function(n, sigma, r2) {
  x1 <- rnorm(n, mean = 7, sd = 1)
  x2 <- rnorm(n, mean = 3, sd = 1)
  beta <- sqrt(r2 / (2 * (1 - r2)))
  y <- beta * x1 + beta * x2 + rnorm(n, sd = sigma)
  
  model1 <- lm(y ~ x1)
  model2 <- lm(y ~ x2)
  reset1 <- resettest(model1, power = 2:3, type = "regressor")
  reset2 <- resettest(model2, power = 2:3, type = "regressor")
  
  return(c(x1_sig = reset1$p.value < 0.05, x2_sig = reset2$p.value < 0.05))
}

power_test <- function(n, sigma, r2, reps = 500) {
  hasil <- replicate(reps, simulasi_once(n, sigma, r2))
  power_x1 <- mean(hasil[1, ])
  power_x2 <- mean(hasil[2, ])
  return(c(power_x1 = power_x1, power_x2 = power_x2))
}

5 Simulasi 125 Kombinasi

kombinasi <- expand.grid(n = n_values, sigma = sigma_values, r2 = r2_values)

set.seed(123)
hasil_simulasi <- pmap_dfr(kombinasi, function(n, sigma, r2) {
  power <- power_test(n, sigma, r2)
  tibble(n = n, sigma = sigma, r2 = r2,
         power_x1 = power[1], power_x2 = power[2])
})

print(hasil_simulasi,n=125)
## # A tibble: 125 × 5
##         n sigma    r2 power_x1 power_x2
##     <dbl> <dbl> <dbl>    <dbl>    <dbl>
##   1    20     1   0.1    0.034    0.054
##   2    40     1   0.1    0.072    0.05 
##   3    60     1   0.1    0.044    0.036
##   4    80     1   0.1    0.05     0.078
##   5   100     1   0.1    0.056    0.048
##   6    20     5   0.1    0.036    0.038
##   7    40     5   0.1    0.048    0.07 
##   8    60     5   0.1    0.04     0.046
##   9    80     5   0.1    0.05     0.05 
##  10   100     5   0.1    0.048    0.048
##  11    20    10   0.1    0.036    0.056
##  12    40    10   0.1    0.056    0.044
##  13    60    10   0.1    0.06     0.064
##  14    80    10   0.1    0.066    0.066
##  15   100    10   0.1    0.042    0.044
##  16    20    15   0.1    0.038    0.05 
##  17    40    15   0.1    0.044    0.048
##  18    60    15   0.1    0.06     0.036
##  19    80    15   0.1    0.052    0.046
##  20   100    15   0.1    0.038    0.046
##  21    20    20   0.1    0.068    0.044
##  22    40    20   0.1    0.05     0.036
##  23    60    20   0.1    0.064    0.05 
##  24    80    20   0.1    0.042    0.046
##  25   100    20   0.1    0.052    0.046
##  26    20     1   0.3    0.054    0.038
##  27    40     1   0.3    0.05     0.048
##  28    60     1   0.3    0.048    0.042
##  29    80     1   0.3    0.038    0.048
##  30   100     1   0.3    0.06     0.062
##  31    20     5   0.3    0.044    0.046
##  32    40     5   0.3    0.032    0.038
##  33    60     5   0.3    0.056    0.036
##  34    80     5   0.3    0.046    0.046
##  35   100     5   0.3    0.052    0.058
##  36    20    10   0.3    0.05     0.046
##  37    40    10   0.3    0.05     0.046
##  38    60    10   0.3    0.052    0.056
##  39    80    10   0.3    0.038    0.066
##  40   100    10   0.3    0.048    0.034
##  41    20    15   0.3    0.066    0.052
##  42    40    15   0.3    0.046    0.064
##  43    60    15   0.3    0.038    0.054
##  44    80    15   0.3    0.054    0.064
##  45   100    15   0.3    0.068    0.052
##  46    20    20   0.3    0.052    0.064
##  47    40    20   0.3    0.052    0.042
##  48    60    20   0.3    0.06     0.06 
##  49    80    20   0.3    0.068    0.062
##  50   100    20   0.3    0.038    0.062
##  51    20     1   0.5    0.066    0.05 
##  52    40     1   0.5    0.048    0.068
##  53    60     1   0.5    0.064    0.052
##  54    80     1   0.5    0.054    0.052
##  55   100     1   0.5    0.042    0.06 
##  56    20     5   0.5    0.054    0.046
##  57    40     5   0.5    0.048    0.052
##  58    60     5   0.5    0.064    0.052
##  59    80     5   0.5    0.044    0.052
##  60   100     5   0.5    0.052    0.044
##  61    20    10   0.5    0.036    0.052
##  62    40    10   0.5    0.048    0.066
##  63    60    10   0.5    0.05     0.046
##  64    80    10   0.5    0.038    0.03 
##  65   100    10   0.5    0.054    0.044
##  66    20    15   0.5    0.04     0.044
##  67    40    15   0.5    0.054    0.048
##  68    60    15   0.5    0.048    0.07 
##  69    80    15   0.5    0.042    0.054
##  70   100    15   0.5    0.056    0.058
##  71    20    20   0.5    0.062    0.038
##  72    40    20   0.5    0.06     0.036
##  73    60    20   0.5    0.034    0.064
##  74    80    20   0.5    0.042    0.06 
##  75   100    20   0.5    0.066    0.052
##  76    20     1   0.7    0.062    0.04 
##  77    40     1   0.7    0.058    0.048
##  78    60     1   0.7    0.046    0.06 
##  79    80     1   0.7    0.046    0.05 
##  80   100     1   0.7    0.052    0.06 
##  81    20     5   0.7    0.046    0.056
##  82    40     5   0.7    0.036    0.044
##  83    60     5   0.7    0.054    0.074
##  84    80     5   0.7    0.048    0.048
##  85   100     5   0.7    0.05     0.06 
##  86    20    10   0.7    0.054    0.078
##  87    40    10   0.7    0.044    0.052
##  88    60    10   0.7    0.056    0.066
##  89    80    10   0.7    0.04     0.046
##  90   100    10   0.7    0.058    0.056
##  91    20    15   0.7    0.044    0.052
##  92    40    15   0.7    0.054    0.044
##  93    60    15   0.7    0.054    0.054
##  94    80    15   0.7    0.056    0.036
##  95   100    15   0.7    0.054    0.052
##  96    20    20   0.7    0.044    0.054
##  97    40    20   0.7    0.05     0.06 
##  98    60    20   0.7    0.042    0.052
##  99    80    20   0.7    0.05     0.05 
## 100   100    20   0.7    0.064    0.044
## 101    20     1   0.9    0.042    0.058
## 102    40     1   0.9    0.06     0.064
## 103    60     1   0.9    0.056    0.038
## 104    80     1   0.9    0.06     0.048
## 105   100     1   0.9    0.066    0.062
## 106    20     5   0.9    0.054    0.038
## 107    40     5   0.9    0.044    0.04 
## 108    60     5   0.9    0.05     0.064
## 109    80     5   0.9    0.06     0.04 
## 110   100     5   0.9    0.04     0.034
## 111    20    10   0.9    0.048    0.042
## 112    40    10   0.9    0.054    0.054
## 113    60    10   0.9    0.062    0.052
## 114    80    10   0.9    0.052    0.056
## 115   100    10   0.9    0.046    0.044
## 116    20    15   0.9    0.034    0.046
## 117    40    15   0.9    0.066    0.04 
## 118    60    15   0.9    0.06     0.052
## 119    80    15   0.9    0.056    0.044
## 120   100    15   0.9    0.028    0.038
## 121    20    20   0.9    0.054    0.066
## 122    40    20   0.9    0.05     0.04 
## 123    60    20   0.9    0.05     0.048
## 124    80    20   0.9    0.052    0.05 
## 125   100    20   0.9    0.05     0.056

6 Visualisasi Heatmap

hasil_long <- hasil_simulasi %>%
  pivot_longer(cols = c(power_x1, power_x2),
               names_to = "variabel",
               values_to = "power")

ggplot(hasil_long, aes(x = sigma, y = r2, fill = power)) +
  geom_tile(color = "white") +
  facet_grid(variabel ~ n) +
  scale_fill_gradient(low = "white", high = "red") +
  labs(title = "Heatmap Power Uji Ramsey RESET",
       x = "Keragaman Galat (σ)",
       y = "Keeratan Hubungan (R²)",
       fill = "Power") +
  theme_minimal()

7 Visualisasi 3d Interaktif

plot_ly(hasil_long, 
        x = ~sigma, 
        y = ~r2, 
        z = ~power, 
        color = ~variabel, 
        type = "scatter3d", 
        mode = "markers") %>%
  layout(title = "Visualisasi 3D Power Uji Ramsey RESET",
         scene = list(xaxis = list(title = "σ"),
                      yaxis = list(title = "R²"),
                      zaxis = list(title = "Power")))
## Warning in RColorBrewer::brewer.pal(N, "Set2"): minimal value for n is 3, returning requested palette with 3 different levels
## Warning in RColorBrewer::brewer.pal(N, "Set2"): minimal value for n is 3, returning requested palette with 3 different levels

8 Uji Kombinasi Manual dan Visualisasi Permukaan (Opsional)

uji_kombinasi_manual_final <- function(n, sigma, r2, reps = 500, seed = 123) {
  set.seed(seed)
  beta <- sqrt(r2 / (2 * (1 - r2)))
  
  simulasi_once <- function() {
    x1 <- rnorm(n)
    x2 <- rnorm(n)
    y <- beta * x1 + beta * x2 + rnorm(n, sd = sigma)
    model1 <- lm(y ~ x1)
    model2 <- lm(y ~ x2)
    reset1 <- resettest(model1, power = 2:3, type = "regressor")
    reset2 <- resettest(model2, power = 2:3, type = "regressor")
    c(x1_sig = reset1$p.value < 0.05, x2_sig = reset2$p.value < 0.05)
  }
  
  hasil <- replicate(reps, simulasi_once())
  power_x1 <- mean(hasil[1, ])
  power_x2 <- mean(hasil[2, ])
  
  cat("===== Hasil Power Uji Ramsey RESET (Manual) =====\n")
  cat("Ukuran Sampel (n):", n, "\n")
  cat("Standar Deviasi Galat (σ):", sigma, "\n")
  cat("Keeratan Hubungan (R²):", r2, "\n")
  cat("--------------------------------------------------\n")
  cat("Power X1:", round(power_x1, 3), "\n")
  cat("Power X2:", round(power_x2, 3), "\n")
  
  x1 <- rnorm(n)
  x2 <- rnorm(n)
  y <- beta * x1 + beta * x2 + rnorm(n, sd = sigma)
  model <- lm(y ~ x1 + x2)
  
  x1_seq <- seq(-2, 2, length.out = 40)
  x2_seq <- seq(-2, 2, length.out = 40)
  grid <- expand.grid(x1 = x1_seq, x2 = x2_seq)
  grid$y_hat <- predict(model, newdata = grid)
  
  z_matrix <- matrix(grid$y_hat, nrow = length(x1_seq), ncol = length(x2_seq))
  
  persp(x = x1_seq, y = x2_seq, z = z_matrix,
        xlab = "x1", ylab = "x2", zlab = "ŷ",
        main = paste0("Permukaan Prediksi (n=", n, ", σ=", sigma, ", R²=", r2, ")"),
        theta = 30, phi = 20, col = "skyblue", shade = 0.5, border = NA)
}
uji_kombinasi_manual_final(n = 20, sigma = 1, r2 = 0.1)
## ===== Hasil Power Uji Ramsey RESET (Manual) =====
## Ukuran Sampel (n): 20 
## Standar Deviasi Galat (σ): 1 
## Keeratan Hubungan (R²): 0.1 
## --------------------------------------------------
## Power X1: 0.034 
## Power X2: 0.054

9 Kesimpulan

Simulasi ini menunjukkan bahwa: