# 1102計量經濟學（二）期末報告《匯率預測》#######
# Packages ########
library(forecast)
## Warning: package 'forecast' was built under R version 4.0.5
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(tseries)
library(readr)
library(readxl)
library(ggplot2)
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.0.5
##
## 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(pander)
## Warning: package 'pander' was built under R version 4.0.5
# Import Data #######
df <- read_xlsx("~/Desktop/Econometrics/data_第12組.xlsx")
# Plot #######
plot(x = df$month, y = df$E, type = "l",
main = "Nominal Exchange Rate Trend (2011.1-2022.3)",
xlab = "Year", ylab = " Nominal Exchange Rate")
plot(x = df$month, y = df$E_real, type = "l",
main = "Real Exchange Rate Trend (2011.1-2022.3)",
xlab = "Year", ylab = "Real Exchange Rate")
# Real Exchange Rate ######
# Augmented Dickey-Fuller Test (Unit Root Test) #######
df$E_real %>%
na.omit %>%
adf.test %>%
pander
Augmented Dickey-Fuller Test: .
| -2.264 |
5 |
0.4668 |
stationary |
# First-Order Difference
df <- df %>%
mutate(dE_real = c(NA, diff(E_real)))
df$dE_real %>%
na.omit %>%
adf.test %>%
pander
## Warning in adf.test(.): p-value smaller than printed p-value
Augmented Dickey-Fuller Test: .
| -4.374 |
5 |
0.01 * * |
stationary |
# Box-Pierce and Ljung-Box Tests (White Noise Test) #######
df$dE_real %>%
na.omit %>%
Box.test(type = "Ljung-Box") %>%
pander
Box-Ljung test: .
| 5.767 |
1 |
0.01633 * |
# Find Best ARIMA #######
md_r <- auto.arima(df$E_real)
auto.arima(df$E_real, trace = TRUE)
##
## ARIMA(2,1,2) with drift : 87.08205
## ARIMA(0,1,0) with drift : 86.2018
## ARIMA(1,1,0) with drift : 82.1572
## ARIMA(0,1,1) with drift : 82.62584
## ARIMA(0,1,0) : 84.85514
## ARIMA(2,1,0) with drift : 84.27802
## ARIMA(1,1,1) with drift : 84.28061
## ARIMA(2,1,1) with drift : Inf
## ARIMA(1,1,0) : 80.57237
## ARIMA(2,1,0) : 82.66264
## ARIMA(1,1,1) : 82.66414
## ARIMA(0,1,1) : 81.0765
## ARIMA(2,1,1) : Inf
##
## Best model: ARIMA(1,1,0)
## Series: df$E_real
## ARIMA(1,1,0)
##
## Coefficients:
## ar1
## 0.2219
## s.e. 0.0870
##
## sigma^2 = 0.1044: log likelihood = -38.24
## AIC=80.48 AICc=80.57 BIC=86.28
accuracy(md_r)
## ME RMSE MAE MPE MAPE MASE
## Training set 0.02040839 0.3206395 0.2488294 0.06229534 0.7948412 0.9520995
## ACF1
## Training set -0.002464701
# Forecasting ######
forecast(md_r, 10)
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 136 32.65785 32.24386 33.07185 32.02470 33.29100
## 137 32.70500 32.05132 33.35868 31.70528 33.70471
## 138 32.71546 31.87627 33.55465 31.43203 33.99889
## 139 32.71778 31.72482 33.71074 31.19918 34.23639
## 140 32.71830 31.59190 33.84469 30.99562 34.44097
## 141 32.71841 31.47270 33.96413 30.81325 34.62357
## 142 32.71844 31.36385 34.07302 30.64678 34.79009
## 143 32.71844 31.26311 34.17377 30.49270 34.94418
## 144 32.71844 31.16890 34.26799 30.34862 35.08827
## 145 32.71844 31.08009 34.35680 30.21280 35.22409
plot(forecast(md_r, 10),
xlab = 'month', ylab = 'real exchange rate',
main='Real Exchange Rate Forecast from ARIMA(1,1,0)',
sub = "t = 10")
# Nominal Exchange Rate #########
# Augmented Dickey-Fuller Test######
df$E %>%
na.omit %>%
adf.test %>%
pander
Augmented Dickey-Fuller Test: .
| -1.691 |
5 |
0.7052 |
stationary |
df <- df %>%
mutate(dE = c(NA, diff(E)))
df$dE %>%
na.omit %>%
adf.test %>%
pander
## Warning in adf.test(.): p-value smaller than printed p-value
Augmented Dickey-Fuller Test: .
| -4.93 |
5 |
0.01 * * |
stationary |
# Box-Pierce and Ljung-Box Tests#######
df$dE %>%
na.omit %>%
Box.test(type = "Ljung-Box") %>%
pander
Box-Ljung test: .
| 22.69 |
1 |
1.908e-06 * * * |
md_n <- auto.arima(df$E)
auto.arima(df$E, trace = T)
##
## ARIMA(2,1,2) with drift : Inf
## ARIMA(0,1,0) with drift : 49.37992
## ARIMA(1,1,0) with drift : 26.48812
## ARIMA(0,1,1) with drift : 25.3912
## ARIMA(0,1,0) : 47.44342
## ARIMA(1,1,1) with drift : 26.80114
## ARIMA(0,1,2) with drift : 26.95079
## ARIMA(1,1,2) with drift : Inf
## ARIMA(0,1,1) : 23.36598
## ARIMA(1,1,1) : 24.72466
## ARIMA(0,1,2) : 24.87963
## ARIMA(1,1,0) : 24.4267
## ARIMA(1,1,2) : Inf
##
## Best model: ARIMA(0,1,1)
## Series: df$E
## ARIMA(0,1,1)
##
## Coefficients:
## ma1
## 0.4444
## s.e. 0.0771
##
## sigma^2 = 0.06801: log likelihood = -9.64
## AIC=23.27 AICc=23.37 BIC=29.07
# Forecasting ######
forecast(md_n, 10)
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 136 28.64073 28.30652 28.97494 28.12960 29.15187
## 137 28.64073 28.05359 29.22788 27.74277 29.53869
## 138 28.64073 27.88061 29.40085 27.47823 29.80323
## 139 28.64073 27.74028 29.54119 27.26361 30.01786
## 140 28.64073 27.61904 29.66242 27.07819 30.20328
## 141 28.64073 27.51073 29.77073 26.91255 30.36892
## 142 28.64073 27.41194 29.86953 26.76145 30.52001
## 143 28.64073 27.32051 29.96095 26.62163 30.65983
## 144 28.64073 27.23502 30.04644 26.49089 30.79058
## 145 28.64073 27.15444 30.12702 26.36765 30.91381
plot(forecast(md_n, 10),
xlab = 'month', ylab = 'real exchange rate',
main='Nominal Exchange Rate Forecast from ARIMA(0, 1, 1)',
sub = "t = 10")
# ARIMA(1, 1, 1)#######
# real exchange rate
md_r_111 <- arima(df$E_real, order=c(1, 1, 1))
plot(forecast(md_r_111, 10), xlab = 'month', ylab = 'real exchange rate',
main = ' Alternative Real Exchange Rate Forecast from ARIMA(1,1,1)',
sub = "t = 10")
# nominal exchange rate
md_n_111 <- arima(df$E, order = c(1, 1, 1))
plot(forecast(md_n_111, 10), xlab = 'month', ylab = 'nominal exchange rate',
main=' Alternative Nominal Exchange Rate Forecast from ARIMA(1,1,1)',
sub = "t = 10")
# Output #######
md_r_110 <- arima(df$E_real, order = c(1, 1, 0))
md_n_011 <- arima(df$E, order = c(0, 1, 1))
stargazer::stargazer(md_r_110, md_r_111 , md_n_011, md_n_111,
dep.var.labels = rep(c("exchange rate", "real exchange rate"), each = 2),
type = "html", title = "ARIMA model", out = "result.doc")
##
## <table style="text-align:center"><caption><strong>ARIMA model</strong></caption>
## <tr><td colspan="5" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"></td><td colspan="4"><em>Dependent variable:</em></td></tr>
## <tr><td></td><td colspan="4" style="border-bottom: 1px solid black"></td></tr>
## <tr><td style="text-align:left"></td><td>exchange rate</td><td>exchange rate</td><td>real exchange rate</td><td>real exchange rate</td></tr>
## <tr><td style="text-align:left"></td><td>(1)</td><td>(2)</td><td>(3)</td><td>(4)</td></tr>
## <tr><td colspan="5" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">ar1</td><td>0.222<sup>**</sup></td><td>0.212</td><td></td><td>0.193</td></tr>
## <tr><td style="text-align:left"></td><td>(0.087)</td><td>(0.278)</td><td></td><td>(0.213)</td></tr>
## <tr><td style="text-align:left"></td><td></td><td></td><td></td><td></td></tr>
## <tr><td style="text-align:left">ma1</td><td></td><td>0.010</td><td>0.444<sup>***</sup></td><td>0.277</td></tr>
## <tr><td style="text-align:left"></td><td></td><td>(0.276)</td><td>(0.077)</td><td>(0.213)</td></tr>
## <tr><td style="text-align:left"></td><td></td><td></td><td></td><td></td></tr>
## <tr><td colspan="5" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">Observations</td><td>134</td><td>134</td><td>134</td><td>134</td></tr>
## <tr><td style="text-align:left">Log Likelihood</td><td>-38.240</td><td>-38.240</td><td>-9.637</td><td>-9.270</td></tr>
## <tr><td style="text-align:left">sigma<sup>2</sup></td><td>0.104</td><td>0.104</td><td>0.067</td><td>0.067</td></tr>
## <tr><td style="text-align:left">Akaike Inf. Crit.</td><td>80.481</td><td>82.480</td><td>23.274</td><td>24.540</td></tr>
## <tr><td colspan="5" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"><em>Note:</em></td><td colspan="4" style="text-align:right"><sup>*</sup>p<0.1; <sup>**</sup>p<0.05; <sup>***</sup>p<0.01</td></tr>
## </table>
# 4 Plots in one ######
par(mfrow = c(2, 2))
plot(forecast(md_r, 10),
xlab = 'month', ylab = 'real exchange rate',
main='Real Exchange Rate Forecast from ARIMA(1,1,0)',
sub = "t = 10")
plot(forecast(md_r_111, 10), xlab = 'month', ylab = 'real exchange rate',
main = ' Alternative Real Exchange Rate Forecast from ARIMA(1,1,1)',
sub = "t = 10")
plot(forecast(md_n, 10),
xlab = 'month', ylab = 'real exchange rate',
main='Nominal Exchange Rate Forecast from ARIMA(0, 1, 1)',
sub = "t = 10")
plot(forecast(md_n_111, 10), xlab = 'month', ylab = 'nominal exchange rate',
main=' Alternative Nominal Exchange Rate Forecast from ARIMA(1,1,1)',
sub = "t = 10")
par(mfrow = c(1, 1))