IBM Stock prices Forecasting:Univariate Analysis_Eviews

Department of Economics,University of Pune,India
nishantup@gmail.com

This is a univariate time series forecasting exercise done entirely in Eviews.

1 Preliminary Data Exporation

Let’s import the file as shown below:

File->Open–>Foreign Data as workfile

Data description:

IBMclosingprices.xlsis a data set of closing daily prices of IBM stock from 4th January, 2000 to 20th August, 2002 thereby having a total no of observations equal to 739 including NAs.The closing prices are given uptill 20/08/2002 (687 obs).There are 52 NAs obs from 31/08/2002 which will be further used for forecasting.

There are two things to be done either before importing or just after the data is imported before any manipulation/analysis is done:

• As these are daily closing prices and there are days which do not appear may be due to holidays etc., we have to create a continuous trend variable either before importing the data or in eviews itself (We will make the time periods equally spaced)..I made the changes in excel itsef by replacing the dates column obs by nos starting from 1 till 739 (made it continuous with equally spaced intervals).This column will be my time trend variable.

In Eviews it can be done as:

Generate series as: \(t = @trend--1\)

• Second thing is to tackle NAs by either removing the NAs in excel and then importing the series with 687 obs or after importing restructuring it.Normally i hate to touch original data set provided. So lets do that in Eviews after the data is imported.

click Proc tab on workfile–>Structure/Resize current page–>unstructured/undated–>data range…change it to 687

2 Stationarity Assessment

Let’s Identify whether the given time series is stationary or not.Broadly there are two ways:

• Graphical method
• various Statistical tests for autocorelation

2.1 Graphical Method

Graph the IBM time series:

2.1.1 Method 1: Stationarity identification from time series plot

Visually inspecting the series we can say that:

• The series in mean non-stationary as the mean is changing (i.e. mean not parallel to the X axis). We can see a perceptible downward trend.
• The series is variance non-stationary as we can see different dispersion widths between two successive peaks or troughs.
• Auto covariance stationarity can only be gauged using correlogram and related statistical tests for autocorrelation significance.

Before going further, let’s see how the log transformed series looks in comparison to the original series:

First generate the log transformed series

Quick–>Generate Series–> lnclose=log(close)

Lets plot both curves close & lnclose by selecting both series using cntrl and opening them as group–>then plot a line graph by stacking them in single graph

We can see that the log transformed series did not smoothen the curve and it is almost exactly superimposed (with slight displacement) on the original series i.e.log series still remains mean & variance Non stationary. But henceforth, we will still use the lnclose variable as our dependent var.

Reasons to take log transformation:

• Used to handle exponential growth of a series
• Used to stabilize the variability

Values must all be positive before the log is taken

• If not all values are positive a positive constant can be added to every data point

2.1.2 Method 2: Stationarity assessment using Correlogram (ACF/PACF plots)

Open lnclose series object–>correlogram–>level–>172 lags–>OK–>save to disk as png

Note:No of lags is normally 1/3rd to 1/4th of total no of obs.Here we are taking 1/4th of 687 which is \(\small \approx\) 172 lags.If suppose we see that our total obs are on lesser side, then we could go for 1/3rd criteria.

From the ACF plot we can see the gradual exponential decay and it reaches zero only around 114th lag. Corresponding auto correlation values are also very high for near lags (If the autocorrelation dies out slowly this indicates that the process is non-stationary).

As the ACF is decaying gradually and PACF drops to near zero after lag 1 (1 significant lag before dropping to zero), it could be an AR(p) model.

Representative time series (Wold decomposition) identification from ACF/PACF plots:

2.2 Non Graphical Methods:Statistical tests

Following methods relate to testing the significance of the auto correlation values.

2.2.1 Method 3: Barlett test for autocorrelation

Sampling distribution of autocorrelations

\[H_0:\rho_k=0\] \[H_a:\rho_k\ne 0\]

Let’s construct the 95% CI for the representative white noise series for similar sample size as our given IBM series is and check for each auto correlation values, whether they fall in the interval or not.

\(0 ± 1.96 *\sqrt{\frac{1}{n}}\)

For our univariate series, the n=687

0 ± 1.96* 0.03815238 = ± 0.0748

i.e. CI–>(-0.0748, 0.0748)

As we can see that ACF values for upto lags 100 (values are in correlogram plot) fall outside the representative white noise CI, hence we can conclude that majority of the sample auto correlations are significantly different from zero, hence the given series is NON stationary. As we can check, even at lag 20 the estimated \(\rho_{20}=0.675\) is statistically significant at the 5% level.

2.2.2 Method 4: Box-Pierce test for autocorrelation (Q statistics)

Returning to the IBM stock prices series, the value of the Q statistic up to lag 172 is about 18444. The probability of obtaining such a Q value under the null hypothesis that the sum of 172 squared estimated autocorrelation coefficients is zero is practically zero, as the last column of that figures shows. Therefore, the conclusion is that the IBM closing stock prices series is non-stationary, therefore reinforcing our hunch from visual inspection that series may be non-stationary.

2.2.3 Method 5: Unit Root test for autocorrelation :DF/ADF tests

Unit Root tests

Note:EViews refers to both the Dickey-Fuller and the Augmented Dickey-Fuller tests as ADF tests. You will face two practical issues in performing the ADF test.

• First, you will have to specify the number of lagged first difference terms to add to the test regression (selecting zero yields the DF test; choosing numbers greater than zero generates ADF tests). The usual (though not particularly useful) advice is to include lags sufficient to remove any serial correlation in the residuals.
• Second, EViews asks you whether to include other exogenous variables in the test regression. You have the choice of including a constant, a constant and a linear time trend, or neither in the test regression.

Openlncloseobject–>View–>Unit Root test–>ADF–>level–>Intercept and Trend

???

For the model, the estimated \(\tau\) value is -2.109594 (prob>0.05), which in absolute value is below even the 10 percent critical value of -3.130503. Since, in absolute terms, the former is smaller than the latter, our conclusion is that the IBM time series is NON stationary (as we fail to reject the null of series having the UNIT ROOT).

Note: The corresponding t statistic of LNCLOSE(-1) is same as \(\tau\) value but the p value is different which is indicative of non stationarity which will not be the case for stationary series i.e. both values will be same.

3 Adjusting for Non stationarity

3.1 Detrending

How this concept works, let’s see:

First regress the depending var on time component

\(lnclose=a_1 + a_2 t + \epsilon_t\)

rearranging, we get

\(\epsilon_t= lnclose - a_1 - a_2 t\) , which is the detrended series (TSP) as we are removing the intercept as well the time component ‘t’ which is the trend.

So in Eviews, We are regressing lnclose on time variable with an intercept

**lnclose=C(1)+C(2)*time**

Let’s see the residual graph ‘\(\epsilon_t\)’ to see whether the deterministic trend has been removed (detrending)

View–>Actual Fitted residual–>Actual Fiited residual Graph

We can still see some trend leftover especially on the right hand side of the plot, which is stochastic trend. But on the whole we have somewhat achieved mean stationary process as the residual series is oscillating around the zero mean. But we still have variance non stationarity as the dispersion between successive troughs or peaks are different.

3.2 Differencing

Let’s now create a first differenced series using

Quick–>Generate Series–>Dlnclose=d(lnclose)

Graphing Dlnclose:

We seem to have achieved both mean as well as variance stationarity after first differencing.

Let’s go through the correlogram to assess the covariance stationarity of the differenced series.

Analysing the ACF and PACF plots we can see a distinct sine/cosine type curve of decay.We can see that majority of bars in ACF & PACF plot are within the 5% critical band implying stationary process.Let’s confirm this through serial correlation significance tests.

Let’s employ the statistical tests:

3.2.1 Barlett test

\[H_0:\rho_k=0\] \[H_a:\rho_k\ne 0\]

Let’s construct the 95% CI for the representative white noise series for similar sample size as our given IBM series is and check for each auto correlation values whether they fall in the interval or not.

\(0 ± 1.96 *\sqrt{\frac{1}{n}}\)

For our univariate series, the n=687

0 ± 1.96* 0.03815238 = ± 0.0748

i.e. CI–>(-0.0748, 0.0748)

As we can see that ACF values for lags 4,18,22,25,123,140 fall outside the representative white noise CI, hence we can conclude some of sample auto correlations are significantly different from zero.Since we have 166 obs (more than 95% of the sample autocorrelations falling in the two standard error band(white noise CI) hence we can say we have a stationary series.

3.2.2 Box Pierce Q statistic test

Correlogram of Differenced IBM stock prices series Dlnclose, the value of the Q statistic up to lag 172 is about 165.68 which is statistically insignificant. Therefore, we fail to reject the null that all the population \(\rho_k\) upto lags 172 are zeros.Hence Dlnclose is stationary .

3.2.3 Dickey-Fuller(DF) and the Augmented Dickey-Fuller(ADF) Unit Root tests

Open Dlnclose object–>View–>Unit Root test–>ADF–>level–>None

From the ADF test output,the estimated \(\tau\) value is -27.69210 (prob<0.00001), which in absolute value is way way higher than the 1 percent critical value of -2.568342. Since, in absolute terms, the former is considerably higher than the latter, our conclusion is that the first differenced IBM time series is stationary (as we reject the null in favour of alternate that the series is NOT having a UNIT ROOT).The corresponding t statistic of DLNCLOSE(-1) is same as \(\tau\) value and also the p values are same which is indicative of stationary series.

So we have achieved a stationary series: Dlnclose We will be using this staionarised series for identifying the appropriate time series specification and finally for forecasting.

4 Identifying the appropriate model specification:Box-Jenkins methodology

Once having found the stationary form of the series, we now look at the ACF and PACF plots to determine the form of the ARIMA(p,d,q) model for the stationary series.

4.1 Introduction to Box & Jenkins procedure

Given an observed time series, the aim of the Box and Jenkins (1970) procedure is to build an ARIMA model.

The Box and Jenkins (hereafter, BJ) procedure is an empirical approach that necessarily:-

• identifies the order of the model (among the infinite set of ARIMA(p,d,q) models identifies the one that is more adapt to the data at hand),
• estimates the parameters
• assesses the goodness of fit of the model

Box and Jenkins procedure’s steps:

• Preliminary analysis: create conditions such that the data at hand can be considered as the realization of a gaussian, stationary, invertible stochastic process.
• identification: specify the orders p, d, q of the ARIMA model so that it is clear the number of parameters to estimate. Recognizing the behaviour of empirical autocorrelation functions plays an extremely important role.
• Estimate: efficient, consistent, sufficient estimate of the parameters of the ARIMA model (maximum likelihood estimator).
• Diagnostics: check if the model is a good one using tests on the parameters and residuals of the model. Note that also when the model is rejected, still this is a very useful step to obtain information to improve the model.
• Usage of the model: if the model passes the diagnostics step, then it can be used to interpret a phenomenon, forecast…

This procedure is a iteration-interaction procedure

The Box-Jenkins model assumes that the time series is stationary. Box and Jenkins recommend differencing non-stationary series one or more times to achieve stationarity. Doing so produces an ARIMA model, with the “I” standing for “Integrated”.

4.1.1 I: Model Identification

The following table summarizes how one can use the sample autocorrelation function for model identification:

Our aim is to find the appropriate ARMA(p,d,q) where p,q is the order of AR & MA processes respectively & ‘d’ here denotes the order of Integration or the level of differencing.In our case since we achieved stationarity after differencing once hence d=1 in our case.

So our model would be of the form ARIMA(p,1,q).We have to just identify what p and q orders would be.

In our differenced series Dlnclose (which is stationary), we found from the correlogram, significant lags at 4,18,22,25,123,140.These lags may have Unit Roots.

Let’s start with the ARIMA(p,1,0) model with significant lags as covariates:

Quick–>Estimate Equation

d(lnclose) c ar(4) ar(18) ar(22) ar(25) ar(123) ar(140) –>ARall

Now we do the ARIMA(0,1,q) model:

d(lnclose) c ma(4) ma(18) ma(22) ma(25) ma(123) ma(140) –>MAall

Observing the model outputs, we can see that ar(18) coefficient is insignificant in ARall.Also AIC and SIC is lower for ARallvis a vis MAall (lower the AIC and SIC values, better the model).This tells us that overall AR model is better than the MA model.

Something about negative AIC’s:

Lets remove ar(18) and remodel ARall

d(lnclose) c ar(4) ar(22) ar(25) ar(123) ar(140) –>ARless18

We get all significant ar coefficients and also there is slight improvement in AIC/SIC & therefore ARless18 is a better model.

Our main intention is to what extent we can improve our Auto Regressive model.Our focus is more on AR model rather than MA stems from the fact that the AR terms are simply lagged Y terms. They hold a lot of information about the series. The MA terms are calculated from the model and may hold less information.

In general, one AR term will often explain more than one MA term. By the law of parsimony, we want to have the least number of terms in our model. Therefore, we look at the AR terms first.

Let’s run thro’ the correlogram again of Dlnclose. and start checking for each significant lag one by one.

First we analyse lag 140:

We notice that at lag 140 we see a MA process with spike at lag140 in ACF and cutsoff at 141st lag.The PACF shows a exponential decay.This looks like a MA process.

Lets remodel the ARless18by replacing the ar(140) with ma(140)

d(lnclose) c ar(4) ar(22) ar(25) ar(123) ma(140) –>ARless18MA140

We get a highly significant ma(140)coefficient (p<0.0001).This model has the lowest AIC/SIC among the previous models.

Lets see if we can achieve better model than ARless18MA140 model.

Now lets go for lag 123:

Since we are getting damped oscillations (sine wave structure), it is suggestive of a mixed AR & MA process.Hence we check for both AR 123 as well as MA 123.So all we now require to do is to incorporate for MA 123 as earlier models has already taken into account AR 123 lag.

Lets remodel the ARless18MA140by replacing the ar(123) with ma(123)

d(lnclose) c ar(4) ar(22) ar(25) ma(123) ma(140) –>ARless18MA123MA140

We get a highly significant ma(123)coefficient (p<0.0001).However the performance has come down as compared to ARless18MA140 model.

Now lets see for lag 25:

Since we are getting exponential decay for both ACF & PACF, it is suggestive of a mixed AR & MA process.Hence we check for both AR 25 as well as MA 25

Lets remodel the ARless18MA140by replacing the ar(25) with ma(25)

d(lnclose) c ar(4) ar(22) ma(25) ar(123) ma(140) –>ARless18MA25MA140

We get an insignificant ma(25)coefficient .However the performance has improved over ARless18MA123MA140 model but it is slightly lower than ARless18MA140 model.

Now lets see for lag 22:

Since we are getting exponential decay for both ACF & PACFat lag 22, it is suggestive of a mixed AR & MA process.Hence we check for both AR 22 as well as MA 22

Lets remodel the ARless18MA140by replacing the ar(25) with ma(25)

d(lnclose) c ar(4) ma(22) ar(25) ar(123) ma(140) –>ARless18MA22MA140

We get an insignificant ma(22)coefficient .However the performance has improved over ARless18MA25MA140 model but it is still slightly lower than ARless18MA140 model.

Now lets see for lag 4:

Since we are getting sinusoidal wave structure for both ACF & PACF, it is suggestive of a mixed AR & MA process.Hence we check for both AR 4 as well as MA 4.

Lets remodel the ARless18MA140by replacing the ar(4) term with ma(4)

d(lnclose) c ma(4) ar(22) ar(25) ar(123) ma(140) –>ARless18MA4MA140

We still get an insignificant ma(4)coefficient .Also the performance has reduced over ARless18MA22MA140 model and it is still slightly lower than ARless18MA140 model.

AIC/SIC summary of all the models

From the summary, we can conclude based on AIC/SIC that the best model is ARless18MA140 model.Henceforth we will utilise this model for further analysis and forecasting.

4.1.2 II: Residual Diagnostics

1.Lets see if we have normality of residuals:

We can see that we have got a near normal distribution for the residuals for the final selected ARless18MA140 model

2.Lets check for serial correlation:Breusch-Godfrey Lagrangian Multiplier(LM) test for serial correlation

The Durbin-Watson statistic can be difficult to interpret, we will perform a more general Breusch-Godfrey test for serial correlation in the residuals.This test is an alternative to the Q-Statistic for testing for serial correlation. It is available for residuals from OLS, and the original regression may include autoregressive (AR) terms.Unlike the Durbin-Watson Test, the Breusch-Godfrey Test may be used to test for serial correlation beyond the first order, and is valid in the presence of lagged dependent variables

open eqarless18ma140 object–>Residual diagnostics–>Serial Correlation LM test

The statistic labeled **Obs*R-squared is the LM test statistic for the null hypothesis of no serial correlation**. The probability value of 0.2532 >0.05 indicates NO serial correlation in the residuals.Hence we fail to reject the null of no serial correlation among the residuals.

3.Performing a test for Heteroskedasticity in EViews

From the Breusch-Pagan-Godfrey (BPG) Heteroskedasticity Test output we get the p value=0.2842 which means the at 5% significance level we do not reject the null of no heteroscedasticity.

4.1.2.1 Checking the stationarity of residuals

Lets first create a separate residual object resma140 by extracting the residuals from the ARless18MA140 model.

Quick–>Generate Series–>resma140=resid

4.1.2.1.1 Graphically

open resma140 object–>View tab–>Graph

open ARless18MA140 object–>View tab–>Actual Fitted Residual–>Actual fitted Residual graph

We can see that it is mean as well as variance stationary.So far so good.

4.1.2.1.2 Correlogram

Correlogram of resma140

From the plot we can see that majority of bars in ACF & PACF plot are withing the 5% critical band implying stationary process.Let’s confirm this through serial correlation significance tests.

4.1.2.1.3 Box pierce Q statistic

ForARless18MA140 model, the value of the Q statistic up to lag 172 is about 169.72 which is statistically insignificant. Therefore, we fail to reject the null that all the population \(\rho_k\) upto lags 172 are zeros.Hence ARless18MA140 is stationary .

4.1.2.1.4 ADF test

From the ADF test output,the estimated \(\tau\) value is -24.2455 (prob<0.00001), which in absolute value is way way higher than the 1 percent critical value of -2.5691. Since, in absolute terms, the former is considerably higher than the latter, our conclusion is that the first ARless18MA140 series is stationary (as we reject the null in favour of alternate that the series is NOT having the UNIT ROOT).We can also notice that t statistic is having the same statistic value as well as SAME p value which is an indicator of stationarity being achieved.

4.1.2.2 Check for Model Stability

So we have narrowed down to an ARMA model which we think is a good model for forecating.Lets check for the stability of this model wherein we check if the roots of AR and MA components lie within the UNIT ROOT CIRCLE.

open eqarless18ma140 object–>ARMA Structure–>Roots–>Table

From the ARMA structure table output we can see that all AR roots are inside the unit circle implying covariance stationarity.Also all MA roots lie inside the unit circle implying that ARMA structure is invertible.

5 Forecasting

This is a famous quote by Lao Tzu.

Forecasting is a process of estimation of an unknown event or parameter such as demand for a product. Now, it is commonly used to refer time series data, where time series is sequence of data points measured at successive interval. For example, these intervals could be weekly or monthly or quarterly. So, people may like to forecast quarterly demand or monthly demand or even weekly demand.

Forecasting method. You have a choice between Dynamic and Static forecast methods.

• Dynamic calculates dynamic, multi-step forecasts starting from the first period in the forecast sample. Here, previously forecasted values for the lagged dependent variables are used in forming forecasts of the current value (out of sample forecasting)
• Static calculates a sequence of one-step ahead forecasts, using the actual, rather than forecasted values for lagged dependent variables, if available (within sample forecasting).

Forecast Evaluation

So if we define

• forecast error as the difference between our forecast and the actual value, then the bias is just the average of those forecast errors.

• Standard error or standard forecast error is defined according to this formula. And you can think of it as a standard deviation for forecast errors.
• Mean squared error is just an average of squared forecast errors, and the root mean square error the root of this measure.
• Mean absolute error and mean absolute percentage error define to this formula. So you see that you take the average of absolute forecast errors or absolute value of forecast errors relative to the actual observations.

So you see that this is just the BIAS squared. So it tells us how far the mean of the forecast is from the mean of actual series. So you want this number to be as small as possible.

Variance proportion measures how far the variation of the forecast from variation of the actual series. And again, if your forecast is good, this number should be as small as possible.

And covariance proportion measures the remaining unsystematic forecasting error. So you want this number to be as close to 1 as possible, because biased proportion, plus variance proportion, plus covariance proportion sum up to 1.And so what properties it should satisfy is that the BIAS in variance proportion should be small, and most of the mean squared error proportion should be concentrated for covariance forecast.

5.1 Static Forecasting

Static forecast is a within sample forecast

To forecast lnclose from eqarless18ma140 model (We are using the final model and using that model for forecasting on the original series), push the Forecast button on the equation toolbar, or select Proc/Forecast

Evaluating the static forecast, we can see that Theil Inequality coefficient is way less than 1 implying a good forecast.Also the forecast is “good”, as the bias and variance proportions are small with most of the bias concentrated in the covariance proportions.

5.2 Dynamic Forecasting

Since dynamic forecasting is for out of sample method i.e. we are going to use the eqarless18ma140 model to forecast lnclose from 688 to 739 i.e. 52 out of sample obs from the available 687 obs.

Hence we have to expand the sample range from 1 to 687 to 739 i.e. add 52 more obs which will contain the forecasted values.

In order to do so we have to change/specify the sample/estimate range in our eqarless18ma140 model.

open eqarless18ma140 object–>proc–>specify/estimate–>Equation specification window opens up—>change the sample from 1 687 to 1 739

You will notice that sample range in the workfile has changed from 687 obs to 739 obs.

open eqarless18ma140 object–>forecast tab

Evaluating the dynamic forecast, we can see that Theil Inequality coefficient is way less than 1 implying a good forecast.Also the forecast is still “good”, as the bias and variance proportions are small with most of the bias still concentrated in the covariance proportions.But if you comapre with the static forecast, we see that bias proportion has increased and covariance proportion has reduced showing that static forecast is always a better forecast when compared to dynamic forecast.

Fan charts

Fan charts originated from the Bank of England in its Inflation Report of 1997. A fan chart consists of a mode forecast (not mean!) of the variable, which might be absent, as well as confidence bands that depict the uncertainty of the forecast over the time horizon. Fan charts shownn below are darkest closest to the mode forecast and get wider and lighter over the time horizon.The fan charts in EViews are not explicitly implemented, but we can use existing tools(codes) to construct standard errors and confidence intervals for the forecast.

so this is our fan chart. So you see that it plots lnclose up to 687th observation.It plots fan charts after that, and actual lnclose, which this data is not used for estimation.

It’s just to see how good our forecast is. So you see that the lightest green band corresponds to the 90% confidence interval. The next confidence interval corresponds to 60% confidence, so it’s slightly darker area. And the one in the middle corresponds to 30%. And the line in the very middle, it corresponds to the mode.

This writeup has been done in r markdown which gives me more flexibility but also interactivity and customisation.

This document can be downloaded from:

IBM Stock prices Forecasting:Univariate Analysis_Eviews