Tugas Individu STA1341 MPDW

Library yang digunakan

library(dplyr)
library(readxl)
library(forecast)
library(TTR)
library(tseries)
library(lmtest) #uji-Durbin Watson
library(orcutt) #Cochrane-Orcutt
library(HoRM)#Hildreth Lu

Input Data

#Membuka data IPM Jawa Barat Tahun 2010-2021
datareg <- read_excel("D:/IPB/Semester 5/STA1341 - Metode Peramalan Deret Waktu/Tugas Individu MPDW.xlsx", sheet = "Jabar")
datareg

## # A tibble: 12 x 2
##    Tahun IPMJabar
##    <dbl>    <dbl>
##  1  2010     66.2
##  2  2011     66.7
##  3  2012     67.3
##  4  2013     68.2
##  5  2014     68.8
##  6  2015     69.5
##  7  2016     70.0
##  8  2017     70.7
##  9  2018     71.3
## 10  2019     72.0
## 11  2020     72.1
## 12  2021     72.4

Eksplorasi Data

x <- datareg$Tahun
y <- datareg$IPMJabar

#diagram pencar identifikasi model
plot(x,y,pch = 20, col = "blue", main = "Scatter Plot X vs Y",
     ylab = "Nilai Peubah Y", xlab = "Nilai Peubah X")

Diagram di atas menunjukkan bahwa antara peubah x dan peubah y memiliki hubungan linear yang kuat, dapat dilihat berdasarkan titik-titik amatan yang hampir dapat membentuk garis lurus sempurna.

cor(x,y) #melihat korelasi x dan y

## [1] 0.9939809

Pernyataan sebelumnya dibuktikan dengan korelasi antar kedua peubah yang hampir mendekati nilai 1.

Model Regresi Deret Waktu

#model regresi
model <- lm(IPMJabar~Tahun, data = datareg)
summary(model)

## 
## Call:
## lm(formula = IPMJabar ~ Tahun, data = datareg)
## 
## 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 ***
## Tahun        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-11

Model regresi deret waktu nilai IPM Provinsi Jawa Barat dengan nilai IPM sebagai peubah respon (y) dan tahun sebagai peubah penjelas (x) memiliki persamaan model y = -1146 - 0.6032x.

Deteksi Autokorelasi

Secara Eksploratif

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 = "steelblue", lwd = 2)

plot(fit, resi1, col = "steelblue", pch = 20, xlab = "Sisaan", 
     ylab = "Fitted Values", main = "Sisaan vs Fitted Values")
abline(a = 0, b = 0, lwd = 2)

hist(resi1, col = "steelblue")

plot(seq(1,12,1), resi1, col = "steelblue", pch = 20, 
     xlab = "Sisaan", ylab = "Order", main = "Sisaan vs Order")
lines(seq(1,12,1), resi1, col = "red")
abline(a = 0, b = 0, lwd = 2)

Secara eksploratif, plot sisaan vs order menunjukkan adanya autokorelasi pada model regresi ini karena terdapat pola menaik pada ujung kanan dan kembali menurun di ujung kiri.

2. ACF dan PACF Plot

#ACF dan PACF identifikasi autokorelasi
acf(resi1)

pacf(resi1)

Plot ACF dan PACF tidak menunjukkan adanya autokorelasi karena tidak adanya garis yang garis hitam yang melewati garis biru (kecuali pada lag 0 yang pasti memiliki skor ACF sebesar 1).

Uji Statistik

Durbin-Watson Test

H0: Tidak ada bukti autokorelasi H1: Terdapat autokorelasi

lmtest::dwtest(model, alternative = 'two.sided')

## 
##  Durbin-Watson test
## 
## data:  model
## DW = 0.79521, p-value = 0.003536
## alternative hypothesis: true autocorrelation is not 0

Nilai p-value dari hasil uji Durbin-Watson kurang dari 0.05 (Tolak H0) atau cukup bukti untuk menunjukkan bahwa terdapat autokorelasi dalam model ini pada taraf 5%. Lalu, nilai DW menunjukkan angka sebesar 0.79521 yang berada di bawah nilai dL pada n = 12 (taraf 5%) yang sebesar 0.971, artinya terdapat autokorelasi positif pada model ini.

Penanganan Autokorelasi

1. Cochrane-Orcutt

# Penanganan Autokorelasi
## Cochrane-Orcutt
modelco <- orcutt::cochrane.orcutt(model, convergence = 6, max.iter = 1000)
modelco 

## Cochrane-orcutt estimation for first order autocorrelation 
##  
## Call:
## lm(formula = IPMJabar ~ Tahun, data = datareg)
## 
##  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)       Tahun 
##   -633.0361      0.3495

Menggunakan metode Cochrane-orcutt, sebanyak 748 iterasi dilakukan dan mendapatkan rho optimum sebesar 0.8797.

### 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.45

y[-1]

##  [1] 66.67 67.32 68.25 68.80 69.50 70.05 70.69 71.30 72.03 72.09 72.45

y[-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.9383

### 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.05592

Model baru yang didapatkan pada metode Cochrane-orcutt dengan persamaan y.trans = -76.133 + 0.3495 x.trans.

### DW Test Model Baru
lmtest::dwtest(modelcorho,alternative = 'two.sided')

## 
##  Durbin-Watson test
## 
## data:  modelcorho
## DW = 2.0818, p-value = 0.8059
## alternative hypothesis: true autocorrelation is not 0

Uji Durbin-Watson pada model yang autokorelasinya sudah ditangani dengan Cochrane-Orcutt menunjukkan kenaikan hasil p-value menjadi lebih besar dari taraf 5% (Tak tolak H0) dan nilai DW yang mendekati nilai 2, sehingga dapat dikatakan tidak terbukti ada autokorelasi.

2. 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, 5) 

##     rho     SSE
## 1  0.10 0.55246
## 2  0.20 0.50710
## 3  0.30 0.46796
## 4  0.40 0.43503
## 5  0.50 0.40832
## 6  0.60 0.38782
## 7  0.70 0.37353
## 8  0.80 0.36546
## 9  0.90 0.36361
## 10 0.91 0.36376
## 11 0.92 0.36398
## 12 0.93 0.36426
## 13 0.94 0.36460
## 14 0.95 0.36501
## 15 0.96 0.36547
## 16 0.97 0.36600
## 17 0.98 0.36659
## 18 0.99 0.36725

min(tab$SSE) #0.3636

## [1] 0.3636054

Didapatkan SSE minimum sebesar 0.36361 pada rho 0.90

#grafik rho dan SSE
plot(tab$SSE ~ tab$rho , type = "l")
abline(v = tab[tab$SSE==min(tab$SSE),"rho"], lty = 3)

Plot menunjukkan bahwa rho optimal berada di kisaran 0.85 sampai 0.95.

r <- seq(0.85, 0.95, by= 0.01)
tab <- data.frame("rho" = r, "SSE" = sapply(r, function(i){deviance(hildreth.lu.func(i, model))}))
round(tab, 5) #0.36348 memiliki SSE Terkecil

##     rho     SSE
## 1  0.85 0.36376
## 2  0.86 0.36360
## 3  0.87 0.36351
## 4  0.88 0.36348
## 5  0.89 0.36351
## 6  0.90 0.36361
## 7  0.91 0.36376
## 8  0.92 0.36398
## 9  0.93 0.36426
## 10 0.94 0.36460
## 11 0.95 0.36501

#grafik SSE optimum
plot(tab$SSE ~ tab$rho , type = "l")
abline(v = tab[tab$SSE==min(tab$SSE),"rho"], lty = 3)

Grafik di atas menunjukkan bahwa rho optimum berada di angka 0.88.

# Model terbaik
modelhl <- hildreth.lu.func(0.88, model)
summary(modelhl)

## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.35859 -0.08506 -0.00473  0.11187  0.26568 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept) -75.8130    38.7695  -1.955   0.0822 .
## x             0.3489     0.1597   2.185   0.0567 .
## ---
## 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.3466, Adjusted R-squared:  0.274 
## F-statistic: 4.774 on 1 and 9 DF,  p-value: 0.05672

# Deteksi autokorelasi
lmtest::dwtest(modelhl)

## 
##  Durbin-Watson test
## 
## data:  modelhl
## DW = 2.0823, p-value = 0.4033
## alternative hypothesis: true autocorrelation is greater than 0

Uji Durbin-Watson pada model yang autokorelasinya sudah ditangani dengan Hildreth-Lu menunjukkan kenaikan hasil p-value menjadi lebih besar dari taraf 5% (Tak tolak H0) dan nilai DW yang mendekati nilai 2, sehingga dapat dikatakan tidak terbukti ada autokorelasi.

Transformasi balik dilakukan untuk mendapatkan persamaan model baru dengan persamaan sebagai berikut.

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

## y = -631.7754+0.3488788x