Tugas STA1341 : Pendeteksian dan Penanganan Autokorelasi pada IPM Provinsi Jawa Barat Tahun 2010 - 2021
Indeks Pembangunan Manusia (IPM) merupakan salah satu cara untuk mengukur taraf kualitas fisik dan non fisik penduduk. Kualitas fisik tercermin dari angka harapan hidup. Sedangkan kualitas non fisik (intelektualitas) melalui lamanya rata-rata penduduk bersekolah dan angka melek huruf dan mempertimbangkan kemampuan ekonomi masyarakat yang tercermin dari nilai purcashing power parity index (ppp). Indeks Pembangunan Manusia memberi informasi tentang tingkat pembangunan sumber daya manusia di suatu wilayah sehingga dapat digunakan sebagai tolok ukur angka kesejahteraan suatu wilayah. Indeks Pembangunan Manusia merupakan indikator gabungan dari tiga indikator, yaitu kesehatan, pendidikan, dan ekonomi. Pada penelitian ini akan dilakukan pendeteksian autokorelasi pada Indeks Pembangunan Manusia Provinsi Jawa Barat tahun 2010-2021.
Kemudian, diperlukan beberapa package dalam menjalankan regresi time series dan melakukan pendeteksian autokorelasi serta penanganannya. Adapun packages yang digunakan sebagai berikut :
readxl
dplyr
forecast
TTR
tseries
lmtest
lawstat
orcutt
HoRM
Memanggil Data
Data yang digunakan merupakan data IPM Provinsi Jawa Barat tahun 2010 sampai dengan 2021 yang diperoleh dari situs resmi Badan Pusat Statistik (BPS). Data terdiri atas dua peubah, yaitu IPM sebagai peubah respon (Y) dan Tahun sebagai peubah penjelas (X).
jabar <- read_excel("C:/Users/LENOVO/Documents/MK SEM 5/STA1341 - MPDW/Praktikum/3/Tugas Individu MPDW.xlsx", sheet = 2)
knitr::kable(jabar, align = "c")| TAHUN | IPM |
|---|---|
| 2010 | 66.15 |
| 2011 | 66.67 |
| 2012 | 67.32 |
| 2013 | 68.25 |
| 2014 | 68.80 |
| 2015 | 69.50 |
| 2016 | 70.05 |
| 2017 | 70.69 |
| 2018 | 71.30 |
| 2019 | 72.03 |
| 2020 | 72.09 |
| 2021 | 72.45 |
Eksplorasi Data
x <- jabar$TAHUN
y <- jabar$IPM
p <- ggplot(jabar, aes(x=TAHUN, y=IPM)) +
geom_line(lwd=1.2,col="#006600")
p +labs(x="Year",y = "Indeks Pembangunan",
title="Time Series Plot Indeks Pembangunan Manusia Jawa Barat ",
subtitle = "Periode 2010 - 2021",
caption = "sumber:bps.go.id")+
theme_bw()+
theme(
plot.title = element_text(size = 14L,
face = "bold",
hjust = 0.5),
plot.subtitle = element_text(size = 11L,
face = "plain",
hjust = 0.5)
)+geom_point(size=2) +geom_text(aes(label=paste(IPM,"%")),vjust=-0.8,size=3)
Berdasarkan time series plot di atas, dapat dilihat adanya hubungan positif antara periode waktu (tahun) dengan nilai IPM Jawa Barat yang ditunjukkan oleh titik-titik dan nilai pada plot yang menghasilkan hubungan linear positif. Dengan demikian, dapat disimpulkan bahwa semakin meningkatnya tahun maka nilai Indeks Pembangunan Manusia Provinsi Jawa Barat semakin meningkat pula.
Korelasi Peubah
cor(x,y)## [1] 0.9939809Pada peubah X dan Y didapatkan nilai korelasi sebesar 0.9939809 (mendekati 1) yang menunjukkan hubungan linier positif yang kuat pada kedua peubah.
Model Regresi
#model regresi
model <- lm(y~x, data = jabar)
summary(model)##
## Call:
## lm(formula = y ~ x, data = jabar)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.4760 -0.1888 0.1183 0.1785 0.3104
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.146e+03 4.237e+01 -27.05 1.10e-10 ***
## x 6.032e-01 2.102e-02 28.69 6.16e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2514 on 10 degrees of freedom
## Multiple R-squared: 0.988, Adjusted R-squared: 0.9868
## F-statistic: 823.2 on 1 and 10 DF, p-value: 6.159e-11Pada model regresi di atas dapat diketahui bahwa hasil uji-F memiliki p-value sebesar 6.159e-11 < α = 0.05 sehingga model tersebut layak untuk digunakan. Selain itu, uji-t menunjukkan bahwa seluruh peubah signifikan dalam model. Hal tersebut ditunjukkan oleh p-value yang bernilai sangat kecil. Koefisien determinasi bernilai sangat besar, yaitu 0.988. Artinya, keragaman nilai IPM yang dapat dijelaskan oleh peubah tahun sebesar 98.8%, sisanya dijelaskan oleh faktor lain yang tidak digunakan dalam penelitian. Dapat disimpulkan bahwa model regresi tersebut merupakan model yang baik.
Pendeteksi Autokorelasi
Autokorelasi merupakan kondisi di mana terdapat hubungan antara sisaan satu pengamatan dengan pengamatan lain pada model regresi. Model regresi yang baik merupakan model yang bebas dari autokorelasi.
Grafik
1. Residual Plot
#sisaan dan fitted value
resi1 <- residuals(model)
fit <- predict(model)
#Diagnostik dengan eksploratif
par(mfrow = c(2,2))
qqnorm(resi1)
qqline(resi1, col = "#006600", lwd = 2)
plot(fit, resi1, col = "#006600", pch = 20, ylab = "Sisaan",
xlab = "Fitted Values", main = "Sisaan vs Fitted Values")
abline(a = 0, b = 0, lwd = 2)
hist(resi1, col = "#006600")
plot(seq(1,12,1), resi1, col = "#006600", pch = 20,
ylab = "Sisaan", xlab = "Order", main = "Sisaan vs Order")
lines(seq(1,12,1), resi1, col = "#003366")
abline(a = 0, b = 0, lwd = 2)
Terlihat bahwa titik-titik pada Q-Q Plot tidak menyebar mengikuti garis lurus yang menandakan sisaan tidak menyebar normal. Hasil plot sisaan vs fitted value menunjukkan lebar pita sisaan relatif sama sehingga ragam sisaan homogen. Hasil plot sisaan vs order memperlihatkan tebaran tidak membentuk pola tertentu yang mengindikasikan terpenuhinya asumsi non-autokorelasi.
2. ACF dan PACF Plot
#ACF dan PACF identifikasi autokorelasi
par(mfrow = c(2,1))
acf(resi1)
pacf(resi1)
Berdasarkan plot ACF dan PACF, data tidak melewati batas garis berwarna biru sehingga dapat disimpulkan tidak adanya autokorelasi.
Deteksi Autokorelasi dengan Statistik Uji
1. Durbin Watson Test
H0: tidak ada autokorelasi || H1: terdapat autokorelasi
lmtest::dwtest(model, alternative = 'two.sided') #ada autokorelasi##
## Durbin-Watson test
##
## data: model
## DW = 0.79521, p-value = 0.003536
## alternative hypothesis: true autocorrelation is not 0p-value (0.003536) < 0.05 : Tolak H0 sehingga cukup bukti untuk menyatakan bahwa terdapat autokorelasi pada taraf 5%
2. Breusch-Godfrey Test
H0: tidak ada autokorelasi || H1: terdapat autokorelasi
lmtest::bgtest(y ~ x, data=jabar, order=1) #Perform Breusch-Godfrey test for first-order serial correlation##
## Breusch-Godfrey test for serial correlation of order up to 1
##
## data: y ~ x
## LM test = 3.9403, df = 1, p-value = 0.04714p-value (0.04714) < 0.05 : Tolak H0 sehingga cukup bukti untuk menyatakan bahwa terdapat autokorelasi pada taraf 5%
3. Run’s Test
H0: tidak ada autokorelasi || H1: terdapat autokorelasi
lawstat::runs.test(resid(model), alternative = 'two.sided')##
## Runs Test - Two sided
##
## data: resid(model)
## Standardized Runs Statistic = -1.2111, p-value = 0.2259p-value (0.2259) > 0.05 : Tak tolak H0 sehingga tidak cukup bukti untuk menyatakan bahwa terdapat autokorelasi pada taraf 5%
Penanganan Autokorelasi
Berdasarkan hasil uji formal yang telah dilakukan, dapat disimpulkan bahwa terdapat autokorelasi pada data IPM Provinsi Jawa Barat sehingga perlu dilakukan penanganan terhadap kondisi tersebut. Penanganan autokorelasi dilakukan dengan menggunakan dua metode, yaitu Cochrane-Orcutt dan Hildreth-Lu.
1. Cochrane-Orcutt
#Interactive method using to solve first order autocorrelation problems
modelco <- orcutt::cochrane.orcutt(model,convergence=6,max.iter = 1000)
modelco ## Cochrane-orcutt estimation for first order autocorrelation
##
## Call:
## lm(formula = y ~ x, data = jabar)
##
## number of interaction: 748
## rho 0.879734
##
## Durbin-Watson statistic
## (original): 0.79521 , p-value: 1.768e-03
## (transformed): 2.08182 , p-value: 4.03e-01
##
## coefficients:
## (Intercept) x
## -633.0361 0.3495Didapatkan nilai rho optimum sebesar 0.879734
#rho optimum
rho <- modelco$rho
y## [1] 66.15 66.67 67.32 68.25 68.80 69.50 70.05 70.69 71.30 72.03 72.09 72.45y[-1]## [1] 66.67 67.32 68.25 68.80 69.50 70.05 70.69 71.30 72.03 72.09 72.45y[-12]## [1] 66.15 66.67 67.32 68.25 68.80 69.50 70.05 70.69 71.30 72.03 72.09#transformasi terhadap y dan x
(y.trans <- y[-1]-y[-12]*rho)## [1] 8.475627 8.668166 9.026339 8.758187 8.974333 8.908520 9.064666 9.111637
## [9] 9.304999 8.722794 9.030010(x.trans <- x[-1]-x[-12]*rho)## [1] 242.7356 242.8559 242.9761 243.0964 243.2167 243.3369 243.4572 243.5775
## [9] 243.6977 243.8180 243.9383Model Baru Hasil Transformasi
#model baru
modelcorho <- lm(y.trans~x.trans)
summary(modelcorho)##
## Call:
## lm(formula = y.trans ~ x.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.35855 -0.08514 -0.00469 0.11189 0.26569
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -76.1330 38.7692 -1.964 0.0812 .
## x.trans 0.3495 0.1593 2.194 0.0559 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.201 on 9 degrees of freedom
## Multiple R-squared: 0.3484, Adjusted R-squared: 0.276
## F-statistic: 4.812 on 1 and 9 DF, p-value: 0.05592Deteksi Autokorelasi Cochrane-Orcutt
lmtest::dwtest(modelcorho,alternative = 'two.sided')##
## Durbin-Watson test
##
## data: modelcorho
## DW = 2.0818, p-value = 0.8059
## alternative hypothesis: true autocorrelation is not 0bgtest(modelcorho)##
## Breusch-Godfrey test for serial correlation of order up to 1
##
## data: modelcorho
## LM test = 0.17522, df = 1, p-value = 0.6755lawstat::runs.test(resid(modelcorho), alternative = 'two.sided')##
## Runs Test - Two sided
##
## data: resid(modelcorho)
## Standardized Runs Statistic = -0.93314, p-value = 0.3507H0: tidak ada autokorelasi || H1: terdapat autokorelasi
Berdasarkan ketiga uji formal Durbin-Watson Test, Breusch-Godfrey Test dan Run Test yang telah dilakukan, diperoleh p-value > 0.05 : Tak Tolak H0 sehingga belum cukup bukti untuk menyatakan bahwa terdapat autokorelasi pada taraf 5%
Transformasi Balik
cat("IPM = ", coef(modelcorho)[1]/(1-rho), "+", coef(modelcorho)[2]," Tahun", sep = "") #persamaan regresi setelah di transformasi ke persamaan awal## IPM = -633.0361+0.3494998 Tahun2. 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))
}
#mencari rho yang meminimumkan SSE (iteratif)
r <- c(seq(0.1,0.8, by= 0.1), seq(0.9,0.99, by= 0.01))
tab <- data.frame("rho" = r, "SSE" = sapply(r, function(i){deviance(hildreth.lu.func(i, model))}))
round(tab, 4) #0,90 memiliki SSE Terkecil## rho SSE
## 1 0.10 0.5525
## 2 0.20 0.5071
## 3 0.30 0.4680
## 4 0.40 0.4350
## 5 0.50 0.4083
## 6 0.60 0.3878
## 7 0.70 0.3735
## 8 0.80 0.3655
## 9 0.90 0.3636
## 10 0.91 0.3638
## 11 0.92 0.3640
## 12 0.93 0.3643
## 13 0.94 0.3646
## 14 0.95 0.3650
## 15 0.96 0.3655
## 16 0.97 0.3660
## 17 0.98 0.3666
## 18 0.99 0.3672#grafik rho dan SSE
plot(tab$SSE ~ tab$rho , type = "l")
abline(v = tab[tab$SSE==min(tab$SSE),"rho"], lty = 3)
Diperoleh nilai rho optimum berada di sekitar 0.90
#rho optimal di sekitar 0.9
r <- seq(0.8, 1, by= 0.01)
tab <- data.frame("rho" = r, "SSE" = sapply(r, function(i){deviance(hildreth.lu.func(i, model))}))
round(tab, 4) #0.87 memiliki SSE Terkecil## rho SSE
## 1 0.80 0.3655
## 2 0.81 0.3650
## 3 0.82 0.3646
## 4 0.83 0.3643
## 5 0.84 0.3640
## 6 0.85 0.3638
## 7 0.86 0.3636
## 8 0.87 0.3635
## 9 0.88 0.3635
## 10 0.89 0.3635
## 11 0.90 0.3636
## 12 0.91 0.3638
## 13 0.92 0.3640
## 14 0.93 0.3643
## 15 0.94 0.3646
## 16 0.95 0.3650
## 17 0.96 0.3655
## 18 0.97 0.3660
## 19 0.98 0.3666
## 20 0.99 0.3672
## 21 1.00 0.4924#grafik SSE optimum
plot(tab$SSE ~ tab$rho , type = "l")
abline(v = tab[tab$SSE==min(tab$SSE),"rho"], lty = 3)
Diperoleh nilai rho optimum adalah 0.87
Model Baru Hasil Transformasi
modelhl <- hildreth.lu.func(0.87, model)
summary(modelhl)##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.35696 -0.08791 -0.00323 0.11264 0.26630
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -87.8191 38.7588 -2.266 0.0497 *
## x 0.3704 0.1474 2.513 0.0331 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.201 on 9 degrees of freedom
## Multiple R-squared: 0.4124, Adjusted R-squared: 0.3471
## F-statistic: 6.316 on 1 and 9 DF, p-value: 0.03314Deteksi Autokorelasi Hildreth-Lu
lmtest::dwtest(modelhl)##
## Durbin-Watson test
##
## data: modelhl
## DW = 2.063, p-value = 0.3902
## alternative hypothesis: true autocorrelation is greater than 0bgtest(modelhl)##
## Breusch-Godfrey test for serial correlation of order up to 1
##
## data: modelhl
## LM test = 0.15368, df = 1, p-value = 0.695lawstat::runs.test(resid(modelhl), alternative = 'two.sided')##
## Runs Test - Two sided
##
## data: resid(modelhl)
## Standardized Runs Statistic = -0.93314, p-value = 0.3507H0: tidak ada autokorelasi || H1: terdapat autokorelasi
Berdasarkan ketiga uji formal Durbin-Watson Test, Breusch-Godfrey Test dan Run Test yang telah dilakukan, diperoleh p-value > 0.05 : Tak Tolak H0 sehingga belum cukup bukti untuk menyatakan bahwa terdapat autokorelasi pada taraf 5%
Transformasi Balik
cat("y = ", coef(modelhl)[1]/(1-0.87), "+", coef(modelhl)[2],"x", sep = "") #persamaan regresi setelah di transformasi ke persamaan awal## y = -675.5317+0.3704406xBerdasarkan dua metode yang telah digunakan, keduanya memiliki model regresi hasil transformasi balik cenderung sama. Hal ini mengisyaratkan bahwa kedua metode menghasilkan dugaan model regresi linear yang serupa, yaitu b0 = -675.53 dan b1 = 0.37. Dengan demikian, dapat disimpulkan bahwa kenaikan satu tahun meningkatkan dugaan rataan besar IPM Provinsi Jawa Barat sebesar 0.37.
Daftar Pustaka
Chatterjee S, Hadi AS. 2012. Regression Analysis by Example. New Jersey (US): John Wiley and Son Inc.
Prasetyoningrum AK, Sukmawati US. Analisis pengaruh indeks pembangunan manusia, pertumbuhan ekonomi dan pengagguran terhadap kemiskinan di Indonesia. Jurnal Ekonomi Syariah. 6(2): 217 - 240.