Data dan Package

Package

library(dplyr)
## Warning: package 'dplyr' was built under R version 4.1.2
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(readxl)
## Warning: package 'readxl' was built under R version 4.1.2
library(forecast)
## Warning: package 'forecast' was built under R version 4.1.3
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(TTR)
## Warning: package 'TTR' was built under R version 4.1.3
library(tseries)
## Warning: package 'tseries' was built under R version 4.1.3
library(lmtest) 
## Warning: package 'lmtest' was built under R version 4.1.2
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.1.2
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(orcutt) 
## Warning: package 'orcutt' was built under R version 4.1.3
library(HoRM)
## Warning: package 'HoRM' was built under R version 4.1.3
library(lawstat)
## Warning: package 'lawstat' was built under R version 4.1.3
## 
## Attaching package: 'lawstat'
## The following object is masked from 'package:tseries':
## 
##     runs.test
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.1.2

Deklarasi Data IPM Sumatera Barat 2010-2021

IPM <- readxl::read_xlsx("C:/Users/HP/Documents/COLLAGE/SMT 5/Metode Peramalan Deret Waktu/Tugas Individu MPDW.xlsx", sheet=2)
knitr::kable(IPM,align = "c")
Tahun IPM
2010 67.25
2011 67.81
2012 68.36
2013 68.91
2014 69.36
2015 69.98
2016 70.73
2017 71.24
2018 71.73
2019 72.39
2020 72.38
2021 72.65

Indeks Pembangunan Manusia

Indeks Pembangunan Manusia (IPM) merupakan suatu indeks komposit yang mencakup tiga bidang pembangunan manusia yang mendasar yaitu usia hidup, pengetahuan, dan standar hidup layak. IPM digunakan untuk mengklarifikasikan apakah sebuah negara adalah negara maju, negara berkembang atau negara terbelakang dan juga untuk mengukur pengaruh dari kebijaksanaan ekonomi terhadap kualitas hidup.

Eksplorasi Data

Nilai Korelasi

x<-IPM$Tahun
y<-IPM$IPM
cor(x,y)
## [1] 0.9932977

Plot Deret Waktu

a <- ggplot(IPM, aes(x=Tahun, y=IPM)) +
  geom_line(lwd=1.2,col="orange")
a +labs(x="Tahun",y = "IPM",
         title="Plot Deret Waktu IPM Sumatera Barat",
        subtitle = "Periode 2010 - 2021")+
  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)

Model Regresi

model <- lm(y~x, data = IPM)
summary(model)
## 
## Call:
## lm(formula = y ~ x, data = IPM)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.45231 -0.09554 -0.03215  0.20099  0.33126 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -981.4216    38.6982  -25.36 2.08e-10 ***
## x              0.5218     0.0192   27.18 1.05e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2296 on 10 degrees of freedom
## Multiple R-squared:  0.9866, Adjusted R-squared:  0.9853 
## F-statistic: 738.5 on 1 and 10 DF,  p-value: 1.053e-10

Deteksi Autokorelasi

1. Residual Plot

#sisaan dan fitted value
res <- residuals(model)
fit <- predict(model)

#Diagnostik dengan eksploratif
par(mfrow = c(2,2))

qqnorm(res)
qqline(res, col = "Orange", lwd = 2)

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

hist(res, col = "#FFB58A")

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

Plot Sisaan vs Order terlihat membentuk suatu pola tertentu, maka asumsi kebebasan tidak terpenuhi. Berdasarkan hasil uji ekspolrasi, terdeteksi adanya gejala autokolerasi.

2. ACF dan PACF Plot

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

Terlihat bahwa plot ACF dan PACF tidak ada yang signifikan. Artinya, sisaan saling bebas.

Uji Statistik

1. Durbin Watson Test

H0: Tidak ada autokorelasi H1: Ada autokorelasi

lmtest::dwtest(model, alternative = 'two.sided') 
## 
##  Durbin-Watson test
## 
## data:  model
## DW = 0.82777, p-value = 0.004736
## alternative hypothesis: true autocorrelation is not 0

P-value < 5% : Tolak H0 sehingga cukup bukti untuk menyatakan bahwa terdapat autokorelasi (positif) pada taraf 5%

2. Breusch-Godfrey Test

H0: Tidak ada autokorelasi H1: Ada autokorelasi

lmtest::bgtest(y ~ x, data=IPM, order=1) 
## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  y ~ x
## LM test = 3.7277, df = 1, p-value = 0.05352

P-value > 5% : Tak tolak H0 sehingga tidak cukup bukti untuk menyatakan bahwa ada autokorelasi pada taraf 5% (tidak terdapat autokorelasi)

Dikarenakan pada uji Durbin-Watson terdeteksi adanya Autokorelasi, sehingga perlu adanya penanganan menggunakan Cochran Orcutt dan Hildreth-lu.

Penanganan Aurokorelasi

1. Cochrane-Orcutt

modelco <- orcutt::cochrane.orcutt(model)
## Warning in orcutt::cochrane.orcutt(model): Did not converge
modelco
## Cochrane-orcutt estimation for first order autocorrelation 
##  
## Call:
## lm(formula = y ~ x, data = IPM)
## 
##  number of interaction: 100
##  rho 0.871476
## 
## Durbin-Watson statistic 
## (original):    0.82777 , p-value: 2.368e-03
## (transformed): NA , p-value: NA
##  
##  coefficients: 
## [1] NA
rho <- modelco$rho
rho
## [1] 0.8714765
y
##  [1] 67.25 67.81 68.36 68.91 69.36 69.98 70.73 71.24 71.73 72.39 72.38 72.65
y[-1]
##  [1] 67.81 68.36 68.91 69.36 69.98 70.73 71.24 71.73 72.39 72.38 72.65
y[-12]
##  [1] 67.25 67.81 68.36 68.91 69.36 69.98 70.73 71.24 71.73 72.39 72.38
(y.trans <- y[-1]-y[-12]*rho)
##  [1] 9.203205 9.265179 9.335866 9.306554 9.534390 9.744075 9.600467 9.646014
##  [9] 9.878991 9.293816 9.572531
(x.trans <- x[-1]-x[-12]*rho)
##  [1] 259.3322 259.4608 259.5893 259.7178 259.8463 259.9749 260.1034 260.2319
##  [9] 260.3604 260.4889 260.6175
#model baru
modelcorho <- lm(y.trans~x.trans)
summary(modelcorho)
## 
## Call:
## lm(formula = y.trans ~ x.trans)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.35303 -0.09636 -0.03508  0.08131  0.27156 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept) -70.2368    36.3099  -1.934   0.0851 .
## x.trans       0.3067     0.1397   2.196   0.0557 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1883 on 9 degrees of freedom
## Multiple R-squared:  0.3488, Adjusted R-squared:  0.2765 
## F-statistic: 4.821 on 1 and 9 DF,  p-value: 0.05573
lmtest::dwtest(modelcorho,alternative = 'two.sided') 
## 
##  Durbin-Watson test
## 
## data:  modelcorho
## DW = 1.8467, p-value = 0.5078
## alternative hypothesis: true autocorrelation is not 0
bgtest(modelcorho)
## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  modelcorho
## LM test = 0.022598, df = 1, p-value = 0.8805

P-value > 0.05, maka tak tolak H0. Artinya, tidak cukup bukti untuk menyatakan bahwa terdapat autokorelasi pada taraf 5% (Tidak terdapat autokorelasi). Disimpulkan bahwa hasil transformasi sudah memenuhi asumsi bebas autokorelasi.

Transformasi Balik

cat("IPM = ", coef(modelcorho)[1]/(1-rho), "+", coef(modelcorho)[2]," Tahun", sep = "")
## IPM = -546.4902+0.3066682 Tahun

2. 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.999, by= 0.001))
tab <- data.frame("rho" = r, "SSE" = sapply(r, function(i){deviance(hildreth.lu.func(i, model))}))
tab$rho[which.min(tab$SSE)]
## [1] 0.896
#grafik rho dan SSE
plot(tab$SSE ~ tab$rho , type = "l")
abline(v = tab[tab$SSE==min(tab$SSE),"rho"], lty = 3)

#Model terbaik
modelhl <- hildreth.lu.func(0.896, model)
summary(modelhl)
## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.35773 -0.09587 -0.03475  0.08096  0.26965 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept) -44.8943    36.3347  -1.236    0.248
## x             0.2501     0.1726   1.449    0.181
## 
## Residual standard error: 0.1882 on 9 degrees of freedom
## Multiple R-squared:  0.1893, Adjusted R-squared:  0.09918 
## F-statistic: 2.101 on 1 and 9 DF,  p-value: 0.1811
lmtest::dwtest(modelhl,alternative = 'two.sided') 
## 
##  Durbin-Watson test
## 
## data:  modelhl
## DW = 1.8879, p-value = 0.556
## alternative hypothesis: true autocorrelation is not 0
bgtest(modelhl)
## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  modelhl
## LM test = 0.0082115, df = 1, p-value = 0.9278

P-value > 0.05, maka tak tolak H0. Artinya, tidak cukup bukti untuk menyatakan bahwa terdapat autokorelasi pada taraf 5% (Tidak terdapat autokorelasi). Disimpulkan bahwa hasil transformasi sudah memenuhi asumsi bebas autokorelasi.

Transformasi Balik

cat("IPM = ", coef(modelhl)[1]/(1-0.896), "+", coef(modelhl)[2]," Tahun", sep = "") #persamaan regresi setelah di transformasi ke persamaan awal
## IPM = -431.6757+0.2501259 Tahun