Urgensi Penanganan Autokorelasi

Asumsi Gauss-Markov dari model regresi linear OLS (Ordinary Least Squares/Metode Kuadrat Terkecil):

  1. Sisaan menyebar normal,
  2. Sisaan homoskedastis (ragam tidak berubah-ubah untuk tiap amatan), dan
  3. Sisaan saling bebas (tidak ada autokorelasi).

Apabila asumsi-asumsi di atas dapat terpenuhi, maka metode OLS akan menghasilkan model BLUE (Best Linear Unbiased Estimator). Best berarti ragam model akan minimum dan pada kasus autokorelasi ragam OLS akan tidak minimum. Oleh karena itu, penanganan autokorelasi menjadi penting untuk dilakukan.

Data

Data yang digunakan adalah Indeks Pembangunan Manusia Metode Baru di Provinsi Bali dari tahun 2010 sampai dengan 2021. Nilai minimum dari IPM adalah 0 dan nilai maksimumnya adalah 100. Nilai 70 - 79.9 dikategorikan sebagai “tinggi”, sementara 80 ke atas adalah sangat tinggi.

Menggunakan data-data IPM di waktu lampau, kita dapat menerapkan regresi linier dan memprediksi kapan Provinsi Bali dapat memiliki IPM 80 atau lebih.

library(readxl)
df <- read_excel("IPM Bali.xlsx")
str(df)
## tibble [12 x 2] (S3: tbl_df/tbl/data.frame)
##  $ Tahun   : num [1:12] 2010 2011 2012 2013 2014 ...
##  $ IPM Bali: chr [1:12] "70.10" "70.87" "71.62" "72.09" ...

Praproses Data

Data yang kita miliki masih belum bersih. Kolom kedua (nilai IPM) masih bertipe “character”. Kolom perlu diubah menjadi data “numeric” agar dapat dimasukkan ke dalam regresi linier.

colnames(df) <- c("tahun", "ipm")
df <- as.data.frame(lapply(df, as.numeric))
str(df)
## 'data.frame':    12 obs. of  2 variables:
##  $ tahun: num  2010 2011 2012 2013 2014 ...
##  $ ipm  : num  70.1 70.9 71.6 72.1 72.5 ...

Eksplorasi Data

Data menunjukkan tren naik. Selama 11 tahun, IPM di Provinsi Bali mengalami kenaikan sebesar 5.6 poin. Artinya, kenaikan rata-rata IPM Bali adalah 0.51 per tahun. Jika tren tersebut terus berlanjut, maka Bali akan memiliki IPM senilai 80 di tahun 2030 (sekitar 8.45 tahun setelah 2021). Apakah regresi linier akan menghasilkan prediksi yang sama? Mari kita cek.

Regresi Linier Sederhana

fit1 <- lm(ipm ~ tahun, data = df)
summary(fit1)
## 
## Call:
## lm(formula = ipm ~ tahun, data = df)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.4800 -0.1125  0.0800  0.1725  0.2500 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -974.75000   41.22291  -23.65 4.15e-10 ***
## tahun          0.52000    0.02045   25.42 2.03e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2446 on 10 degrees of freedom
## Multiple R-squared:  0.9848, Adjusted R-squared:  0.9832 
## F-statistic: 646.4 on 1 and 10 DF,  p-value: 2.033e-10

Model hasil regresi linier sederhana adalah IPM = -974.75 + 0.52*(Tahun). Artinya, setiap perubahan tahun akan diiringi oleh perubahan sebesar 0.52 rataan poin IPM. Nilai tersebut tidak jauh berbeda dari hasil perhitungan perubahan rata-rata di atas.

futureYear <- 2022:2030
pred1 <- cbind(futureYear, predict(fit1, data.frame(tahun = futureYear)))
colnames(pred1) <- c('tahun', 'ipm')
pred1
##   tahun   ipm
## 1  2022 76.69
## 2  2023 77.21
## 3  2024 77.73
## 4  2025 78.25
## 5  2026 78.77
## 6  2027 79.29
## 7  2028 79.81
## 8  2029 80.33
## 9  2030 80.85

Prediksi regresi linier ialah Bali akan menyentuh IPM sebesar 80 pada tahun 2028 - 2029. Hasil prediksi regresi linier rupanya sedikit lebih optimistis dibandingkan metode prediksi primitif menggunakan perubahan rata-rata. Walaupun begitu, tidak serta-merta hasil tersebut dapat kita terima karena masih diperlukan uji asumsi sebelum regresi tersebut dapat dikatakan sebagai valid.

car::durbinWatsonTest(fit1)
##  lag Autocorrelation D-W Statistic p-value
##    1       0.2673019     0.8754597   0.002
##  Alternative hypothesis: rho != 0

Berdasarkan eksplorasi visual, dapat dilihat bahwa nilai-nilai autokorelasi tidak ada yang melewati batas toleransi sehingga tidak ada yang signifikan berbeda dari nol. Akan tetapi, uji formal Durbin-Watson justru menunjukkan p-value yang berada di bawah taraf nyata alpha = 0.05. Oleh karena itu, regresi linier sederhana di atas mengandung autokorelasi yang cukup signifikan.

Penyesuaian Cochrane–Orcutt

fit2 <- orcutt::cochrane.orcutt(fit1, max.iter = 1000)
summary(fit2)
## Call:
## lm(formula = ipm ~ tahun, data = df)
## 
##                Estimate  Std. Error t value  Pr(>|t|)    
## (Intercept) -682.272588  127.751754  -5.341 0.0004681 ***
## tahun          0.375121    0.063282   5.928 0.0002213 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1757 on 9 degrees of freedom
## Multiple R-squared:  0.7961 ,  Adjusted R-squared:  0.7734
## F-statistic: 35.1 on 1 and 9 DF,  p-value: < 2.213e-04
## 
## Durbin-Watson statistic 
## (original):    0.87546 , p-value: 3.527e-03
## (transformed): 2.13019 , p-value: 4.363e-01
futureYear2 <- 2022:2035
pred2 <- cbind(futureYear2, coefficients(fit2)[[1]] + coefficients(fit2)[[2]]*
                   futureYear2)
colnames(pred2) <- c('tahun', 'ipm')
pred2
##       tahun      ipm
##  [1,]  2022 76.22225
##  [2,]  2023 76.59737
##  [3,]  2024 76.97249
##  [4,]  2025 77.34761
##  [5,]  2026 77.72273
##  [6,]  2027 78.09785
##  [7,]  2028 78.47297
##  [8,]  2029 78.84809
##  [9,]  2030 79.22322
## [10,]  2031 79.59834
## [11,]  2032 79.97346
## [12,]  2033 80.34858
## [13,]  2034 80.72370
## [14,]  2035 81.09882
plot(df, type = "l", main = "Cochrane-Orcutt Prediction", xlim = c(2010, 2035),
     ylim = c(68, 82))
lines(rbind(df[dim(df)[1], ], pred2), col = "blue", lty = "dashed")
abline(v = 2032.5, h = 80, lty = "dotted")

Hasil regresi yang telah disesuaikan menggunakan metode Cochrane-Orcutt memiliki hasil yang berbeda dari sebelum-sebelumnya. Regresi ini justru memprediksikan bahwa IPM Bali baru akan mencapai poin 80 di tahun 2032 atau 2033. Regresi ini juga berhasil menangani autokorelasi, ditunjukkan oleh plot ACF, plot PACf, dan uji Durbin-Watson di bawah.

fit2$DW[3:4]
##        DW           
## 2.1301908 0.4363058

Walaupun plot ACF menunjukkan pola siklikal dan PACF ada yang melawati batas, uji Durbin-Watson menunjukkan p-value sebesar 0.44 atau jauh berada di atas alpha = 0.05 sehingga autokorelasi tidak ada.

Penyesuaian Hildreth-Lu

Metode Hildreth-Lu akan menggunakan rho dari metode Cochrane-OOrcutt yang sebelumnya.

fit3 <- HoRM::hildreth.lu(df$ipm, df$tahun, fit2$rho)
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
summary(fit3)
## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.18524 -0.15403 -0.03408  0.13616  0.24765 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -180.56944   33.81063  -5.341 0.000468 ***
## x              0.37512    0.06328   5.928 0.000221 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1757 on 9 degrees of freedom
## Multiple R-squared:  0.7961, Adjusted R-squared:  0.7734 
## F-statistic: 35.14 on 1 and 9 DF,  p-value: 0.0002213
pred3 <- cbind(futureYear2, coef(fit3)[[1]]/(1 - fit2$rho) + coef(fit3)[[2]]*
                   futureYear2)
colnames(pred3) <- c('tahun', 'ipm')
pred3
##       tahun      ipm
##  [1,]  2022 76.22225
##  [2,]  2023 76.59737
##  [3,]  2024 76.97249
##  [4,]  2025 77.34761
##  [5,]  2026 77.72273
##  [6,]  2027 78.09785
##  [7,]  2028 78.47297
##  [8,]  2029 78.84809
##  [9,]  2030 79.22322
## [10,]  2031 79.59834
## [11,]  2032 79.97346
## [12,]  2033 80.34858
## [13,]  2034 80.72370
## [14,]  2035 81.09882
plot(df, type = "l", main = "Hildreth-Lu Prediction", xlim = c(2010, 2035),
     ylim = c(68, 82))
lines(rbind(df[dim(df)[1], ], pred3), col = "blue", lty = "dashed")
abline(v = 2032.5, h = 80, lty = "dotted")

Metode Hildreth-Lu rupanya menghasilkan prediksi yang tidak jauh berbeda dengan metode Cochrane-Orcutt di atas. Hasil prediksi juga sama, yaitu Bali akan mencapai poin ipm 80 pada tahun 2032 atau 2033. Bahkan hasil uji Durbin-Watson pun memiliki kesimpulan yang sama, yaitu sisaan tidak mengalami autokorelasi.

car::durbinWatsonTest(fit3)
##  lag Autocorrelation D-W Statistic p-value
##    1      -0.1240508      2.130191   0.904
##  Alternative hypothesis: rho != 0

Uji Asumsi Kondisi Gauss-Markov

Walaupun telah berhasil melalui uji autokorelasi, model regresi hasil penanganan juga harus berhasil melalui asumsi-asumsi lainnya. Asumsi yang akan diuji adalah normalitas menggunakan Shapiro-Wilkinson dan homoskedastisitas menggunakan Breusch-Pagan.

Uji Normalitas Sisaan Regresi Hasil Penanganan

shapiro.test(fit2$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  fit2$residuals
## W = 0.93952, p-value = 0.4919
shapiro.test(fit3$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  fit3$residuals
## W = 0.89104, p-value = 0.1433

Model hasil penyesuaian Cochrane-Orcutt adalah fit2 dan Hildreth-Lu adalah fit3. Kedua model memiliki sisaan yang menyebar normal karena p-value > alpha. Akan tetapi, model Cochrane-Orcutt cenderung lebih menyerupai sebaran normal dibandingkan model Hildreth-Lu.

Uji Homoskedastisitas Sisaan Regresi Hasil Penanganan

lmtest::bptest(fit2)
## 
##  studentized Breusch-Pagan test
## 
## data:  fit2
## BP = 2.8908, df = 1, p-value = 0.08909
lmtest::bptest(fit3)
## 
##  studentized Breusch-Pagan test
## 
## data:  fit3
## BP = 2.8908, df = 1, p-value = 0.08909

Secara uji formal, kedua model memiliki sisaan yang homoskedastis. Akan tetapi, secara eksploratif, kita dapat melihat bahwa sisaan model Hildreth-Lu terkesan cenderung lebih acak terhadap nilai amatannya dibandingkan Cochrane-Orcutt.

par(mfrow = c(1, 2))
plot(fit2$fitted.values, fit2$residuals, ylim = c(-3, 3)*sd(fit2$residuals),
     main = "Sisaan Model Cochrane-Orcutt\ntiap Nilai Amatan")
abline(h = 0, col = "blue", lty = "dotted")

plot(fit3$fitted.values, fit3$residuals, ylim = c(-3, 3)*sd(fit3$residuals),
     main = "Sisaan Model Hildreth-Lu\ntiap Nilai Amatan")
abline(h = 0, col = "blue", lty = "dotted")

Simpulan

Regresi linier yang mengalami autokorelasi dapat ditangani menggunakan metode Cochrane-Orcutt dan Hildreth-Lu. Kedua metode saling berhubungan dan hasil prediksinya kurang lebih sama, yaitu Provinsi Bali kemungkinan besar akan mencapai IPM sebesar 80 pada tahun 2032 atau 2033.