g13dj:: Time Series Analysis (NAG Toolbox)

nag_tsa_multi_varma_forecast (g13dj) computes forecasts of a multivariate time series. It is assumed that a vector ARMA model has already been fitted to the appropriately differenced/transformed time series using nag_tsa_multi_varma_estimate (g13dd). The standard deviations of the forecast errors are also returned. A reference vector is set up so that, should future series values become available, the forecasts and their standard errors may be updated by calling nag_tsa_multi_varma_update (g13dk).

Syntax

Description

Let the vector

Z_{t} = {(z_{1 t}, z_{2 t}, \dots, z_{k t})}^{T}

, for

t = 1, 2, \dots, n

, denote a

k

-dimensional time series for which forecasts of

Z_{n + 1}, Z_{n + 2}, \dots, Z_{n + l_{\max}}

are required. Let

W_{t} = {(w_{1 t}, w_{2 t}, \dots, w_{k t})}^{T}

be defined as follows:

w_{i t} = δ_{i} (B) z_{i t}^{*}, i = 1, 2, \dots, k,

where

δ_{i} (B)

is the differencing operator applied to the

i

th series and where

z_{i t}^{*}

is equal to either

z_{i t}

\sqrt{z_{i t}}

\log_{e} (z_{i t})

depending on whether or not a transformation was required to stabilize the variance before fitting the model.

If the order of differencing required for the

i

th series is

d_{i}

, then the differencing operator for the

i

th series is defined by

δ_{i} (B) = 1 - δ_{i 1} B - δ_{i 2} B^{2} - \dots - δ_{i d_{i}} B^{d_{i}}

where

B

is the backward shift operator; that is,

B Z_{t} = Z_{t - 1}

. The differencing parameters

δ_{i j}

, for

i = 1, 2, \dots, k

and

j = 1, 2, \dots, d_{i}

, must be supplied by you. If the

i

th series does not require differencing, then

d_{i} = 0

W_{t}

is assumed to follow a multivariate ARMA model of the form:

W_{t} - μ = ϕ_{1} (W_{t - 1} - μ) + ϕ_{2} (W_{t - 2} - μ) + \dots + ϕ_{p} (W_{t - p} - μ) + ε_{t} - θ_{1} ε_{t - 1} - \dots - θ_{q} ε_{t - q},

(1)

where

ε_{t} = {(ε_{1 t}, ε_{2 t}, \dots, ε_{k t})}^{T}

, for

t = 1, 2, \dots, n

, is a vector of

k

residual series assumed to be Normally distributed with zero mean and positive definite covariance matrix

Σ

. The components of

ε_{t}

are assumed to be uncorrelated at non-simultaneous lags. The

ϕ_{i}

and

θ_{j}

are

k

k

matrices of parameters. The matrices

ϕ_{i}

, for

i = 1, 2, \dots, p

, are the autoregressive (AR) parameter matrices, and the matrices

θ_{i}

, for

i = 1, 2, \dots, q

, the moving average (MA) parameter matrices. The parameters in the model are thus the

p

(

k

k

)

ϕ

-matrices, the

q

(

k

k

)

θ

-matrices, the mean vector

μ

and the residual error covariance matrix

Σ

. The ARMA model (1) must be both stationary and invertible; see nag_tsa_uni_arma_roots (g13dx) for a method of checking these conditions.

Forecasts are computed using a multivariate version of the procedure described in Box and Jenkins (1976). If

{\hat{Z}}_{n}^{*} (l)

denotes the forecast of

Z_{n + l}^{*}

, then

{\hat{Z}}_{n}^{*} (l)

is taken to be that linear function of

Z_{n}^{*}, Z_{n - 1}^{*}, \dots

which minimizes the elements of

E \{e_{n} (l) e_{n}^{'} (l)\}

where

e_{n} (l) = Z_{n + l}^{*} - {\hat{Z}}_{n}^{*} (l)

is the forecast error.

{\hat{Z}}_{n}^{*} (l)

is referred to as the linear minimum mean square error forecast of

Z_{n + l}^{*}

The linear predictor which minimizes the mean square error may be expressed as

{\hat{Z}}_{n}^{*} (l) = τ + ψ_{l} ε_{n} + ψ_{l + 1} ε_{n - 1} + ψ_{l + 2} ε_{n - 2} + \dots .

The forecast error at

t

for lead

l

is then

e_{n} (l) = Z_{n + l}^{*} - {\hat{Z}}_{n}^{*} (l) = ε_{n + l} + ψ_{1} ε_{n + l - 1} + ψ_{2} ε_{n + l - 2} + \dots + ψ_{l - 1} ε_{n + 1} .

Let

d = \max (d_{i})

, for

i = 1, 2, \dots, k

. Unless

q = 0

the function requires estimates of

ε_{t}

, for

t = d + 1, \dots, n

, which are obtainable from nag_tsa_multi_varma_estimate (g13dd). The terms

ε_{t}

are assumed to be zero, for

t = n + 1, \dots, n + l_{\max}

. You may use nag_tsa_multi_varma_update (g13dk) to update these

l_{\max}

forecasts should further observations,

Z_{n + 1}, Z_{n + 2}, \dots

, become available. Note that when

l_{\max}

or more further observations are available then nag_tsa_multi_varma_forecast (g13dj) must be used to produce new forecasts for

Z_{n + l_{\max} + 1}, Z_{n + l_{\max} + 2}, \dots

, should they be required.

References

Parameters

Compulsory Input Parameters

Optional Input Parameters

Output Parameters

Error Indicators and Warnings

Accuracy

Further Comments

The same differencing operator does not have to be applied to all the series. For example, suppose we have

k = 2

, and wish to apply the second order differencing operator

\nabla^{2}

to the first series and the first-order differencing operator

\nabla

to the second series:

\begin{matrix} w_{1 t} = \nabla^{2} z_{1 t} = {(1 - B)}^{2} z_{1 t} = (1 - 2 B + B^{2}) Z_{1 t},   and \\ w_{2 t} = \nabla z_{2 t} = (1 - B) z_{2 t} . \end{matrix}

Then

d_{1} = 2, d_{2} = 1

d = \max (d_{1}, d_{2}) = 2

, and

delta = [\begin{array}{r} δ_{11} & δ_{12} \\ δ_{21} \end{array}] = [\begin{array}{r} 2 & - 1 \\ 1 \end{array}] .

Note: although differencing may already have been applied prior to the model fitting stage, the differencing parameters supplied in delta are part of the model definition and are still required by this function to produce the forecasts.

Example

function g13dj_example


fprintf('g13dj example results\n\n');

% Series
z = [-1.490 -1.620  5.200  6.230  6.210  5.860  4.090  3.180  ...
      2.620  1.490  1.170  0.850 -0.350  0.240  2.440  2.580  ...
      2.040  0.400  2.260  3.340  5.090  5.000  4.780  4.110  ...
      3.450  1.650  1.290  4.090  6.320  7.500  3.890  1.580  ...
      5.210  5.250  4.930  7.380  5.870  5.810  9.680  9.070  ...
      7.290  7.840  7.550  7.320  7.970  7.760  7.000  8.350;
      7.340  6.350  6.960  8.540  6.620  4.970  4.550  4.810  ...
      4.750  4.760 10.880 10.010 11.620 10.360  6.400  6.240  ...
      7.930  4.040  3.730  5.600  5.350  6.810  8.270  7.680  ...
      6.650  6.080 10.250  9.140 17.750 13.300  9.630  6.800  ...
      4.080  5.060  4.940  6.650  7.940 10.760 11.890  5.850  ...
      9.010  7.500 10.020 10.380  8.150  8.370 10.730 12.140];
[k,n] = size(z);

% Difference /transform series
tr    = {'N'; 'N'};
id    = [int64(0);0];
delta = [0;  0];
[w, nd, ifail] = g13dl( ...
			z, tr, id, delta);

% VARMA info
ip    = int64(1);
iq    = int64(0);
mean_p = true;

% Initial parameter estimates and free parameter flags
par    = zeros(6, 1);
parhld = [false;  false;  true;  false;  false;  false];

% Exact likelihood
exact  = true;
% control parameters
iprint = int64(-1);
cgetol = 0.0001;
ishow  = int64(0);

qq = [0, 0; 0, 0];

% Fit VARMA
[par, qq, ~, ~, v, ~, ~, ifail] = ...
  g13dd( ...
	 ip, iq, mean_p, par, qq, w, parhld, exact, iprint, cgetol, ...
	 ishow, 'n', nd);

% Perform forecast
lmax = int64(5);
lref = int64(150);
mean_p = 'M';
[qq, predz, sefz, ref, ifail] = ...
  g13dj( ...
	 z, tr, id, delta, ip, iq, mean_p, par, qq, v, lmax, lref);

% Display results
fprintf(' Forecast Summary Table\n');
fprintf(' ----------------------\n\n');
fprintf(' Forecast origin is set at t = %4d\n\n', n);
loop = lmax/5;
if mod(lmax,5)~=0
  loop = loop + 1;
end
for j = 1:loop
  i2 = (j-1)*5;
  l2 = min(i2+5,lmax);
  fprintf('Lead Time %14s',' ');
  fprintf('%7d', [i2+1:l2]);
  fprintf('\n\n');
  for i = 1:k
    fprintf('Series %d : Forecast       ', i);
    fprintf('%7.2f',predz(i,i2+1:l2));
    fprintf('\n%8s : Standard Error ', ' ');
    fprintf('%7.2f',sefz(i,i2+1:l2));
    fprintf('\n');
  end
end

On entry,	$k < 1$ ,
or	$n < 3$ ,
or	$kmax < k$ ,
or	$id (i) < 0$ for some $i = 1, 2, \dots, k$ ,
or	$id (i) \geq n - \max (ip, iq)$ for some $i = 1, 2, \dots, k$ ,
or	$ip < 0$ ,
or	$iq < 0$ ,
or	$mean_p \neq'M'$ or $'Z'$ ,
or	$lpar < (ip + iq) \times k \times k + k$ , and $mean_p ='M'$ ,
or	$lpar < (ip + iq) \times k \times k$ and $mean_p ='Z'$ ,
or	$n \times k \leq (ip + iq) \times k \times k + k + k (k + 1) / 2$ , and $mean_p ='M'$ ,
or	$n \times k \leq (ip + iq) \times k \times k + k (k + 1) / 2$ and $mean_p ='Z'$ ,
or	$lmax < 1$ ,
or	$lref < (lmax - 1) \times k \times k + 2 \times k \times lmax + k$ ,
or	lwork is too small,
or	liwork is too small.

NAG Toolbox: nag_tsa_multi_varma_forecast (g13dj)

▸▿ Contents

Purpose

Syntax

Description

References

Parameters

Compulsory Input Parameters

Optional Input Parameters

Output Parameters

Error Indicators and Warnings

Accuracy

Further Comments

Example