This dossier is a continuation of previous work done by our group in an attempt to sufficiently predict the movement of the ASX 200 for the following month using a sample of a number of key economic and financial variables. However, this piece will only be focusing on one macroeconomic variable, the cash rate.
To maintain coherence of the writing, the aim of this work is to provide an example of evaluating the tools available associated with time series forecasting to help sufficiently predict macroeconomic data for the future as well as an overview of a selection assessment of the available methodologies.
There are extensive research and rich literature that support the empirical evidence of the relationship between stock market indexes and macroeconomic activities. Cheung & Ng (1998) have worked on international evidence of the explanatory power of aggregate macroeconomic variables with five national stock market indexes. Moreover, Asprem (1989) had found an inverse relationship between interest rates and stock prices indices in ten different European countries. To globalise the assumption, Gupta, Chevalier, and Sayekt (2001) has found a causality between interest rate and stock market prices in an emerging market.
The author of this article personally believes that when the interest rates are low, people are taking more risks by investing in the stock market because the return in their safe deposits are less favourable. In addition to that, other factor believed by the author is, since the cost of money is low, companies listed in stock markets have lower expenses on their financial liabilities which means higher profits, as well as the ability to borrow more money to expand the business, both of which increases investors confidence which in turn increases the stock market indices in general.
In the previous work, the group has concluded that using delayed and freely available inputs of economic and financial data do not provide the informational advantage to successfully consistently predict the returns of the ASX200.
There are three main recommendations that was introduced, in this piece we are going to focus on the third option by continuing the work and focus the analyses by: Reducing the scope of the task (focus on a more specific input) Incorporate more advanced modelling techniques
To summarise, we are going to focus on macroeconomic data, the cash rate, and using time series analysis and modelling to forecast this macroeconomic data. This will provide an example of a successful or unsuccessful use of time series modelling to be used to predict the ASX200 returns in future studies.
The source of the data will follow the previous work, focusing on the Reserve Bank of Australia’s cash rate announcement, with time period between January 2005 to June 2019. Then the data will be put into forward filling exercises to allow for monthly data and avoid missing values.
The granularity of the most commonly used in time series model. Using ARIMA, the seasonality of monthly data are not difficult to assess as it contains 12 data points per year. Though, I doubt there would be a seasonality in the cash rate datasets that we have. Conversely, the number of data points spans thirteen-and-a-half years, accumulating 174 data points. Lastly, we have zero missing values in our data.
From the visual features of the datasets, we can see extract key information below * An upward trend from 2005 until around second half of 2008 * Another upward trend from mid first-half of 2009 to the end of 2011 * A downward trend from end of 2011until the end of the dataset (mid 2019)
plot(og_cr, main = "RBA Cash Rate")
Let’s ensure that our datasets are in time series object.
# convert to ts object
cr_ts = ts(og_cr, start = 2005, frequency = 12)
cr_ts
## Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
## 2005 5.25 5.25 5.50 5.50 5.50 5.50 5.50 5.50 5.50 5.50 5.50 5.50
## 2006 5.50 5.50 5.50 5.50 5.75 5.75 5.75 6.00 6.00 6.00 6.25 6.25
## 2007 6.25 6.25 6.25 6.25 6.25 6.25 6.25 6.50 6.50 6.50 6.75 6.75
## 2008 6.75 7.00 7.25 7.25 7.25 7.25 7.25 7.25 7.00 6.00 5.25 4.25
## 2009 4.25 3.25 3.25 3.00 3.00 3.00 3.00 3.00 3.00 3.25 3.50 3.75
## 2010 3.75 3.75 4.00 4.25 4.50 4.50 4.50 4.50 4.50 4.50 4.75 4.75
## 2011 4.75 4.75 4.75 4.75 4.75 4.75 4.75 4.75 4.75 4.75 4.50 4.25
## 2012 4.25 4.25 4.25 4.25 3.75 3.50 3.50 3.50 3.50 3.25 3.25 3.00
## 2013 3.00 3.00 3.00 3.00 2.75 2.75 2.75 2.50 2.50 2.50 2.50 2.50
## 2014 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50 2.50
## 2015 2.50 2.25 2.25 2.25 2.00 2.00 2.00 2.00 2.00 2.00 2.00 2.00
## 2016 2.00 2.00 2.00 2.00 1.75 1.75 1.75 1.50 1.50 1.50 1.50 1.50
## 2017 1.50 1.50 1.50 1.50 1.50 1.50 1.50 1.50 1.50 1.50 1.50 1.50
## 2018 1.50 1.50 1.50 1.50 1.50 1.50 1.50 1.50 1.50 1.50 1.50 1.50
## 2019 1.50 1.50 1.50 1.50 1.50 1.25
Here is a seasonal monthly plot of the RBA cash rate.
# seasonal plot
ggseasonplot(cr_ts, main = "Monthly Seasonal Plot of RBA Cash Rate")
As suggested earlier, there is no coherent seasonality in the cash rate data. While in some years, the rates are higher than when its started, other years like in recent years were flat. Some even fell lower than it started. This does not conform to the definition of seasonal pattern prescribed by (Hyndman & Athanasopoulos 2018, chap. 2.3) which states that a seasonal pattern occurs when a time series is affected by seasonal factors such as the time of the year or the day of the week. Furthermore, they stressed that seasonality is always of a fixed and known frequency.
Adopting the distinction of seasonality above, there is no sign of consistent seasonality pattern throughout the datasets. Though in nine out of fourteen years, the first four months of the year goes flat. As this is just a majority finding, not a fixed pattern, seasonality is still ruled out in this dataset.
In order to conform to the foundation of statistical inference in time series analysis, we need to ensure stationarity in our data. The stationarity of our data is also required to fit an ARIMA model for time series forecasting.
# KPPS test of stationarity
cr_ts %>% ur.kpss %>% summary()
##
## #######################
## # KPSS Unit Root Test #
## #######################
##
## Test is of type: mu with 4 lags.
##
## Value of test-statistic is: 3.0587
##
## Critical value for a significance level of:
## 10pct 5pct 2.5pct 1pct
## critical values 0.347 0.463 0.574 0.739
As shown above, the test-statistic value is so much higher than the 1pct critical value. This also implies that the time-series object of our datasets fails the test and the data is not stationary. Next, we can apply differencing to our data to determine if the result will make our data a stationary series.
# Adding differencing for stationarity series
cr_ts %>% diff() %>% ur.kpss %>% summary()
##
## #######################
## # KPSS Unit Root Test #
## #######################
##
## Test is of type: mu with 4 lags.
##
## Value of test-statistic is: 0.0732
##
## Critical value for a significance level of:
## 10pct 5pct 2.5pct 1pct
## critical values 0.347 0.463 0.574 0.739
From the above, we can see that the test-statistic became smaller. Though it is not much lower than the 1pct of the critical values, it is still within the range of what we would expect for stationarity data. So we can conclude that the differenced data are stationary.
To delve even further, let’s find out the number of first differences carried out to make our time series data a stationary series.
# Finding out number of first differencing carried out
ndiffs(cr_ts)
## [1] 1
From the above function, we are informed that there is only one first difference data required and carried out.
There is another function we can use to prove that there is no seasonality on our time-series data.
# Determining number of seasonal differencing required
nsdiffs(cr_ts)
## [1] 0
The above shows that there is zero number of seasonal differencing required for our time-series data, cementing out earlier inference that our time-series data has no seasonality.
As expected from the non-stationary time-series dataset, we would be seeing definitive level of autocorrelations in the datasets shown below.
# Correlogram on non-stationary time-series
ggAcf(cr_ts) + ggtitle("ACF of a non-stationary time series") + theme_bw()
We could also detect trend from the above correlogram, shown by the deliberate reduction of ACF values while the lags are increasing. This supports our initial inference of the existence of trend in our dataset.
Furthermore, we can see in the partial autocorrelation plot below that a majority of the correlations are not significant except on lag 1 and 3, which indicates an autoregressive term of order 2.
#PACF
ggPacf(cr_ts) + ggtitle("Partial ACF of a non-stationary time series") + theme_bw()
We have split our data into 75% train and 25% test data. Let us use the benchmark forecasting method of simple tools such as the average, naive, seasonal naive, and drift method.
# Plot analysis
autoplot(window(cr_ts, start=2005)) +
autolayer(fit1, series="Mean", PI=FALSE) +
autolayer(fit2, series="Naïve", PI=FALSE) +
autolayer(fit3, series="Seasonal Naïve", PI=FALSE) +
autolayer(fit4, series="Drift", PI=FALSE) +
xlab("Year") + ylab("Percentage Point") +
ggtitle("42 months Forecast: 2016-2019") +
guides(colour=guide_legend(title = "Forecast"))
From the chart above, we are predicting the values using the benchmark tools for 42 months, 25% of the data, while we also plot the actual data for those 42 months. Intuitively, Both Naïve and Drift method proves as better predictor. Naïve method predicted that the cash rate will be flat which was true for the majority of the 42 months, while Drift method also follows the downward trend of the actual data.
Let us check the accuracy of these methods.
# Check Benchmark Forecasting Accuracy
# Mean
accuracy(fit1, cr_ts_test)
## ME RMSE MAE MPE MAPE MASE
## Training set 8.558839e-17 1.578359 1.372661 -16.12982 38.52908 1.447425
## Test set -2.765129e+00 2.769106 2.765129 -180.72981 180.72981 2.915735
## ACF1 Theil's U
## Training set 0.9823924 NA
## Test set 0.7338497 45.39631
# Naïve
accuracy(fit2, cr_ts_test)
## ME RMSE MAE MPE MAPE
## Training set -0.02462121 0.202961 0.08143939 -0.8463268 1.999091
## Test set -0.45121951 0.474984 0.45121951 -30.1509872 30.150987
## MASE ACF1 Theil's U
## Training set 0.08587509 0.4049244 NA
## Test set 0.47579574 0.7338497 7.971388
# Seasonal Naïve
accuracy(fit3, cr_ts_test)
## ME RMSE MAE MPE MAPE MASE
## Training set -0.3367769 1.3660963 0.9483471 -15.31527 26.70219 1.0000000
## Test set -0.5243902 0.5466082 0.5243902 -34.72416 34.72416 0.5529518
## ACF1 Theil's U
## Training set 0.9698040 NA
## Test set 0.6790777 9.1445
# Drift
accuracy(fit4, cr_ts_test)
## ME RMSE MAE MPE MAPE
## Training set -1.749320e-16 0.2014621 0.1000918 -0.182059 2.486708
## Test set 6.582594e-02 0.2378345 0.1988636 4.501153 13.193477
## MASE ACF1 Theil's U
## Training set 0.1055435 0.4049244 NA
## Test set 0.2096950 0.9500327 3.975859
As mentioned above, both Naïve and Drift has the lowest scores among the four benchmark methods, where the values for both methods are very close across all the accuracy scores in the training set but different in the test set with Naïve method has larger errors than Drift. Moreover, the seasonal method is not very useful in our datasets as we have established earlier that our macroeconomic data is not seasonal.
We can further explore the benchmark forecasting by checking the residuals of these methods, focusing on both Naïve and Drift as shown to have the best accuracy scores.
# Residual diagnostics
# Naïve
checkresiduals(fit2)
##
## Ljung-Box test
##
## data: Residuals from Random walk
## Q* = 121.98, df = 24, p-value = 4.33e-15
##
## Model df: 0. Total lags used: 24
# Drift
checkresiduals(fit4)
##
## Ljung-Box test
##
## data: Residuals from Random walk with drift
## Q* = 121.98, df = 23, p-value = 1.776e-15
##
## Model df: 1. Total lags used: 24
Both the Naïve and Drift have their mean for the residuals quiet close to zero, except in the event of the end of 2008 and early 2009, which corresponds to the GFC where the cash rate slumped in response to the crisis. Moreover, the histogram of both methods appears to be normally distributed if we exclude the small number of outliers, hence we can assume the residual variance as constant. This is also supported by the constant variation in the time plot of the residuals, with the exception of the sudden slump. Though, from the ACF chart, we can see some correlations in the residuals from both models.
These charts indicates that both methods constructs predictions that seems to account for most of the available information. However, we have to bear in mind that since the p-values of both methods are quite small, the residuals are distinguishable from a white noise series.
Let us start the advanced modelling with ARIMA and it’s variations.
### More advanced modelling
# ARIMA
cr_ts_arima = auto.arima(cr_ts_train, parallel = TRUE, stepwise = FALSE, approximation = FALSE)
checkresiduals(cr_ts_arima)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(2,1,0)
## Q* = 19.643, df = 22, p-value = 0.6054
##
## Model df: 2. Total lags used: 24
# ARIMA with Box-Cox transformation
cr_ts_arima_box = auto.arima(cr_ts_train, lambda = "auto", parallel = TRUE, stepwise = FALSE, approximation = FALSE)
checkresiduals(cr_ts_arima_box)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(1,1,2)
## Q* = 18.987, df = 21, p-value = 0.586
##
## Model df: 3. Total lags used: 24
Inside forecast package, we have auto.arima function that we have used above to fit into the model. The function auto.arima relies on three inputs(p,d,q) which represents the lags, the degree of differencing, and the moving average respectively.
Moreover, we have specified the stepwise and approximation parameters to FALSE as we are only analysing just one time series. The parallel parameter is also set for parallel computing.
Using the checkresiduals function, we can see that both auto.arima and auto.arima with Box.Cox transformation has a similar theme on their charts. Both time plots moved away from zero during 2009 due to the sudden flop of the interest rates. Otherwise, the residuals are quite close to zero.
There are few correlations from the ACF graph (bottom left) and the histogram (bottom right) are normally distributed if we exclude the one outlier. Pay attention to the p-values, as it is not significant, we can conclude that the residuals are not distinguishable from the white noise
Let’s check manual ARIMA model
# ARIMA Manual
cr_arima_manual = Arima(cr_ts_train, order=c(10,1,10))
checkresiduals(cr_arima_manual)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(10,1,10)
## Q* = 9.3654, df = 4, p-value = 0.05259
##
## Model df: 20. Total lags used: 24
You can see that the time plots generate the same theme as previous arima functions, as well as the histogram. We can also see improvement in the correlogram plot as there is no correlation in the residuals. The p-value is still not significant, even though it is very close.
Next, let’s use three exponential smoothing models.
# Simple exponential smoothing
cr_ses <- ses(cr_ts_train, h=42)
checkresiduals(cr_ses)
##
## Ljung-Box test
##
## data: Residuals from Simple exponential smoothing
## Q* = 122.88, df = 22, p-value = 5.551e-16
##
## Model df: 2. Total lags used: 24
# Holt's linear trend method - because it also smoothing trend
cr_holt <- holt(cr_ts_train, h=42)
checkresiduals(cr_holt)
##
## Ljung-Box test
##
## data: Residuals from Holt's method
## Q* = 24.808, df = 20, p-value = 0.2089
##
## Model df: 4. Total lags used: 24
# ETS method
cr_ets <- ets(cr_ts_train)
summary(cr_ets)
## ETS(M,Ad,N)
##
## Call:
## ets(y = cr_ts_train)
##
## Smoothing parameters:
## alpha = 0.8728
## beta = 0.2908
## phi = 0.8
##
## Initial states:
## l = 5.3376
## b = -0.1823
##
## sigma: 0.0414
##
## AIC AICc BIC
## 180.0542 180.7209 197.3963
##
## Training set error measures:
## ME RMSE MAE MPE MAPE
## Training set -0.009639467 0.179153 0.09739388 -0.2518063 2.421613
## MASE ACF1
## Training set 0.1026986 0.1894643
checkresiduals(cr_ets)
##
## Ljung-Box test
##
## data: Residuals from ETS(M,Ad,N)
## Q* = 21.736, df = 19, p-value = 0.2976
##
## Model df: 5. Total lags used: 24
Simple exponential smoothing has a lot more correlations in the residuals in the correlogram suggesting that there is information left in the residuals which should be used in computing forecasts. The histogram is not normal at all.
However, the Holt’s linear model shows a very normal distribution on the histogram as this model also smooth trends within the time series, while the ETS model has its left tail a little too long for normal distribution. Both models have lack of correlations on their residuals as shown on their correlogram.
ETS model uses maximum likelihood for its default estimation method not rather than minimum sum of squares.
Let’s see the chart for the estimated states over time.
# Visually check estimated check overtime
autoplot(cr_ets)
The below is showing us the residuals and one-step forecasts errors from the ETS model.
# Residuals and one-step forecasts errors
cbind('Residuals' = residuals(cr_ets), 'Forecast errors' = residuals(cr_ets, type = 'response')) %>%
autoplot(facet=TRUE)
Let us both visually forecasts these advanced models and check their accuracy.
## Forecast charts for advanced models
autoplot(window(cr_ts, start=2005)) +
forecast::autolayer(fc_arima, series="Auto Arima", PI=FALSE) +
forecast::autolayer(fc_box, series="Box-Cox", PI=FALSE) +
forecast::autolayer(fc_man, series="Manual Arima", PI=FALSE) +
forecast::autolayer(fc_ses, series="SES", PI=FALSE) +
forecast::autolayer(fc_holt, series="Holt's Trend", PI=FALSE) +
forecast::autolayer(fc_ets, series="ETS", PI=FALSE) +
xlab("Year") + ylab("Percentage Point") +
ggtitle("42 months Forecast: 2016-2019") +
guides(colour=guide_legend(title = "Forecast"))
The above seemingly shows only two lines but in fact it is manual ARIMA is the only one predicting differently compared to the rest of the models. This could suggest that all the tested models above would do a better job than manual ARIMA in predicting the cash rate.
Let us check the accuracy below.
## Checking accuracy on the advanced models
# Auto ARIMA
accuracy(fc_arima,cr_ts_test)
## ME RMSE MAE MPE MAPE
## Training set -0.009146908 0.1703826 0.09563574 -0.2219174 2.382698
## Test set -0.451219512 0.4749840 0.45121951 -30.1509872 30.150987
## MASE ACF1 Theil's U
## Training set 0.1008447 0.009837688 NA
## Test set 0.4757957 0.733849703 7.971388
# Auto ARIMA with Box-Cox Transformation
accuracy(fc_box,cr_ts_test)
## ME RMSE MAE MPE MAPE
## Training set -0.01900226 0.1780056 0.09743696 -0.4721108 2.39897
## Test set -0.45034049 0.4740897 0.45037668 -30.0928142 30.09462
## MASE ACF1 Theil's U
## Training set 0.102744 0.2322395 NA
## Test set 0.474907 0.7336652 7.956394
# Manual ARIMA
accuracy(fc_man,cr_ts_test)
## ME RMSE MAE MPE MAPE
## Training set -0.01479111 0.1559287 0.09623257 -0.4811149 2.492096
## Test set -0.65699891 0.6836289 0.65699891 -43.7438498 43.743850
## MASE ACF1 Theil's U
## Training set 0.1014740 -0.008672941 NA
## Test set 0.6927832 0.769139298 11.46714
# Simple Exponential Smoothing
accuracy(fc_ses,cr_ts_test)
## ME RMSE MAE MPE MAPE
## Training set -0.02443895 0.202205 0.0808356 -0.8400628 1.984271
## Test set -0.45121951 0.474984 0.4512195 -30.1509872 30.150987
## MASE ACF1 Theil's U
## Training set 0.08523841 0.4051082 NA
## Test set 0.47579574 0.7338497 7.971388
# Holt's Linear Trend Method
accuracy(fc_holt,cr_ts_test)
## ME RMSE MAE MPE MAPE
## Training set -0.001222241 0.1804123 0.1105777 0.1635016 2.811532
## Test set -0.450641340 0.4743733 0.4506528 -30.1122877 30.112859
## MASE ACF1 Theil's U
## Training set 0.1166004 -0.01627111 NA
## Test set 0.4751981 0.73360284 7.961053
# ETS (Error Trend Seasonal) Model
accuracy(fc_ets,cr_ts_test)
## ME RMSE MAE MPE MAPE
## Training set -0.009639467 0.1791530 0.09739388 -0.2518063 2.421613
## Test set -0.449529288 0.4732865 0.44963611 -30.0394839 30.044825
## MASE ACF1 Theil's U
## Training set 0.1026986 0.1894643 NA
## Test set 0.4741261 0.7335335 7.942943
# Drift
accuracy(fit4, cr_ts_test)
## ME RMSE MAE MPE MAPE
## Training set -1.749320e-16 0.2014621 0.1000918 -0.182059 2.486708
## Test set 6.582594e-02 0.2378345 0.1988636 4.501153 13.193477
## MASE ACF1 Theil's U
## Training set 0.1055435 0.4049244 NA
## Test set 0.2096950 0.9500327 3.975859
We can see that, with reference to Drift accuracy scores, all the errors for the advanced models are very similar and in fact smaller than Drift errors in the training set. However, just like Naïve model, it does not work very well in the test set, suggesting again a overfitting on these advanced modelling.
In this post, we have compared the simple forecasting methods with the more advanced methods using a real macroeconomic data as announced by the RBA, the cash rate. The advanced modelling techniques have not proven themselves to be better than the simple benchmark modelling in the forecasting, though it shows similar performance in the training set.
Therefore, if we want to model forecasting for the ASX200 using macroeconomic data, we should simply use the Drift method, as shown by its superior performance both on training and test set.
However, we need to remember that the time series data being used in this piece is one without seasonality. Perhaps the performance would yield different results with time-series data containing seasonality.
Asprem, M. 1989, ‘Stock prices, asset portfolios and macroeconomic variables in ten European countries’, Journal of Banking & Finance, vol.13, viewed 28 October 2019, https://doi.org/10.1016/0378-4266(89)90032-0.
Cheung, Y. & Ng, L. 1998, ‘International evidence on the stock market and aggregate economic activity’, Journal of Empirical Finance, vol.5, viewed 28 October 2019, https://doi.org/10.1016/S0927-5398(97)00025-X.
Hyndman, R. J. & Athanasopoulos, G. 2018, Forecasting: principles and practice, electronic book, 2nd edition, OTexts:Melbourne, viewed 28 October 2019, https://otexts.org/fpp2/
Gupta, J, Chevalier, A., & Sayekt, F. 2001, Fuzzy sets in Management, Economics and Marketing, electronic book, World Scientific, viewed 28 October 2019, https://doi.org/10.1142/9789812810892_0010.