library(dLagM)
## Warning: package 'dLagM' was built under R version 4.3.3
## Loading required package: nardl
## Warning: package 'nardl' was built under R version 4.3.3
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
## Loading required package: dynlm
## Warning: package 'dynlm' was built under R version 4.3.3
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.3.2
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(dynlm)
library(MLmetrics)
## Warning: package 'MLmetrics' was built under R version 4.3.3
##
## Attaching package: 'MLmetrics'
## The following object is masked from 'package:dLagM':
##
## MAPE
## The following object is masked from 'package:base':
##
## Recall
library(lmtest)
## Warning: package 'lmtest' was built under R version 4.3.2
library(car)
## Warning: package 'car' was built under R version 4.3.3
## Loading required package: carData
Data ini terdiri atas angka harapan hidup sebagai peubah dependen (Y) dan supply pangan sebagai peubah independen (X) Indonesia periode 1970 - 2019.
library(readxl)
## Warning: package 'readxl' was built under R version 4.3.2
data <- read_excel("D:/Sem 5/MPDW/Tugas/Individu/3/AngkaHarapanHidup~SupplyPangan Indonesia.xlsx")
str(data)
## tibble [50 × 4] (S3: tbl_df/tbl/data.frame)
## $ t : num [1:50] 1 2 3 4 5 6 7 8 9 10 ...
## $ Yt : num [1:50] 53 53.6 54.2 54.8 55.4 ...
## $ Y(t-1): num [1:50] NA 53 53.6 54.2 54.8 ...
## $ Xt : num [1:50] 1.15e+08 1.18e+08 1.22e+08 1.25e+08 1.28e+08 ...
data
## # A tibble: 50 × 4
## t Yt `Y(t-1)` Xt
## <dbl> <dbl> <dbl> <dbl>
## 1 1 53.0 NA 115228390
## 2 2 53.6 53.0 118347140
## 3 3 54.2 53.6 121504136
## 4 4 54.8 54.2 124709060
## 5 5 55.4 54.8 127945190
## 6 6 56.0 55.4 131213220
## 7 7 56.5 56.0 134521020
## 8 8 57.1 56.5 137861540
## 9 9 57.6 57.1 141250960
## 10 10 58.2 57.6 144693090
## # ℹ 40 more rows
# Split Data
training <- data[1:40,]
testing <- data[41:50,]
# Data Time Series
training.ts <- ts(training)
testing.ts <- ts(testing)
data.ts <- ts(data)
Model Koyck didasarkan pada asumsi bahwa semakin jauh jarak lag peubah independen dari periode sekarang maka semakin kecil pengaruh peubah lag terhadap peubah dependen.
Koyck mengusulkan suatu metode untuk menduga model dinamis distributed lag dengan mengasumsikan bahwa semua koefisien \(\beta\) mempunyai tanda sama.
Model kyock merupakan jenis paling umum dari model infinite distributed lag dan juga dikenal sebagai geometric lag
\[ y_t=a(1-\lambda)+\beta_0X_t+\beta_1Z_t+\lambda Y_{t-1}+V_t \]
dengan \[V_t=u_t-\lambda u_{t-1}\]
Pemodelan model Koyck dengan R dapat menggunakan
dLagM::koyckDlm() . Fungsi umum dari koyckDlm
adalah sebagai berikut.
koyckDlm(x , y , intercept)
Fungsi koyckDlm() akan menerapkan model lag
terdistribusi dengan transformasi Koyck satu prediktor. Nilai
x dan y tidak perlu sebagai objek time
series (ts). intercept dapat dibuat
TRUE untuk memasukkan intersep ke dalam model.
# Model Koyck
modelkoyck <- koyckDlm(x = training$Xt, y = training$Yt)
summary(modelkoyck)
##
## Call:
## "Y ~ (Intercept) + Y.1 + X.t"
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.92147 -0.07364 0.01528 0.11772 1.28354
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.294e+00 3.791e+00 2.188 0.0352 *
## Y.1 8.180e-01 9.491e-02 8.619 2.81e-10 ***
## X.t 1.873e-08 1.180e-08 1.587 0.1212
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4114 on 36 degrees of freedom
## Multiple R-Squared: 0.992, Adjusted R-squared: 0.9916
## Wald test: 2237 on 2 and 36 DF, p-value: < 2.2e-16
##
## Diagnostic tests:
## NULL
##
## alpha beta phi
## Geometric coefficients: 45.57815 1.873248e-08 0.8180297
AIC(modelkoyck)
## [1] 46.27588
BIC(modelkoyck)
## [1] 52.93012
Dari hasil tersebut, didapat bahwa peubah \(y_{t-1}\) memiliki nilai \(P-Value<0.05\). Hal ini menunjukkan bahwa peubah \(y_{t-1}\) berpengaruh signifikan terhadap \(y\) sementara \(x_t\) tidak berpengaruh signifikan terhadap \(y_t\). Artinya, angka harapan hidup satu periode sebelumnya berpengaruh signifikan terhadap angka harapan hidup periode sekarang. Adapun model keseluruhannya adalah sebagai berikut
\[ \hat{Y_t}=8.294+0.00000001873X_t+0.818Y_{t-1} \]
Berikut adalah hasil peramalan y untuk 10 periode kedepan menggunakan model koyck
forekoyck <- forecast(model = modelkoyck, x = training$Xt, h = 10)
forekoyck
## $forecasts
## [1] 66.48218 64.89520 63.65615 62.70260 61.98319 61.45591 61.08654 60.84696
## [9] 60.71447 60.67057
##
## $call
## forecast.koyckDlm(model = modelkoyck, x = training$Xt, h = 10)
##
## attr(,"class")
## [1] "forecast.koyckDlm" "dLagM"
# Akurasi Data Testing
mapekoyck <- MAPE(forekoyck$forecasts, testing$Yt)
mapekoyck
## [1] 0.1018948
# Akurasi Data Training
GoF(modelkoyck)
## n MAE MPE MAPE sMAPE MASE
## modelkoyck 39 0.1942202 -4.634154e-05 0.003058909 0.003053095 0.404581
## MSE MRAE GMRAE
## modelkoyck 0.1562271 0.5061085 0.250326
Pemodelan model Regression with Distributed Lag dengan R
dapat menggunakan dLagM::dlm() . Fungsi umum dari
dlm adalah sebagai berikut.
dlm(formula , data , x , y , q , remove )
Fungsi dlm() akan menerapkan model lag terdistribusi
dengan satu atau lebih prediktor. Nilai x dan
y tidak perlu sebagai objek time series
(ts). \(q\) adalah integer
yang mewakili panjang lag yang terbatas.
finiteDLMauto(formula = Yt ~ Xt,
data = data.frame(training), q.min = 1, q.max = 10,
model.typ = "dlm", error.type = "AIC", trace = FALSE)
## q - k MASE AIC BIC GMRAE MBRAE R.Adj.Sq Ljung-Box
## 10 10 0.50952 37.61211 55.82768 0.41145 -0.1325 0.98336 0.2917478
Berdasarkan output tersebut, lag optimum didapatkan ketika lag = 10. Selanjutnya dilakukan pemodelan untuk lag = 10.
modeldlm <- dlm(x = training$Xt, y = training$Yt, q = 10)
summary(modeldlm)
##
## Call:
## lm(formula = model.formula, data = design)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.77888 -0.19534 0.01772 0.20152 0.66379
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.119e+01 8.122e+00 2.609 0.0178 *
## x.t 5.366e-06 2.879e-06 1.864 0.0788 .
## x.1 -1.429e-05 7.888e-06 -1.811 0.0868 .
## x.2 1.861e-05 1.157e-05 1.609 0.1250
## x.3 -1.553e-05 1.347e-05 -1.153 0.2639
## x.4 1.173e-05 1.436e-05 0.817 0.4248
## x.5 -6.058e-06 1.497e-05 -0.405 0.6906
## x.6 -2.646e-06 1.613e-05 -0.164 0.8715
## x.7 5.389e-06 1.812e-05 0.297 0.7696
## x.8 -8.511e-06 1.738e-05 -0.490 0.6303
## x.9 1.300e-05 1.141e-05 1.139 0.2696
## x.10 -6.942e-06 3.500e-06 -1.984 0.0627 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3791 on 18 degrees of freedom
## Multiple R-squared: 0.9897, Adjusted R-squared: 0.9834
## F-statistic: 156.8 on 11 and 18 DF, p-value: 1.563e-15
##
## AIC and BIC values for the model:
## AIC BIC
## 1 37.61211 55.82768
AIC(modeldlm)
## [1] 37.61211
BIC(modeldlm)
## [1] 55.82768
Dari hasil tersebut peubah \(x\) secara simultan berpengaruh terhadap peubah \(y_t\) pada taraf nyata 5%. Namun, tidak terdapat peubah yang berpengaruh signifikan secara parsial terhadap taraf nyata 5% terhadap \(y_t\). Adapun keseluruhan model yang terbentuk adalah
\[ \hat{Y_t}=21.119+0.000005366X_t-...-0.000006942X_{t-10} \]
Berikut merupakan hasil peramalan \(y\) untuk 10 periode kedepan.
foredlm <- forecast(model = modeldlm, x = training$Xt, h = 10)
foredlm
## $forecasts
## [1] -621.27843 1219.30294 -1177.75095 822.39679 -688.72014 93.04647
## [7] 434.71199 -258.68342 838.07115 -836.62050
##
## $call
## forecast.dlm(model = modeldlm, x = training$Xt, h = 10)
##
## attr(,"class")
## [1] "forecast.dlm" "dLagM"
# Akurasi Data Testing
mapedlm <- MAPE(foredlm$forecasts, testing$Yt)
mapedlm
## [1] 10.06573
# Akurasi Data Training
GoF(modeldlm)
## n MAE MPE MAPE sMAPE MASE MSE
## modeldlm 30 0.2295306 -2.175937e-05 0.00354587 0.003545361 0.5095214 0.086224
## MRAE GMRAE
## modeldlm 1.0362 0.4114548
Peubah dependen dipengaruhi oleh peubah independen pada waktu sekarang, serta dipengaruhi juga oleh peubah dependen itu sendiri pada satu waktu yang lalu maka model tersebut disebut autoregressive (Gujarati 2004).
Pemodelan Autoregressive dilakukan menggunakan fungsi
dLagM::ardlDlm() . Fungsi tersebut akan menerapkan
autoregressive berordo \((p,q)\) dengan satu prediktor. Fungsi umum
dari ardlDlm() adalah sebagai berikut.
ardlDlm(formula = NULL , data = NULL , x = NULL , y = NULL , p = 1 , q = 1 ,
remove = NULL )
Dengan \(p\) adalah integer yang mewakili panjang lag yang terbatas dan \(q\) adalah integer yang merepresentasikan ordo dari proses autoregressive.
modelardl <- ardlBoundOrders(data = data.frame(data), ic = "AIC",
formula = Yt ~ Xt )
min_p=c()
for(i in 1:5){
min_p[i]=min(modelardl$Stat.table[[i]])
}
q_opt=which(min_p==min(min_p, na.rm = TRUE))
p_opt=which(modelardl$Stat.table[[q_opt]] ==
min(modelardl$Stat.table[[q_opt]], na.rm = TRUE))
data.frame("q_optimum" = q_opt, "p_optimum" = p_opt,
"AIC"=modelardl$min.Stat)
## q_optimum p_optimum AIC
## 1 1 4 -75.75334
Dari tabel di atas, dapat terlihat bahwa nilai AIC terendah didapat
ketika \(p=4\) dan \(q=1\), yaitu sebesar
-75,75334. Artinya, model autoregressive optimum didapat
ketika \(p=4\) dan \(q=1\).
Selanjutnya dapat dilakukan pemodelan dengan nilai \(p\) dan \(q\) optimum.
modelardl <- ardlDlm(formula = Yt ~ Xt, data = training, p = 4, q = 1)
summary(modelardl)
##
## Time series regression with "ts" data:
## Start = 5, End = 40
##
## Call:
## dynlm(formula = as.formula(model.text), data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.24409 -0.16700 -0.01354 0.22831 0.55611
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.074e+00 4.886e+00 1.652 0.109269
## Xt.t 6.805e-06 2.000e-06 3.402 0.001969 **
## Xt.1 -1.883e-05 5.960e-06 -3.159 0.003687 **
## Xt.2 2.038e-05 8.080e-06 2.523 0.017377 *
## Xt.3 -9.657e-06 6.203e-06 -1.557 0.130373
## Xt.4 1.337e-06 2.138e-06 0.625 0.536753
## Yt.1 6.631e-01 1.493e-01 4.442 0.000119 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3585 on 29 degrees of freedom
## Multiple R-squared: 0.9931, Adjusted R-squared: 0.9917
## F-statistic: 700.2 on 6 and 29 DF, p-value: < 2.2e-16
AIC(modelardl)
## [1] 36.52248
BIC(modelardl)
## [1] 49.19063
Dari hasil tersebut, didapat bahwa peubah \(x_t\), \(x_{t-1}\), \(x_{t-2}\), dan \(y_{t-1}\) memiliki nilai \(P-Value<0.05\). Hal ini menunjukkan bahwa peubah \(x_t\), \(x_{t-1}\), \(x_{t-2}\), dan \(y_{t-1}\) berpengaruh signifikan terhadap \(y_t\) sementara \(x_{t-3}\) dan \(x_{t-4}\) tidak berpengaruh signifikan terhadap \(y_t\). Artinya, supply pangan pada periode sekarang, satu periode sebelumnya, dua periode sebelumnya, dan angka harapan hidup satu periode sebelumnya berpengaruh signifikan terhadap angka harapan hidup periode sekarang. Adapun model keseluruhannya adalah sebagai berikut
\[ \hat{Y_t}=8.074+0.000006805X_t-...+0.6631Y_{t-1} \]
foreardl <- forecast(model = modelardl, x = training$Xt, h = 10)
foreardl
## $forecasts
## [1] -807.058041 1037.532878 -365.345969 -51.190781 -15.039224 9.260458
## [7] 25.595721 36.536743 44.162002 49.352246
##
## $call
## forecast.ardlDlm(model = modelardl, x = training$Xt, h = 10)
##
## attr(,"class")
## [1] "forecast.ardlDlm" "dLagM"
# Akurasi Data Testing
mapeardl <- MAPE(foreardl$forecasts, testing$Yt)
mapeardl
## [1] 3.872911
# Akurasi Data Training
GoF(modelardl)
## n MAE MPE MAPE sMAPE MASE
## modelardl 36 0.2320023 -2.698237e-05 0.003615361 0.003612415 0.4954591
## MSE MRAE GMRAE
## modelardl 0.1035389 0.7500387 0.3739476
dynlmPemodelan regresi dengan peubah lag tidak hanya dapat
dilakukan dengan fungsi pada packages dLagM ,
tetapi terdapat packages dynlm yang dapat
digunakan. Fungsi dynlm secara umum adalah sebagai
berikut.
dynlm(formula, data, subset, weights, na.action, method = "qr",
model = TRUE, x = FALSE, y = FALSE, qr = TRUE, singular.ok = TRUE,
contrasts = NULL, offset, start = NULL, end = NULL, ...)
Untuk menentukan formula model yang akan digunakan,
tersedia fungsi tambahan yang memungkinkan spesifikasi dinamika (melalui
d() dan L()) atau pola linier/siklus dengan
mudah (melalui trend(), season(), dan
harmon()). Semua fungsi formula baru mengharuskan
argumennya berupa objek deret waktu (yaitu, "ts" atau
"zoo").
#sama dengan model dlm q=10
cons_lm1 <- dynlm(Yt ~ Xt+L(Xt)+L(Xt,2)+L(Xt,3)+L(Xt,4)+L(Xt,5)+L(Xt,6)+L(Xt,7)+L(Xt,8)+L(Xt,9)+L(Xt,10),data = training.ts)
cons_lm1
##
## Time series regression with "ts" data:
## Start = 11, End = 40
##
## Call:
## dynlm(formula = Yt ~ Xt + L(Xt) + L(Xt, 2) + L(Xt, 3) + L(Xt,
## 4) + L(Xt, 5) + L(Xt, 6) + L(Xt, 7) + L(Xt, 8) + L(Xt, 9) +
## L(Xt, 10), data = training.ts)
##
## Coefficients:
## (Intercept) Xt L(Xt) L(Xt, 2) L(Xt, 3) L(Xt, 4)
## 2.119e+01 5.366e-06 -1.429e-05 1.861e-05 -1.553e-05 1.173e-05
## L(Xt, 5) L(Xt, 6) L(Xt, 7) L(Xt, 8) L(Xt, 9) L(Xt, 10)
## -6.058e-06 -2.646e-06 5.389e-06 -8.511e-06 1.300e-05 -6.942e-06
summary(cons_lm1)
##
## Time series regression with "ts" data:
## Start = 11, End = 40
##
## Call:
## dynlm(formula = Yt ~ Xt + L(Xt) + L(Xt, 2) + L(Xt, 3) + L(Xt,
## 4) + L(Xt, 5) + L(Xt, 6) + L(Xt, 7) + L(Xt, 8) + L(Xt, 9) +
## L(Xt, 10), data = training.ts)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.77888 -0.19534 0.01772 0.20152 0.66379
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.119e+01 8.122e+00 2.609 0.0178 *
## Xt 5.366e-06 2.879e-06 1.864 0.0788 .
## L(Xt) -1.429e-05 7.888e-06 -1.811 0.0868 .
## L(Xt, 2) 1.861e-05 1.157e-05 1.609 0.1250
## L(Xt, 3) -1.553e-05 1.347e-05 -1.153 0.2639
## L(Xt, 4) 1.173e-05 1.436e-05 0.817 0.4248
## L(Xt, 5) -6.058e-06 1.497e-05 -0.405 0.6906
## L(Xt, 6) -2.646e-06 1.613e-05 -0.164 0.8715
## L(Xt, 7) 5.389e-06 1.812e-05 0.297 0.7696
## L(Xt, 8) -8.511e-06 1.738e-05 -0.490 0.6303
## L(Xt, 9) 1.300e-05 1.141e-05 1.139 0.2696
## L(Xt, 10) -6.942e-06 3.500e-06 -1.984 0.0627 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3791 on 18 degrees of freedom
## Multiple R-squared: 0.9897, Adjusted R-squared: 0.9834
## F-statistic: 156.8 on 11 and 18 DF, p-value: 1.563e-15
#sama dengan model ardl p=4 q=1
cons_lm2 <- dynlm(Yt ~ Xt+L(Xt)+L(Xt,2)+L(Xt,3)+L(Xt,4)+L(Yt),data = training.ts)
cons_lm2
##
## Time series regression with "ts" data:
## Start = 5, End = 40
##
## Call:
## dynlm(formula = Yt ~ Xt + L(Xt) + L(Xt, 2) + L(Xt, 3) + L(Xt,
## 4) + L(Yt), data = training.ts)
##
## Coefficients:
## (Intercept) Xt L(Xt) L(Xt, 2) L(Xt, 3) L(Xt, 4)
## 8.074e+00 6.805e-06 -1.883e-05 2.039e-05 -9.657e-06 1.337e-06
## L(Yt)
## 6.631e-01
summary(cons_lm2)
##
## Time series regression with "ts" data:
## Start = 5, End = 40
##
## Call:
## dynlm(formula = Yt ~ Xt + L(Xt) + L(Xt, 2) + L(Xt, 3) + L(Xt,
## 4) + L(Yt), data = training.ts)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.24409 -0.16700 -0.01354 0.22831 0.55611
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.074e+00 4.886e+00 1.652 0.109269
## Xt 6.805e-06 2.000e-06 3.402 0.001969 **
## L(Xt) -1.883e-05 5.960e-06 -3.159 0.003687 **
## L(Xt, 2) 2.038e-05 8.080e-06 2.523 0.017377 *
## L(Xt, 3) -9.657e-06 6.203e-06 -1.557 0.130373
## L(Xt, 4) 1.337e-06 2.138e-06 0.625 0.536753
## L(Yt) 6.631e-01 1.493e-01 4.442 0.000119 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3585 on 29 degrees of freedom
## Multiple R-squared: 0.9931, Adjusted R-squared: 0.9917
## F-statistic: 700.2 on 6 and 29 DF, p-value: < 2.2e-16
deviance(cons_lm1)
## [1] 2.58672
deviance(cons_lm2)
## [1] 3.727399
dwtest(cons_lm1)
##
## Durbin-Watson test
##
## data: cons_lm1
## DW = 1.6142, p-value = 0.03518
## alternative hypothesis: true autocorrelation is greater than 0
dwtest(cons_lm2)
##
## Durbin-Watson test
##
## data: cons_lm2
## DW = 2.5269, p-value = 0.8099
## alternative hypothesis: true autocorrelation is greater than 0
Berdasarkan hasil uji diagnostik model yaitu autokorelasi diperoleh \(P-Value<0.05\) untuk kedua model. Artinya, tidak terjadi autokorelasi pada model dlm dan ardl.
bptest(cons_lm1)
##
## studentized Breusch-Pagan test
##
## data: cons_lm1
## BP = 18.569, df = 11, p-value = 0.06929
bptest(cons_lm2)
##
## studentized Breusch-Pagan test
##
## data: cons_lm2
## BP = 11.614, df = 6, p-value = 0.07114
Berdasarkan hasil uji diagnostik model yaitu heterogenitas diperoleh \(P-Value<0.05\) untuk kedua model. Artinya, tidak terjadi heterogenitas pada model dlm dan ardl.
shapiro.test(residuals(cons_lm1))
##
## Shapiro-Wilk normality test
##
## data: residuals(cons_lm1)
## W = 0.98891, p-value = 0.9845
ks.test(residuals(cons_lm1), "pnorm", mean=mean(residuals(cons_lm1)), sd=sd(residuals(cons_lm1)))
##
## Exact one-sample Kolmogorov-Smirnov test
##
## data: residuals(cons_lm1)
## D = 0.074553, p-value = 0.9918
## alternative hypothesis: two-sided
shapiro.test(residuals(cons_lm2))
##
## Shapiro-Wilk normality test
##
## data: residuals(cons_lm2)
## W = 0.90139, p-value = 0.003724
ks.test(residuals(cons_lm2), "pnorm", mean=mean(residuals(cons_lm2)), sd=sd(residuals(cons_lm2)))
##
## Exact one-sample Kolmogorov-Smirnov test
##
## data: residuals(cons_lm2)
## D = 0.11885, p-value = 0.646
## alternative hypothesis: two-sided
Berdasarkan hasil uji diagnostik model yaitu normalitas dengan menggunakan uji shapiro-wilk diperoleh \(P-Value<0.05\) untuk model dlm tetapi \(P-Value>0.05\) untuk model ardl. Artinya, normalitas terpenuhi pada model dlm tetapi tidak pada model ardl. Namun, berdasarkan hasil uji diagnostik model yaitu normalitas dengan menggunakan uji kolmogorov-smirnov diperoleh \(P-Value<0.05\) untuk kedua model. Artinya, normalitas terpenuhi pada model dlm dan ardl.
akurasi <- matrix(c(mapekoyck, mapedlm, mapeardl))
row.names(akurasi)<- c("Koyck","DLM","Autoregressive")
colnames(akurasi) <- c("MAPE")
akurasi
## MAPE
## Koyck 0.1018948
## DLM 10.0657312
## Autoregressive 3.8729112
Berdasarkan nilai MAPE, model paling optimum didapat pada Model Koyck karena memiliki nilai MAPE yang terkecil.
par(mfrow=c(1,1))
plot(testing$Xt, testing$Yt, type="b", col="black", ylim=c(0,200))
points(testing$Xt, forekoyck$forecasts,col="red")
lines(testing$Xt, forekoyck$forecasts,col="red")
points(testing$Xt, foredlm$forecasts,col="blue")
lines(testing$Xt, foredlm$forecasts,col="blue")
points(testing$Xt, foreardl$forecasts,col="green")
lines(testing$Xt, foreardl$forecasts,col="green")
legend("topleft",c("aktual", "koyck","DLM", "autoregressive"), lty=1, col=c("black","red","blue","green"), cex=0.8)
Berdasarkan plot tersebut, terlihat bahwa plot yang paling mendekati data aktualnya adalah Model koyck, sehingga dapat disimpulkan model terbaik dalam hal ini adalah model regresi koyck.