Anggota Kelompok:
Muhammad Nur Sholihin (2207016011)
Dewi Maharani (2207016018)
Gita Rahayu (2207016019)
Nabila Husnul Listia (2207016038)
Dalam analisis deret waktu (time series), pemilihan metode yang tepat sangat krusial untuk menghasilkan prediksi yang akurat serta interpretasi yang bermakna terhadap fenomena sosial, ekonomi, maupun politik. Salah satu peristiwa penting yang kerap memengaruhi indikator-indikator makroekonomi dan sosial di banyak negara, termasuk Indonesia, adalah pemilihan umum (pemilu). Pemilu sering kali menjadi titik balik dalam kebijakan fiskal, stabilitas politik, hingga sentimen publik, yang pada akhirnya tercermin dalam data deret waktu seperti pertumbuhan ekonomi, inflasi, nilai tukar, dan indeks kepercayaan konsumen.
Dalam konteks tersebut, analisis dampak peristiwa pemilu terhadap variabel-variabel tersebut menjadi relevan, terutama ketika perubahan terjadi tidak hanya pada tahun pemilu, tetapi juga sebelum dan sesudahnya. Oleh karena itu, diperlukan model yang mampu menangkap dampak temporer (transitory) dan struktural dari momen pemilu. Adapun pendekatan yang digunakan dalam analisis ini adalah ARIMAX (Autoregressive Integrated Moving Average with Exogenous Variables) dan Time Series Regression.
Dengan menggunakan variabel dummy yaitu satu tahun sebelum pemilu, tahun ketika pemilu, dan satu tahun setelah pemilu, analisis ini bertujuan untuk mengukur dan membandingkan bagaimana masing-masing metode menangkap dampak peristiwa pemilu terhadap variabel yang diamati. Analisis ini bertujuan untuk melihat perbandingan akurasi prediksi yaitu dengan melihat nilai MAPE dari model TSR dan ARIMAX.
Data yang digunakan pada penelitian ini adalah data Indeks Harga Saham Gabungan bulanan dari bulan Januari 2010 hingga Mei 2025. Kemudian variabel dummy yang digunakan adalah 1 tahun sebelum pemilu, 1 tahun saat pelaksanaan pemilu, dan 1 tahun setelah pemilu. Data penelitian dapat ditampilkan seperti berikut:
data<-read.table(file.choose(),header=TRUE)
data
## Tanggal Terakhir D1 D2 D3 t
## 1 1/8/2010 3081.88 0 0 0 1
## 2 1/9/2010 3501.30 0 0 0 2
## 3 1/10/2010 3635.32 0 0 0 3
## 4 1/11/2010 3531.21 0 0 0 4
## 5 1/12/2010 3703.51 0 0 0 5
## 6 1/1/2011 3409.17 0 0 0 6
## 7 1/2/2011 3470.35 0 0 0 7
## 8 1/3/2011 3678.67 0 0 0 8
## 9 1/4/2011 3819.62 0 0 0 9
## 10 1/5/2011 3836.97 0 0 0 10
## 11 1/6/2011 3888.57 0 0 0 11
## 12 1/7/2011 4130.80 0 0 0 12
## 13 1/8/2011 3841.73 0 0 0 13
## 14 1/9/2011 3549.03 0 0 0 14
## 15 1/10/2011 3790.85 0 0 0 15
## 16 1/11/2011 3715.08 0 0 0 16
## 17 1/12/2011 3821.99 0 0 0 17
## 18 1/1/2012 3941.69 0 0 0 18
## 19 1/2/2012 3985.21 0 0 0 19
## 20 1/3/2012 4121.55 0 0 0 20
## 21 1/4/2012 4180.73 0 0 0 21
## 22 1/5/2012 3832.82 0 0 0 22
## 23 1/6/2012 3955.58 0 0 0 23
## 24 1/7/2012 4142.34 0 0 0 24
## 25 1/8/2012 4060.33 0 0 0 25
## 26 1/9/2012 4262.56 0 0 0 26
## 27 1/10/2012 4350.29 0 0 0 27
## 28 1/11/2012 4276.14 0 0 0 28
## 29 1/12/2012 4316.69 0 0 0 29
## 30 1/1/2013 4453.70 0 0 0 30
## 31 1/2/2013 4795.79 0 0 0 31
## 32 1/3/2013 4940.99 0 0 0 32
## 33 1/4/2013 5034.07 0 0 0 33
## 34 1/5/2013 5068.63 0 0 0 34
## 35 1/6/2013 4818.90 0 0 0 35
## 36 1/7/2013 4610.38 1 0 0 36
## 37 1/8/2013 4195.09 0 0 0 37
## 38 1/9/2013 4316.18 0 0 0 38
## 39 1/10/2013 4510.63 0 0 0 39
## 40 1/11/2013 4256.44 0 0 0 40
## 41 1/12/2013 4274.18 0 0 0 41
## 42 1/1/2014 4418.76 0 0 0 42
## 43 1/2/2014 4620.22 0 0 0 43
## 44 1/3/2014 4768.28 0 0 0 44
## 45 1/4/2014 4840.15 0 0 0 45
## 46 1/5/2014 4893.91 0 0 0 46
## 47 1/6/2014 4878.58 0 0 0 47
## 48 1/7/2014 5088.80 0 1 0 48
## 49 1/8/2014 5136.86 0 0 0 49
## 50 1/9/2014 5137.58 0 0 0 50
## 51 1/10/2014 5089.55 0 0 0 51
## 52 1/11/2014 5149.89 0 0 0 52
## 53 1/12/2014 5226.95 0 0 0 53
## 54 1/1/2015 5289.40 0 0 0 54
## 55 1/2/2015 5450.29 0 0 0 55
## 56 1/3/2015 5518.67 0 0 0 56
## 57 1/4/2015 5086.42 0 0 0 57
## 58 1/5/2015 5216.38 0 0 0 58
## 59 1/6/2015 4910.66 0 0 0 59
## 60 1/7/2015 4802.53 0 0 1 60
## 61 1/8/2015 4509.61 0 0 0 61
## 62 1/9/2015 4223.91 0 0 0 62
## 63 1/10/2015 4455.18 0 0 0 63
## 64 1/11/2015 4446.46 0 0 0 64
## 65 1/12/2015 4593.01 0 0 0 65
## 66 1/1/2016 4615.16 0 0 0 66
## 67 1/2/2016 4770.96 0 0 0 67
## 68 1/3/2016 4845.37 0 0 0 68
## 69 1/4/2016 4838.58 0 0 0 69
## 70 1/5/2016 4796.87 0 0 0 70
## 71 1/6/2016 5016.65 0 0 0 71
## 72 1/7/2016 5215.99 0 0 0 72
## 73 1/8/2016 5386.08 0 0 0 73
## 74 1/9/2016 5364.80 0 0 0 74
## 75 1/10/2016 5422.54 0 0 0 75
## 76 1/11/2016 5148.91 0 0 0 76
## 77 1/12/2016 5296.71 0 0 0 77
## 78 1/1/2017 5294.10 0 0 0 78
## 79 1/2/2017 5386.69 0 0 0 79
## 80 1/3/2017 5568.11 0 0 0 80
## 81 1/4/2017 5685.30 0 0 0 81
## 82 1/5/2017 5738.15 0 0 0 82
## 83 1/6/2017 5829.71 0 0 0 83
## 84 1/7/2017 5840.94 0 0 0 84
## 85 1/8/2017 5864.06 0 0 0 85
## 86 1/9/2017 5900.85 0 0 0 86
## 87 1/10/2017 6005.78 0 0 0 87
## 88 1/11/2017 5952.14 0 0 0 88
## 89 1/12/2017 6355.65 0 0 0 89
## 90 1/1/2018 6605.63 0 0 0 90
## 91 1/2/2018 6597.22 0 0 0 91
## 92 1/3/2018 6188.99 0 0 0 92
## 93 1/4/2018 5994.60 1 0 0 93
## 94 1/5/2018 5983.59 0 0 0 94
## 95 1/6/2018 5799.24 0 0 0 95
## 96 1/7/2018 5936.44 0 0 0 96
## 97 1/8/2018 6018.46 0 0 0 97
## 98 1/9/2018 5976.55 0 0 0 98
## 99 1/10/2018 5831.65 0 0 0 99
## 100 1/11/2018 6056.12 0 0 0 100
## 101 1/12/2018 6194.50 0 0 0 101
## 102 1/1/2019 6532.97 0 0 0 102
## 103 1/2/2019 6443.35 0 0 0 103
## 104 1/3/2019 6468.75 0 0 0 104
## 105 1/4/2019 6455.35 0 1 0 105
## 106 1/5/2019 6209.12 0 0 0 106
## 107 1/6/2019 6358.63 0 0 0 107
## 108 1/7/2019 6390.50 0 0 0 108
## 109 1/8/2019 6328.47 0 0 0 109
## 110 1/9/2019 6169.10 0 0 0 110
## 111 1/10/2019 6228.32 0 0 0 111
## 112 1/11/2019 6011.83 0 0 0 112
## 113 1/12/2019 6299.54 0 0 0 113
## 114 1/1/2020 5940.05 0 0 0 114
## 115 1/2/2020 5452.70 0 0 0 115
## 116 1/3/2020 4538.93 0 0 0 116
## 117 1/4/2020 4716.40 0 0 1 117
## 118 1/5/2020 4753.61 0 0 0 118
## 119 1/6/2020 4905.39 0 0 0 119
## 120 1/7/2020 5149.63 0 0 0 120
## 121 1/8/2020 5238.49 0 0 0 121
## 122 1/9/2020 4870.04 0 0 0 122
## 123 1/10/2020 5128.23 0 0 0 123
## 124 1/11/2020 5612.42 0 0 0 124
## 125 1/12/2020 5979.07 0 0 0 125
## 126 1/1/2021 5862.35 0 0 0 126
## 127 1/2/2021 6241.80 0 0 0 127
## 128 1/3/2021 5985.52 0 0 0 128
## 129 1/4/2021 5995.62 0 0 0 129
## 130 1/5/2021 5947.46 0 0 0 130
## 131 1/6/2021 5985.49 0 0 0 131
## 132 1/7/2021 6070.04 0 0 0 132
## 133 1/8/2021 6150.30 0 0 0 133
## 134 1/9/2021 6286.94 0 0 0 134
## 135 1/10/2021 6591.35 0 0 0 135
## 136 1/11/2021 6533.93 0 0 0 136
## 137 1/12/2021 6581.48 0 0 0 137
## 138 1/1/2022 6631.15 0 0 0 138
## 139 1/2/2022 6888.17 0 0 0 139
## 140 1/3/2022 7071.44 0 0 0 140
## 141 1/4/2022 7228.91 0 0 0 141
## 142 1/5/2022 7148.97 0 0 0 142
## 143 1/6/2022 6911.58 0 0 0 143
## 144 1/7/2022 6951.12 0 0 0 144
## 145 1/8/2022 7178.59 0 0 0 145
## 146 1/9/2022 7040.80 0 0 0 146
## 147 1/10/2022 7098.89 0 0 0 147
## 148 1/11/2022 7081.31 0 0 0 148
## 149 1/12/2022 6850.62 0 0 0 149
## 150 1/1/2023 6839.34 0 0 0 150
## 151 1/2/2023 6843.24 1 0 0 151
## 152 1/3/2023 6805.28 0 0 0 152
## 153 1/4/2023 6915.72 0 0 0 153
## 154 1/5/2023 6633.26 0 0 0 154
## 155 1/6/2023 6661.88 0 0 0 155
## 156 1/7/2023 6931.36 0 0 0 156
## 157 1/8/2023 6953.26 0 0 0 157
## 158 1/9/2023 6939.89 0 0 0 158
## 159 1/10/2023 6752.21 0 0 0 159
## 160 1/11/2023 7080.74 0 0 0 160
## 161 1/12/2023 7272.80 0 0 0 161
## 162 1/1/2024 7207.94 0 0 0 162
## 163 1/2/2024 7316.11 0 1 0 163
## 164 1/3/2024 7288.81 0 0 0 164
## 165 1/4/2024 7234.20 0 0 0 165
## 166 1/5/2024 6970.74 0 0 0 166
## 167 1/6/2024 7063.58 0 0 0 167
## 168 1/7/2024 7255.76 0 0 0 168
## 169 1/8/2024 7670.73 0 0 0 169
## 170 1/9/2024 7527.93 0 0 0 170
## 171 1/10/2024 7574.02 0 0 0 171
## 172 1/11/2024 7114.27 0 0 0 172
## 173 1/12/2024 7079.90 0 0 0 173
## 174 1/1/2025 7109.20 0 0 0 174
## 175 1/2/2025 6270.60 0 0 1 175
## 176 1/3/2025 6510.62 0 0 0 176
## 177 1/4/2025 6766.80 0 0 0 177
## 178 1/5/2025 6839.55 0 0 0 178
Pada penelitian ini, data dibagi menjadi data In Sample dan Out Sample dengan proporsi 90:10, sehingga data in sample merupakan 161 periode pertama, dan data out sample adalah 17 periode berikutnya
Yt=data$Terakhir[1:161]
t=data$t[1:161]
D1=data$D1[1:161]
D2=data$D2[1:161]
D3=data$D3[1:161]
Berdasarkan Data In Sample dapat diperoleh grafik runtun waktu sebagai berikut:
ts.plot(as.ts(Yt),lwd=5,ylab="Zt",col="red")
Berdasarkan grafik runtun waktu, dapat dilihat bahwa data IHSG membentuk pola tren naik, dan akan diselidiki apakah kejadian pemilu memengaruhi nilai IHSG.
Untuk mengetahui apakah nilai IHSG dipengaruhi oleh kejadian pemilu dilakukan pemodelan TSR dengan hasil seperti berikut:
#IMPORT PACKAGE
library(MASS)
library(hms)
library(car)
## Loading required package: carData
library(forecast)
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(tseries)
library(timeSeries)
## Warning: package 'timeSeries' was built under R version 4.4.1
## Loading required package: timeDate
##
## Attaching package: 'timeSeries'
## The following objects are masked from 'package:graphics':
##
## lines, points
library(quadprog)
library(zoo)
##
## Attaching package: 'zoo'
## The following object is masked from 'package:timeSeries':
##
## time<-
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(fracdiff)
library(fUnitRoots)
## Warning: package 'fUnitRoots' was built under R version 4.4.1
library(lmtest)
library(RcmdrMisc)
## Warning: package 'RcmdrMisc' was built under R version 4.4.1
## Loading required package: sandwich
##
## Attaching package: 'RcmdrMisc'
## The following object is masked from 'package:forecast':
##
## CV
#MODEL AWAL
model1=lm(Yt~t+D1+D2+D3)
summary(model1)
##
## Call:
## lm(formula = Yt ~ t + D1 + D2 + D3)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1570.68 -252.33 14.12 305.81 1032.41
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3716.4649 70.7499 52.530 <2e-16 ***
## t 20.6306 0.7555 27.308 <2e-16 ***
## D1 174.0864 259.6565 0.670 0.5036
## D2 477.3697 316.8011 1.507 0.1339
## D3 -782.8074 316.8396 -2.471 0.0146 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 445.1 on 156 degrees of freedom
## Multiple R-squared: 0.8286, Adjusted R-squared: 0.8242
## F-statistic: 188.5 on 4 and 156 DF, p-value: < 2.2e-16
#SELEKSI BACKWARD
model2=lm(Yt~t+D2+D3)
summary(model2)
##
## Call:
## lm(formula = Yt ~ t + D2 + D3)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1574.66 -257.03 14.02 306.94 1028.92
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3718.2865 70.5737 52.687 <2e-16 ***
## t 20.6492 0.7537 27.399 <2e-16 ***
## D2 474.1262 316.2083 1.499 0.136
## D3 -786.2740 316.2415 -2.486 0.014 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 444.3 on 157 degrees of freedom
## Multiple R-squared: 0.8281, Adjusted R-squared: 0.8248
## F-statistic: 252.1 on 3 and 157 DF, p-value: < 2.2e-16
model3=lm(Yt~t+D3)
summary(model3)
##
## Call:
## lm(formula = Yt ~ t + D3)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1580.20 -262.10 18.54 311.83 1023.06
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3725.2226 70.6995 52.691 <2e-16 ***
## t 20.6372 0.7566 27.277 <2e-16 ***
## D3 -792.1466 317.4639 -2.495 0.0136 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 446.1 on 158 degrees of freedom
## Multiple R-squared: 0.8256, Adjusted R-squared: 0.8234
## F-statistic: 374 on 2 and 158 DF, p-value: < 2.2e-16
Berdasarkan output dapat diperoleh bahwa kejadian 1 tahun setelah pemilu berpengaruh terhadap IHSG, dimana pengaruhnya merupakan pengaruh negatif. Selanjutnya akan dilihat nilai akurasi hasil prediksi in sample dan out sample untuk model TSR:
# Prediksi in-sample
y_pred_in <- predict(model1)
# Nilai aktual in-sample
y_actual_in <- Yt
# MAPE
mape_in_tsr <- mean(abs((y_actual_in - y_pred_in) / y_actual_in)) * 100
# Tampilkan hasil
cat("In-sample:\n")
## In-sample:
cat("MAPE:", mape_in_tsr, "%\n")
## MAPE: 6.381513 %
Yt_out = data$Terakhir[162:178]
t_out = data$t[162:178]
D3_out = data$D3[162:178]
data_out = as.data.frame(cbind(t_out,D3_out))
# Prediksi out-of-sample
y_pred_out <- 3725.2226 + 20.6372 * data_out$t_out - 792.1466 * data_out$D3_out
y_pred_out
## [1] 7068.449 7089.086 7109.723 7130.361 7150.998 7171.635 7192.272 7212.909
## [9] 7233.547 7254.184 7274.821 7295.458 7316.095 6544.586 7357.370 7378.007
## [17] 7398.644
# Nilai aktual dari data uji
y_actual_out <- Yt_out
# MAPE
mape_out_tsr <- mean(abs((y_actual_out - y_pred_out) / y_actual_out)) * 100
# Tampilkan hasil
cat("Out-of-sample:\n")
## Out-of-sample:
cat("MAPE:", mape_out_tsr, "%\n")
## MAPE: 4.165688 %
Kemudian dapat diperoleh grafik perbandingan data aktual, prediksi in sample, dan prediksi out sample sebagai berikut:
#Plot Perbandingan
ts.plot(data$Terakhir[1:178], xlab = "Periode", ylab = "IHSG", main = "Grafik Prediksi dan peramalan",ylim=c(2000,13000)) +
points(data$Terakhir[1:178], cex = 1, pch = 1) +
lines(1:161,y_pred_in[1:161], col="red",lwd=2) +
lines(162:178,y_pred_out[1:17], col="blue",lwd=2)
## integer(0)
legend("topleft", legend = c("Data aktual","Prediksi In Sample", "Prediksi Out Sample"),
col = c("black","red", "blue"), lty = 1, pch = 1, lwd = 3, bty = "o", cex = 0.8)
abline(v = 161, col = "black", lty = 2)
abline(v = 178, col = "black", lty = 2)
grid()
Berdasarkan hasil yang diperoleh, didapatkan bahwa pada model TSR memiliki nilai akurasi In Sample dan Out Sample yang tidak jauh berbeda, dan nilai akurasi Out Sample lebih baik dari nilai akurasi Out Sample. Berdasarkan grafik perbandingan, dapat dilihat bahwa data prediksi In Sample, dan prediksi Out Sample mengikuti pola data aktualnya.
res=resid(model3) #Residual dari model TSR
#Identifikasi Apakah Residual dari TSR Sudah WN
acf(res,237,lwd=1,col="blue")
x<-res
pv2<-rep(1,length(x))
for(i in 1:length(x))
pv2[i]=Box.test(x,lag=i,type="Ljung")$p.value
plot(pv2,type="p",main="Ljung-Box TSR",xlab="lag",ylab="p value",ylim=c(0,1))
abline(h=0.01,lty=1,col="red")
Berdasarkan output, dapat diperoleh bahwa residual model TSR belum memiliki sifat white noise, akan dilakukan pemodelan ARIMA pada residual TSR:
ts.plot(res,ylab="Residual TSR",col="red",lwd=3)
acf(res,237,col="blue",lwd=1) #Membuat ACF
pacf(res,237,col="blue",lwd=1) #Membuat PACF
#Stasioneritas dalam Variansi
min(res)
## [1] -1580.204
#Drifting agar residual tidak negatif
res=res+4681
#Transformasi Box-Cox
boxcox(res~1)
p<-powerTransform(res)
p
## Estimated transformation parameter
## res
## 3.270962
Berdasarkan grafik Box-Cox diperoleh bahwa nilai likelihood maksimum tercapai pada saat lambda = 3,270962 dimana nilai tersebut masih jauh dari 1, yang berarti data residual TSR belum stasioner dalam variansi sehingga akan dilakukan transformasi Box-Cox sebagai berikut:
res_trans = res^(3.270962)
head(res_trans)
## 1 2 3 4 5 6
## 614091525450 836928223236 909291177618 829903516437 927085915351 733419669840
boxcox(res_trans~1)
p1 <- powerTransform(res_trans)
p1
## Estimated transformation parameter
## res_trans
## 0.9999993
Diperoleh nilai likelihood maksimum tercapai pada saat lambda = 1 yang berarti data residual TSR sudah stasioner dalam variansi. Selanjutnya akan diperiksa stasioneritas data dalam rata-rata.
#Stasioneritas dalam Rata-Rata ADF
adf.test(res_trans)
##
## Augmented Dickey-Fuller Test
##
## data: res_trans
## Dickey-Fuller = -3.2585, Lag order = 5, p-value = 0.08059
## alternative hypothesis: stationary
acf(res_trans,237,col="blue",lwd=1) #Membuat ACF
Diperoleh bahwa data masih belum stasioner dalam rata-rata karena nilai p-value uji ADF masih > taraf siginifikansi 5%, sehingga akan dilakukan differencing orde 1 sebagai berikut
diff1 = diff(res_trans, differences = 1)
head(diff1)
## 2 3 4 5 6
## 222836697785 72362954382 -79387661181 97182398914 -193666245511
## 7
## 23182952357
adf.test(diff1)
## Warning in adf.test(diff1): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: diff1
## Dickey-Fuller = -5.3351, Lag order = 5, p-value = 0.01
## alternative hypothesis: stationary
acf(diff1,161*(1/4),col="blue",lwd=1) #Membuat ACF
pacf(diff1,161*(1/4),col="blue",lwd=1) #Membuat PACF
Diperoleh dari uji ADF bahwa data sudah stasioner dalam rata-rata akan tetapi grafik ACF dan PACF menunjukkan bahwa tidak ada lag yang signifikan sehingga peneliti kesulitan dalam menentukan orde model ARIMA, oleh karena itu akan dilakukan differencing orde 2 sebagai berikut:
diff2 = diff(res_trans, differences = 2)
head(diff2)
## 3 4 5 6 7
## -150473743403 -151750615564 176570060095 -290848644424 216849197868
## 8
## 90798351116
adf.test(diff2)
## Warning in adf.test(diff2): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: diff2
## Dickey-Fuller = -8.8888, Lag order = 5, p-value = 0.01
## alternative hypothesis: stationary
acf(diff2,161*(1/4),col="blue",lwd=1) #Membuat ACF
pacf(diff2,161*(1/4),col="blue",lwd=1) #Membuat PACF
Berdasarkan grafik ACF dan PACF dapat diperoleh model ARIMA sementara sebagai berikut:
#Pemodelan ARIMA
library(lmtest)
orders <- expand.grid(p = 0:3, q = 0:1)
orders <- orders[!(orders$p == 0 & orders$q == 0), ]
#List menyimpan model
models <- list()
#Looping model dengan tiap kombinasi p dan q
for (i in 1:nrow(orders)) {
p <- orders$p[i]
q <- orders$q[i]
model_name <- paste0("ARIMA_", p, "2", q)
models[[model_name]] <- arima(x = res_trans, order = c(p,2,q))
cat("\nModel", model_name, "\n")
print(models[[model_name]])
print(coeftest(models[[model_name]]))
}
##
## Model ARIMA_120
##
## Call:
## arima(x = res_trans, order = c(p, 2, q))
##
## Coefficients:
## ar1
## -0.5583
## s.e. 0.0654
##
## sigma^2 estimated as 3.845e+22: log likelihood = -4360.08, aic = 8724.16
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.558268 0.065411 -8.5348 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Model ARIMA_220
##
## Call:
## arima(x = res_trans, order = c(p, 2, q))
##
## Coefficients:
## ar1 ar2
## -0.7194 -0.2852
## s.e. 0.0762 0.0765
##
## sigma^2 estimated as 3.532e+22: log likelihood = -4353.43, aic = 8712.87
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.719373 0.076218 -9.4383 < 2.2e-16 ***
## ar2 -0.285241 0.076514 -3.7279 0.000193 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Model ARIMA_320
##
## Call:
## arima(x = res_trans, order = c(p, 2, q))
##
## Coefficients:
## ar1 ar2 ar3
## -0.8108 -0.5181 -0.3197
## s.e. 0.0754 0.0909 0.0754
##
## sigma^2 estimated as 3.168e+22: log likelihood = -4344.95, aic = 8697.9
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.810817 0.075364 -10.7586 < 2.2e-16 ***
## ar2 -0.518123 0.090905 -5.6996 1.201e-08 ***
## ar3 -0.319708 0.075380 -4.2413 2.222e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Model ARIMA_021
##
## Call:
## arima(x = res_trans, order = c(p, 2, q))
##
## Coefficients:
## ma1
## -1.0000
## s.e. 0.0171
##
## sigma^2 estimated as 2.538e+22: log likelihood = -4329.41, aic = 8662.81
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ma1 -0.999999 0.017134 -58.365 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Model ARIMA_121
##
## Call:
## arima(x = res_trans, order = c(p, 2, q))
##
## Coefficients:
## ar1 ma1
## -0.1086 -1.0000
## s.e. 0.0794 0.0176
##
## sigma^2 estimated as 2.505e+22: log likelihood = -4328.48, aic = 8662.95
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.108582 0.079384 -1.3678 0.1714
## ma1 -1.000000 0.017568 -56.9228 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Model ARIMA_221
##
## Call:
## arima(x = res_trans, order = c(p, 2, q))
##
## Coefficients:
## ar1 ar2 ma1
## -0.1072 0.0098 -1.0000
## s.e. 0.0801 0.0803 0.0175
##
## sigma^2 estimated as 2.505e+22: log likelihood = -4328.47, aic = 8664.94
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.1072433 0.0800889 -1.3391 0.1806
## ar2 0.0098076 0.0803133 0.1221 0.9028
## ma1 -0.9999827 0.0175245 -57.0620 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Model ARIMA_321
##
## Call:
## arima(x = res_trans, order = c(p, 2, q))
##
## Coefficients:
## ar1 ar2 ar3 ma1
## -0.1068 -0.0022 -0.0962 -1.000
## s.e. 0.0796 0.0804 0.0798 0.018
##
## sigma^2 estimated as 2.479e+22: log likelihood = -4327.75, aic = 8665.5
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.106794 0.079629 -1.3411 0.1799
## ar2 -0.002172 0.080450 -0.0270 0.9785
## ar3 -0.096189 0.079764 -1.2059 0.2279
## ma1 -1.000000 0.018003 -55.5463 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# (1,2,0) (2,2,0) (3,2,0) (0,2,1)
Diperoleh model yang signifikan adalah sebagai berikut:
#Pemodelan ARIMA
fit1=arima(x=res_trans,order=c(1,2,0))
coeftest(fit1)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.558268 0.065411 -8.5348 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
fit2=arima(x=res_trans,order=c(2,2,0))
coeftest(fit2)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.719373 0.076218 -9.4383 < 2.2e-16 ***
## ar2 -0.285241 0.076514 -3.7279 0.000193 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
fit3=arima(x=res_trans,order=c(3,2,0))
coeftest(fit3)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.810817 0.075364 -10.7586 < 2.2e-16 ***
## ar2 -0.518123 0.090905 -5.6996 1.201e-08 ***
## ar3 -0.319708 0.075380 -4.2413 2.222e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
fit4=arima(x=res_trans,order=c(0,2,1))
coeftest(fit4)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ma1 -0.999999 0.017134 -58.365 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Selanjutnya akan dilakukan pemeriksaan diagostik untuk model ARIMA yang signifikan:
#Pemeriksaan Diagnostik
res1=resid(fit1)
res2=resid(fit2)
res3=resid(fit3)
res4=resid(fit4)
#WN
x<-res1
pv2<-rep(1,length(x))
for(i in 1:length(x))
pv2[i]=Box.test(x,lag=i,type="Ljung")$p.value
plot(pv2,type="p",main="Ljung-Box",xlab="lag",ylab="p value",ylim=c(0,1))
abline(h=0.01,lty=1,col="red")
x<-res2
pv2<-rep(1,length(x))
for(i in 1:length(x))
pv2[i]=Box.test(x,lag=i,type="Ljung")$p.value
plot(pv2,type="p",main="Ljung-Box",xlab="lag",ylab="p value",ylim=c(0,1))
abline(h=0.01,lty=1,col="red")
x<-res3
pv2<-rep(1,length(x))
for(i in 1:length(x))
pv2[i]=Box.test(x,lag=i,type="Ljung")$p.value
plot(pv2,type="p",main="Ljung-Box",xlab="lag",ylab="p value",ylim=c(0,1))
abline(h=0.01,lty=1,col="red")
x<-res4
pv2<-rep(1,length(x))
for(i in 1:length(x))
pv2[i]=Box.test(x,lag=i,type="Ljung")$p.value
plot(pv2,type="p",main="Ljung-Box",xlab="lag",ylab="p value",ylim=c(0,1))
abline(h=0.01,lty=1,col="red")
Diperoleh bahwa hanya model ARIMA (3,2,0) dan ARIMA (0,2,1) yang memenuhi asumsi white noise, selanjutnya akan dilakukan uji residual berdistribusi normal untuk model ARIMA (3,2,0) dan ARIMA (0,2,1)
#Cek Normalitas
library(nortest)
ks.test(res3,"pnorm",mean(res3),sqrt(var(res3)))
##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: res3
## D = 0.068813, p-value = 0.4309
## alternative hypothesis: two-sided
ks.test(res4,"pnorm",mean(res4),sqrt(var(res4)))
##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: res4
## D = 0.13385, p-value = 0.006247
## alternative hypothesis: two-sided
lillie.test(res3)
##
## Lilliefors (Kolmogorov-Smirnov) normality test
##
## data: res3
## D = 0.068813, p-value = 0.06013
lillie.test(res4)
##
## Lilliefors (Kolmogorov-Smirnov) normality test
##
## data: res4
## D = 0.13385, p-value = 1.919e-07
Diperoleh bahwa model ARIMA (3,2,0) memenuhi asumsi independensi residual (white noise) dan asumsi residual berdistribusi normal, sehingga dapat dipilih model ARIMA terbaik adalah model ARIMA (3,2,0). Selanjutnya akan dilakukan pemodelan ARIMAX dengan orde model ARIMA (3,2,0)
#Pemodelan ARIMAX
library(TSA)
## Warning: package 'TSA' was built under R version 4.4.1
## Registered S3 methods overwritten by 'TSA':
## method from
## fitted.Arima forecast
## plot.Arima forecast
##
## Attaching package: 'TSA'
## The following objects are masked from 'package:timeDate':
##
## kurtosis, skewness
## The following objects are masked from 'package:stats':
##
## acf, arima
## The following object is masked from 'package:utils':
##
## tar
xreg<-data.frame(t,D3)
xregmat <- as.matrix(xreg)
xreg_trimmed <- xreg[3:length(Yt), ]
Yt=data$Terakhir[1:161]
modelarimax=arima(diff(Yt, differences = 2),order=c(3,0,0),xreg = xreg_trimmed,include.mean = FALSE)
coeftest(modelarimax)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.609283 0.078886 -7.7235 1.131e-14 ***
## ar2 -0.389731 0.087943 -4.4317 9.351e-06 ***
## ar3 -0.233631 0.078859 -2.9627 0.003050 **
## t -0.032161 0.084669 -0.3798 0.704063
## D3 354.605917 127.206552 2.7876 0.005309 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Berdasarkan pemodelan ARIMAX diperoleh bahwa variabel t tidak berpengaruh secara signifikan terhadap IHSG, maka variabel t akan dieliminasi dan diperoleh:
xreg1<-data.frame(D3)
xreg1_trimmed <- xreg1[3:length(Yt), ]
Yt=data$Terakhir[1:161]
modelarimax=arima(diff(Yt, differences = 2),order=c(3,0,0),xreg=xreg1_trimmed,include.mean = FALSE)
coeftest(modelarimax)
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.609355 0.078880 -7.7251 1.118e-14 ***
## ar2 -0.388897 0.087970 -4.4208 9.834e-06 ***
## ar3 -0.233112 0.078873 -2.9555 0.003121 **
## xreg 345.543969 125.036903 2.7635 0.005718 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Telah diperoleh model ARIMAX dengan seluruh parameternya signifikan, selanjutnya akan dilakukan pemeriksaan diagnostik model ARIMAX:
#WN
x<-resid(modelarimax)
pv2<-rep(1,length(x))
for(i in 1:length(x))
pv2[i]=Box.test(x,lag=i,type="Ljung")$p.value
plot(pv2,type="p",main="Ljung-Box",xlab="lag",ylab="p value",ylim=c(0,1))
abline(h=0.04,lty=1,col="red")
#Normalitas
ks.test(x,"pnorm",mean(x),sqrt(var(x)))
##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: x
## D = 0.074579, p-value = 0.3394
## alternative hypothesis: two-sided
Diperoleh bahwa model ARIMAX (3,2,0) telah memenuhi asumsi independensi residual (White noise) dan asumsi normalitas residual. Selanjutnya akan diperoleh nilai prediksi in sample dan out sample dari model ARIMAX (3,2,0)
#Prediksi In Sample
fitted_diff2 <- fitted(modelarimax) # output dari arima(Yt_diff2, ...)
# Inisialisasi: panjang fitted = length(Yt) - 2
n <- length(fitted_diff2)
Y_hat <- numeric(n + 2) # kita akan rekonstruksi hingga panjang Yt
# Gunakan rumus rekursif untuk kembalikan fitted ke skala Yt
for (t in 3:(n + 2)) {
Y_hat[t] <- 2 * Yt[t - 1] - Yt[t - 2] + fitted_diff2[t - 2]
}
Pred = as.data.frame(Y_hat)
Pred
## Y_hat
## 1 0.000
## 2 0.000
## 3 3876.563
## 4 3871.284
## 5 3634.179
## 6 3866.516
## 7 3347.195
## 8 3431.933
## 9 3767.848
## 10 3861.524
## 11 3921.536
## 12 3983.072
## 13 4272.362
## 14 3794.291
## 15 3420.725
## 16 3832.222
## 17 3625.808
## 18 3816.490
## 19 4056.587
## 20 4027.592
## 21 4227.974
## 22 4268.589
## 23 3741.342
## 24 3967.838
## 25 4201.957
## 26 4007.488
## 27 4381.192
## 28 4459.904
## 29 4278.901
## 30 4376.993
## 31 4525.061
## 32 4948.663
## 33 5103.925
## 34 5187.673
## 35 5205.016
## 36 4777.311
## 37 4500.950
## 38 3956.041
## 39 4181.230
## 40 4499.982
## 41 4122.065
## 42 4283.592
## 43 4484.880
## 44 4674.302
## 45 4897.191
## 46 4965.954
## 47 5000.784
## 48 4930.154
## 49 5192.671
## 50 5212.123
## 51 5177.632
## 52 5127.438
## 53 5174.188
## 54 5263.041
## 55 5328.988
## 56 5552.979
## 57 5608.544
## 58 4972.261
## 59 5220.013
## 60 5114.029
## 61 4822.933
## 62 4488.394
## 63 4040.164
## 64 4411.701
## 65 4381.248
## 66 4617.765
## 67 4708.674
## 68 4857.503
## 69 4946.398
## 70 4881.767
## 71 4826.990
## 72 5109.599
## 73 5334.233
## 74 5520.986
## 75 5476.272
## 76 5513.370
## 77 5091.082
## 78 5298.158
## 79 5296.496
## 80 5381.523
## 81 5693.440
## 82 5784.891
## 83 5834.477
## 84 5937.676
## 85 5901.064
## 86 5902.151
## 87 5943.412
## 88 6061.101
## 89 5965.439
## 90 6526.377
## 91 6808.345
## 92 6699.402
## 93 6160.669
## 94 5885.628
## 95 5870.877
## 96 5599.351
## 97 5902.365
## 98 6049.462
## 99 5956.660
## 100 5810.567
## 101 6124.455
## 102 6265.701
## 103 6696.890
## 104 6556.843
## 105 6543.902
## 106 6520.655
## 107 6093.043
## 108 6366.586
## 109 6394.428
## 110 6277.157
## 111 6132.985
## 112 6214.086
## 113 5901.027
## 114 6336.280
## 115 5843.124
## 116 5177.421
## 117 4431.139
## 118 4635.116
## 119 4685.693
## 120 4868.073
## 121 5325.669
## 122 5359.367
## 123 4819.127
## 124 5218.641
## 125 5821.802
## 126 6183.376
## 127 6033.202
## 128 6534.287
## 129 6036.345
## 130 5974.971
## 131 5979.403
## 132 5931.561
## 133 6106.305
## 134 6194.991
## 135 6380.049
## 136 6772.603
## 137 6618.605
## 138 6666.671
## 139 6723.052
## 140 6993.546
## 141 7218.518
## 142 7382.447
## 143 7240.922
## 144 6868.475
## 145 6938.486
## 146 7220.550
## 147 6987.942
## 148 7135.859
## 149 7118.809
## 150 6733.556
## 151 6794.879
## 152 6802.241
## 153 6735.777
## 154 6948.472
## 155 6542.261
## 156 6619.146
## 157 7024.682
## 158 6959.838
## 159 6988.148
## 160 6742.177
## 161 7170.725
#Prediksi Out Sample
Ytnew=data$Terakhir[162:178]
tnew=data$t[162:178]
D3new=data$D3[162:178]
xregnewbgt=data.frame(D3new)
fore.ARIMA=predict(modelarimax,newYt=Ytnew,newxreg=xregnewbgt,17)$pred
se.fore.ARIMA=predict(modelarimax,newxreg=xregnewbgt,26)$se
## Warning in z[[1L]] + xm: longer object length is not a multiple of shorter
## object length
Pred_out=fore.ARIMA
Pred_out
## Time Series:
## Start = 160
## End = 176
## Frequency = 1
## [1] -76.96012798 -20.36572016 74.15223210 -19.32456268 -12.31457103
## [6] -2.26653484 10.67500028 -2.75274617 -1.94571959 -0.23229762
## [11] 1.53993278 -0.39445599 -0.30436012 345.52385846 0.22257111
## [16] -0.05685418 -0.04722495
Y_hatout <- numeric(17) # kita akan rekonstruksi hingga panjang Yt
Y_hatout
## [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
# Gunakan rumus rekursif untuk kembalikan fitted ke skala Yt
Y_hatout[1] = 2*Yt[161] - Yt[160] + Pred_out[1]
Y_hatout[1]
## [1] 7387.9
Y_hatout[2] = 2*Y_hatout[1] - Yt[161] + Pred_out[2]
Y_hatout[2]
## [1] 7482.634
for (t in 3:17) {
Y_hatout[t] <- 2 * Y_hatout[t - 1] - Y_hatout[t - 2] + Pred_out[t]
}
Y_hatout
## [1] 7387.900 7482.634 7651.520 7801.082 7938.329 8073.310 8218.966
## [8] 8361.869 8502.826 8643.551 8785.816 8927.686 9069.252 9556.342
## [15] 10043.655 10530.911 11018.119
Predd = as.data.frame(Y_hatout)
Predd
## Y_hatout
## 1 7387.900
## 2 7482.634
## 3 7651.520
## 4 7801.082
## 5 7938.329
## 6 8073.310
## 7 8218.966
## 8 8361.869
## 9 8502.826
## 10 8643.551
## 11 8785.816
## 12 8927.686
## 13 9069.252
## 14 9556.342
## 15 10043.655
## 16 10530.911
## 17 11018.119
Kemudian dapat diperoleh nilai akurasi prediksi In Sample dan Out Sample sebagai berikut:
#Akurasi In Sample
actual_values <- Yt[3:161]
mape_in_arimax <- mean(abs((actual_values - Pred$Y_hat[3:161])/actual_values)) * 100
print(mape_in_arimax)
## [1] 3.371419
#Akurasi Out Sample
actual_values <- data$Terakhir[162:178]
mape_out_arimax <- mean(abs((actual_values - Predd$Y_hatout)/actual_values)) * 100
print(mape_out_arimax)
## [1] 23.27459
Kemudian dapat diperoleh grafik perbandingan data aktual dengan prediksi In Sample dan prediksi Out Sample adalah sebagai berikut:
#Plot Perbandingan
ts.plot(data$Terakhir[1:178], xlab = "Periode", ylab = "IHSG", main = "Grafik Prediksi dan peramalan",ylim=c(2000,13000)) +
points(data$Terakhir[1:178], cex = 1, pch = 1) +
lines(3:161,Pred$Y_hat[3:161], col="red",lwd=2) +
lines(162:178,Predd$Y_hatout[1:17], col="blue",lwd=2)
## integer(0)
legend("topleft", legend = c("Data aktual","Prediksi In Sample", "Prediksi Out Sample"),
col = c("black","red", "blue"), lty = 1, pch = 1, lwd = 3, bty = "o", cex = 0.8)
abline(v = 161, col = "black", lty = 2)
abline(v = 178, col = "black", lty = 2)
grid()
Berdasarkan output dapat diperoleh bahwa hasil akurasi prediksi in sample dengan out sample berbeda signifikan, dimana model ARIMAX berarti lebih baik dalam memprediksi data in sample dibandingkan data out sample.
list_model = c("Model TSR","Model ARIMAX")
mape_in_model1 = mape_in_tsr
mape_in_model2 = mape_in_arimax
mape_out_model1 = mape_out_tsr
mape_out_model2 = mape_out_arimax
mape_in_model = c(mape_in_model1,mape_in_model2)
mape_out_model = c(mape_out_model1, mape_out_model2)
pemilihan_model = data.frame(row.names = list_model, "MAPE In Sample" = mape_in_model, "MAPE Out Sample" = mape_out_model)
pemilihan_model
## MAPE.In.Sample MAPE.Out.Sample
## Model TSR 6.381513 4.165688
## Model ARIMAX 3.371419 23.274590
Berdasarkan nilai MAPE yang diperoleh, MAPE In Sample terkecil terdapat pada model ARIMAX, sedangkan MAPE Out Sample terkecil terdapat pada model TSR. Selain itu dapat diketahui bahwa nilai MAPE In Sample dan Out Sample untuk model TSR tidak berbeda signifikan, sedangkan untuk model ARIMAX MAPE In Sample dan Out Sample nya berbeda signifikan.
Berdasarkan hasil penelitian, dapat diperoleh dari model TSR ataupun ARIMAX bahwa kejadian 1 tahun setelah pemilu berpengaruh terhadap nilai IHSG dan termasuk pengaruh negatif. Model TSR maupun ARIMAX memiliki kemampuan yang sangat baik dalam memprediksi data In Sample, tetapi model TSR lebih baik dalam memprediksi data Out Sample dibandingkan model ARIMAX. Sehingga pada penelitian ini, diperoleh model terbaik untuk memprediksi nilai IHSG adalah model TSR.