Penanganan Autokorelasi pada Model Regresi Sederhana Menggunakan Data IPM Provinsi DKI Jakarta

Persiapan Data

Install Library

library("readxl")
library("car")

## Loading required package: carData

library("lmtest")

## Loading required package: zoo

## 
## Attaching package: 'zoo'

## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

library("orcutt")
library("HoRM")

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

Membaca Data

df <- read_excel("D:/IPM DKI Jakarta.xlsx")
df

## # A tibble: 12 x 2
##    Tahun   IPM
##    <dbl> <dbl>
##  1  2010  76.3
##  2  2011  77.0
##  3  2012  77.5
##  4  2013  78.1
##  5  2014  78.4
##  6  2015  79.0
##  7  2016  79.6
##  8  2017  80.1
##  9  2018  80.5
## 10  2019  80.8
## 11  2020  80.8
## 12  2021  81.1

Data tersebut terdiri dari kolom Tahun dan kolom IPM. Rentang tahun yang digunakan mulai dari tahun 2010 hingga tahun 2019.

Analisis Eksplorasi Data

Scater Plot IPM

plot(df$Tahun,df$IPM, main = "Scatter Plot IPM Provinsi DKI Jakarta Tahun 2010-2019",xlab= "Tahun", ylab = "IPM")

Koefesien Korelasi Antara Tahun dengan IPM

cor(df)

##           Tahun       IPM
## Tahun 1.0000000 0.9879386
## IPM   0.9879386 1.0000000

Dari grafik scatter plot diatas dapat kita lihat bahwa IPM memiliki hubungan linear dengan Tahun. Hal tersebut dapat dibuktikan oleh nilai korelasi yang dihasilkan kedua peubah tersebut yaitu sebesar 0.9879 .

Model Awal Regrei Linear Sederhana

Pemodelan Regresi Linear Sederhana

model <- lm(IPM ~ Tahun, data = df)
summary(model)

## 
## Call:
## lm(formula = IPM ~ Tahun, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.42154 -0.16017  0.05061  0.16141  0.30594 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -816.54150   44.39144  -18.39 4.86e-09 ***
## Tahun          0.44437    0.02202   20.18 1.97e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2634 on 10 degrees of freedom
## Multiple R-squared:  0.976,  Adjusted R-squared:  0.9736 
## F-statistic: 407.1 on 1 and 10 DF,  p-value: 1.97e-09

Persamaan model regresi awal yaitu: \(y = -816.54150 + 0.44437(Tahun)\) Didapatkan \(R^2\) sebesar 0.97 yang artinya hanya 97% keragaman dari IPM bisa dijelaskan oleh peubah Tahun.

Uji Asumsi Gauss-Markov Model Awal

Berikut ini merupakan asumsi Gauss-Markov yang harus dipenuhi ketika menggunakan model regresi OLS:

• Residual menyebar normal
• Ragam homogen
• Residual saling bebas (tidak ada autokorelasi)

Uji normalitas

shapiro.test(model$residuals)

## 
##  Shapiro-Wilk normality test
## 
## data:  model$residuals
## W = 0.91712, p-value = 0.263

Hipotesis

H0 : Residual berdistribusi normal

H1 : Residual tidak berdistribusi normal

• Keputusan : Karena p-Value(0.263) > alpha(0.05), maka terima H0
• Kesimpulan : Dengan tingkat keyakinan 95%, kita yakin bahwa residual berdistribusi normal

Uji Homoskedastisitas

bptest(model)

## 
##  studentized Breusch-Pagan test
## 
## data:  model
## BP = 3.0615, df = 1, p-value = 0.08017

Hipotesis

H0 : Ragam homogen

H1 : Ragam tidak homogen

• Keputusan : Karena p-Value(0.08017) > alpha(0.05), maka terima H0
• Kesimpulan : Dengan tingkat keyakinan 95%, kita yakin bahwa ragam homogen

Uji Autokorelasi

residual <- residuals(model)

par(mfrow=c(1,2))
acf(residual)
pacf(residual )

dwtest(model)

## 
##  Durbin-Watson test
## 
## data:  model
## DW = 0.53015, p-value = 6.42e-05
## alternative hypothesis: true autocorrelation is greater than 0

Hipotesis

H0 : Tidak ada autokorelasi

H1 : Ada Autokoerelasi

• Keputusan : Karena p-Value(6.42e-05) < alpha(0.05), maka tolak H0
• Kesimpulan : Dengan tingkat keyakinan 95%, kita yakin bahwa ada autokorelasi

Penanganan Gejala Autokorelasi

Mencari Nilai Rho Terbaik Metode 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))
}


r <- c(seq(0.85,0.99, by= 0.01), 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)

##     rho    SSE
## 1  0.85 0.1775
## 2  0.86 0.1771
## 3  0.87 0.1768
## 4  0.88 0.1765
## 5  0.89 0.1764
## 6  0.90 0.1763
## 7  0.91 0.1763
## 8  0.92 0.1764
## 9  0.93 0.1766
## 10 0.94 0.1769
## 11 0.95 0.1773
## 12 0.96 0.1778
## 13 0.97 0.1783
## 14 0.98 0.1790
## 15 0.99 0.1797
## 16 0.90 0.1763
## 17 0.91 0.1763
## 18 0.92 0.1764
## 19 0.93 0.1766
## 20 0.94 0.1769
## 21 0.95 0.1773
## 22 0.96 0.1778
## 23 0.97 0.1783
## 24 0.98 0.1790
## 25 0.99 0.1797

Dapat dilihat dari tabel diatas didapatkan nilai rho dengan error terkecil yaitu rho sebesar 0.90

Model Baru dengan Metode Hildreth-Lu (rho=0.90)

model_hl <- hildreth.lu(df$IPM, df$Tahun, 0.90)
summary(model_hl)

## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.262945 -0.031786  0.002045  0.082941  0.182273 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.92123   27.02291  -0.108    0.916
## x            0.05555    0.13345   0.416    0.687
## 
## Residual standard error: 0.14 on 9 degrees of freedom
## Multiple R-squared:  0.01889,    Adjusted R-squared:  -0.09013 
## F-statistic: 0.1733 on 1 and 9 DF,  p-value: 0.687

Penggunaan rho sebesar 0.90 menghasilkan model regresi dengan p-value > alpha(0.05) sehingga model tersebut kurang significant. Maka dilakukan pemodelan lain menggunakan rho sebesar 0.85

Model Baru dengan Metode Hildreth-Lu (rho=0.85)

model_hl <- hildreth.lu(df$IPM, df$Tahun, 0.85)
summary(model_hl)

## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.263600 -0.035007  0.001823  0.094207  0.186591 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)  
## (Intercept) -46.20020   27.07096  -1.707   0.1221  
## x             0.19282    0.08927   2.160   0.0591 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1404 on 9 degrees of freedom
## Multiple R-squared:  0.3414, Adjusted R-squared:  0.2682 
## F-statistic: 4.665 on 1 and 9 DF,  p-value: 0.05907

Model baru yang dihasilkan sudah significant karena memiliki nilai p-value < alpha(0.05)

Uji Asumsi Gauss-Markov Model Baru

Uji normalitas

shapiro.test(model_hl$residuals)

## 
##  Shapiro-Wilk normality test
## 
## data:  model_hl$residuals
## W = 0.93467, p-value = 0.46

Hipotesis

H0 : Residual berdistribusi normal

H1 : Residual tidak berdistribusi normal

• Keputusan : Karena p-Value(0.46) > alpha(0.05), maka terima H0
• Kesimpulan : Dengan tingkat keyakinan 95%, kita yakin bahwa residual berdistribusi normal

Uji Homoskedastisitas

bptest(model_hl)

## 
##  studentized Breusch-Pagan test
## 
## data:  model_hl
## BP = 0.89753, df = 1, p-value = 0.3434

Hipotesis

H0 : Ragam homogen

H1 : Ragam tidak homogen

• Keputusan : Karena p-Value(0.3434) > alpha(0.05), maka terima H0
• Kesimpulan : Dengan tingkat keyakinan 95%, kita yakin bahwa ragam homogen

Uji Autokorelasi

residual <- residuals(model_hl)

par(mfrow=c(1,2))
acf(residual)
pacf(residual)

dwtest(model_hl)

## 
##  Durbin-Watson test
## 
## data:  model_hl
## DW = 1.7859, p-value = 0.2203
## alternative hypothesis: true autocorrelation is greater than 0

Hipotesis

H0 : Tidak ada autokorelasi

H1 : Ada Autokoerelasi

• Keputusan : Karena p-Value(0.2203) > alpha(0.05), maka terima H0
• Kesimpulan : Dengan tingkat keyakinan 95%, kita yakin bahwa tidak ada autokorelasi

Kesimpulan

Didapatkan persamaan regresi setelah penanganan gejala autokorelasi yaitu:

cat("y = ", coef(model_hl)[1]/(1-0.85), "+", coef(model_hl)[2],"x", sep = "")

## y = -308.0014+0.1928182x