円ドルダービー5月末提出分
本田力也
2021/5/12
Preparation
全コードは、右上のcodeをクリックすることで表示される。必要に応じて参照されたい。
Packages
Import data
# import exchange rate data -------------------------------------------------------------
# with quantmod
## but only access to the past 180 days
## so this is for a short run analysis
getFX("USD/JPY")
getFX("USD/EUR")
getFX("USD/CNY")
getFX("USD/GBP")
# with European Central Bank
## but all currencies are against euro
## so, this is for a long run analysis
url <- "http://www.ecb.europa.eu/stats/eurofxref/eurofxref-hist.zip"
download.file(url, "eurofxref-hist.zip")
rates <- read_csv(unz("eurofxref-hist.zip", "eurofxref-hist.csv"), na = "N/A")
rates
rates <- rates %>%
mutate(USDJPY = JPY / USD, # JPY against USD
USDEUR = 1 / USD, # EUR against USD
USDCNY = CNY / USD, # CNY (yuan) against USD
USDGBP = GBP / USD, # GBP (Pond sterling) against USD
)
# create xts objects
USDJPY_long <- as.xts(rates$USDJPY, order.by = rates$Date)
# Bank of Japan -----------------------------------------------------------
# url : https://www.stat-search.boj.or.jp/ssi/cgi-bin/famecgi2?cgi=$nme_a000&lstSelection=FM08
# full dat; "https://www.stat-search.boj.or.jp/ssi/html/nme_R031.15138.20210513000306.02.csv"
rates_boj <- read_csv("C:/Users/rikit/Documents/R/Musashi_University/2021_first_semester/exchage_forcast/Data/nme_R031.1110.20210513135227.02.csv") %>%
`colnames<-`(c("Date", "Open", "High", "Low", "Close", "Volume")) %>%
na.omit() %>%
read.zoo() %>%
as.xts()
# import economic data ----------------------------------------------------
# interest rate
# countrycode
names <- c("Japan", "USA", "China", "Great Britain")
iso3c_names <- countrycode::countrycode(sourcevar = names,
origin = "country.name",
destination = "iso3c")
iso3c_names
# get interest rate data
start_year <- 1999
end_year <- 2021
interest_rate <- WDI(country = iso3c_names,
indicator = "FR.INR.RINR", # see the more indicator from "https://data.worldbank.org/indicator"
start = start_year,
end = end_year) %>%
as_tibble() %>%
`colnames<-`(c("iso2c", "country", "real_interest_rate", "year"))
# get inflation (GDP deflator)
inflation_rate <- WDI(country = iso3c_names,
indicator = "NY.GDP.DEFL.KD.ZG", # see the more indicator from "https://data.worldbank.org/indicator"
start = start_year,
end = end_year) %>%
as_tibble() %>%
`colnames<-`(c("iso2c", "country", "inflation_rate", "year"))
# get GDP growth rate
gdp_growth_rate <- WDI(country = iso3c_names,
indicator = "NY.GDP.MKTP.KD.ZG", # see the more indicator from "https://data.worldbank.org/indicator"
start = start_year,
end = end_year) %>%
as_tibble() %>%
`colnames<-`(c("iso2c", "country", "gdp_growth_rate", "year"))
# get unemployment rate
unemployment_rate <- WDI(country = iso3c_names,
indicator = "SL.UEM.TOTL.ZS", # see the more indicator from "https://data.worldbank.org/indicator"
start = start_year,
end = end_year) %>%
as_tibble() %>%
`colnames<-`(c("iso2c", "country", "unemployment_rate", "year"))
# get industry, value added
industry_value_added <- WDI(country = iso3c_names,
indicator = "NV.IND.TOTL.ZS", # see the more indicator from "https://data.worldbank.org/indicator"
start = start_year,
end = end_year) %>%
as_tibble() %>%
`colnames<-`(c("iso2c", "country", "industry_value_added", "year"))
economic_df <-
full_join(interest_rate, inflation_rate) %>%
full_join(gdp_growth_rate) %>%
full_join(unemployment_rate) %>%
full_join(industry_value_added) %>%
select(year, everything())
economic_dfBasic Visualization
a Candle Chart
candleChart(rates_boj, major.ticks = "months", subset = "2020-12::",
theme = "white", type = "candles", TA = "addBBands();addVo();addCCI()")
currencies against USD
# USD...
g1 <- rates %>%
select(Date, starts_with("USD")) %>%
pivot_longer(cols = -Date) %>%
filter(name != "USD") %>%
ggplot(aes(Date, value, colour = name)) +
geom_line() +
facet_wrap(~name, ncol = 1, scales = "free_y")
# ggplotly(g1)
g1
Modeling
hypothesis1:Correlation between currencies
まず、通貨ごとの相関関係を確認してみる。長期（1999年～2020年）までの相関関係と、短期（過去180日間）のデータを 使って相関関係を図示化する。
短期的見ると、Euroと正の相関関係と、GBPとの負の相関関係がみられる。その為、ユーロ圏の経済とイギリス経済を予測することで、円ドルへの影響を予測することができるかもしれない。
in the long run
mrates_long <- rates[,c(45:47,44)] %>%
na.omit() %>%
cor() %>%
round(digits = 2)
corrplot(mrates_long,
type = "lower", # 下半分を表示
method = "shade", # 四角いパネルで表示
shade.col = NA, # パネルに政府を示す斜線を加えない
tl.col = "black", # textlabelの色を変える
tl.srt = 45, # 上部テキストラベルの角度を変更
addCoef.col = "black", # 相関係数を黒で表示
cl.pos = "n", # 凡例を表示させない
order = "original") # "hclust", "AOE", "FPC"
in the short run
mrates_short <- tibble(
USDEUR = USDEUR,
USDCNY = USDCNY,
USDGBP = USDGBP,
USDJPY = USDJPY
) %>%
cor() %>%
round(digits = 2)
corrplot(mrates_short,
type = "lower", # 下半分を表示
method = "shade", # 四角いパネルで表示
shade.col = NA, # パネルに政府を示す斜線を加えない
tl.col = "black", # textlabelの色を変える
tl.srt = 45, # 上部テキストラベルの角度を変更
addCoef.col = "black", # 相関係数を黒で表示
cl.pos = "n", # 凡例を表示させない
order = "original") # "hclust", "AOE", "FPC"
Forecasting with ACF&P
ACF
自己相関関係：機能の値は、今日の値に影響している。そして、一昨日の値は、昨日の値に影響しているとする考え方。
過去に遡るほど、現在との相関関係が弱くなる。
USDJPY_return <- periodReturn(USDJPY, period = "daily")
acf(USDJPY_return, plot = F, lag.max = 40) # F or T##
## Autocorrelations of series 'USDJPY_return', by lag
##
## 0 1 2 3 4 5 6 7 8 9 10
## 1.000 0.357 0.013 0.061 0.063 0.021 -0.091 -0.036 0.019 0.015 -0.042
## 11 12 13 14 15 16 17 18 19 20 21
## -0.028 -0.064 -0.024 0.001 -0.038 0.098 0.039 0.030 0.072 0.040 0.023
## 22 23 24 25 26 27 28 29 30 31 32
## -0.040 -0.133 -0.056 0.007 0.042 0.125 0.169 0.152 0.017 -0.026 0.069
## 33 34 35 36 37 38 39 40
## 0.001 -0.084 -0.001 0.033 -0.018 -0.021 -0.054 -0.052
PAFC
偏自己相関係数：昨日の影響を除去した、昨日と今日の値の感あ系を調べる方法
先ほどとは異なり、一日前の影響のみが95％有意となっている。つまり、偏自己相関係数を基に30日後を予測するのは至難の技であることが分かる。
##
## Partial autocorrelations of series 'USDJPY_return', by lag
##
## 1 2 3 4 5 6 7 8 9 10 11
## 0.357 -0.132 0.120 -0.004 0.006 -0.115 0.047 -0.002 0.021 -0.055 0.017
## 12 13 14 15 16 17 18 19 20 21 22
## -0.095 0.047 -0.015 -0.021 0.141 -0.071 0.063 0.030 0.004 -0.001 -0.044
## 23 24 25 26 27 28 29 30 31 32 33
## -0.133 0.047 -0.007 0.091 0.103 0.134 0.040 -0.077 0.016 0.044 -0.062
## 34 35 36 37 38 39 40
## -0.045 0.033 -0.013 -0.022 0.024 -0.007 -0.032
Ljung-Box検定
Ling-Box検定は、自己相関関係を有しているかどうかを調べる。帰無仮説は、自己相関関係を有していない。
今回は、p-valueが0.05以下である為、自己相関関係があるという対立仮説を採用する。
##
## Box-Ljung test
##
## data: USDJPY_return
## X-squared = 23.196, df = 1, p-value = 1.463e-06ARモデル
推定アルゴリズムは、ユールウォーカー法を使用している。
##
## Call:
## ar(x = USDJPY)
##
## Coefficients:
## 1 2
## 1.0994 -0.1097
##
## Order selected 2 sigma^2 estimated as 0.1114usdjpy_pr <- predict(usdjpy_ar, n.ahead = 20)
usdjpy_se1 <- usdjpy_pr$pred + 2 * usdjpy_pr$se
usdjpy_se2 <- usdjpy_pr$pred - 2 * usdjpy_pr$se
ts.plot(usdjpy_pr$pred, usdjpy_se1, usdjpy_se2,
gpars=list(lt=c(1,2,3,3),col=c(1,2,4,4)))
重回帰分析
まず、多重共線性の確認をする。変数間で異常な相関関係がないかを調べる。
# rates と economic_dfをfull_joinする
df <- rates %>%
mutate(year = year(Date)) %>%
group_by(year) %>%
summarise(mean_USDJPY = mean(USDJPY)) %>%
full_join(economic_df %>% filter(country == "Japan")) %>%
mutate(return_USDJPY = )
# 多重共線性の確認
df %>%
select(mean_USDJPY, real_interest_rate, inflation_rate,
gdp_growth_rate, unemployment_rate, industry_value_added) %>%
na.omit() %>%
cor()## mean_USDJPY real_interest_rate inflation_rate
## mean_USDJPY 1.00000000 -0.1078971 0.2055002
## real_interest_rate -0.10789706 1.0000000 -0.9804042
## inflation_rate 0.20550021 -0.9804042 1.0000000
## gdp_growth_rate 0.07549084 0.0507250 -0.1057339
## unemployment_rate -0.07375917 0.7622289 -0.6997876
## industry_value_added 0.68778812 0.4200985 -0.3065225
## gdp_growth_rate unemployment_rate industry_value_added
## mean_USDJPY 0.07549084 -0.07375917 0.6877881
## real_interest_rate 0.05072500 0.76222891 0.4200985
## inflation_rate -0.10573392 -0.69978762 -0.3065225
## gdp_growth_rate 1.00000000 -0.14727598 0.2071851
## unemployment_rate -0.14727598 1.00000000 0.2757046
## industry_value_added 0.20718510 0.27570456 1.0000000以上が、が行列になるが、分かりにくいので図示する。
corrplot(df %>%
select(mean_USDJPY, real_interest_rate, inflation_rate,
gdp_growth_rate, unemployment_rate, industry_value_added) %>%
na.omit() %>%
cor(),
type = "lower", # 下半分を表示
method = "shade", # 四角いパネルで表示
shade.col = NA, # パネルに政府を示す斜線を加えない
tl.col = "black", # textlabelの色を変える
tl.srt = 45, # 上部テキストラベルの角度を変更
addCoef.col = "black", # 相関係数を黒で表示
cl.pos = "n", # 凡例を表示させない
order = "original") # "hclust", "AOE", "FPC"
以上が、変数の事前分析になる。
以下では、AICを用いて、重回帰分析を行った。
最終的に、以下の式が定義された。
\[ y = \hat{\alpha} + \hat{\beta_1}X1 + \hat{\beta_2}X2\] \[ USDJPY = \hat{\alpha} + \hat{\beta_1}RealInterestRate + \hat{\beta_2}Industry Value Added \]
# 重回帰分析
## mean_USDJPY
mod1 <- lm(mean_USDJPY ~ real_interest_rate + inflation_rate + gdp_growth_rate +
unemployment_rate + industry_value_added, data = df %>% na.omit())
summary(mod1)##
## Call:
## lm(formula = mean_USDJPY ~ real_interest_rate + inflation_rate +
## gdp_growth_rate + unemployment_rate + industry_value_added,
## data = df %>% na.omit())
##
## Residuals:
## Min 1Q Median 3Q Max
## -14.2077 -4.6153 0.4899 5.0397 14.1418
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -88.9963 49.0339 -1.815 0.09266 .
## real_interest_rate -8.8292 14.7378 -0.599 0.55941
## inflation_rate -3.7157 16.3010 -0.228 0.82324
## gdp_growth_rate -0.6331 1.2939 -0.489 0.63279
## unemployment_rate 2.0726 5.1843 0.400 0.69581
## industry_value_added 6.9821 1.9136 3.649 0.00295 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.996 on 13 degrees of freedom
## Multiple R-squared: 0.6758, Adjusted R-squared: 0.5511
## F-statistic: 5.42 on 5 and 13 DF, p-value: 0.00654## Start: AIC=88.27
## mean_USDJPY ~ real_interest_rate + inflation_rate + gdp_growth_rate +
## unemployment_rate + industry_value_added
##
## Df Sum of Sq RSS AIC
## - inflation_rate 1 4.20 1056.2 86.343
## - unemployment_rate 1 12.93 1065.0 86.499
## - gdp_growth_rate 1 19.37 1071.4 86.614
## - real_interest_rate 1 29.05 1081.1 86.784
## <none> 1052.0 88.267
## - industry_value_added 1 1077.34 2129.4 99.664
##
## Step: AIC=86.34
## mean_USDJPY ~ real_interest_rate + gdp_growth_rate + unemployment_rate +
## industry_value_added
##
## Df Sum of Sq RSS AIC
## - unemployment_rate 1 8.95 1065.2 84.503
## - gdp_growth_rate 1 15.18 1071.4 84.614
## <none> 1056.3 86.343
## - real_interest_rate 1 354.27 1410.5 89.838
## - industry_value_added 1 2125.07 3181.3 105.292
##
## Step: AIC=84.5
## mean_USDJPY ~ real_interest_rate + gdp_growth_rate + industry_value_added
##
## Df Sum of Sq RSS AIC
## - gdp_growth_rate 1 24.25 1089.5 82.931
## <none> 1065.2 84.503
## - real_interest_rate 1 629.59 1694.8 91.327
## - industry_value_added 1 2120.82 3186.0 103.320
##
## Step: AIC=82.93
## mean_USDJPY ~ real_interest_rate + industry_value_added
##
## Df Sum of Sq RSS AIC
## <none> 1089.5 82.931
## - real_interest_rate 1 620.55 1710.0 89.496
## - industry_value_added 1 2117.89 3207.3 101.446##
## Call:
## lm(formula = mean_USDJPY ~ real_interest_rate + industry_value_added,
## data = df %>% na.omit())
##
## Coefficients:
## (Intercept) real_interest_rate industry_value_added
## -73.459 -4.826 6.508# AICの結果
mod1_1 <- lm(formula = mean_USDJPY ~ real_interest_rate + industry_value_added,
data = df %>% na.omit())
summary(mod1_1)##
## Call:
## lm(formula = mean_USDJPY ~ real_interest_rate + industry_value_added,
## data = df %>% na.omit())
##
## Residuals:
## Min 1Q Median 3Q Max
## -15.1107 -4.8170 0.8384 5.0178 16.1088
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -73.459 32.969 -2.228 0.04057 *
## real_interest_rate -4.826 1.599 -3.019 0.00815 **
## industry_value_added 6.508 1.167 5.577 4.17e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.252 on 16 degrees of freedom
## Multiple R-squared: 0.6643, Adjusted R-squared: 0.6223
## F-statistic: 15.83 on 2 and 16 DF, p-value: 0.0001614how to predict real interest rate?
For a central bank, there are three possible ways to control the real interest rate.
- Monetary policies
- Reserve Ration
- Open market operations
- bank rate operation
Fisher equation:
\[ i + \pi^{e} = r\]