0% found this document useful (0 votes)

15 views249 pages

Univariate Time Series Forecasting Guide

This document discusses univariate time series models for economic forecasting, focusing on the specification, estimation, and forecasting of ARMA(p,q) models. It covers concepts such as stationarity, autocorrelation, and partial autocorrelation, providing examples and empirical applications. The chapter emphasizes the importance of past behavior in predicting future developments in time series data.

Uploaded by

saidul

We take content rights seriously. If you suspect this is your content, claim it here.

Available Formats

Download as PDF, TXT or read online on Scribd

0% found this document useful (0 votes)

15 views249 pages

Univariate Time Series Forecasting Guide

Uploaded by

saidul

We take content rights seriously. If you suspect this is your content, claim it here.

Available Formats

Download as PDF, TXT or read online on Scribd

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 1 / 248

Applied Economic Forecasting using Time Series

Methods

Eric Ghysels and Massimiliano Marcellino

Companion Slides - Chapter 5 Univariate Time Series Models

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 2 / 248

The Baseline Linear Regression Model - Outline
Overview
Representation
Model specification
Estimation
Unit Root Tests
Diagnostic Checking
Forecasting, Known Parameters
Forecasting, Estimated Parameters
Multi-steps (or Direct) Estimation
Permanent-transitory Decomposition
Exponential Smoothing
Seasonality
Examples With Simulated Data
Empirical examples
Concluding Remarks
c Eric Ghysels & Massimiliano Marcellino 2018 Edition 3 / 248
Overview

Box and Jenkins popularized the use of univariate time series

models for forecasting.
The key idea is to exploit the past behavior of a time series to
forecast its future development, which requires the future to be
rather similar to the past.
In this chapter we will consider the specification, estimation,
diagnostic checking, and use for forecasting of ARMA(p,q) models.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 4 / 248

A time series process is strictly stationary if

F {yt , . . . , yt+T } = F {yt+k , . . . , yt+T+k } , ∀ t, T, k. (1)

where F(·) indicates the joint density for a segment of length T of

the process y.
A time series process is weakly stationary if

E (yt ) = E (yt+k ) , ∀ t, k, (2)

Var (yt ) = Var (yt+k ) , ∀ t, k,
Cov (yt , yt−m ) = Cov (yt+k , yt−m+k ) , ∀ t, m, k.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 6 / 248

Stationarity

A strictly stationary process is also a weakly stationary process

provided F(·) has finite first and second moments.
A weakly stationary process is not in general strictly stationary.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 7 / 248

Representation

A weakly stationary process can be represented as:

yt = εt + c1 εt−1 + c2 εt−2 + . . . (3)

X∞ ∞
X
= ci εt−i = ci Li εt
i=0 i=0
= c (L) εt ,

L is the lag operator: Lεt = εt−1 and Li εt = εt−i

c0 = 1
the error process εt is uncorrelated across time and has a constant
variance, hence: εt ∼ WN(0, σ 2 ), meaning white noise with mean zero
and variance σ 2 .
The infinite moving average representation in (3) is known as the
Wold decomposition.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 8 / 248

Normalization

For this chapter we assume that the overall mean of the process yt is
zero.
This my not be true of all time series. However, if we denote the
unconditional mean by µ, all the analysis in this chapter goes
through with yt replaced by (yt − µ) whenever we have a non-zero
mean process.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 9 / 248

Representation

A problem with the model in (3) is that it has an infinite number of

parameters. However, in general, we can approximate c (L) via a
ratio of two finite polynomials, namely:

ψ (L)
c (L) =
φ (L)

ψ (L) = 1 - ψ1 L - . . . - ψq Lq
φ (L) = 1 - φ1 L - . . . - φp Lp

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 10 / 248

Representation

Because ofP weak stationarity φ (L) is invertible, i.e.,

φ(z) = 0 = pj=0 φj zj has all the roots outside the unit circle. Thus we
can rewrite (3) as
φ (L) yt = ψ (L) εt ,
or equivalently as

yt = φ1 yt−1 + . . . + φp yt−p + εt − ψ1 εt−1 − . . . − ψq εt−q .

This is a moving average autoregressive process of order p and

q, henceforth called ARMA(p,q).

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 11 / 248

Autocorrelation

Two useful tools to study ARMA processes are the autocorrelation

(AC) and partial autocorrelation (PAC) functions.
The autocovariance function reports the covariance of yt with its own
lags:

Cov (yt , yt−1 ) = γ (1) ,

Cov (yt , yt−2 ) = γ (2) ,
.. ..
. . (4)
Cov (yt , yt−k ) = γ (k) .

Then, the AC is defined as:

Cov (yt , yt−k ) γ (k)
AC(k) = p p = .
Var(yt ) Var(yt−k ) γ (0)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 12 / 248

Partial Autocorrelation

The kth value of the PAC measures the correlation between yt and yt−k ,
conditional on yt−1 , . . . , yt−k+1 .
The elements of the PAC can be considered as specific coefficients
in regression equations. In particular, they are

PAC (1) : coefficient of yt−1 in the regression of yt on yt−1 ,

PAC (2) : coefficient of yt−2 in the regression of yt on yt−1 , yt−2 ,
.. .. ..
. . .
PAC (k) : coefficient of yt−k in the regression of yt on yt−1 , . . . , yt−k .

Note that inserting a deterministic component (like a constant) in the

model changes the expected value of yt , while both the AC and PAC
remain the same.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 13 / 248

Representation

Before studying the characteristics of the class of ARMA processes, we

analyze the pure AR and MA processes seperately.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 14 / 248

Autoregressive processes

We assume that c (L) in (3) is invertible, i.e., c(z) = 0 = ∞ j

P
j=0 cj z has all
the roots outside the unit circle. Then we can rewrite (3) as an AR(∞):
∞
X
yt = φj yt−j +εt , εt ∼ WN(0, σ 2 ). (5)
j=1

Weak stationarity implies that the effect of yt−j onto yt fades as j

becomes large.
In practice we can reasonably approximate (5) with an AR(p):

yt = φ1 yt−1 + . . . + φp yt−p + εt . (6)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 15 / 248

Stationarity

If we view the AR(p) process in (6) as an approximation of (3), then it

is weakly stationary by definition.
If instead (6) is the data generating process (DGP), then for yt to be
weakly stationary the roots of φ (z) = 1 - φ1 z - . . . - φp zp = 0 are all
larger than one in absolute value.
Weak stationarity also guarantees that φ(L) can be inverted, namely,
yt can be represented as
1
yt = εt , (7)
φ (L)

which is an MA(∞) representation.

Adding a deterministic component in (6) does not create additional
complications.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 16 / 248

AR(1) Example

As an example, consider the following AR(1) process:

yt = φ1 yt−1 + εt (8)
φ (L) = 1 − φ1 L.

Weak stationarity requires |φ1 | < 1, since

1
φ (z) = 1 − φ1 z = 0 → z = .
φ1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 17 / 248

AR(1) Example

The MA (∞) of (8) representation is

yt = εt + φ1 yt−1 (9)
= εt + φ1 εt−1 + φ21 yt−2
= εt + φ1 εt−1 + φ21 εt−2+ φ31 yt−3
..
.
∞
X 1
= φi1 εt−i = εt .
φ(L)
i=0

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 18 / 248

AR(1) Example

We can use the MA(∞) representation to compute the mean and

variance of of yt .

E (yt ) = 0
∞ ∞
!
X X
φi1 εt−i Var φi1 εt−i

Var (yt ) = Var =
i=0 i=0
∞
X σ2
= φ2i
1 Var (εt−i ) =
i=0
1 − φ21

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 19 / 248

AR(1) Example

The autocovariance function appearing in (5) for the AR(1) is as

follows:

γ (1) = Cov (φ1 yt−1 + εt , yt−1 ) = φ1 Var (yt ) = φ1 γ (0) ,

γ (2) = Cov (φ1 yt−1 + εt , yt−2 ) = φ1 γ (1) = φ21 γ (0) ,
..
.
γ (k) = Cov (φ1 yt−1 + εt , yt−k ) = φ1 γ (k − 1) = φk1 γ (0) .

Therefore, for j = 1, 2, . . . , the AC is defined as:

Cov (yt , yt−j ) γ (j)

AC(j) = = = φj1 .
σy σy−j γ (0)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 20 / 248

AR(1) Example

Finally, for the elements of the PAC we have:

PAC (1) = φ1 ,
PAC (j) = 0, j > 1.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 21 / 248

AR(2) Example

Let us consider now the case of an AR(2) process:

yt = φ1 yt−1 + φ2 yt−2 + εt . (10)

To derive the weak stationarity conditions let us consider the

solutions of
φ (z) = 1 − φ1 z − φ2 z2 = 0,
which are q
−φ1 ± φ21 + 4φ2
z1,2 = .
2φ2

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 22 / 248

AR(2) Example

Therefore, to have |z1 | > 1 and |z2 | > 1, we need:

φ1 + φ2 < 1,
φ2 − φ1 < 1,
−φ2 < 1.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 23 / 248

AR(2) Example

The autocovariance function for an AR(2) process is given by:

γ (0) = φ1 γ (1) + φ2 γ (2) + σε2 ,

γ (1) = φ1 γ (0) + φ2 γ (1) ,
γ (2) = φ1 γ (1) + φ2 γ (0) ,
..
.
γ (k) = φ1 γ (k − 1) + φ2 γ (k − 2) .

The general element j of the autocovariance function can be

obtained as follows. First, one multiplies both sides of (10) with yt−j .
Second, one takes expectations on both sides and uses the fact that
yt−j and εt are uncorrelated and that γ(i) = γ(−i).

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 24 / 248

AR(2) Example

The γ’s can be solved for in terms of the AR parameters using a

system of equation for the first three lags of the autocovariance
function:
(1 − φ1 ) σε2
γ (0) = h i,
(1 + φ2 ) (1 − φ2 )2 − φ21
φ1 γ (0)
γ (1) = ,
1 − φ2
γ (2) = φ1 γ (1) + φ2 γ (0) ,

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 25 / 248

AR(2) Example
Thus, the autocorrelation function for an AR(2) process is:
φ1
AC(1) =
(1 − φ2 )
φ21
AC(2) = φ2 + ,
(1 − φ2 )
..
.
AC(k) = φ1 AC(k − 1) + φ2 AC(k − 2).

These are known as the Yule-Walker equations, and they can be

used to obtain estimators of φ1 and φ2 .
If we substitute AC(1) and AC(2) with their estimated counterparts,
the first two equations can be solved for φ̂1 and φ̂2 .
There are more efficient estimators than those based on the
Yule-Walker equations, but the latter can provide a simple initial
estimate.
c Eric Ghysels & Massimiliano Marcellino 2018 Edition 26 / 248
AR(2) Example

It can be easily shown that the first two lags of the PAC are different
from zero, while PAC(j) = 0 for j > 2.
Using methods similar to those seen in the examples considered so
far, we can calculate the AC and PAC for any AR (p) process.
The equations which define the first p lags of the AC function can
also be used to obtain the initial (Yule-Walker) estimators of the
parameters.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 27 / 248

MA Process

The qth-order moving average process, or MA(q), is defined as

εt ∼ WN 0, σε2 .

yt = εt − ψ1 εt−1 − . . . − ψq εt−q = ψ (L) εt , (11)

It can easily be shown that an MA(q) process is always weakly

stationary.
Its first two moments are

E (yt ) = 0,
1 + ψ12 + . . . + ψq2 σε2 = γ (0) ,

Var (yt ) =

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 28 / 248

MA Process

The autocovariance function is as follows:

γ(k) = Cov (yt , yt−k )

= Cov (εt − ψ1 εt−1 − . . . − ψq εt−q , εt−k − ψ1 εt−k−1 − . . . − ψq εt−k−q )

(−ψk + ψk+1 ψ1 + . . . + ψq ψq−k ) σε2 k = 1, . . . , q
= .
0 k>q

Dividing γ(k) by γ(0) yields the autocorrelation function.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 29 / 248

MA Process

Another relevant property for an MA process is invertibility, that is, the

possibility to represent an MA process as an AR(∞). This requires
that all the roots of ψ (z) = 0 are larger than one in absolute value.
When the MA process is invertible, we can write:
1
yt = ε t .
ψ (L)

The AR(∞) representation is useful to derive the PAC for an MA

process, which will coincide with that of the AR(∞) process.
The PAC elements decay (possibly non-monotonically) towards zero,
but are always different from zero except in the limit.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 30 / 248

MA(1) Example

As an example, let us consider the MA(1) process

yt = εt − ψ1 εt−1 . (12)

Its first two moments are

E (yt ) = 0,
1 + ψ12 σ2 ,

Var (yt ) =

The autocovariance function is

γ(1) = −ψ1 σ2

γ(k) = 0, k > 1,

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 31 / 248

MA(1) Example

For the AC we have AC(1) = −ψ1 /(1 + ψ12 ) and AC(k) = 0, k > 1.
The condition for invertibility is |ψ1 | < 1, with associated AR(∞)
representation

1 + ψ1 L + ψ12 L2 + ψ13 L3 + . . . yt = εt .

We can see that the PAC of the MA(1) process declines

exponentially.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 32 / 248

Remarks on Pure AR and MA Process

Note the different shapes of the AC and PAC for AR and MA

processes.
For an AR(p) process the AC decays but is always different from
zero, except in the limit.
For a MA(q) process the AC is only different from zero up to q lags.

The PAC of an AR(p) process is only different from zero up to p lags

The PAC of an MA(q) process decays but is always different from
zero, except in the limit.
These considerations suggest that the estimated AC and PAC could be
used to determine whether the underlying process is of the AR or MA
type, and what is the order of the lags. We will revisit this issue latter.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 33 / 248

ARMA Process

An ARMA(p, q) process is defined as

φ (L) yt = ψ (L) εt . (13)

The weak stationarity condition is, as for the AR processes,

φ (z) = 0 → |zi | > 1, i = 1, . . . , p,

where zi are the roots of the AR polynomial.

Likewise, the invertibility condition is, as for the MA processes,

ψ (z) = 0 → |zi | > 1, i = 1, . . . , q.

where zi are the roots of the MA polynomial.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 34 / 248

ARMA Process

For a stationary ARMA, we can write the infinite MA representation

as:
yt = φ−1 (L) ψ (L) εt = c (L) εt
Therefore, the expected value and variance are:

E (yt ) = 0 (14)
∞
X
Var (yt ) = σε2 c2i . (15)
i=0

The AC and the PAC are like those of an MA (∞) or AR (∞) , so both
decline exponentially.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 35 / 248

ARMA (1,1) Example
As an example, let us derive the autocovariance structure for an
ARMA(1,1) process:

yt = φ1 yt−1 + εt − ψ1 εt−1 . (16)

We can write γ0 and γ1 as

γ0 = Var (yt ) = E (φ1 yt−1 + εt − ψ1 εt−1 )2 =

= φ21 γ0 + σε2 + ψ12 σε2 − 2φ1 ψ1 E (yt−1 εt−1 )
= φ21 γ0 + σε2 + ψ12 σε2 − 2φ1 ψ1 σε2 ,

and

γ1 = Cov (yt , yt−1 ) = E (yt−1 (φ1 yt−1 + εt − ψ1 εt−1 ))

(1 − φ1 ψ1 ) (φ1 − ψ1 ) 2
= φ1 γ0 − ψ1 σε2 = σε .
1 − φ21

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 36 / 248

ARMA (1,1) Example

We can write γ2 as:

γ2 = Cov (yt , yt−2 ) = E (yt−2 (φ1 yt−1 + εt − ψ1 εt−1 )) = φ1 γ1 .

In general γk can be written as:

γk = Cov (yt , yt−k ) = φ1 γk−1 ,

The the kth autocovariance function is

AC(k) = γk /γ0 .

Therefore, for k > 2 the shape of the AC is similar to that of an AR(1)

process.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 37 / 248

ARMA (1,1) Example

It can be easily shown that the MA component plays a similar role for
the PAC, in the sense that for k > 2 the shape of the PAC is similar to
that of an MA(1) process.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 38 / 248

ARMA (1,1) Example

Finally, to obtain a first estimate of the ARMA(1, 1) parameters φ1 , ψ1 ,

σε2 we can use the Yule-Walker approach.
Form a system with the equations for γ0 , γ1 , γ2 .
Using sample counterparts γ̂0 , γ̂1 , γ̂2 , one obtains parameter
estimates by solving the system for φ̂1 , ψ̂1 , σ̂ε2 .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 39 / 248

Integrated Process

An integrated process yt is a non stationary process such that

(1 − L)d yt is stationary, where d is the order of integration. The
process is typically labeled I(d).
These are also called unit root processes, since integration is
associated with roots of the AR polynomial φ(L) exactly equal to one,
namely, one or more values for z in φ(z) = 0 are equal to one.
In the weakly stationary case, all the roots must be larger than one
in absolute value.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 40 / 248

Random Walk

The most common integrated process is the Random Walk (RW):

yt = yt−1 + εt , (17)

for which d = 1, since

(1 − L)yt = ∆yt = εt .

A RW can also be written as

yt = εt + εt−1 + εt−2 + . . . , (18)

so that the effects of a shock do not decay over time, contrary to the
case of a weakly stationary process (compare (18) with (9)).

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 41 / 248

Random Walk

From (18) we see that

E(yt ) = 0.
The variance is not properly defined:

Var(yt ) = Var(εt + εt−1 + εt−2 + . . .) → ∞.

The AC is also not properly defined, but the persistence of the

effects of the shocks is such that if we computed empirically the AC,
its elements would not decay as in the weakly stationary case but
would remain close to one at all lags.
For the PAC, we can easily see that PAC(1) = 1, while PAC(j) = 0 for
j > 1.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 42 / 248

From (18) we see that

E(yt ) = 0,

while the variance is not properly defined:

Var(yt ) = Var(εt + εt−1 + εt−2 + . . .) → ∞.

The AC is also not properly defined, but the persistence of the effects
of the shocks is such that if we computed empirically the AC, its
elements would not decay as in the weakly stationary case but would
remain close to one at all lags. For the PAC, we can easily see that
PAC(1) = 1, while PAC(j) = 0 for j > 1.
c Eric Ghysels & Massimiliano Marcellino 2018 Edition 43 / 248
Random Walk with a Drift

Inserting a deterministic component into an integrated process can

change its features substantially. As an example, let us consider the
RW with drift:
yt = µ + yt−1 + εt . (19)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 43 / 248

Random Walk with a Drift

Repeated substitution yields:

yt = µ + εt + µ + εt−1 + µ + εt−2 + . . .

so that

E (yt ) = E(µ + εt + µ + εt−1 + µ + εt−2 + . . .) → ∞,

Var (yt ) = Var(µ + εt + µ + εt−1 + µ + εt−2 + . . .) → ∞,

and the AC and PAC are as in the RW case.

Also in this case first differencing eliminates the non-stationarity,
since
(1 − L) yt = ∆yt = µ + εt .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 44 / 248

Final Notes

Recall that throughout the chapter we assumed for convenience that

the mean of yt is zero. It is clear that from the above discussion that
this assumption can be maintained for integrated processes once we
subtract µ from (1 − L) yt .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 45 / 248

ARIMA Process

An ARIMA(p, d, q) process is

φ (L) ∆d yt = ψ (L) εt . (20)

with ∆d ≡ (1 − L)d , whereas φ (L) and ψ (L) are polynomials in the lag
operator of order, respectively, p and q, while yt is I(d).
If we define xt = ∆d yt , xt is an ARMA(p, q).
Hence, the only additional complication with respect to the ARMA
case is the determination of the order of integration d.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 46 / 248

In order to specify an ARIMA model, we need to determine d, p and q,

namely, the order of integration and the lag length of the AR and MA
components. Three main approaches are available, based respectively
on the AC/PAC, testing, and information criteria.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 48 / 248

AC/PAC Based Specification

Using the available data on yt , we can easily estimate the AC and PAC
for different values of k. These values are typically graphed, with
increasing values of k on the right axis.
If the estimated values of both AC and PAC decay when k increases,
there is evidence that yt is weakly stationary and we can set d = 0.
If the AC declines very slowly and PAC(1) is close to one, then there
is evidence for a unit root. Therefore, we difference yt once, and
repeat the analysis with ∆yt . If necessary, we further difference yt ,
otherwise we set d = 1 and move to determine p and q.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 49 / 248

AC/PAC Based Specification

Evidence of pure AR(p)

PAC presents some peaks and then it is close to zero.
AC declines exponentially.
p is equal to the number of peaks (coefficients statistically different from
zero) in the PAC.
Evidence of pure MA(q)
AC present peaks and then it is close to zero.
PAC declines exponentially.
q is equal to the number of peaks in the AC.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 50 / 248

AC/PAC Based Specification

Often, we have both AR and MA components (an ARMA process).

The AC and PAC can provide an idea on the order of the AR and MA
components.
Identification of p and q from AC/PAC is more complex in the ARMA
case.
One possibility is to make a guess on p and q, estimate the ARMA
model, control whether the resulting residuals are WN, and if not go
back and try with a higher value for p and/or q.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 51 / 248

Testing Based Specification

Another diagnostic tool is testing for autocorrelations in the residuals.

If we have a correctly specified an ARMA(p, q) model, then we should
expect that the estimated errors, ε̂t , are temporally uncorrelated. Two
tests are:
Ljung-Box Q test (Ljung and Box (1978))
Box-Pierce test (Box and Pierce (1970))
For both tests, the null hypothesis is that the correlations in the
population from which the sample is taken are zero

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 52 / 248

Testing Based Specification

The Ljung-Box Q test is defined as:

h
X ρ̂2k
QLB = T (T + 2)
T −k
k=1

The Box-Pierce test is defined as:

h
X
QBP = T ρ̂2k ,
k=1

ρ̂k is the sample autocorrelation at lag k,.

h is the number of lags being tested.
For both, under the null Q follows a χ2(h) .
The degrees of freedom should be set to h − p − q.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 53 / 248

Testing Based Specification

Many regression-based model selection procedures exist for

ARIMA(p,d,0) type models, i.e., models without MA terms, based on
formal testing procedures, such as the Wald or LR statistics.
Testing the lag length of the MA component is more difficult.
Its presence prevents the use of OLS estimation.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 54 / 248

Testing for ARCH

We can also test for ARCH (autoregressive conditional

heteroskedasticity).
The test involves regressing ε2t onto a constant and h lags ε2t−1 . . . ,
ε2t−h .
Using a straightforward derivation of the LM test leads to the TR2 test
statistic, where the R2 pertains to the aforementioned regression.
Under the null hypothesis that there is no ARCH, the test statistic is
asymptotically distributed as chi-square distribution with h degrees
of freedom.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 55 / 248

Specification with Information Criteria

The third option for ARIMA model specification is information criteria

(IC), which combine goodness of fit with a penalty function related to
the number of model parameters.
Start with an ARIMA(pMAX ,d, qMAX ) wher it is assumed that pMAX > p0
and qMAX > q0 , where p0 and q0 are the true orders.
Compute the IC for all combinations of p, q, and possibly d, selecting
the model with the lowest IC.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 56 / 248

Specification with Information Criteria
A technical condition related to the penalty function guarantees the
procedure provides consistent results asymptotically, in the sense of
selecting the true orders p0 and q0 with probability approaching one
when the sample size T diverges.
This condition is satisfied by the BIC criterion

σε2 ) + (p + q) log(T)/T
BIC = log(b

but not by the AIC criterion

σε2 ) + 2(p + q)/T

AIC = log(b

In finite samples the relative performance of the BIC and AIC is not
uniquely defined.
The BIC generally leads to more parsimonious specifications than
the AIC

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 57 / 248

Concluding Remarks on Specification

The model specification methods can be applied jointly.

It may be that several model specifications perform similarly, in which
case we can proceed to forecasting with several specifications.
In general, lower orders of p and q are preferred.
Often, pure AR approximations are used in forecasting exercises.
The choice of d is not a major practical issue in practice.
If d0 = 1 and we set d = 0, then an AR root will be close to one.
If d0 = 0 and we set d = 1, then an MA root will be close to one,
canceling the effected of the imposed unit root.
Precise specification may be more relevant when considering
structural models.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 58 / 248

Once the order of the ARIMA(p,d,q) is specified, we need to estimate

the parameters in φ (L) and ψ (L) , where

φ (L) ∆d yt = ψ (L) εt .

If we define wt = ∆d yt , then

φ (L) wt = ψ (L) εt

and again
εt = ψ (L)−1 φ (L) wt .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 60 / 248

Estimation

Suppose the objective function to be minimized is the usual sum of

squared errors: X
St = ε2t .

If there is an MA component, St is non-linear in the parameters.

For example, for an MA(1) process we have

yt − ψ1 yt−1 − ψ12 yt−2 − . . . = εt .

We cannot find an analytical expression for the parameter estimators.

Non-linear least squares (NLS) estimators need to be used.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 61 / 248

Estimation

Typically, the estimator used for ARMA models is not NLS, but MLE.
iid
Make the additional assumptions that εt ∼ N(0, σ 2 ) for all t.
Collect all parameters to be estimated into a single vector θ.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 62 / 248

Estimation

An additional complication is dealing with the starting values. For

example, for an ARMA(1,1) how do we treat y0 and ε0 which are in the
specification for y1 ? Solution following Box and Jenkins (1976, p.211):
Compute the likelihood of an ARMA(p,q) conditional on y1 , . . . , yp ,
equal to the first p observations and ε1 = . . . = εMax(p,q) = 0.
The resulting sample log likelihood for an ARMA(p,q) model is
therefore:
T −p T
L(yp+1 , . . . , yT ; θ) = − log (2π) − log (σ 2 ) (1)
2 2
T 2
X ε t
−
2σ 2
t=p+1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 63 / 248

Estimation

Under mild conditions,

√ the MLE estimator is consistent and, when
multiplied by T, it has an asymptotically normal distribution,
centered at the true parameter values denoted by θ0 .
There is no analytical expression for the estimator in the presence of
an MA component.
For convergence of a numerical optimization method, we need a
reasonable initial value. One solution is to use the Yule-Walker
equations for an initial value.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 64 / 248

Consider the following model:

yt = Tt + zt
Tt = ν0 + ν1 t
zt = ρzt−1 + εt (2)

with εt ∼ WN(0, σ 2 ), and Tt is a deterministic linear trend.

If ρ < 1 then yt is I(0) about the deterministic trend Tt , sometimes
referred to as trend-stationary.
If ρ = 1 and ν1 = 0, then zt is a random walk and yt is I(1) with drift.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 66 / 248

Unit Root Tests

Let us start with ν0 = ν1 = 0 and hence Tt is zero, yielding:

yt = ρyt−1 + εt , (3)

or
∆yt = (ρ − 1) yt−1 + εt . (4)

We are interested in whether there is a unit root, i.e. ρ = 1. Hence

we are interested in the following hypothesis:

H0 : ρ=1 → yt is I(1) without drift

H1 : |ρ| < 1 → yt is I(0) with mean zero

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 67 / 248

Unit Root Tests

Under the unit root null, the usual sample moments of yt used to
compute the t−test do not converge to fixed constants.
Instead, Dickey and Fuller (1979) and Phillips (1987) showed that
the sample moments of yt converge to random functions of Brownian
motions on the unit [Link] particular:
T
d R1
X
T −3/2 yt−1 →σ 0 W(τ )dτ
t=1
T
d R1
X
T −2 y2t−1 → σ2 0 W(τ )2 dτ
t=1
T
d R1
X
T −1 yt−1 εt → σ 2 0 W(τ )dW(τ )
t=1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 68 / 248

Unit Root Tests

Using these results Phillips showed that under the unit root null:
R1
d 0 W(τ )dW(τ )
T(ρ̂ − 1) → R1
2
0 W(τ ) dτ
R1
d W(τ )dW(τ )
t − test(ρ = 1) → R01
( 0 W(τ )2 dτ )1/2

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 69 / 248

Unit Root Tests

√
The convergence rate of ρ̂ is not the standard T but instead T,
typically referred to as ‘super-consistency.
Neither ρ̂ nor t − test(ρ = 1) are asymptotically normally distributed.
The limiting distribution of t − test(ρ = 1) is called the Dickey-Fuller
(DF) (Dicky and Fuller (1979)) distribution, which does not have a
closed form representation. Consequently,
The p-values of the distribution must be computed by numerical
approximation or by simulation.
Critical values are available in tables or are automatically generated by
econometric software packages.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 70 / 248

Unit Root Tests

The T(ρ̂ − 1) and t − test(ρ = 1) statistics are called respectively

Dickey-Fuller (DF) normalized bias test and DF t-test.
The critical values of the DF distribution are generally larger than the
standard t distribution ones, such that using the standard critical
values will lead to rejecting the null hypothesis of a unit root too
often.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 71 / 248

Unit Root Tests with Deterministic Components

Including a deterministic trend to the model yields a different

distribution of the test statistics. The two most common trend cases
are a constant only and a constant and time trend

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 72 / 248

Unit Root Tests with Constant Only

The test regression for a constant only is

∆yt = ν0 + (ρ − 1) yt−1 + εt (5)

The hypotheses to be tested is:

H0 : ρ = 1 and ν0 = 0 → yt is I(1) without drift

H1 : |ρ| < 1 and ν0 6= 0 → yt is I(0) with non-zero mean

Under H0 the asymptotic distributions of the normalized biased and t

test statistics are influenced by the presence but not the coefficient
value of the intercept in the DF test regression.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 73 / 248

Unit Root Tests with Constant and Time Trend

The test regression for a constant and time trend is:

Deltayt = ν0 + ν1 t + (ρ − 1) yt−1 + εt (6)

The hypotheses to be tested are

H0 : ρ = 1 and ν1 = 0 → yt is I(1) with drift

H1 : |ρ| < 1 and ν1 6= 0 → yt is I(0) with deterministic trend

Under H0 the asymptotic distributions of the normalized biased and t

test statistics are again influenced by the presence but not the
parameter values of the intercept and trend slope coefficients in the
DF test regression.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 74 / 248

Unit Root Tests with p > 1

If the AR process is of order p , with p > 1, we can rewrite the process

as
∆yt = γyt−1 + ρ1 ∆yt−1 + . . . + ρp−1 ∆yt−p+1 + εt , (7)
where γ = φ (1) , and the ρ coefficients are related to the original φ
coefficients.
If one of the roots of φ (z) = 0 is z = 1, it follows that φ (1) = γ = 0
and therefore the coefficient of yt−1 in (7) will be equal to zero.
We can apply the same test as before, namely a t-test for γ = 0, or
normalized bias test, both known as Augmented Dickey Fuller (ADF)
tests.
The ADF test has the same distribution as the DF test.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 75 / 248

Unit Root Tests with p > 1

We can run test regressions:

∆yt = Dt β + γyt−1 + ρ1 ∆yt−1 + . . . + ρp−1 ∆yt−p+1 + εt (8)

where Dt is a vector of deterministic terms (constant, trend).

The specification of the deterministic terms depends on the
assumed behavior under the alternative hypothesis of trend
stationarity.
The ADF t-statistic and normalized bias statistic are based on the
least squares estimates of the above regression and have the same
limiting distribution as with p = 1.
In finite samples, the size and power of the ADF test are affected by
the number of lags included in (7), so that the determination of p is
relevant.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 76 / 248

Unit Root Tests

Other unit root tests include:

Phillips-Perron (PP) unit root tests (Phillips and Perron (1988)).
Elliot, Rothenburg, and Stock (1996) proposed a modification of the DF
test statistic based on the generalized least squares (GLS) principle,
known as the DF-GLS.

The size and power of unit root tests are also affected by the
presence of deterministic breaks in the model parameters.
More complex statistical procedures to allow for breaks when testing for
unit roots, see e.g.m Perron (1989) and Stock (1994).

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 77 / 248

Unit Root Tests with More than 1 Unit Root

There are also procedures to test for the presence of more than one
unit root, e.g., d = 2.
The simplest approach is first to test whether ∆yt has a unit root (so
that d = 2) and then, if the hypothesis is rejected, whether yt has a
unit root (so that d = 1).
At each step we can use the DF or ADF tests that we have
described.
One should use proper critical values that control for the sequential
applications of the procedure.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 78 / 248

Unit Root Tests

For forecasting, it is worth mentioning that if yt is I(d), so that

∆d yt = wt , (9)

then
yt = Σd wt , (10)
where, as an example,
t
X
Σwt = wi , (11)
i=−∞
t
X j
X
Σ2 wt = wj . (12)
i=−∞ i=−∞

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 79 / 248

Conclusions

There are many alternative procedures available for testing for unit
roots (or for stationarity), notably where the null is trend-stationarity
see e.g., Kwiatkowski, Philips, Schmidt, and Shin (1992).
It is difficult in general to outperform the ADF procedure, which has
the advantage of being easy to implement and understand.
For further details, there are a number of excellent surveys one can
read, including Campbell and Perron (1991), Stock (1994), and
Phillips and Xiao (1998).

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 80 / 248

We noted in Section 3 that the Ljung-Box and Box-Pierce statistics can

be used to test whether estimated residuals are white noise.
More informally, we can verify whether the estimated values of the
residual AC and PAC lie within the approximate asymptotic 95%
confidence
√ bands (around the null hypothesis of zero), which are
±1.96/ T.
We can use also the other tests for no correlation in the residuals
that we have seen in the context of the linear regression model, as
well as those for homoskedasticity and parameter stability, which all
have an asymptotic justification.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 82 / 248

Diagnostic Checking

If the model fails to pass the diagnostic checks, we need to

re-specify it, for example by increasing the AR or MA orders.
If the diagnostic checks are passed but some of the model
coefficients are not statistically significant, we could also assess
whether the diagnostic checks remain satisfactory for a more
parsimonious specification.
After determining the ARIMA model order, estimating its parameters,
and verifying the underlying assumptions are not rejected, we can
use it for forecasting purposes.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 83 / 248

As we have seen in the case of the linear regression model, the

optimal forecast of yT+h in the MSFE sense is

byT+h = E (yT+h |yT , yT−1 , . . . , y1 ) . (13)

We consider optimal linear forecasts for ARIMA(p, d, q) models,

which coincide with E (yT+h |yT , yT−1 , . . . , y1 ) if we assume that {εt } is
normal.
We also assume for the moment that the ARIMA parameters are
known.
We first discuss the general case, then present some examples, and
finally make a set of additional comments related to ARIMA
forecasts.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 85 / 248

General Formula

To calculate the optimal linear forecast for an ARIMA(p,d,q), we can

proceed as follows. We start by defining ∆d yt = wt , so that wt is
ARMA(p,q):

wT = φ1 wT−1 + . . . + φp wT−p + εT − ψ1 εT−1 − . . . − ψq εT−q .

From

wT+1 = φ1 wT + . . . + φp wT−p+1 + εT+1 − ψ1 εT − . . . − ψq εT−q+1 ,

it follows that

b T+1 = E (wT+1 |IT ) = φ1 wT + . . . + φp wT−p+1 − ψ1 εT − . . . − ψq εT−q+1 .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 86 / 248

General Formula

Similarly,

b T+2 = E (wT+2 |IT ) = φ1 w

w b T+1 + . . .
+φp wT−p+2 − ψ2 εT − . . . − ψq εT−q+2
... (14)
b T+h = E (wT+h |IT ) = φ1 w
w b T+h−1 + . . .
b T−p+h − ψh εT − . . . − ψq εT−q+h
+φp w

b T−j = wT−j if j ≤ 0 and there is no MA component for h > q.

where w

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 87 / 248

General Formula

We obtain the forecast of yT+h summing appropriately those for wT+j ,

j = 1, . . . , h. For example, for d = 1, we have

byT+h = yT + w
b T+1 + . . . + w
b T+h . (15)

To clarify the general derivation of the optimal ARIMA forecasts, we

now discuss a few examples.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 88 / 248

AR(1) Example

We consider a number of specific examples which are relatively simple

to work with beginning with the AR(1) process:

yt = φyt−1 + εt . (16)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 89 / 248

AR(1) Example

The formula in (14) simplifies to:

byT+1 = φyT ,
byT+2 = φbyT+1 = φ2 yT ,
...
byT+k = φk yT .

Since

yT+1 = φyT + εT+1 ,

yT+2 = φ2 yT + εT+2 + φεT+1 ,
...
yT+h = φh yT + εT+h + φεT+h−1 + . . . + φh−1 εT+1 .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 90 / 248

AR(1) Example

The forecast errors are

eT+1 = εT+1 ,
eT+2 = εT+2 + φεT+1 ,
...
eT+h = εT+h + φεT+h−1 + . . . + φh−1 εT+1 ,

and their variances are

Var (eT+1 ) = σε2 ,

1 + φ2 σε2 ,

Var (eT+2 ) =
...
1 + φ2 + . . . + φ2k−2 σε2 .

Var (eT+k ) =

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 91 / 248

AR(1) Example

Moreover, we have

lim byT+h = 0 = E (yt )

h→∞
1
lim Var (byT+h ) = σ 2 = Var (yt )
h→∞ 1 − φ2 ε

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 92 / 248

MA(1) Example

Now Let us consider the case of an MA(1) process,

yt = εt − ψ1 εt−1 (17)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 93 / 248

MA(1) Example

From the general formula in (14), we have

byT+1 = −ψ1 εT ,
byT+2 = 0,
...
byT+h = 0.

Since

yT+1 = εT+1− ψ1 εT ,
yT+2 = εT+2− ψ1 εT+1 ,
...
yT+h = εT+h− ψ1 εT+h−1 .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 94 / 248

MA(1) Example

The forecast errors are

eT+1 = εT+1
eT+2 = εT+2 − ψ1 εT+1
...
eT+k = εT+k− ψ1 εT+k−1 ,

with variances

Var (eT+1 ) = σε2

1 + ψ12 σε2

Var (eT+2 ) =
...
1 + ψ12 σε2

Var (eT+h ) =

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 95 / 248

MA(1) Example

Again, we have

lim byT+h = 0 = E (yt )

h→∞
lim Var (beT+h ) = 1 + ψ12 σε2 = Var (yt ) .

h→∞

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 96 / 248

ARMA(p,q) Example

Let us now consider an alternative derivation of the optimal forecast for

an ARMA(p,q), and show that it is equivalent to that in (14). We can
write the model as

a(L)yt = b(L)εt
p
X q
X
j
a(L) = 1 − aj L , b(L) = bj Lj , b0 = 1
j=1 j=0

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 97 / 248

ARMA(p,q) Example

Rewriting as an MA(∞) representation we have

yt = c(L)εt ,

with a(z)c(z) = b(z), namely,

p
X
ck − aj ck−j = bk , (18)
j=1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 98 / 248

ARMA(p,q) Example

We have that
∞
X
ŷT+h = ci+h εT−i .
i=0

The above equation implies that

p
X ∞
X p
X
ŷT+h − aj ŷT+h−j = (ci+h εT−i − aj ci+h−j εT−i ). (19)
j=1 i=0 j=1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 99 / 248

ARMA(p,q) Example

If we now use (18) in (19), we have

p
X ∞
X
ŷT+h − aj ŷT+h−j = bj+h εT−j , (20)
j=1 j=0

with bj+h ≡ 0, j + h > q, which is indeed equivalent to what we would

obtain with (14).
With this formula we can generate optimal forecasts for all h, p and q.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 100 / 248

Random Walk Example

Now consider the random walk

yt = yt−1 + εt , (21)

The optimal forecast is

byT+h = yT ,
for any h, and

eT+h = εT+1 + εT+2 + . . . + εT+h. .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 101 / 248

Random Walk Example

The variance of the forecast error is

Var (eT+h ) = hσε2 .

From these expressions, it follows that

lim byT+k = yT ,
k→∞
lim Var (byT+k ) = ∞.
k→∞

More generally, the presence of a unit root in the AR component

implies that the variance of the forecast error grows linearly over
time, while in the stationary it converges to the unconditional
variance of the variable.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 102 / 248

Additional Comments

What follows are a few interesting implications of what we have

derived:
Using the MA(∞) representation, it can be shown that the forecast
error is:
h−1
X
eT+h = yT+h − ŷT+h = cj εT+h−j . (22)
j=0

Therefore, even when using an optimal forecast, the h-steps ahead

forecast error is serially correlated. In particular, it is an MA(h-1)
process.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 103 / 248

Additional Comments

Moreover

E(eT+h ) = 0,
h−1
X
Var(eT+h ) = σε2 c2j ,
j=0

lim Var(eT+h ) = Var(yt ), Var(eT+h+1 ) − Var(eT+h ) = σε2 c2h ≥ 0.

h→∞

Thus the forecast error variance increases monotonically with the

forecast horizon. This result is no longer necessarily true if the
parameters are estimated.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 104 / 248

Additional Comments

If the error ε is Gaussian, so is the forecast error. In particular

y − ŷT+h
pT+h ∼ N(0, 1).
Var(eT+h )

We can use this result to construct (1 − α)% interval forecasts as:

p p
ŷT+h − cα/2 Var(eT+h ); ŷT+h + cα/2 Var(eT+h )

where cα/2 are critical values from the standard normal distribution.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 105 / 248

Additional Comments

From (22), for h = 1 we have

eT+1 = εT+1 , (23)

which can also be read as

εT+1 = yT+1 − ŷT+1, (24)

which provides an interpretation of the errors in the MA(∞ )

representation of a weakly stationary process.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 106 / 248

Additional Comments

Consider ŷT+h and ŷT+h+k , i.e., forecasts of yT+h and yT+h+k made in
period T. From (22) it can be easily shown that:
h−1
X
E(eT+h eT+h+k ) = σε2 cj cj+k ,
j=0

so that the forecast errors for different horizons are in general

correlated.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 107 / 248

Additional Comments

From (22) and the fact that εt is white noise, and considering the
predictor ŷT+h as an estimator (hence random), it follows that

Cov(ŷT+h , eT+h ) = 0.

Therefore,
Var(yT+h ) = Var(ŷT+h ) + Var(eT+h )
and
Var(yT+h ) ≥ Var(ŷT+h ).
Hence, the forecast is always less volatile than the actual realized
value.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 108 / 248

We now move on to the case where parameters are estimated.

If we use consistent parameter estimators, the formulas we have
derived for the optimal forecasts remain valid.
However, there is an increase in the variance of the forecast error
due to the estimation uncertainty.
To illustrate this we focus on only a few examples.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 110 / 248

AR(1) with Drift Example

The first case we assess is that of a stationary AR(1) with drift:

yt = µ + ayt−1 + εt . (25)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 111 / 248

AR(1) with Drift Example
If the parameters µ and a have to be estimated, by µ̂ and â, the
forecast error for h = 1 is:

eT+1 = εT+1 + (µ − µ̂) + (a − â)yT = εT+1 + (Θ − Θ̂)0 xT

1 µ − µ̂
xT = Θ − Θ̂ =
yT a − â

and
Var(eT+1 ) = σε2 + xT0 Var(Θ̂)xT , (26)
where
PT −1
µ̂ T yt−1
Var(Θ̂) = Var = σε2 E PT t=1
PT 2
â t=1 yt−1 t=1 yt−1
σε2 + µ2 (1 + a)(1 − a)−1

∼ − µ(1 + a)
= T −1 , (27)
−µ(1 + a) (1 − a2 )

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 112 / 248

AR(1) with Drift Example
The h-steps ahead prediction is:

(1 − âh )
ŷT+h = µ̂ + âh yT .
(1 − â)
Therefore,
h−1
X h−1
X
i i h h
êT+h = (µa − µ̂â ) + (a − â )yT + ai εT+h−i
i=0 i=0

" h−1 #2
− a2h ) (1
X
Var(êT+h ) = σε2
+E i i
(µa − µ̂â )
(1 − a2 )
i=0
" h−1 #
X
+Var (a − âh ) y2T + 2E
h
(µai − µ̂âi )(ah − âh ) yT .

i=0

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 113 / 248

Unit Root Example

Consider a unit root with a drift by setting a = 1 in equation (25).

We have already seen that

ŷT+h = µh + yT
h−1
X
eT+h = εT+h−i
i=0

and Var(eT+h ) = hσε2 increases with the forecast horizon h.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 114 / 248

AR(1) with Drift Example

Using estimated parameters from an AR(1) with drift, without

imposing the unit root, the forecast error becomes
h−1
X
eT+h = (µ − µ̂)h + (1 − âh )yT + εT+h−i . (28)
i=0

It can be shown that the OLS estimator of a converges at a rate of

T 3/2 , rather than T 1/2 as in the stationary case, so that its variance
can be neglected.
The variance of µ̂ decreases instead with T, so that the key
determinant of the forecast uncertainty is the estimation of µ, the
“local trend,” combined with the cumulated future errors.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 115 / 248

The idea of multi-steps (or direct) estimation is to estimate the

parameters that will be used in forecasting by minimizing the same
loss function as in the forecast period.
As an example, let us consider the AR(1):

yt = ayt−1 + εt

so that
h−1
X
h
yT+h = a yT + ai εT+h−i
i=0

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 117 / 248

Standard or Iterated Forecast

The standard forecast is:

ah yT ,
byT+h = b

where
T PT
X 2 yt yt−1
a = argmin
b (yt − ayt−1 ) = Pt=1
T 2
and (29)
a
t=1 t=1 yt−1
h h
E(yT+h − byT+h ) = (a − E(b
a ))yT . (30)

The forecast byT+h is also called “iterated” as it can be derived by

replacing the unknown future values of y with their forecasts for
T + 1, . . . , T + h − 1.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 118 / 248

Direct Forecast

The alternative forecast is

eyT+h = e
ah yT

where
T PT
X 2 yt yt−h
ah = argmin
e (yt − ah yt−h ) = Pt=h
T 2
ah
t=1 t=h yt−h

and
E(yT+h − eyT+h ) = (ah − E(e
ah ))yT .

The forecast eyT+h is labeled “direct” since it is derived from a model

where the target variable yT+h is directly related to the available
information set in period T.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 119 / 248

Comparison

How do the forecasts byT+h and eyT+h compare?

The relative performance in terms of bias and efficiency depends on
the bias and efficiency of the alternative estimators of ah , b
ah and e
ah .
In the presence of correct model specification, both estimators of ah
ah is more efficient than e
are consistent, but b ah since it coincides with
the maximum likelihood estimator.
In the presence of model mis-specification the ranking can change.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 120 / 248

Example

As an example, assume that the DGP is an MA(1):

yt = εt + ψεt−1 (31)
εt ∼ WN(0, σε2 )

Suppose that the chosen model for yt is the AR(1):

yt = ρyt−1 + vt , (32)
vt ∼ WN(0, σv2 )

For illustration, we will compare standard and direct estimation

based forecasts, assuming h = 2 and using the MSFE as a
comparison criterion.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 121 / 248

Example

Standard estimation yields

T T
!−1
X X
ρ̂ = yt yt−1 y2t−1
t=1 t=1

And to a first approximation

ψ
E(ρ̂) ∼
= = ρ.
(1 + ψ 2 )

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 122 / 248

Example

Then
ŷT+2 = ρ̂2 yT , E(ŷT+2 ) ∼
= ρ2 yT
It can be shown that

\ = E (yT+2 − ρ̂2 yT )2 |yT

MSFE
∼
= (1 + ψ 2 )σ 2 + Var(ρ̂2 ) + ρ4 y2

ε T

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 123 / 248

Example

In the case of direct estimation:

ρ̃2 =
PT PT 2 −1 ∼
t=2 yt yt−2 ( t=2 yt−2 ) =0

ỹT+2 = ρ̃2 yT ∼
=0

^ = E[(yT+2 − ỹT+2 )|yT ] ∼

MSFE = (1 + ψ 2 )σε2 + Var(ρ̃2 )y2T

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 124 / 248

Example

It can be shown that for particular values of the parameters it is

possible that
^ ≤ MSFE.
MSFE \

A necessary condition for this is that the AR(1) model is

mis-specified. Otherwise, ρ̂2 is the ML estimator of ρ2 and the
associated forecast cannot be beaten.
The relative performance of the iterated forecasts improves with the
forecast horizon.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 125 / 248

It is sometimes of interest to decompose a process yt into two

components:
Permanent- captures the long-run, trend-like behavior of yt .
Transitory- measures short term deviations from the trend.

There is no unique way to achieve a permanent-transitory. We

consider the to most common approaches:
The Beveridge and Nelson (1981) (BN) decomposition features a
permanent component which behaves as a random walk.
The Hodrick and Prescott (1997) (HP) has potentially more
complex dynamics for the permanent component.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 127 / 248

Beveridge and Nelson Decomposition

Recall that a weakly stationary process can be written as an MA(∞),

and that if yt ∼ I (d) , then ∆d yt is weakly stationary.
Assuming d = 1, we have:
iid
εt ∼ 0, σε2 .

∆yt = µ + c (L) εt , (33)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 128 / 248

Beveridge and Nelson Decomposition

Let us consider another polynomial in L, defined as

d (L) = c (L) − c (1) . (34)

Since d (1) = 0, 1 is a root of d (L) , which can therefore been

rewritten as
d (L) = c̃ (L) (1 − L) . (35)
Combining equations (34) and (35) we obtain

(L) = c̃ (L) (1 − L) + c (1) and

∆yt = µ + c̃ (L) ∆εt + c (1) εt .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 129 / 248

Beveridge and Nelson Decomposition

To obtain a representation for yt , we need to integrate both sides:

t
X
yt = µt + c (1) εj + c̃ (L) εt .
| {z }
j=1
|{z } cycle
trend (transitory component)
(permanent component) (CC)
(PC)

It also follows that the permanent component is a random walk with

drift:
PCt = PCt−1 + µ + c (1) εt . (36)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 130 / 248

Beveridge and Nelson Decomposition

Moreover, the variance of the trend innovation is c (1)2 σε2 , which is

larger (smaller) than the innovation in yt if c (1) is larger (smaller)
than one.
The innovation in the cyclical component is c̃ (0) εt . Since
c̃ (L) = c(L)−c(1)
1−L , then c̃ (0) = c (0) − c (1) = 1 − c (1) . Therefore, the
innovation in the cyclical component is (1 − c (1)) εt .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 131 / 248

BN Decomposition Example

As an example, let us derive the BN decomposition for an

ARIMA(1,1,1) model,

∆xt = φ∆xt−1 + εt + ψεt−1 .

From the MA representation for ∆xt , we have:

1 + ψL 1+ψ
c (L) = , c (1) = ,
1 − φL 1−φ
c (L) − c (1) φ+ψ
c̃ (L) = =− ;
1−L (1 − ψ) (1 − ψL)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 132 / 248

BN Decomposition Example

It follows that the BN decomposition is

t
1 + ψX φ+ψ
yt = PC + CC = εj − εt . (37)
1−φ (1 − ψ) (1 − ψL)
j=1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 133 / 248

The Hodrick-Prescott Filter

The Hodrick-Prescott filter is an alternative way to compute the

permanent component. The permanent component is obtained as

XT T−1 h
X i
min (yt − PCt )2 + λ (PCt+1 − PCt )2 + (PCt − PCt−1 )2 . (38)
PC
t=1 t=2
| {z }
Variance of CC

The bigger is λ, the smoother is the trend.

In practice the following values are used:
λ = 100 for annual data
λ = 1600 for quarterly data
λ = 14400 for monthly data
Note that if λ = 0, it is PCt = yt .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 134 / 248

Exponential smoothing(ES) is a method to produce short-term

forecasts quickly and with sufficient accuracy.
ES decomposes a time series into a “level” component and an
unpredictable residual component.
Once the level at the end of the estimation sample is obtained, say
yLT , it is used as a forecast for yT+h , h > 1.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 136 / 248

Exponential Smoothing

If yt is an i.i.d. process with a non-zero mean, we could estimate yLT

as the sample mean of y.
If instead yt is persistent, then the more recent observations should
receive a greater weight. Hence, we could use
T−1
X
yLT = α(1 − α)t yT−t ,
t=1

with 0 < α < 1 and

eyT+h = yLT . (39)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 137 / 248

Exponential Smoothing

Since
T−1
X
(1 − α)yLT = α(1 − α)t+1 yT−t−1 ,
t=1

we have
yLT = αyT + (1 − α)yLT−1 ,
with the starting condition yL1 = y1 .
The larger α the larger the weight on the most recent observations.
Note that in the limiting case where α = 1, it is yLT = yT and the ES
forecasts coincides with that from a RW model for y.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 138 / 248

Exponential Smoothing

A more elaborate model underlies the Holt-Winters procedure

yLt = ayt + (1 − a)(yLt−1 + Tt−1 ), (40)

Tt = c(yLt − yLt−1 ) + (1 − c)Tt−1 ,

with 0 < a < 1, and starting conditions T2 = y2 − y1 and yL2 = y2 . In

this case, we have
eyT+h = yLT + hTT . (41)
The coefficients α and a control the smoothness of yLt .
In practice, the smoothing coefficients, and c, are selected by
minimizing the in-sample MSFE:
T−1
X
(yt − eyt )2 .
t=3

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 139 / 248

Exponential Smoothing

When can the ES forecasts be considered optimal?

Focusing on the more general case in (41) for which we have:

et = yt − eyt|t−1| = yt − yLt−1 − Tt−1 . (42)

From (40) and (42) we have:

yLt − yLt−1 = a(yt − yLt−1 ) + (1 − a)Tt−1 = Tt−1 + aet , (43)

Combining (40) and (43) we have:

Tt − Tt−1 = c(yLt − yLt−1 ) − cTt−1 = caet (44)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 140 / 248

Exponential Smoothing

Using (43) we have

(yLt − yLt−1 ) − (yLt−1 − yLt−2 ) = Tt−1 − Tt−2 + aet − aet−1

= caet−1 + aet − aet−1

so that
(1 − L)2 yLt = a [1 − (1 − c)L] et .
From (44) we have

(1 − L)2 Tt = ca(1 − L)et .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 141 / 248

Exponential Smoothing

Putting together the expressions for (1 − L)2 yLt and (1 − L)2 Tt we

obtain
(1 − L)2 (yLt + Tt ) = a(1 + c)et − aet−1
From (42) et = yt − yLt−1 − Tt−1 , we can write

(1 − L)2 yt = (1 − L2 )et + a(1 + c)et−1 − aet−2 .

Therefore, in conclusion, the Holt-Winter procedure forecast in (41)

is optimal (in the MSFE sense) if yt is an ARIMA(0,2,2).
Similarly, it can be shown that the forecast in (39) is optimal when yt
is an ARIMA(0,1,1).

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 142 / 248

Seasonality is a systematic but possibly stochastic intra-year

variation in the behavior of an economic variable, related, e.g., to the
weather or the calendar and their impact on economic decisions.
e.g. sales are higher in December, electricity consumption peaks in the
summer
A common approach in the analysis and forecasting of economic
time series is to work with seasonally adjusted variables.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 144 / 248

Seasonality

However, seasonal adjustment procedures can spuriously alter the

dynamic behavior of a variable, in particular when the seasonal
pattern changes over time.
Ghysels, Granger, and Siklos (1996) show that seasonal adjusted
white noise produces a series that has predictable patterns.
We will discuss two alternative approaches in this context:
Deterministic seasonality
Stochastic seasonality

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 145 / 248

Deterministic Seasonality

Assume that the time series xt can be decomposed into a deterministic

seasonal component, st , and an ARMA component, yt :

xt = st + yt . (45)

Under the additional assumption that t is measured in months and

there is a deterministic monthly seasonality, we can write
12
X
st = γi Dit ,
i=1

where Dit are seasonal dummy variables taking values 1 in each

month i and 0 otherwise, i = 1, . . . , 12.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 146 / 248

Deterministic Seasonality

The seasonal coefficients γi can be estimated in the ARMA model

for xt as previously discussed.
Alternative, we could seasonally adjust the variable xt by subtracting
from each of its values est , where
12
X
est = γ
ei Dit , γ
ei = xi ,
i=1

and xi is the sample mean of all month-i observations of x.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 147 / 248

Deterministic Seasonality

We can handle a slow evolution in the shape or amplitude of the

seasonal component, by using weighted averages of the month-i
observations of x for the construction of γ
ei .
For example, also based on the discussion of exponential smoothing
in the previous section we could use

γ γi (j − 1) + (1 − α)xj,i ,
ei (j) = αe

where α ∈ [0, 1] , j is measured in years, and xj,i indicates the value

of x in month i of year j.
A similar dummy variable approach can be used for different
patterns of seasonality (e.g., quarterly), and also to handle other
data irregularities such as working days effects and moving festivals.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 148 / 248

Stochastic Seasonality

If the seasonal pattern could be stochastic rather than deterministic,

we could extend the ARMA(p,q) specification to take into explicit
account the possibility of seasonal dynamics.
We could use a model such as

φ1 (L) φ2 (Ls ) yt = ψ1 (L) ψ2 (Ls ) εt , (46)

where Ls is the seasonal lag operator, so that

φ2 (Ls ) = 1 − φ21 Ls − . . . − φ2p2 Lsps ,

ψ2 (Ls ) = 1 + ψ21 Ls + . . . + ψ2ψ2 Lsqs ,

.
This model is typically known as a seasonal ARMA model and it can
properly represent several seasonal time series.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 149 / 248

Stochastic Seasonality

One example is the seasonal random walk model:

yt = yt−s + εt ,

year-to-year differences are white noise.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 150 / 248

Stochastic Seasonality

Another example is the case where p = ps = q = qs =1 and s =12:

(1 − φ11 L)(1 − φ21 L12 )yt = (1 + ψ11 L)(1 + ψ21 L12 )εt ,

(1 − φ11 L − φ21 L12 + φ11 φ21 L13 )yt = (1 + ψ11 L + ψ21 L12 + ψ11 ψ21 L13 )εt .

This can also be considered as an ARMA(13,13), though with a set

of zero restrictions on the coefficients.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 151 / 248

Stochastic Seasonality

Since a seasonal ARMA is just a constrained version of a general

ARMA model, the same tools for specification, estimation and
diagnostic testing that we have seen for the ARMA case can be
applied in this context.
Finally, some of the roots of the seasonal AR polynomial φ2 (Ls )
could be equal to one as in the seasonal random walk model. In the
general case, these are known as seasonal unit roots.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 152 / 248

We use standard Box–Jenkins procedure to specify and estimate

several types of ARMA models. Once the appropriate ARMA
specification is found, it is used for an out-of-sample forecasting
exercise.
The data set contains 600 simulated observations.
In order to avoid dependence on the starting values, the first 100
observations are discarded.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 154 / 248

Stationary ARMA Process

We consider three stationary ARMA process:

y1 , with an AR(2) DGP
y2 , with an MA(2) DGP
y3 , with an ARMA(2,2) DGP

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 155 / 248

Y1 Y2
4 4

2 2

0 0

-2 -2

-4 -4
150 200 250 300 350 400 450 500 550 600 150 200 250 300 350 400 450 500 550 600

Y3
6

-2

-4
150 200 250 300 350 400 450 500 550 600

Figure 1: Stationary ARMA process y1 through y3 , sample 101 - 600

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 156 / 248

Stationary ARMA Process

Note from Figure 1 that the dynamic behavior of the variables looks
similar although they have been generated from different DGPs.
Figures 2 - 4 show the (partial) correlograms of the three series for
the entire sample under consideration.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 157 / 248

Correlogram for y1

Autocorrelation Partial Correlation AC PAC

1 0.370 0.370
2 0.350 0.246
3 0.149 -0.049
4 0.115 -0.002
5 0.054 -0.001
6 0.040 0.003
7 0.029 0.010
8 -0.024 -0.054
9 -0.058 -0.058
10 -0.124 -0.089
11 -0.099 -0.011
12 -0.154 -0.074

Figure 2: AC and PAC functions for y1 entire sample (101 -600)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 158 / 248

Correlogram for y2

Autocorrelation Partial Correlation AC PAC

1 0.145 0.145
2 -0.335 -0.364
3 -0.022 0.114
4 -0.014 -0.184
5 -0.029 0.054
6 0.009 -0.074
7 0.055 0.089
8 0.004 -0.060
9 0.003 0.079
10 -0.010 -0.067
11 -0.056 -0.001
12 0.005 -0.010

Figure 3: AC and PAC functions for y2 entire sample (101 - 600)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 159 / 248

Correlogram for y3

Autocorrelation Partial Correlation AC PAC

1 0.227 0.227
2 -0.234 -0.301
3 0.002 0.161
4 0.001 -0.140
5 -0.001 0.097
6 0.052 -0.013
7 0.079 0.105
8 0.020 -0.029
9 0.032 0.097
10 -0.014 -0.082
11 -0.072 -0.001
12 0.008 -0.003

Figure 4: AC and PAC functions for y3 entire sample (101 - 600)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 160 / 248

Stationary ARMA Process

Figure 2, y1 has an exponentially decaying autocorrelation function

and significant peaks at lags one and two in the PAC, strongly
supporting an AR(2) specification.
Figure 3 has an oscillating PAC function while the AC displays a
significant peak at lag two, suggesting that y2 could in fact follow a
MA(2) process.
From Figure 4 it is harder to make a proposal for a model for y3
based on the reported AC and PAC functions. However, the
presence of significant peaks in both of them suggests a combined
ARMA specification.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 161 / 248

Stationary ARMA Process

Let’s use BIC to select an appropriate ARIMA specification and

compare with the AC and PAC analysis.
Below is the BIC for various specifications for y1 .
AR / MA 0
0 3.042114
1 2.909087
2 2.860312
3 2.871777
4 2.885768
5 2.896720
6 2.910012
Table 1: The ARIMA selection with BIC for y1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 162 / 248

BIC for y2

AR / MA 0 1 2 3
0 3.038216 2.965412 2.817112 2.820006
1 3.024794 2.871028 2.819972 2.823934
2 2.886401 2.836022 2.823942 2.827933
3 2.877178 2.830599 2.827901 2.823809
Table 2: The ARIMA selection with BIC for y2

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 163 / 248

BIC for y3

AR / MA 0 1 2 3
0 3.187754 3.059214 2.981900 2.985839
1 3.142723 2.995391 2.985840 2.988778
2 3.051139 2.987988 2.981839 2.991453
3 3.028714 2.989539 2.993529 2.987754
Table 3: The ARIMA selection with BIC for y3

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 164 / 248

Stationary ARMA Process

we see that the procedure always selects the true DGPs when the
entire sample is taken as a reference, as shown in Tables 1 - 3 for
BIC.
Next, we will estimate the specifications selected by BIC.
Finally, we will illustrate the forecasting performance of the selected
models.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 165 / 248

Estimation: Full sample analysis

Tables 4 through 6 report the estimation results of the three models

suggested by the ARIMA selection procedure;
y1 is modeled as an AR(2)
y2 is modeled as an MA(2)
y3 is modeled as an ARMA(2,2)

Note that the selected AR and MA terms are always significant and
Durbin-Watson statistics are always pretty close to the value of 2,
signaling almost with noise residuals.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 166 / 248

Estimation results for y1
Dep Var: Y1

Variable Coefficient Std. Error t-Statistic Prob.

C 0.042014 0.093569 0.449014 0.653

AR(1) 0.279133 0.043545 6.410269 0.000
AR(2) 0.245699 0.043530 5.644410 0.000

R-squared 0.189038 Mean dep var 0.03953

Adjusted R-squared 0.185774 S.D. dep var 1.10175
S.E. of regression 0.994163 Akaike IC 2.83215
Sum squared resid 491.2152 Schwarz IC 2.85743
Log likelihood -705.0378 Hannan-Quinn 2.84207
F-statistic 57.92603 DW stat 1.96986
Prob(F-statistic) 0.000000

Table 4: Estimation results for y1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 167 / 248

Estimation results for y2

Dep Var: Y2
Variable Coefficient Std. Error t-Statistic Prob.

C -0.009993 0.038195 -0.261621 0.793

MA(1) 0.281535 0.040840 6.893551 0.000
MA(2) -0.414625 0.040851 -10.14959 0.000

R-squared 0.209487 Mean dep var -0.01033

Adjusted R-squared 0.206306 S.D. dep var 1.10425
S.E. of regression 0.983775 Akaike IC 2.81114
Sum squared resid 481.0032 Schwarz IC 2.83643
Log likelihood -699.7857 Hannan-Quinn 2.82106
F-statistic 65.85277 DW stat 1.96942
Prob(F-statistic) 0.000000

Table 5: Estimation results for y2

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 168 / 248

Estimation results for y3

Dep Var: Y3
Variable Coefficient Std. Error t-Statistic Prob.

C 0.011417 0.054939 0.207816 0.835

AR(1) -0.193760 0.283357 -0.683803 0.494
AR(2) -0.076268 0.132683 -0.574815 0.565
MA(1) 0.594037 0.284018 2.091550 0.037
MA(2) -0.135818 0.241172 -0.563158 0.573

R-squared 0.198335 Mean dep var 0.01040

Adjusted R-squared 0.191857 S.D. dep var 1.18998
S.E. of regression 1.069757 Akaike IC 2.98269
Sum squared resid 566.4681 Schwarz IC 3.02483
Log likelihood -740.672 Hannan-Quinn 2.99922
F-statistic 30.61621 DW stat 1.99818
Prob(F-statistic) 0.000000
Table 6: Estimation results for y3

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 169 / 248

Diagnostics

Diagnostics appear in Tables 7 though 9, showing some routine

statistics for the residuals coming from the estimated models for y1 ,
y2 and y3 .

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 170 / 248

Diagnostic tests on the AR(2) model for y1

Heteroskedasticity Test: ARCH

F-statistic 0.009042 F(1,497) 0.924

Obs*R-squared 0.009078 Chi-Sq(1) 0.924
Breusch-Godfrey Serial
Correlation LM Test:
F-statistic 0.628107 F(2,495) 0.534
Obs*R-squared 1.265691 Chi-Sq(2) 0.531
Heteroskedasticity Test: White

F-statistic 0.776695 F(5,494) 0.566

Obs*R-squared 3.899984 Chi-Sq(5) 0.563
Scaled explained SS 3.132991 Chi-Sq(5) 0.679
Jarque-Bera normality Test: 3.516553 0.172
Table 7: Diagnostic tests on the AR(2) model for y1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 171 / 248

Diagnostic tests on the MA(2) model for y2

Heteroskedasticity Test: ARCH

F-statistic 0.169720 F(1,497) 0.680

Obs*R-squared 0.170345 Chi-Sq(1) 0.679
Breusch-Godfrey Serial
Correlation LM Test:
F-statistic 0.281619 F(2,495) 0.754
Obs*R-squared 0.568031 Chi-Sq(2) 0.752
Heteroskedasticity Test: White

F-statistic 0.776314 F(9,490) 0.638

Obs*R-squared 7.029190 Chi-Sq(9) 0.634
Scaled explained SS 7.486218 Chi-Sq(9) 0.586
Jarque-Bera normality Test: 0.945034 Prob. 0.623
Table 8: Diagnostic tests on the MA(2) model for y2

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 172 / 248

Diagnostic tests on the ARMA(2,2) model for y3

Heteroskedasticity Test: ARCH

F-statistic 0.087968 F(1,497) 0.766

Obs*R-squared 0.088306 Chi-Sq.(1) 0.766
Breusch-Godfrey Serial
Correlation LM Test:
F-statistic 3.380311 F(2,495) 0.034
Obs*R-squared 6.736898 Chi-Sq.(2) 0.034
Heteroskedasticity Test: White

F-statistic 1.475130 F(9,490) 0.154

Obs*R-squared 13.18975 Chi-Sq.(9) 0.154
Scaled explained SS 13.83836 Chi-Sq.(9) 0.128
Jarque-Bera normality Test: 0.504597 Prob. 0.777
Table 9: Diagnostic tests on the ARMA(2,2) model for y3

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 173 / 248

Residuals

For all the three models, the actual vs fitted and residuals (see
Figures 5 through 9) show a good fit for the selected models.
The residuals display frequent sign changes, although there might
be a few outliers.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 174 / 248

Actual vs fitted residuals

3
-2
2

1 -4

-1

-2

-3
150 200 250 300 350 400 450 500 550 600

Residual Actual Fitted

Figure 5: Actual vs fitted residuals

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 175 / 248

Correlogram of the residuals for the model for y1

Autocorrelation Partial Correlation AC PAC

1 0.012 0.012
2 0.019 0.019
3 -0.040 -0.041
4 -0.011 -0.011
5 -0.018 -0.016
6 0.007 0.007
7 0.038 0.038
8 -0.001 -0.004
9 -0.021 -0.023
10 -0.078 -0.075
11 -0.035 -0.032
12 -0.124 -0.122

Figure 6: Correlogram of the residuals for the AR(2) model for y1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 176 / 248

Actual vs fitted residuals for y2

0
4

-2
2

-4
0

-2

-4
150 200 250 300 350 400 450 500 550 600

Residual Actual Fitted

Figure 7: Actual vs fitted residuals for y2

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 177 / 248

Correlogram of the residuals for the model for y2

Autocorrelation Partial Correlation AC PAC

1 0.014 0.014
2 -0.009 -0.009
3 -0.023 -0.023
4 -0.011 -0.011
5 -0.022 -0.022
6 0.002 0.002
7 0.044 0.043
8 0.001 -0.002
9 -0.010 -0.010
10 0.005 0.007
11 -0.080 -0.080
12 0.005 0.009

Figure 8: Correlogram of the residuals for the MA(2) model for y2

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 178 / 248

Actual vs fitted residuals for y3

0
4
-2

2 -4

-2

-4
150 200 250 300 350 400 450 500 550 600

Residual Actual Fitted

Figure 9: Actual vs fitted residuals for y3

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 179 / 248

Correlogram of the residuals for the model for y3

Autocorrelation Partial Correlation AC PAC

1 -0.000 -0.000
2 -0.000 -0.000
3 -0.025 -0.025
4 -0.004 -0.004
5 0.011 0.011
6 0.030 0.029
7 0.083 0.083
8 0.002 0.002
9 0.035 0.037
10 -0.002 0.002
11 -0.072 -0.073
12 0.013 0.011

Figure 10: Correlogram of the residuals for the ARMA(2,2) model for y3

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 180 / 248

Model Misspecification

Having a simulated dataset allows us to answer also another

interesting question: What if the researcher “gets it wrong” and
estimates the following models:
AR(1) for y1
MA(1) for y2
ARMA(1,1) for y3
For y1 estimating an AR(1) yields a reasonable good fit with only a
slight decrease in the

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 181 / 248

Model Misspecification

For y1 estimating an AR(1) yields a reasonable good fit with only a

slight decrease in the in the adjusted R2 to 14% but no major signs
of misspecification (see next slide).
Nevertheless, the correlogram of the residuals of the AR(1) question
features a relevant peak exactly at lag 2, as the AR(2) term was
omitted in the estimation.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 182 / 248

Estimation results for the mis-specified model for y1

Variable Coefficient Std. Error t-Statistic Prob.

C 0.025639 0.045844 0.559268 0.576

AR(1) 0.370487 0.041659 8.893349 0.000

R-squared 0.137052 Mean dep var 0.03953

Adjusted R-squared 0.135319 S.D. dep var 1.10175
S.E. of regression 1.024503 Akaike IC 2.89028
Sum squared resid 522.7037 Schwarz IC 2.90714
Log likelihood -720.570 Hannan-Quinn 2.89689
F-statistic 79.09166 DW stat 2.17887
Prob(F-statistic) 0.000000
Table 10: Estimation results for the mis-specified model for y1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 183 / 248

Correlograms for Mis-specified Model for y1

Autocorrelation Partial Correlation AC PAC

1 -0.091 -0.091
2 0.238 0.231
3 -0.004 0.036
4 0.065 0.013
5 0.003 0.002
6 0.018 0.000
7 0.031 0.031
8 -0.020 -0.022
9 -0.012 -0.031
10 -0.096 -0.098
11 -0.013 -0.021
12 -0.140 -0.105

Figure 11: Correlograms for mis-specified model for y1

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 184 / 248

Correlograms for Mis-specified Model for y2

Autocorrelation Partial Correlation AC PAC

1 -0.179 -0.179
2 -0.259 -0.300
3 0.112 -0.003
4 -0.062 -0.130
5 0.011 0.003
6 -0.025 -0.084
7 0.067 0.071
8 -0.027 -0.041
9 0.006 0.051
10 0.018 -0.006
11 -0.071 -0.038
12 0.038 0.005

Figure 12: Correlograms for mis-specified model for y2

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 185 / 248

Correlograms for Mis-specified Model for y3

Autocorrelation Partial Correlation AC PAC

1 -0.036 -0.036
2 -0.101 -0.102
3 -0.026 -0.034
4 0.022 0.009
5 -0.018 -0.023
6 0.044 0.046
7 0.061 0.063
8 0.006 0.019
9 0.025 0.043
10 0.005 0.012
11 -0.078 -0.072
12 0.029 0.026

Figure 13: Correlograms for mis-specified model for y3

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 186 / 248

Subsample Analysis

We now will apply the ARIMA selection procedure to two smaller

samples: observations 101 to 200 and 201 to 300
Subsample 101-200 Subsample 201-300
AIC BIC AIC BIC
y1 AR(2) AR(1) AR(2) AR(2)
y2 MA(2) MA(2) MA(2) MA(10)
y3 ARMA(2,2) ARMA(2,2) ARMA(8,8) ARMA(3,3)

As we can see, in shorter samples the criteria can select models

rather different from the DGPs.
This is not necessarily a bad outcome, as long as the residuals of
the models satisfy the usual assumptions and the forecasting
performance remains satisfactory.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 187 / 248

Forecasting

For illustrative purposes here we chose to estimate the three models

from observation 100 to 499 and to forecast the last 100 observations
of the sample, although 100 periods ahead is quite a large forecast
horizon.
The static forecast is simply the one-step ahead forecast.
The dynamic forecasting method calculates multi-steps ahead
forecasts, from 1- up to 100-period ahead.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 188 / 248

Actual values and forecasts for y1

Figure 14: Actual values and static vs dynamic forecasts for y1 in the sample 501 - 600

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 189 / 248

Actual values and forecasts for y2

Figure 15: Actual values and static vs dynamic forecasts for y2 in the sample 501 - 600

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 190 / 248

Actual values and forecasts for y3

Figure 16: Actual values and static vs dynamic forecasts for y3 in the sample 501 - 600

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 191 / 248

The aim of this empirical application is to formulate an ARMA model

suitable to describe the dynamics of the effective federal funds rate r in
the United States.
Monthly dataset spanning from January 1985 until the end of 2012
from the Federal Reserve Economic Data (FRED) website.
We will first analyze the time frame going from 1985 to 2002, using
the period from 2003 to 2006 for an out-of-sample forecasting
exercise prior the global financial crisis.
Then, we will enlarge the estimation sample to the start of the crisis,
i.e., in August 2007, and we will see whether our model is able to
forecast the steep fall that the federal funds rate displayed from 2008
until the end of the sample.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 193 / 248

US Federal Funds Rate 1985 - 2012
R
10

0
86 88 90 92 94 96 98 00 02 04 06 08 10 12

Figure 17: The effective US federal funds rate 1985 - 2012

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 194 / 248

Modeling and forecasting the US federal funds rate

The Augmented Dickey-Fuller test can help us investigate the

dynamic properties of the series more thoroughly.
Based on the results in the following slide, we cannot reject the null
of a unit root at conventional significance levels.
However, when we perform the same test on the differenced variable
Dr instead, there is strong evidence in favor of stationarity,
suggesting that a specification in differences for r is preferred.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 195 / 248

The Augmented Dickey-Fuller test

ADF Test for r t-Statistic Prob.*

-0.912779 0.782
Test critical values: 1% level -3.460884
5% level -2.874868
10% level -2.573951

ADF Test for Dr t-Statistic Prob.*

-8.810261 0.000
Test critical values: 1% level -2.575813
5% level -1.942317
10% level -1.615712
Table 11: The Augmented Dickey-Fuller test

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 196 / 248

Lag Specification

BIC suggests a ARMA(10,1) for Dr (see next slide).

However, AIC suggest ARMA(5,2) for Dr.
We started specifying an ARMA(10,2) model for Dr, and then we
noticed that some AR lags could be conveniently eliminated. At the
end of this general-to-specific specification procedure we ended up
with a parsimonious ARMA(5,2) model reported in Table 13.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 197 / 248

ARIMA selection with BIC for Dr

AR / MA 0 1 2
0 -0.276287 -0.514755 -0.544955
1 -0.596637 -0.611789 -0.627400
2 -0.600159 -0.622702 -0.622388
3 -0.615346 -0.621617 -0.616286
4 -0.611976 -0.615451 -0.610219
5 -0.619612 -0.615223 -0.621251
6 -0.616315 -0.612577 -0.621058
7 -0.613225 -0.607978 -0.621052
8 -0.611476 -0.607423 -0.611349
9 -0.614132 -0.610542 -0.609949
10 -0.615026 -0.650287 -0.644754
11 -0.615257 -0.645099 -0.641540
12 -0.609648 -0.640053 -0.638794
Table 12: ARIMA selection with BIC for Dr

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 198 / 248

Parameter estimates ARMA(5,2) model for Dr

Variable Coefficient Std. Error t-Statistic Prob.

AR(1) 0.169824 0.063523 2.673423 0.008

AR(2) -0.799708 0.043098 -18.55554 0.000
AR(3) 0.515466 0.076887 6.704226 0.000
AR(5) 0.205128 0.066405 3.089048 0.002
MA(1) 0.262565 0.013437 19.54028 0.000
MA(2) 0.976636 0.009565 102.1080 0.000

R-squared 0.282454 Mean dep var -0.02995

Adjusted R-squared 0.264867 S.D. dep var 0.22866
S.E. of regression 0.196060 Akaike IC -0.39263
Sum squared resid 7.841698 Schwarz IC -0.29700
Log likelihood 47.22639 Hannan-Quinn -0.35397
DW stat 2.004751
Table 13: Parameter estimates ARMA(5,2) model for Dr

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 199 / 248

Correlogram of residuals from ARMA(5,2) for Dr

Autocorrelation Partial Correlation AC PAC

1 -0.016 -0.016
2 0.005 0.005
3 -0.014 -0.014
4 0.064 0.064
5 0.040 0.042
6 0.035 0.036
7 0.000 0.003
8 -0.044 -0.047
9 0.049 0.044
10 0.020 0.016
11 0.011 0.007
12 0.021 0.027

Figure 18: Correlogram of residuals from ARMA(5,2) for Dr

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 200 / 248

Diagnostic Tests

Breusch-Godfrey
Serial Correlation LM Test:
F-statistic 0.090390 F(2,202) 0.913
Obs*R-squared 0.187771 Chi-Sq(2) 0.910
Heteroskedasticity Test: ARCH
F-statistic 35.80092 F(1,207) 0.000
Obs*R-squared 30.81699 Chi-Sq(1) 0.000
Heteroskedasticity Test: White
F-statistic 4.645695 F(21,188) 0.000
Obs*R-squared 71.74515 Chi-Sq(21) 0.000
Table 14: Diagnostic tests on the residuals of ARMA(5,2) model for Dr

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 201 / 248

Forecasting

We now want to produce one-step ahead forecasts for r in the period

2003 - 2006.
The first option is to construct forecasts for Dr and then cumulate
them with the starting value for r, namely,

brT+1 = rT + Dr
c T+1 .

When we use this method, we add to the forecast series the suffix
“differences,” in order to clarify that the forecast for rT+1 has been
obtained indirectly from the model specified in differences.
As an alternative, we can specify a model for the level of the interest
rate r and produce a forecast for rT+1 directly.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 202 / 248

Forecasting

Which model is suitable to represent r in levels?

The unit root test indicated that there was evidence for
non-stationarity, and the correlogram showed a very persistent AR
type of behavior.
It turns out that an AR(11) specification, with many non-significant
intermediate lags removed, is able to capture almost all the
variability of the series.
From The estimation output in Table 15, we note the very high R2 ,
which is not surprising for a specification for the levels of an I(1)
variable.
The diagnostic tests on the residuals reported in Table 16 confirm
the likely correctly specified dynamics.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 203 / 248

Parameter estimates restricted AR(11) model for r

Variable Coefficient Std. Error t-Statistic Prob.

AR(1) 1.410725 0.065970 21.38418 0.000

AR(2) -0.385824 0.071123 -5.424708 0.000
AR(11) -0.028306 0.011827 -2.393426 0.017

R-squared 0.989310 Mean dep var 5.51282

Adjusted R-squared 0.989204 S.D. dep var 1.92831
S.E. of regression 0.200362 Akaike IC -0.36285
Sum squared resid 8.109257 Schwarz IC -0.31422
Log likelihood 40.19299 Hannan-Quinn -0.34318
DW stat 2.004220
Table 15: Parameter estimates restricted AR(11) model for r

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 204 / 248

Diagnostic tests on the residuals of AR(11) model for r

Breusch-Godfrey
Serial Correlation LM Test:
F-statistic 0.660035 F(2,200) 0.518
Obs*R-squared 1.344200 Chi-Sq(2) 0.510
Heteroskedasticity Test: ARCH
F-statistic 26.24542 F(1,202) 0.000
Obs*R-squared 23.45749 Chi-Sq(1) 0.000
Heteroskedasticity Test: White
F-statistic 11.92362 F(6,198) 0.000
Obs*R-squared 54.41108 Chi-Sq(6) 0.000
Table 16: Diagnostic tests on the residuals of AR(11) model for r

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 205 / 248

Forecasts

We can now compare the forecasting performance of the models

specified in differences and in levels.
Figure 19 shows dynamic or multi-steps ahead forecasts.
Figure 20 shows static or one-step ahead forecasts.
Table 17 complements the analysis by showing the standard
forecast evaluation criteria, the RMSFE and the MAFE, for the
one-step ahead forecasts.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 206 / 248

Dynamic or Multi-steps Ahead Forecasts

Figure 19: The actual values of r against the forecasted series: dynamic or multi-steps ahead
forecasts for the period 2003 - 2006

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 207 / 248

Static or One-step Ahead Forecast

Figure 20: The actual values of r against the forecasted series: static or one-step ahead forecast,
for the period 2003 - 2006.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 208 / 248

Forecast Evaluations for the One-step Ahead
Forecasts

RMSFE MAFE
Levels 0.0868 0.0637
Differences 0.1075 0.0834
Table 17: Forecast evaluations for the one-step ahead forecasts

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 209 / 248

The Crisis Period

Will the one-step predictive capacity of the above-selected models be

any different if we enlarge our to include the financial crisis period?
Tables 18 and 20 clearly show that the models’ fit remains quite
similar also using the enlarged estimation sample.
Figures 21, 22, and Table 19 report detailed forecast evaluation
criteria.
We also look at 2-step ahead forecast in Figure 23.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 210 / 248

Parameter estimates ARMA(5,2) for Dr with enlarged
sample

Variable Coefficient Std. Error t-Statistic Prob.

AR(1) 0.309924 0.179822 1.723511 0.086

AR(2) -0.529328 0.201804 -2.622974 0.009
AR(3) 0.385577 0.096964 3.976501 0.000
AR(5) 0.227198 0.056587 4.015022 0.000
MA(1) 0.099835 0.180357 0.553539 0.580
MA(2) 0.613386 0.175845 3.488216 0.000

R-squared 0.290327 Mean dep var -0.00856

Adjusted R-squared 0.276627 S.D. dep var 0.21258
S.E. of regression 0.180808 Akaike IC -0.56037
Sum squared resid 8.467142 Schwarz IC -0.47932
Log likelihood 80.24991 Hannan-Quinn -0.52781
DW stat 2.007645
Table 18: Parameter estimates ARMA(5,2) for Dr with enlarged sample
c Eric Ghysels & Massimiliano Marcellino 2018 Edition 211 / 248
Forecast Evaluations for the One-step Ahead
Forecasts

RMSFE MAFE
Levels 0.1491 0.0694
Differences 0.1552 0.0703
Table 19: Forecast evaluations for the one-step ahead forecasts

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 212 / 248

Parameter estimates AR(11) model for r with enlarged
sample

Variable Coefficient Std. Error t-Statistic Prob.

AR(1) 1.417957 0.058500 24.23864 0.000

AR(2) -0.388422 0.063316 -6.134685 0.000
AR(11) -0.032273 0.010137 -3.183665 0.001

R-squared 0.992941 Mean dep var 4.97984

Adjusted R-squared 0.992886 S.D. dep var 2.15119
S.E. of regression 0.181445 Akaike IC -0.56425
Sum squared resid 8.461069 Schwarz IC -0.52316
Log likelihood 76.35278 Hannan-Quinn -0.54773
DW stat 2.009223
Table 20: Parameter estimates AR(11) model for r with enlarged sample

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 213 / 248

Actual vs 1-step Forecast

Figure 21: Actual vs 1-step forecasted r series, period 2007 - 2012

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 214 / 248

Actual vs h-step Ahead Forecast

Figure 22: Actual vs h-step forecasted r series, period 2007 - 2012

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 215 / 248

Actual vs 2-steps Ahead Forecast
6

-1
III IV I II III IV I II III IV I II III IV I II III IV I II III IV
2007 2008 2009 2010 2011 2012

RF_2STEPS
RF_2STEPS_DIRECT
R

Figure 23: Actual vs 2-steps ahead forecasted r series, period 2007 - 2012

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 216 / 248

Forecast evaluation

One-step Two-steps ahead

ahead iterated direct
RMSFE 0.1552 0.3160 0.1813
MAFE 0.0703 0.1461 0.1813
Table 21: Forecast evaluation

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 217 / 248

Change in Real Private Inventories

As a second empirical example, we consider the quarterly time series

of the change in real private inventories (rcpi).
Data for the United States, for the period 1985 - 2012, downloaded
again from the FRED database.
The series seems stationary around its mean while its variance
seems to increase towards the end of the sample.
We begin our analysis by inspecting the AC and PAC.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 218 / 248

rcpi in the US for 2007-2012

RCPI
1.2

0.8

0.4

0.0

-0.4

-0.8

-1.2

-1.6

-2.0
86 88 90 92 94 96 98 00 02 04 06 08 10 12

Figure 24: The series of real changes in private inventories for the US for 2007 - 2012

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 219 / 248

Correlogram of rcpi

Autocorrelation Partial Correlation AC PAC

1 0.385 0.385
2 -0.083 -0.272
3 -0.087 0.075
4 -0.320 -0.421
5 -0.152 0.263
6 -0.003 -0.303
7 -0.002 0.290
8 0.220 -0.093
9 0.089 0.016
10 -0.186 -0.272
11 -0.206 -0.021
12 -0.165 -0.009
13 -0.048 -0.054
14 0.039 -0.035
15 0.020 -0.112

Figure 25: Correlogram of rcpi

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 220 / 248

Change in Real Private Inventories

The ADF confirms the time series is stationary.

AIC- and BIC based results both suggest an ARMA(3,2).
Estimation and diagnostics follow.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 221 / 248

Augmented Dickey-Fuller Test for rcpi

t-Statistic Prob.*

Augmented Dickey-Fuller test statistic -3.9513 0.002

Test critical values: 1% level -3.4970
5% level -2.8906
10% level -2.5823
Table 22: Augmented Dickey-Fuller test for rcpi

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 222 / 248

Parameter estimates ARMA(3,2) model for rcpi

Variable Coefficient Std. Error t-Statistic Prob.

C 0.311431 0.094479 3.296302 0.001

AR(1) 0.596686 0.113719 5.247034 0.000
AR(2) -0.712288 0.098247 -7.249984 0.000
AR(3) 0.451393 0.106507 4.238135 0.000
MA(1) -0.154118 0.034152 -4.512663 0.000
MA(2) 0.967220 0.020041 48.26132 0.000

R-squared 0.433316 Mean dep var 0.31330

Adjusted R-squared 0.388342 S.D. dep var 0.36738
S.E. of regression 0.287329 Akaike IC 0.42656
Sum squared resid 5.201136 Schwarz IC 0.62083
Log likelihood -8.71634 Hannan-Quinn 0.50363
F-statistic 9.634632 DW stat 2.01732
Prob(F-statistic) 0.000001
Table 23: Parameter estimates ARMA(3,2) model for rcpi

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 223 / 248

Actual vs fitted and residuals from ARMA(3,2) for rcpi,
1986 - 2002 1.5

1.0

0.5

0.0
0.8

-0.5
0.4

-1.0
0.0

-0.4

-0.8

-1.2
1986 1988 1990 1992 1994 1996 1998 2000 2002

Residual Actual Fitted

Figure 26: Actual vs fitted and residuals from ARMA(3,2) for rcpi, 1986 - 2002

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 224 / 248

Correlogram of the residuals from ARMA(3,2)

Autocorrelation Partial Correlation AC PAC

1 -0.013 -0.013
2 0.029 0.029
3 0.054 0.055
4 -0.017 -0.017
5 0.015 0.012
6 -0.071 -0.073
7 0.054 0.053
8 -0.215 -0.214
9 0.064 0.073
10 0.165 0.175
11 0.005 0.031
12 -0.031 -0.070

Figure 27: Correlogram of the residuals from ARMA(3,2)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 225 / 248

Diagnostic tests on the residuals from ARMA(3,2)

Breusch-Godfrey
Serial Correlation LM Test:
F-statistic 0.054017 F(2,61) 0.947
Obs*R-squared 0.121987 Chi-Sq(2) 0.940

Heteroskedasticity Test: ARCH

F-statistic 0.625783 F(1,66) 0.431
Obs*R-squared 0.638691 Chi-Sq(1) 0.424

Heteroskedasticity Test: White

F-statistic 0.443385 F(27,41) 0.985
Obs*R-squared 15.59381 Chi-Sq(27) 0.960
Scaled explained SS 19.89035 Chi-Sq(27) 0.835
Table 24: Diagnostic tests on the residuals from ARMA(3,2)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 226 / 248

Change in Real Private Inventories

We will compare both one-step ahead/static, multi-steps

ahead/dynamic forecasts and two steps ahead forecast.
The quality of these forecasts is not very high.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 227 / 248

The one-step ahead and the h-steps ahead forecast
series against the actuals for rcpi, 1986 - 2002
1.0

0.8

0.6

0.4

0.2

0.0

-0.2
I II III IV I II III IV I II III IV I II III IV
2003 2004 2005 2006

RCPI
RCPIF_1STEP
RCPIF_HSTEPS

Figure 28: The one-step ahead and the h-steps ahead forecast series against the actuals for rcpi,
1986 - 2002

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 228 / 248

The two-steps ahead forecasts against the actuals for
rcpi, 1986 - 2002
1.4

1.2

1.0

0.8

0.6

0.4

0.2

0.0

-0.2
I II III IV I II III IV I II III IV I II III IV
2003 2004 2005 2006

RCPI
RCPI_2STEPS
RCPI_2STEPS_DIRECT

Figure 29: The two-steps ahead forecasts against the actuals for rcpi, 1986 - 2002

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 229 / 248

Forecast evaluation

1-step 2-steps 2-steps 1-to-x

ahead ahead, iterated ahead, direct steps ahead
RMSFE 0.3110 0.5800 0.3152 0.3362
MAFE 0.2563 0.4393 0.2822 0.2961
Table 25: Forecast evaluation

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 230 / 248

The Crisis Period

We now assess whether an ARMA model can produce reliable

forecasts also during the global financial crisis period.
We now include the crisis period in the estimation sample, which
spans from 1985Q1 to 2009Q4, and we produce forecasts for
2010Q1 - 2012Q4.
We use a dummy variable equal to 1 in the period 2007Q4 - 2009Q4
and equal to 0 elsewhere.
The selected most parsimonious specification is reported in Table
26.
A seventh-order AR term had to be included, as otherwise the
correlogram indicated significant serial autocorrelation. The term is
only mildly significant if multiplied by the dummy, but it substantially
improves the overall model fit.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 231 / 248

Estimation results for the AR(7) model for rcpi
Coefficient Std. Error t-Statistic Prob.

C(1) 0.144059 0.048902 2.945845 0.004

C(2) -0.567565 0.410887 -1.381317 0.170
C(3) 0.465738 0.109646 4.247637 0.000
C(4) 2.409279 0.546987 4.404637 0.000
C(5) -2.893790 0.726124 -3.985256 0.000
C(6) 0.134667 0.109553 1.229241 0.222
C(7) 1.155596 0.855702 1.350465 0.180

R-squared 0.689825 Mean dep var 0.24323

Adjusted R-squared 0.668185 S.D. dep var 0.52727
S.E. of regression 0.303729 Akaike IC 0.52692
Sum squared resid 7.933607 Schwarz IC 0.71754
Log likelihood -17.50192 Hannan-Quinn 0.60389
F-statistic 31.87716 DW stat 2.12359
Prob(F-statistic) 0.000000
Table 26: Estimation results for the AR(7) model for rcpi in the estimation period containing the
global financial crisis.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 232 / 248

The Crisis Period

This specification yields an adjusted R2 = 0.67, substantially higher

than the one for the pre-crisis period (which was 0.38), with a Durbin
Watson statistic of 2.12.
Figures 30 and 31, and Table 27 report the diagnostic tests on the
residuals that overall support the model specification.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 233 / 248

The residuals of the estimated AR(7) model for rcpi

0.8 -1

0.4 -2

0.0 -3

-0.4

-0.8

-1.2
86 88 90 92 94 96 98 00 02 04 06 08

Residual Actual Fitted

Figure 30: The residuals of the estimated AR(7) model for rcpi in the extended sample including
the crisis

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 234 / 248

Correlogram of the residuals of the AR(7)

Autocorrelation Partial Correlation AC PAC

1 -0.070 -0.070
2 0.024 0.019
3 0.112 0.115
4 -0.089 -0.075
5 -0.106 -0.126
6 0.053 0.031
7 0.134 0.174
8 -0.138 -0.110
9 0.066 0.004
10 0.105 0.093
11 -0.132 -0.064
12 -0.130 -0.173

Figure 31: Correlogram of the residuals of the AR(7) model for rcpi

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 235 / 248

Diagnostic checking on the residuals from AR(7)

Breusch-Godfrey
Serial Correlation LM Test:
F-statistic 3.429180 F(2,84) 0.037
Obs*R-squared 7.020020 Chi-Sq(2) 0.029

Heteroskedasticity Test: ARCH

F-statistic 0.025863 F(1,90) 0.872
Obs*R-squared 0.026431 Chi-Sq(1) 0.870

Heteroskedasticity Test: White

F-statistic 0.188135 F(14,78) 0.999
Obs*R-squared 3.037834 Chi-Sq(14) 0.999
Scaled explained SS 3.384554 Chi-Sq(14) 0.998
Table 27: Diagnostic checking on the residuals from AR(7) model

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 236 / 248

The Crisis Period

It is interesting to note how the dummy variable affects the estimated

parameter values, since it is as if we had two separate equations

Before the crisis: rcpit = 0.14 + 0.46rcpit−1 + 0.13rcpit−2 + ε1t

(47a)
During the crisis: rcpit = −0.43 + 2.88rcpit−1 − 2.76rcpit−2
+1.15rcpit−7 + ε2t (47b)

In particular, the persistence (as measured by the sum of the AR

coefficients) increases substantially during the crisis period.
Computing forecasts with this specification amounts to exploiting the
information of both models (47a) and (47b) combined into a single
specification.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 237 / 248

Forecasts from AR(7) model for rcpi against the
actuals 2010 - 2012
1.0

0.8

0.6

0.4

0.2

0.0

-0.2

-0.4
I II III IV I II III IV I II III IV
2010 2011 2012

RCPI RCPIF_1STEP
RCPIF_2STEPS RCPIF_HSTEPS

Figure 32: Forecasts from AR(7) model for rcpi against the actuals 2010 - 2012

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 238 / 248

The Crisis Period

The overall impression from these results is that the forecast

performance remains overall stable after the crisis, once the model
is appropriately modified to allow for parameter changes during the
crisis.
It is interesting to consider what would happen if we did not include
the dummies for the financial crisis, and just produced forecasts with
a simple AR(2) model.
The estimation of the model without augmentation for the dummy
produces the following equation:

rcpit = 0.048 + 0.7rcpit−1 + 0.07rcpit−2 + ε1t . (48)

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 239 / 248

The Crisis Period

The associated residuals, graphed in Figure 33, are large and

negative during the global financial crisis.
The proceeding slides show the forecasts obtained through model
(48) against those we obtained through the richer specification
(47a)-(47b), both at one-, two-, and one- to h-steps ahead.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 240 / 248

Residuals of the mis-specified model

1.0 0

0.5 -1

0.0 -2

-0.5

-1.0
86 88 90 92 94 96 98 00 02 04 06 08

Residual Actual Fitted

Figure 33: Residuals of the mis-specified model which does not contain the dummy variable
c accounting
Eric Ghysels &for the crisisMarcellino
Massimiliano 2018 Edition 241 / 248
1-step

Figure 34: 1-step forecasts not accounting for the crisis (no dummy) against those obtained
through the richer specification containing the dummy, against the actuals

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 242 / 248

h-steps

Figure 35: h-steps forecasts not accounting for the crisis (no dummy) against those obtained
through the richer specification containing the dummy, against the actuals

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 243 / 248

2-steps

Figure 36: 2-steps forecasts not accounting for the crisis (no dummy) against those obtained
through the richer specification containing the dummy, against the actuals

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 244 / 248

2-steps direct

Figure 37: 2-steps direct forecasts not accounting for the crisis (no dummy) against those
obtained through the richer specification containing the dummy, against the actuals

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 245 / 248

The models we have considered are overall rather simple, since they
basically explain the variable of interest using its past only. Hence, a
natural question is whether these univariate time series methods are
useful in practice. The answer is yes, for the following reasons.
From a theoretical perspective, any weakly stationary process (or
integrated process after proper differencing) can be written as an
MA(∞), and under mild conditions the latter can be approximated by
an ARMA model.
The high persistence of several economic variables suggests that
forecasting methods that exploit the present and the past behavior of
the variable to predict the future can perform well.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 247 / 248

Concluding Remarks

Forecasts from these models can provide a benchmark for

comparison of more elaborate forecasts or can be combined with
forecasts from other models to assess whether a lower loss function
can be obtained.
Forecast failure of these models can provide an indication of the type
of information that is missing in more elaborate models, e.g., it can
suggest that a substantial amount of dynamics should be included.
More practically, empirically these forecasting methods tend to
perform well for a variety of variables and across different evaluation
samples, when used for short-term forecasting.

c Eric Ghysels & Massimiliano Marcellino 2018 Edition 248 / 248

1 Uts 260303 195852
No ratings yet
1 Uts 260303 195852
105 pages
Univariate Time Series Analysis Overview
No ratings yet
Univariate Time Series Analysis Overview
75 pages
Time Series Analysis: Stationarity & ARMA
No ratings yet
Time Series Analysis: Stationarity & ARMA
6 pages
Ec 221 L 5 Slnewx
No ratings yet
Ec 221 L 5 Slnewx
34 pages
Stationary Univariate Time Series
No ratings yet
Stationary Univariate Time Series
28 pages
ST304 Time Series Models Overview
No ratings yet
ST304 Time Series Models Overview
20 pages
CH 3. Introduction To Time Series Regression and Forecasting 2025
No ratings yet
CH 3. Introduction To Time Series Regression and Forecasting 2025
67 pages
Time Series Analysis in Econometrics
100% (1)
Time Series Analysis in Econometrics
246 pages
Linear Time
No ratings yet
Linear Time
11 pages
Understanding AR, MA, and ARMA Models
No ratings yet
Understanding AR, MA, and ARMA Models
37 pages
Univariate Time Series Models
No ratings yet
Univariate Time Series Models
13 pages
Introduction to AR, MA, and ARMA Models
No ratings yet
Introduction to AR, MA, and ARMA Models
15 pages
ARMA Models and Stationarity Analysis
No ratings yet
ARMA Models and Stationarity Analysis
20 pages
Understanding ARMA Processes in Time Series
No ratings yet
Understanding ARMA Processes in Time Series
28 pages
Module1_TimeSeries_Notes
No ratings yet
Module1_TimeSeries_Notes
5 pages
AR(1) Autoregressive Model Overview
No ratings yet
AR(1) Autoregressive Model Overview
12 pages
Understanding Autoregressive Models
No ratings yet
Understanding Autoregressive Models
18 pages
Understanding ARMA and ARIMA Models
No ratings yet
Understanding ARMA and ARIMA Models
46 pages
Chapter 5
No ratings yet
Chapter 5
7 pages
Box Jenkins Time Series Modeling Steps
No ratings yet
Box Jenkins Time Series Modeling Steps
68 pages
Chapter 4
No ratings yet
Chapter 4
96 pages
Stochastic Time Series Models Explained
No ratings yet
Stochastic Time Series Models Explained
96 pages
Time Series Models Exam Guide
No ratings yet
Time Series Models Exam Guide
13 pages
001 Time Series Reader
No ratings yet
001 Time Series Reader
8 pages
Understanding AR(p) Models in Time Series
No ratings yet
Understanding AR(p) Models in Time Series
12 pages
Univariate Time Series Analysis Guide
No ratings yet
Univariate Time Series Analysis Guide
37 pages
ARMA Process Covariance Functions
No ratings yet
ARMA Process Covariance Functions
9 pages
Autoregressive Models for Time Series
No ratings yet
Autoregressive Models for Time Series
24 pages
AR Model Stationarity and Roots Analysis
No ratings yet
AR Model Stationarity and Roots Analysis
40 pages
Weak Stationarity in Time Series Models
No ratings yet
Weak Stationarity in Time Series Models
32 pages
Time Series Analysis: AR Models Explained
No ratings yet
Time Series Analysis: AR Models Explained
32 pages
MATH5845_Lecture_Slies_2026
No ratings yet
MATH5845_Lecture_Slies_2026
59 pages
Lec 6 Autoregressive (AR) Models and Moving Average (MA) Models
No ratings yet
Lec 6 Autoregressive (AR) Models and Moving Average (MA) Models
41 pages
Univariate Time Series Models Overview
No ratings yet
Univariate Time Series Models Overview
54 pages
Stationarity and Autocorrelation in Time Series
No ratings yet
Stationarity and Autocorrelation in Time Series
24 pages
Understanding ARIMA Time Series Models
No ratings yet
Understanding ARIMA Time Series Models
25 pages
ARMA Process and Wold's Theorem Overview
No ratings yet
ARMA Process and Wold's Theorem Overview
13 pages
Advanced Econometrics Overview
No ratings yet
Advanced Econometrics Overview
48 pages
Examples of Stationary Time Series
No ratings yet
Examples of Stationary Time Series
21 pages
Box Jenkins Technique: ARMA & ARIMA Models
No ratings yet
Box Jenkins Technique: ARMA & ARIMA Models
25 pages
4 Arma
No ratings yet
4 Arma
35 pages
ARIMA Models in Time Series Analysis
100% (1)
ARIMA Models in Time Series Analysis
50 pages
Time Series Analysis Exercises
100% (1)
Time Series Analysis Exercises
30 pages
AoFD 2025 L9
No ratings yet
AoFD 2025 L9
46 pages
AR(1) Model: Stationarity and Forecasting
No ratings yet
AR(1) Model: Stationarity and Forecasting
43 pages
Linear Stochastic Models Overview
No ratings yet
Linear Stochastic Models Overview
10 pages
ARMA Models: Time Series Analysis Guide
No ratings yet
ARMA Models: Time Series Analysis Guide
21 pages
Intro 1
No ratings yet
Intro 1
36 pages
Kalman Filter for Nonstationary Time Series
No ratings yet
Kalman Filter for Nonstationary Time Series
33 pages
Stationarity and Cointegration Explained
No ratings yet
Stationarity and Cointegration Explained
52 pages
Univariate Time Series Models Overview
100% (2)
Univariate Time Series Models Overview
72 pages
Introduction to ARMA Models Explained
No ratings yet
Introduction to ARMA Models Explained
15 pages
Time Series Analysis Fundamentals
No ratings yet
Time Series Analysis Fundamentals
38 pages
Durbin-Watson and Assumption Violation
No ratings yet
Durbin-Watson and Assumption Violation
185 pages
Lect 5 S2020
No ratings yet
Lect 5 S2020
12 pages
Advanced Time-Series Models Overview
No ratings yet
Advanced Time-Series Models Overview
110 pages
Univariate Time Series Modeling Guide
No ratings yet
Univariate Time Series Modeling Guide
31 pages
Univariate Time Series Analysis
No ratings yet
Univariate Time Series Analysis
21 pages
Time SeriesClass Notes
No ratings yet
Time SeriesClass Notes
18 pages
EWMA Methods for Covariate Shift Detection
No ratings yet
EWMA Methods for Covariate Shift Detection
11 pages
Inflation and Unemployment's Impact on Sierra Leone's Growth
No ratings yet
Inflation and Unemployment's Impact on Sierra Leone's Growth
19 pages
Stochastic Hydrology Overview
No ratings yet
Stochastic Hydrology Overview
24 pages
KTU M.Tech AI & Data Science Syllabus
No ratings yet
KTU M.Tech AI & Data Science Syllabus
73 pages
ECE II B.Tech Course Structure 2015-16
No ratings yet
ECE II B.Tech Course Structure 2015-16
16 pages
Understanding Random Signals and Processes
No ratings yet
Understanding Random Signals and Processes
34 pages
Key Probability Concepts and Definitions
No ratings yet
Key Probability Concepts and Definitions
7 pages
Streamflow Measurement Protocols in Amazon
No ratings yet
Streamflow Measurement Protocols in Amazon
27 pages
Time Series Analysis Techniques Explained
No ratings yet
Time Series Analysis Techniques Explained
7 pages
Understanding Unit Root Processes
No ratings yet
Understanding Unit Root Processes
33 pages
Time Series Analysis: Stationarity & Models
No ratings yet
Time Series Analysis: Stationarity & Models
30 pages
ARMA Model Overview and Analysis
No ratings yet
ARMA Model Overview and Analysis
39 pages
Price Discovery in India's Gold Market
No ratings yet
Price Discovery in India's Gold Market
9 pages
O-GARCH Volatility Analysis in Finance
No ratings yet
O-GARCH Volatility Analysis in Finance
12 pages
Time Series Analysis and ARIMA Model Guide
No ratings yet
Time Series Analysis and ARIMA Model Guide
11 pages
Random Processes and Sequences: Chapter Three
No ratings yet
Random Processes and Sequences: Chapter Three
52 pages
Triple Barrier Method in Finance Analysis
No ratings yet
Triple Barrier Method in Finance Analysis
94 pages
Spurious Regression in Time Series Analysis
No ratings yet
Spurious Regression in Time Series Analysis
3 pages
Sustainability 14 09763
No ratings yet
Sustainability 14 09763
15 pages
Stochastic Modelling Overview
No ratings yet
Stochastic Modelling Overview
130 pages
KITS Singapur B.Tech ECE Syllabus R-23
No ratings yet
KITS Singapur B.Tech ECE Syllabus R-23
38 pages
Random Variables and Stochastic Processes
No ratings yet
Random Variables and Stochastic Processes
2 pages
EDA Techniques for Time Series Analysis
No ratings yet
EDA Techniques for Time Series Analysis
20 pages
AR(2) Process and Stationarity Conditions
No ratings yet
AR(2) Process and Stationarity Conditions
18 pages
Quantitative Trading Resource Roadmap
100% (2)
Quantitative Trading Resource Roadmap
21 pages
Bitcoin and Ethereum: Causality and Volatility Analysis
No ratings yet
Bitcoin and Ethereum: Causality and Volatility Analysis
40 pages
Automated Data Selection for Surface Waves
No ratings yet
Automated Data Selection for Surface Waves
18 pages
Recursive Macroeconomic Theory 4th Edition Lars Ljungqvist Ebook Testbank Solutions Full Size Version
100% (4)
Recursive Macroeconomic Theory 4th Edition Lars Ljungqvist Ebook Testbank Solutions Full Size Version
116 pages
Introduction to Time Series Modeling
No ratings yet
Introduction to Time Series Modeling
5 pages
Python for Machine Learning Basics
No ratings yet
Python for Machine Learning Basics
11 pages