#install.packages("https://cran.r-project.org/src/contrib/Archive/orcutt/orcutt_2.1.tar.gz", repos = NULL, type = "source")Penanganan Autokorelasi Deret Waktu
Dalam analisis regresi, salah satu asumsi penting adalah bahwa komponen galat (error term) bersifat independen atau tidak saling berkorelasi antar waktu. Namun, pada data deret waktu (time series), asumsi ini sering kali tidak terpenuhi karena adanya autokorelasi diri (autocorrelation) — yaitu kondisi ketika error pada periode sekarang \(e_{t}\) berkorelasi dengan error pada periode sebelumnya (\(e_{t-1}\)).
Autokorelasi dapat menyebabkan:
Estimasi parameter regresi tetap unbiased, tetapi tidak efisien (varian besar)
Standar error menjadi bias, sehingga uji t dan uji F bisa memberikan kesimpulan yang keliru
Untuk mengatasi hal ini, beberapa metode koreksi dikembangkan agar model regresi tetap valid dan efisien. Dua pendekatan umum yang sering digunakan adalah metode Cochrane–Orcutt dan metode Hildreth–Lu.
1 Packages
#install.packages("orcutt")
#install.packages("HoRM")library(dplyr)Warning: package 'dplyr' was built under R version 4.5.1
Attaching package: 'dplyr'The following objects are masked from 'package:stats':
filter, lagThe following objects are masked from 'package:base':
intersect, setdiff, setequal, unionlibrary(TTR)Warning: package 'TTR' was built under R version 4.5.1library(forecast)Warning: package 'forecast' was built under R version 4.5.1Registered S3 method overwritten by 'quantmod':
method from
as.zoo.data.frame zoo library(lmtest) #digunakan untuk uji formal pendeteksian autokorelasiWarning: package 'lmtest' was built under R version 4.5.1Loading required package: zoo
Attaching package: 'zoo'The following objects are masked from 'package:base':
as.Date, as.Date.numericlibrary(orcutt) #untuk membuat model regresi Cochrane-Orcutt
library(HoRM) #untuk membuat model regresi Hildreth-LuWarning: package 'HoRM' was built under R version 4.5.12 Input Data
Data yang digunakan dalam kesempatan kali ini adalah data IPM Provinsi Gorontalo periode tahun 2010-2021.
tahun <- c(2010,2011,2012,2013,2014,2015,2016,2017,2018,2019,2020,2021)
ipm <- c(62.65,63.48,64.16,64.70,65.17,65.86,66.29,67.01,67.71,68.49,68.68,69.00)
dtGorontalo <- cbind(tahun,ipm)
dtGorontalo <- as.data.frame(dtGorontalo)
dtGorontalo3 Eksplorasi Data
Sebelum melakukan regresi, akan diperlihatkan plot time-series dari IPM Provinsi Gorontalo Periode 2010-2021
#Membentuk objek time series
data.ts<-ts(dtGorontalo$ipm)
data.tsTime Series:
Start = 1
End = 12
Frequency = 1
[1] 62.65 63.48 64.16 64.70 65.17 65.86 66.29 67.01 67.71 68.49 68.68 69.00#Membuat plot time series
ts.plot(data.ts, xlab="Time Period ", ylab="IPM", main= "Time Series Plot of IPM")
points(data.ts)
Selanjutnya akan dilakukan ramalan dan pemulusan dengan metode DMA dan DES karena terlihat pada plot di atas menunjukkan adanya trend.
dt.sma <- SMA(data.ts, n=3)
dma <- SMA(dt.sma, n = 3)
At <- 2*dt.sma - dma
Bt <- 2/(3-1)*(dt.sma - dma)
dt.dma<- At+Bt
dt.ramal<- c(NA, dt.dma)
t = 1:5
f = c()
for (i in t) {
f[i] = At[length(At)] + Bt[length(Bt)]*(i)
}dt.gab <- cbind(aktual = c(data.ts,rep(NA,5)),
pemulusan1 = c(dt.sma,rep(NA,5)),
pemulusan2 = c(dt.dma, rep(NA,5)),
At = c(At, rep(NA,5)),
Bt = c(Bt,rep(NA,5)),
ramalan = c(dt.ramal, f[-1]))
dt.gab aktual pemulusan1 pemulusan2 At Bt ramalan
[1,] 62.65 NA NA NA NA NA
[2,] 63.48 NA NA NA NA NA
[3,] 64.16 63.43000 NA NA NA NA
[4,] 64.70 64.11333 NA NA NA NA
[5,] 65.17 64.67667 65.88333 65.28000 0.6033333 NA
[6,] 65.86 65.24333 66.37444 65.80889 0.5655556 65.88333
[7,] 66.29 65.77333 66.85778 66.31556 0.5422222 66.37444
[8,] 67.01 66.38667 67.55778 66.97222 0.5855556 66.85778
[9,] 67.71 67.00333 68.23444 67.61889 0.6155556 67.55778
[10,] 68.49 67.73667 69.12556 68.43111 0.6944444 68.23444
[11,] 68.68 68.29333 69.52444 68.90889 0.6155556 69.12556
[12,] 69.00 68.72333 69.66778 69.19556 0.4722222 69.52444
[13,] NA NA NA NA NA 69.66778
[14,] NA NA NA NA NA 70.14000
[15,] NA NA NA NA NA 70.61222
[16,] NA NA NA NA NA 71.08444
[17,] NA NA NA NA NA 71.55667#Plot time series
ts.plot(dt.gab[,1], xlab="Time Period ", ylab="IPM",
main= "DMA N=3 Data IPM", ylim=c(62,75))
points(dt.gab[,1])
points(dt.gab[,3])
points(dt.gab[,6])
lines(dt.gab[,3],col="green",lwd=2)
lines(dt.gab[,6],col="red",lwd=2)
legend("topleft",c("data aktual","data pemulusan","data peramalan"),
lty=8, col=c("black","green","red"), cex=0.8)
Selanjutnya akan dilihat keakuratan dari metode DMA
#Menghitung nilai keakuratan
error.dma = data.ts-dt.ramal[1:length(data.ts)]
SSE.dma = sum(error.dma[6:length(data.ts)]^2)
MSE.dma = mean(error.dma[6:length(data.ts)]^2)
MAPE.dma = mean(abs((error.dma[6:length(data.ts)]/data.ts[6:length(data.ts)])*100))
akurasi.dma <- matrix(c(SSE.dma, MSE.dma, MAPE.dma))
row.names(akurasi.dma)<- c("SSE", "MSE", "MAPE")
colnames(akurasi.dma) <- c("Akurasi m = 3")
akurasi.dma Akurasi m = 3
SSE 0.59288889
MSE 0.08469841
MAPE 0.34238965Selanjutnya akan digunakan metode Double Exponential Smoothing dengan cara sebagai berikut.
Pertama akan data akan dibagi menjadi data training dan data testing.
#membagi training dan testing
training<-dtGorontalo[1:10,2]
testing<-dtGorontalo[11:12,2]
#data time series
training.ts<-ts(training)
testing.ts<-ts(testing,start=11)
#eksplorasi data
plot(data.ts, col="red",main="Plot semua data")
points(data.ts)
plot(training.ts, col="blue",main="Plot data training")
points(training.ts)
Selanjutnya akan dilakukan pemulusan dengan DES, kali ini langsung dicari lambda dan gamma optimum sebagai berikut. Nilai lambda dan gamma optimum dapat dilihat pada smoothing parameters alpha untuk nilai lambda dan beta untuk nilai gamma.
#Lamda dan gamma optimum
des.opt<- HoltWinters(training.ts, gamma = FALSE)
des.optHolt-Winters exponential smoothing with trend and without seasonal component.
Call:
HoltWinters(x = training.ts, gamma = FALSE)
Smoothing parameters:
alpha: 0.7826887
beta : 1
gamma: FALSE
Coefficients:
[,1]
a 68.4727317
b 0.7823472plot(des.opt)
legend("topleft", c("Data Aktual", "Peramalan"), col = c("black", "red"),
lty = c(1,1))
#ramalan
ramalandesopt<- forecast(des.opt, h=5)
ramalandesopt Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
11 69.25508 69.04282 69.46734 68.93045 69.57971
12 70.03743 69.64314 70.43171 69.43442 70.64043
13 70.81977 70.18427 71.45528 69.84785 71.79170
14 71.60212 70.68262 72.52163 70.19586 73.00838
15 72.38447 71.14531 73.62363 70.48934 74.27960Selanjutnya akan dicari akurasi dari metode DES.
ssedes.train<-des.opt$SSE
msedes.train<-ssedes.train/length(training.ts)
sisaandes<-ramalandesopt$residuals
head(sisaandes)Time Series:
Start = 1
End = 6
Frequency = 1
[1] NA NA -0.1500000 -0.2051934 -0.1265850 0.2095741mapedes.train <- sum(abs(sisaandes[3:length(training.ts)]/training.ts[3:length(training.ts)])*100)/length(training.ts)
akurasides.opt <- matrix(c(ssedes.train,msedes.train,mapedes.train))
row.names(akurasides.opt)<- c("SSE", "MSE", "MAPE")
colnames(akurasides.opt) <- c("Akurasi lamda dan gamma optimum")
akurasides.opt Akurasi lamda dan gamma optimum
SSE 0.19249578
MSE 0.01924958
MAPE 0.17991246#Akurasi data testing
selisihdesopt<-ramalandesopt$mean-testing.ts
selisihdesoptTime Series:
Start = 11
End = 12
Frequency = 1
[1] 0.575079 1.037426SSEtestingdesopt<-sum(selisihdesopt^2)
SSEtestingdesopt<-SSEtestingdesopt/length(testing.ts)
MAPEtestingdesopt<-sum(abs(selisihdesopt/testing.ts)*100)/length(testing.ts)
akurasiDesTesting <- matrix(c(SSEtestingdesopt,SSEtestingdesopt,MAPEtestingdesopt))
row.names(akurasiDesTesting)<- c("SSE", "MSE", "MAPE")
colnames(akurasiDesTesting) <- c("Akurasi lamda dan gamma optimum")
akurasiDesTesting Akurasi lamda dan gamma optimum
SSE 0.7034845
MSE 0.7034845
MAPE 1.1704237Setelah didapatkan nilai akurasi untuk metode DMA dan DES, selanjutnya akan dibandingkan keakuratan antar metode keduanya.
cbind(akurasi.dma, akurasides.opt) Akurasi m = 3 Akurasi lamda dan gamma optimum
SSE 0.59288889 0.19249578
MSE 0.08469841 0.01924958
MAPE 0.34238965 0.17991246Berdasarkan perbandingan akurasi tersebut, terlihat nilai SSE, MSE, dan MAPE metode DES lebih kecil dibandingkan dengan metode DMA. Oleh karena itu, metode peramalan dan pemulusan yang terbaik antara keduanya adalah dengan metode DES.
Setelah melakukan peramalan, data yang telah dimasukkan kemudian dieksplorasi. Eksplorasi pertama yang dilakukan adalah dengan menggunakan scatter plot.
#Eksplorasi Data
#Pembuatan Scatter Plot
plot(tahun,ipm, pch = 20, col = "blue",
main = "Scatter Plot Tahun vs Nilai IPM",
xlab = "Tahun",
ylab = "Nilai IPM")
#Menampilkan Nilai Korelasi
cor(tahun,ipm)[1] 0.9966848Berdasarkan scatter plot di atas, terlihat adanya hubungan / korelasi positif antara peubah tahun dengan nilai IPM, terlihat titik-titik pada plot yang naik ke arah kanan atas. Hal tersebut juga diperkuat dengan hasil perhitungan aplikasi R di mana didapatkan nilai korelasi sebesar \(0.9966848\).
Setalah mengetahui adanya hubungan antar dua peubah, maka model regresi dapat ditentukan.
4 Model Regresi Linear: IPM terhadap Tahun
#Pembuatan Model Regresi
#model regresi
model<- lm(ipm~tahun, data = dtGorontalo)
summary(model)
Call:
lm(formula = ipm ~ tahun, data = dtGorontalo)
Residuals:
Min 1Q Median 3Q Max
-0.33000 -0.07295 0.02591 0.08000 0.33455
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -1.118e+03 3.055e+01 -36.58 5.55e-12 ***
tahun 5.873e-01 1.516e-02 38.74 3.14e-12 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.1813 on 10 degrees of freedom
Multiple R-squared: 0.9934, Adjusted R-squared: 0.9927
F-statistic: 1501 on 1 and 10 DF, p-value: 3.136e-12Model yang dihasilkan adalah \[y_i=-1118+0.5873x_t\] Berdasarkan ringkasan model dapat diketahui bahwa hasil uji F memiliki p-value < \(\alpha\) (5%). Artinya, minimal terdapat satu variabel yang berpengaruh nyata terhadap model. Hasil uji-t parsial kedua parameter regresi, yaitu intersep dan koefisien regresi juga menunjukkan hal yang sama, yaitu memiliki p-value < \(\alpha\) (5%) sehingga nyata dalam taraf 5%. Selanjutnya dapat dilihat juga nilai \(R^2=0.9934\). Artinya, sebesar 99.34% keragaman nilai IPM dapat dijelaskan oleh peubah tahun. Hasil ini menunjukkan hasil yang bagus, seolah mendapatkan hasil terbaik. Namun, kita perlu melakukan uji terhadap sisaannya seperti berikut ini.
#sisaan dan fitted value
sisaan<- residuals(model)
fitValue<- predict(model)
#Diagnostik dengan eksploratif
par(mfrow = c(2,2))
qqnorm(sisaan)
qqline(sisaan, col = "steelblue", lwd = 2)
plot(fitValue, sisaan, col = "steelblue", pch = 20, xlab = "Sisaan", ylab = "Fitted Values", main = "Sisaan vs Fitted Values")
abline(a = 0, b = 0, lwd = 2)
hist(sisaan, col = "steelblue")
plot(seq(1,12,1), sisaan, col = "steelblue", pch = 20, xlab = "Sisaan", ylab = "Order", main = "Sisaan vs Order")
lines(seq(1,12,1), sisaan, col = "red")
abline(a = 0, b = 0, lwd = 2)
Dua plot di samping kiri digunakan untuk melihat apakah sisaan menyebar normal. Normal Q-Q Plot di atas menunjukkan bahwa sisaan cenderung menyebar normal, tetapi histogram dari sisaan tidak menunjukkan demikian. Selanjutnya, dua plot di samping kanan digunakan untuk melihat autokorelasi. Plot Sisaan vs Fitted Value dan Plot Sisaan vs Order menunjukkan adanya pola pada sisaan. Untuk lebih lanjut akan digunakan uji formal melihat normalitas sisaan dan plot ACF dan PACF untuk melihat apakah ada autokorelasi atau tidak.
#Melihat Sisaan Menyebar Normal/Tidak
#H0: sisaan mengikuti sebaran normal
#H1: sisaan tidak mengikuti sebaran normal
shapiro.test(sisaan)
Shapiro-Wilk normality test
data: sisaan
W = 0.97479, p-value = 0.954ks.test(sisaan, "pnorm", mean=mean(sisaan), sd=sd(sisaan))
Exact one-sample Kolmogorov-Smirnov test
data: sisaan
D = 0.13564, p-value = 0.9589
alternative hypothesis: two-sidedBerdasarkan uji formal Saphiro-Wilk dan Kolmogorov-Smirnov didapatkan nilai p-value > \(\alpha\) (5%). Artinya, cukup bukti untuk menyatakan sisaan berdistribusi normal.
#ACF dan PACF identifikasi autokorelasi
par(mfrow = c(1,2))
acf(sisaan)
pacf(sisaan)
Berdasarkan plot ACF dan PACF, terlihat semua dalam rentang batas dan tidak ada lag yang signifikan. Namun, untuk lebih memastikan akan dilakukan uji formal dengan uji Durbin Watson.
#Deteksi autokorelasi dengan uji-Durbin Watson
#H0: tidak ada autokorelasi
#H1: ada autokorelasi
dwtest(model)
Durbin-Watson test
data: model
DW = 1.2644, p-value = 0.03741
alternative hypothesis: true autocorrelation is greater than 0Berdasarkan hasil DW Test, didapatkan nilai \(DW = 1.2644\) dan p-value = \(0.03741\). Berdasarkan tabel Durbin-Watson diperoleh nilai \(DL = 0.9771\) dan \(DU = 1.331\). Nilai DW masih berada di antara nilai DL dan DU. Artinya, berada di daerah inkonklusif, tidak dapat dikatakan berada di daerah autokorelasi positif maupun bebas dari autokorelasi. Namun, dengan nilai p-value < 0.05 dapat disimpulkan bahwa tolak H0, cukup bukti mengatakan adanya autokorelasi. Oleh karena itu, diperlukan penangan autokorelasi. Penanganan yang akan digunakan menggunakan dua metode, yaitu Cochrane-Orcutt dan Hildret-Lu.
5 Penanganan Autokorelasi
5.1 Metode Cochrane-Orcutt
Penanganan metode Cochrane-Orcutt dapat dilakukan dengan bantuan packages Orcutt pada aplikasi R maupun secara manual. Berikut ini ditampilkan cara menggunakan bantuan library packages Orcutt.
#Penanganan Autokorelasi Cochrane-Orcutt
modelCO<-cochrane.orcutt(model)
modelCOCochrane-orcutt estimation for first order autocorrelation
Call:
lm(formula = ipm ~ tahun, data = dtGorontalo)
number of interaction: 30
rho 0.340921
Durbin-Watson statistic
(original): 1.26437 , p-value: 3.741e-02
(transformed): 1.49014 , p-value: 9.403e-02
coefficients:
(Intercept) tahun
-1062.042646 0.559755 Hasil keluaran model setelah dilakukan penanganan adalah sebagai berikut. \[y_i=-1062.042646+0.559755x_t\] Hasil juga menunjukkan bahwa nilai DW dan p-value meningkat menjadi \(1.49014\) dan \(0.09403\). Nilai DW sudah berada pada rentang DU < DW < 4-DU atau \(1.331 < DW < 2.669\). Hal tersebut juga didukung dengan nilai p-value > 0.05, artinya belum cukup bukti menyatakan bahwa sisaan terdapat autokorelasi pada taraf nyata 5%. Untuk nilai \(ρ ̂\) optimum yang digunakan adalah \(0.3409214\). Nilai tersebut dapat diketahui dengan syntax berikut.
#Rho optimum
rho<- modelCO$rho
rho[1] 0.3409214Selanjutnya akan dilakukan transformasi secara manual dengan syntax berikut ini.
ipm [1] 62.65 63.48 64.16 64.70 65.17 65.86 66.29 67.01 67.71 68.49 68.68 69.00ipm[-1] [1] 63.48 64.16 64.70 65.17 65.86 66.29 67.01 67.71 68.49 68.68 69.00#Transformasi Manual
ipm.trans<- ipm[-1]-ipm[-12]*rho
tahun.trans<- tahun[-1]-tahun[-12]*rho
modelCOmanual<- lm(ipm.trans~tahun.trans)
summary(modelCOmanual)
Call:
lm(formula = ipm.trans ~ tahun.trans)
Residuals:
Min 1Q Median 3Q Max
-0.22770 -0.11619 -0.00273 0.05763 0.33083
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -699.96956 31.38081 -22.31 3.46e-09 ***
tahun.trans 0.55975 0.02361 23.71 2.02e-09 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.1632 on 9 degrees of freedom
Multiple R-squared: 0.9842, Adjusted R-squared: 0.9825
F-statistic: 562 on 1 and 9 DF, p-value: 2.017e-09Hasil model transformasi bukan merupakan model sesungguhnya. Koefisien regresi masih perlu dicari kembali mengikuti \(β_0^*=β_0+ρ ̂β_0\) dan \(β_1^*=β_1\).
#Mencari Penduga Koefisien Regresi setelah Transformasi ke Persamaan Awal
b0bintang <- modelCOmanual$coefficients[-2]
b0 <- b0bintang/(1-rho)
b1 <- modelCOmanual$coefficients[-1]
b0(Intercept)
-1062.043 b1tahun.trans
0.5597546 Hasil perhitungan koefisien regresi tersebut akan menghasilkan hasil yang sama dengan model yang dihasilkan menggunakan packages.
5.2 Metode Hildreth-Lu
Penanganan kedua adalah menggunakan metode Hildreth-Lu. Metode ini akan mencari nilai SSE terkecil dan dapat dicari secara manual maupun menggunakan packages. Jika menggunakan packages, gunakan library packages HORM.
library(HoRM)
# ambil y dan x (satu prediktor saja)
y <- dtGorontalo$ipm
x <- dtGorontalo$tahun
# grid search rho 0–0.9
r <- seq(0, 0.9, by = 0.1)
tab <- data.frame(
rho = r,
SSE = sapply(r, function(rr) deviance(hildreth.lu(y, x, rr)))
)
tab # tabel rho vs SSErho_hat <- tab$rho[which.min(tab$SSE)] # rho terbaik
fit_hl <- hildreth.lu(y, x, rho_hat) # model final
summary(fit_hl)
Call:
lm(formula = y ~ x)
Residuals:
Min 1Q Median 3Q Max
-0.23518 -0.10992 -0.00927 0.05327 0.33320
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -747.18687 31.39851 -23.80 1.95e-09 ***
x 0.56242 0.02224 25.28 1.14e-09 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.1633 on 9 degrees of freedom
Multiple R-squared: 0.9861, Adjusted R-squared: 0.9846
F-statistic: 639.2 on 1 and 9 DF, p-value: 1.139e-09#Penanganan Autokorelasi Hildreth lu
# Hildreth-Lu
hildreth.lu.func<- function(r, model){
x <- model.matrix(model)[,-1]
y <- model.response(model.frame(model))
n <- length(y)
t <- 2:n
y <- y[t]-r*y[t-1]
x <- x[t]-r*x[t-1]
return(lm(y~x))
}
#Pencariab rho yang meminimumkan SSE
r <- c(seq(0.1,0.9, by= 0.1))
tab <- data.frame("rho" = r, "SSE" = sapply(r, function(i){deviance(hildreth.lu.func(i, model))}))
round(tab, 4)Pertama-tama akan dicari di mana kira-kira \(ρ\) yang menghasilkan SSE minimum. Pada hasil di atas terlihat \(ρ\) minimum ketika 0.3. Namun, hasil tersebut masih kurang teliti sehingga akan dicari kembali \(ρ\) yang lebih optimum dengan ketelitian yang lebih. Jika sebelumnya jarak antar \(ρ\) yang dicari adalah 0.1, kali ini jarak antar \(ρ\) adalah 0.001 dan dilakukan pada selang 0.2 sampai dengan 0.5.
#Rho optimal di sekitar 0.4
rOpt <- seq(0.2,0.5, by= 0.001)
tabOpt <- data.frame("rho" = rOpt, "SSE" = sapply(rOpt, function(i){deviance(hildreth.lu.func(i, model))}))
head(tabOpt[order(tabOpt$SSE),])#Grafik SSE optimum
par(mfrow = c(1,1))
plot(tab$SSE ~ tab$rho , type = "l", xlab = "Rho", ylab = "SSE")
abline(v = tabOpt[tabOpt$SSE==min(tabOpt$SSE),"rho"], lty = 2, col="red",lwd=2)
text(x=0.341, y=0.2397500, labels = "rho=0.341", cex = 0.8)
Perhitungan yang dilakukan aplikasi R menunjukkan bahwa nilai \(ρ\) optimum, yaitu saat SSE terkecil terdapat pada nilai \(ρ=0.341\). Hal tersebut juga ditunjukkan pada plot. Selanjutnya, model dapat didapatkan dengan mengevaluasi nilai \(ρ\) ke dalam fungsi hildreth.lu.func, serta dilanjutkan dengan pengujian autokorelasi dengan uji Durbin-Watson. Namun, setelah pengecekan tersebut tidak lupa koefisien regresi tersebut digunakan untuk transformasi balik. Persamaan hasil transformasi itulah yang menjadi persamaan sesungguhnya.
#Model terbaik
modelHL <- hildreth.lu.func(0.341, model)
summary(modelHL)
Call:
lm(formula = y ~ x)
Residuals:
Min 1Q Median 3Q Max
-0.22769 -0.11621 -0.00271 0.05763 0.33083
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -699.87889 31.38081 -22.3 3.47e-09 ***
x 0.55975 0.02361 23.7 2.02e-09 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.1632 on 9 degrees of freedom
Multiple R-squared: 0.9842, Adjusted R-squared: 0.9825
F-statistic: 561.9 on 1 and 9 DF, p-value: 2.02e-09#Transformasi Balik
cat("y = ", coef(modelHL)[1]/(1-0.341), "+", coef(modelHL)[2],"x", sep = "")y = -1062.032+0.5597492xSetelah dilakukan tranformasi balik, didapatkan model dengan metode Hildreth-Lu sebagai berikut. \[y_i=-1062.032+0.5597492x_t\]
#Deteksi autokorelasi
dwtest(modelHL)
Durbin-Watson test
data: modelHL
DW = 1.4902, p-value = 0.09404
alternative hypothesis: true autocorrelation is greater than 0Hasil uji Durbin-Watson juga menunjukkan bawah nilai DW sebesar \(1.4092\) berada pada selang daerah tidak ada autokorelasi, yaitu pada rentang DU < DW < 4-DU atau \(1.331 < DW < 2.669\). Hal tersebut juga didukung oleh p-value sebesar \(0.09404\), di mana p-value > \(\alpha\)=5%. Artinya tak tolak \(H_0\) atau belum cukup bukti menyatakan bahwa ada autokorelasi dalam data nilai IPM dengan metode Hildreth-Lu pada taraf nyata 5%.
Terakhir, akan dibandingkan nilai SSE dari ketiga metode (metode awal, metode Cochrane-Orcutt, dan Hildreth-Lu).
#Perbandingan
sseModelawal <- anova(model)$`Sum Sq`[-1]
sseModelCO <- anova(modelCOmanual)$`Sum Sq`[-1]
sseModelHL <- anova(modelHL)$`Sum Sq`[-1]
mseModelawal <- sseModelawal/length(ipm)
mseModelCO <- sseModelCO/length(ipm)
mseModelHL <- sseModelHL/length(ipm)
akurasi <- matrix(c(sseModelawal,sseModelCO,sseModelHL,
mseModelawal,mseModelCO,mseModelHL),nrow=2,ncol=3,byrow = T)
colnames(akurasi) <- c("Model Awal", "Model Cochrane-Orcutt", "Model Hildreth-Lu")
row.names(akurasi) <- c("SSE","MSE")
akurasi Model Awal Model Cochrane-Orcutt Model Hildreth-Lu
SSE 0.32863636 0.23974997 0.23974997
MSE 0.02738636 0.01997916 0.01997916Berdasarkan hasil tersebut dapat diketahui bahwa hasil penanganan autokorelasi dengan metode Cochrane-Orcutt dan Hildreth-Lu memiliki SSE yang sama, sebesar \(0.23975\) dan lebih baik dibandingkan model awal ketika autokorelasi masih terjadi, yaitu sebesar \(0.3286364\).
6 Simpulan
Autokorelasi yang terdapat pada data IPM terjadi akibat adanya korelasi di antara unsur penyusunnya. Indikator IPM yang erat hubungannya dengan perekonomian sangat rawan menjadi penyebab adanya autokorelasi. Adanya autokorelasi menyebabkan model regresi kurang baik karena akan meingkatkan galatnya. Autokorelasi dapat dideteksi secara eksploratif melalui plot sisaan, ACF, dan PACF, serta dengan uji formal Durbin-Watson. Namun, autokorelasi tersebut dapat ditangani dengan metode Cochrane-Orcutt dan Hildreth-Lu. Kedua metode menghasilkan nilai SSE yang sama, artinya keduanya baik untuk digunakan.
7 Daftar Pustaka
Aprianto A, Debataraja NN, Imro’ah N. 2020. Metode cochrane-orcutt untuk mengatasi autokorelasi pada estimasi parameter ordinary least squares. Bimaster : Buletin Ilmiah Matematika, Statistika dan Terapannya. 9(1):95–102. doi:10.26418/bbimst.v9i1.38590.
BPS. 2021a. Indeks Pembangunan Manusia 2020. Jakarta (ID): Badan Pusat Statistik.
BPS. 2021b. Indeks Pembangunan Manusia (IPM) 2021. Berita Resmi Statistik., siap terbit.