This example shows time series forecasting of Euro-AUD exchange rates with the with ARIMA and STL models. The data used are historical currency exchange rates from January 1999 to June 2014 provided by the European Central Bank.
Download data from the European Central Bank at http://www.ecb.europa.eu/stats/exchange/eurofxref/html/index.en.html.
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"), header = T)
rates[1:5, ]
## Date USD JPY BGN CYP CZK DKK EEK GBP HUF LTL
## 1 2015-04-17 1.0814 128.45 1.9558 N/A 27.486 7.4608 N/A 0.7199 301.36 N/A
## 2 2015-04-16 1.0711 127.64 1.9558 N/A 27.503 7.4637 N/A 0.7189 300.95 N/A
## 3 2015-04-15 1.0579 126.52 1.9558 N/A 27.415 7.4671 N/A 0.7168 297.85 N/A
## 4 2015-04-14 1.0564 126.67 1.9558 N/A 27.345 7.4691 N/A 0.7217 296.64 N/A
## 5 2015-04-13 1.0552 127.20 1.9558 N/A 27.392 7.4715 N/A 0.7215 296.85 N/A
## LVL MTL PLN ROL RON SEK SIT SKK CHF ISK NOK HRK
## 1 N/A N/A 4.0187 N/A 4.418 9.2921 N/A N/A 1.0297 N/A 8.3970 7.569
## 2 N/A N/A 4.0310 N/A 4.4158 9.2606 N/A N/A 1.0327 N/A 8.3725 7.5748
## 3 N/A N/A 4.0178 N/A 4.4197 9.3115 N/A N/A 1.0321 N/A 8.4150 7.5798
## 4 N/A N/A 4.0105 N/A 4.4118 9.3235 N/A N/A 1.0344 N/A 8.4820 7.583
## 5 N/A N/A 4.0136 N/A 4.4156 9.3309 N/A N/A 1.0373 N/A 8.5480 7.568
## RUB TRL TRY AUD BRL CAD CNY HKD IDR INR
## 1 54.9086 N/A 2.9037 1.3842 3.2562 1.3183 6.7013 8.3820 13838.67 67.474
## 2 53.4111 N/A 2.8899 1.3790 3.2312 1.3180 6.6377 8.3031 13758.43 66.6873
## 3 53.611 N/A 2.8711 1.3945 3.2661 1.3277 6.5644 8.2004 13686.9 65.9826
## 4 54.935 N/A 2.8284 1.3936 3.2857 1.3292 6.5633 8.1872 13710.09 65.905
## 5 55.2605 N/A 2.7985 1.3925 3.2694 1.3324 6.5552 8.1779 13723.03 65.8671
## KRW MXN MYR NZD PHP SGD THB ZAR ILS X
## 1 1167.11 16.4205 3.911 1.4025 47.819 1.4522 34.992 12.9192 4.239 NA
## 2 1164.56 16.3075 3.9046 1.4038 47.562 1.4502 34.709 12.9081 4.2206 NA
## 3 1159.50 16.2123 3.9206 1.4104 47.086 1.4402 34.327 12.8567 4.2095 NA
## 4 1156.94 16.2263 3.9092 1.4169 47.117 1.4400 34.291 12.8100 4.2108 NA
## 5 1163.30 16.1398 3.9172 1.4187 47.2 1.4491 34.376 12.8258 4.2239 NA
str(rates$Date)
## Factor w/ 4171 levels "1999-01-04","1999-01-05",..: 4171 4170 4169 4168 4167 4166 4165 4164 4163 4162 ...
## convert into date format
rates$Date <- as.Date(rates$Date, "%Y-%m-%d")
str(rates$Date)
## Date[1:4171], format: "2015-04-17" "2015-04-16" "2015-04-15" "2015-04-14" ...
range(rates$Date)
## [1] "1999-01-04" "2015-04-17"
rates <- rates[order(rates$Date), ]
## plot time series
plot(rates$Date, rates$AUD, type = "l")
The code below shows that there are no data for weekends or public holidays.
head(rates$Date, 20)
## [1] "1999-01-04" "1999-01-05" "1999-01-06" "1999-01-07" "1999-01-08"
## [6] "1999-01-11" "1999-01-12" "1999-01-13" "1999-01-14" "1999-01-15"
## [11] "1999-01-18" "1999-01-19" "1999-01-20" "1999-01-21" "1999-01-22"
## [16] "1999-01-25" "1999-01-26" "1999-01-27" "1999-01-28" "1999-01-29"
years <- format(rates$Date, "%Y")
tab <- table(years)
tab
## years
## 1999 2000 2001 2002 2003 2004 2005 2006 2007 2008 2009 2010 2011 2012 2013
## 259 255 254 255 255 259 257 255 255 256 256 258 257 256 255
## 2014 2015
## 255 74
## number of days per year after removing 2014
mean(tab[1:(length(tab) - 1)])
## [1] 256.0625
Based on above result, there are about 256 values per year, so the windows size is set to 256 in time series analysis in section 5. Another way is to fill weekends and public holidays with values in the previous populated days.
source("forecast.R") ## see code file in section 5
result.arima <- forecastArima(rates, n.ahead = 90)
source("plotForecastResult.R") ## see code file in section 5
plotForecastResult(result.arima, title = "Exchange rate forecasting with ARIMA")
result.stl <- forecastStl(rates, n.ahead = 90)
plotForecastResult(result.stl, title = "Exchange rate forecasting with STL")
## exchange rate in 2014
result <- subset(result.stl, date >= "2014-01-01" )
plotForecastResult(result, title = "Exchange rate forecasting with STL (2014)")
Below are two source files used in section 3 and 4.
It provides functions for forecasting with ARIMA and STL.
> library(forecast)
> forecastStl <- function(x, n.ahead = 30) {
+ myTs <- ts(x$AUD, start = 1, frequency = 256)
+ fit.stl <- stl(myTs, s.window = 256)
+ sts <- fit.stl$time.series
+ trend <- sts[, "trend"]
+ fore <- forecast(fit.stl, h = n.ahead, level = 95)
+ plot(fore)
+ pred <- fore$mean
+ upper <- fore$upper
+ lower <- fore$lower
+ output <- data.frame(actual = c(x$AUD, rep(NA, n.ahead)),
+ trend = c(trend, rep(NA, n.ahead)), pred = c(rep(NA,
+ nrow(x)), pred), lower = c(rep(NA, nrow(x)), lower),
+ upper = c(rep(NA, nrow(x)), upper), date = c(x$Date,
+ max(x$Date) + (1:n.ahead)))
+ return(output)
+ }
> forecastArima <- function(x, n.ahead = 30) {
+ myTs <- ts(x$AUD, start = 1, frequency = 256)
+ fit.arima <- arima(myTs, order = c(0, 0, 1))
+ fore <- forecast(fit.arima, h = n.ahead)
+ plot(fore)
+ upper <- fore$upper[, "95%"]
+ lower <- fore$lower[, "95%"]
+ trend <- as.numeric(fore$fitted)
+ pred <- as.numeric(fore$mean)
+ output <- data.frame(actual = c(x$AUD, rep(NA, n.ahead)),
+ trend = c(trend, rep(NA, n.ahead)), pred = c(rep(NA,
+ nrow(x)), pred), lower = c(rep(NA, nrow(x)), lower),
+ upper = c(rep(NA, nrow(x)), upper), date = c(x$Date,
+ max(x$Date) + (1:n.ahead)))
+ return(output)
+ }
It provides a function for ploting time series forecasting result, incl. trend, forecast and bounds.
> plotForecastResult <- function(x, title = NULL) {
+ x <- x[order(x$date), ]
+ max.val <- max(c(x$actual, x$upper), na.rm = T)
+ min.val <- min(c(x$actual, x$lower), na.rm = T)
+ plot(x$date, x$actual, type = "l", col = "grey", main = title,
+ xlab = "Time", ylab = "Exchange Rate", xlim = range(x$date),
+ ylim = c(min.val, max.val))
+ grid()
+ lines(x$date, x$trend, col = "yellowgreen")
+ lines(x$date, x$pred, col = "green")
+ lines(x$date, x$lower, col = "blue")
+ lines(x$date, x$upper, col = "blue")
+ legend("bottomleft", col = c("grey", "yellowgreen", "green",
+ "blue"), lty = 1, c("Actual", "Trend", "Forecast", "Lower/Upper Bound"))
+ }