NAG Library Routine Document
g13baf (multi_filter_arima)
1
Purpose
g13baf filters a time series by an ARIMA model.
2
Specification
Fortran Interface
Subroutine g13baf ( 
y, ny, mr, nmr, par, npar, cy, wa, nwa, b, nb, ifail) 
Integer, Intent (In)  ::  ny, mr(nmr), nmr, npar, nwa, nb  Integer, Intent (Inout)  ::  ifail  Real (Kind=nag_wp), Intent (In)  ::  y(ny), par(npar), cy  Real (Kind=nag_wp), Intent (Out)  ::  wa(nwa), b(nb) 

C Header Interface
#include nagmk26.h
void 
g13baf_ (const double y[], const Integer *ny, const Integer mr[], const Integer *nmr, const double par[], const Integer *npar, const double *cy, double wa[], const Integer *nwa, double b[], const Integer *nb, Integer *ifail) 

3
Description
From a given series ${y}_{1},{y}_{2},\dots ,{y}_{n}$, a new series ${b}_{1},{b}_{2},\dots ,{b}_{n}$ is calculated using a supplied (filtering) ARIMA model. This model will be one which has previously been fitted to a series ${x}_{t}$ with residuals ${a}_{t}$. The equations defining ${b}_{t}$ in terms of ${y}_{t}$ are very similar to those by which ${a}_{t}$ is obtained from ${x}_{t}$. The only dissimilarity is that no constant correction is applied after differencing. This is because the series ${y}_{t}$ is generally distinct from the series ${x}_{t}$ with which the model is associated, though ${y}_{t}$ may be related to ${x}_{t}$. Whilst it is appropriate to apply the ARIMA model to ${y}_{t}$ so as to preserve the same relationship between ${b}_{t}$ and ${a}_{t}$ as exists between ${y}_{t}$ and ${x}_{t}$, the constant term in the ARIMA model is inappropriate for ${y}_{t}$. The consequence is that ${b}_{t}$ will not necessarily have zero mean.
The equations are precisely:
the appropriate differencing of
${y}_{t}$; both the seasonal and nonseasonal inverted autoregressive operations are then applied,
followed by the inverted moving average operations
Because the filtered series value
${b}_{t}$ depends on present and past values
${y}_{t},{y}_{t1},\dots \text{}$, there is a problem arising from ignorance of
${y}_{0},{y}_{1},\dots \text{}$ which particularly affects calculation of the early values
${b}_{1},{b}_{2},\dots \text{}$, causing ‘transient errors’. The routine allows two possibilities.
(i) 
The equations (1), (2) and (3) are applied from successively later time points so that all terms on their righthand sides are known, with ${v}_{t}$ being defined for $t=\left(1+d+s\times D+s\times P\right),\dots ,n$. Equations (4) and (5) are then applied over the same range, taking any values on the righthand side associated with previous time points to be zero.
This procedure may still however result in unacceptably large transient errors in early values of ${b}_{t}$. 
(ii) 
The unknown values ${y}_{0},{y}_{1},\dots \text{}$ are estimated by backforecasting. This requires that an ARIMA model distinct from that which has been supplied for filtering, should have been previously fitted to ${y}_{t}$. 
For efficiency, you are asked to supply both this ARIMA model for
${y}_{t}$ and a limited number of backforecasts which are prefixed to the known values of
${y}_{t}$. Within the routine further backforecasts of
${y}_{t}$, and the series
${w}_{t}$,
${u}_{t}$,
${v}_{t}$ in
(1),
(2) and
(3) are then easily calculated, and a set of linear equations solved for backforecasts of
${z}_{t},{b}_{t}$ for use in
(4) and
(5) in the case that
$q+Q>0$.
Even if the best model for
${y}_{t}$ is not available, a very approximate guess such as
or
can help to reduce the transients substantially.
The backforecasts which need to be prefixed to
${y}_{t}$ are of length
${Q}_{y}^{\prime}={q}_{y}+{s}_{y}\times {Q}_{y}$, where
${q}_{y}$ and
${Q}_{y}$ are the nonseasonal and seasonal moving average orders and
${s}_{y}$ the seasonal period for the ARIMA model of
${y}_{t}$. Thus you need not carry out the backforecasting exercise if
${Q}_{y}^{\prime}=0$. Otherwise, the series
${y}_{1},{y}_{2},\dots ,{y}_{n}$ should be reversed to obtain
${y}_{n},{y}_{n1},\dots ,{y}_{1}$ and
g13ajf
should be used to forecast
${Q}_{y}^{\prime}$ values,
${\hat{y}}_{0},\dots ,{\hat{y}}_{1{Q}_{y}^{\prime}}$. The ARIMA model used is that fitted to
${y}_{t}$ (as a forward series) except that, if
${d}_{y}+{D}_{y}$ is odd, the constant should be changed in sign (to allow, for example, for the fact that a forward upward trend is a reversed downward trend). The ARIMA model for
${y}_{t}$ supplied to the filtering routine must however have the appropriate constant for the forward series.
The series ${\hat{y}}_{1{Q}_{y}^{\prime}},\dots ,{\hat{y}}_{0},{y}_{1},\dots ,{y}_{n}$ is then supplied to the routine, and a corresponding set of values returned for ${b}_{t}$.
4
References
Box G E P and Jenkins G M (1976) Time Series Analysis: Forecasting and Control (Revised Edition) Holden–Day
5
Arguments
 1: $\mathbf{y}\left({\mathbf{ny}}\right)$ – Real (Kind=nag_wp) arrayInput

On entry: the
${Q}_{y}^{\prime}$ backforecasts, starting with backforecast at time
$1{Q}_{y}^{\prime}$ to backforecast at time
$0$, followed by the time series starting at time
$1$, where
${Q}_{y}^{\prime}={\mathbf{mr}}\left(10\right)+{\mathbf{mr}}\left(13\right)\times {\mathbf{mr}}\left(14\right)$. If there are no backforecasts, either because the ARIMA model for the time series is not known, or because it is known but has no moving average terms, then the time series starts at the beginning of
y.
 2: $\mathbf{ny}$ – IntegerInput

On entry: the total number of backforecasts and time series data points in array
y.
Constraint:
${\mathbf{ny}}\ge \mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1+{Q}_{y}^{\prime},{\mathbf{npar}}\right)$.
 3: $\mathbf{mr}\left({\mathbf{nmr}}\right)$ – Integer arrayInput

On entry: the orders vector for the filtering model, followed by the orders vector for the ARIMA model for the time series if the latter is known. The orders appear in the standard sequence
$\left(p,d,q,P,D,Q,s\right)$ as given in the
G13 Chapter Introduction. If the ARIMA model for the time series is supplied, the routine will assume that the first
${Q}_{y}^{\prime}$ values of the array
y are backforecasts.
Constraints:
the filtering model is restricted in the following ways:
 ${\mathbf{mr}}\left(1\right)+{\mathbf{mr}}\left(3\right)+{\mathbf{mr}}\left(4\right)+{\mathbf{mr}}\left(6\right)>0$, i.e., filtering by a model which contains only differencing terms is not permitted;
 ${\mathbf{mr}}\left(\mathit{k}\right)\ge 0$, for $\mathit{k}=1,2,\dots ,7$;
 if ${\mathbf{mr}}\left(7\right)=0$, ${\mathbf{mr}}\left(4\right)+{\mathbf{mr}}\left(5\right)+{\mathbf{mr}}\left(6\right)=0$;
 if ${\mathbf{mr}}\left(7\right)\ne 0$, ${\mathbf{mr}}\left(4\right)+{\mathbf{mr}}\left(5\right)+{\mathbf{mr}}\left(6\right)\ne 0$;
 ${\mathbf{mr}}\left(7\right)\ne 1$.
the ARIMA model for the time series is restricted in the following ways:
 ${\mathbf{mr}}\left(\mathit{k}\right)\ge 0$, for $\mathit{k}=8,9,\dots ,14$;
 if ${\mathbf{mr}}\left(14\right)=0$, ${\mathbf{mr}}\left(11\right)+{\mathbf{mr}}\left(12\right)+{\mathbf{mr}}\left(13\right)=0$;
 if ${\mathbf{mr}}\left(14\right)\ne 0$, ${\mathbf{mr}}\left(11\right)+{\mathbf{mr}}\left(12\right)+{\mathbf{mr}}\left(13\right)\ne 0$;
 ${\mathbf{mr}}\left(14\right)\ne 1$.
 4: $\mathbf{nmr}$ – IntegerInput

On entry: the number of values specified in the array
mr. It takes the value
$7$ if no ARIMA model for the time series is supplied but otherwise it takes the value
$14$. Thus
nmr acts as an indicator as to whether backforecasting can be carried out.
Constraint:
${\mathbf{nmr}}=7$ or $14$.
 5: $\mathbf{par}\left({\mathbf{npar}}\right)$ – Real (Kind=nag_wp) arrayInput

On entry: the parameters of the filtering model, followed by the parameters of the ARIMA model for the time series, if supplied. Within each model the parameters are in the standard order of nonseasonal AR and MA followed by seasonal AR and MA.
 6: $\mathbf{npar}$ – IntegerInput

On entry: the total number of parameters held in array
par.
Constraints:
 if ${\mathbf{nmr}}=7$, ${\mathbf{npar}}={\mathbf{mr}}\left(1\right)+{\mathbf{mr}}\left(3\right)+{\mathbf{mr}}\left(4\right)+{\mathbf{mr}}\left(6\right)$;
 if ${\mathbf{nmr}}=14$, ${\mathbf{npar}}={\mathbf{mr}}\left(1\right)+{\mathbf{mr}}\left(3\right)+{\mathbf{mr}}\left(4\right)+{\mathbf{mr}}\left(6\right)+\phantom{\rule{0ex}{0ex}}{\mathbf{mr}}\left(8\right)+{\mathbf{mr}}\left(10\right)+{\mathbf{mr}}\left(11\right)+{\mathbf{mr}}\left(13\right)$.
Note: the first constraint (i.e.,
${\mathbf{mr}}\left(1\right)+{\mathbf{mr}}\left(3\right)+{\mathbf{mr}}\left(4\right)+{\mathbf{mr}}\left(6\right)>0$) on the orders of the filtering model, in argument
mr, ensures that
${\mathbf{npar}}>0$.
 7: $\mathbf{cy}$ – Real (Kind=nag_wp)Input

On entry: if the ARIMA model is known (i.e.,
${\mathbf{nmr}}=14$),
cy must specify the constant term of the ARIMA model for the time series. If this model is not known (i.e.,
${\mathbf{nmr}}=7$),
cy is not used.
 8: $\mathbf{wa}\left({\mathbf{nwa}}\right)$ – Real (Kind=nag_wp) arrayWorkspace
 9: $\mathbf{nwa}$ – IntegerInput

On entry: the dimension of the array
wa as declared in the (sub)program from which
g13baf is called. Workspace is only required if the ARIMA model for the time series is known.
Constraints:
let $K={\mathbf{mr}}\left(3\right)+{\mathbf{mr}}\left(6\right)\times {\mathbf{mr}}\left(7\right)+{\mathbf{mr}}\left(8\right)+{\mathbf{mr}}\left(9\right)+\left({\mathbf{mr}}\left(11\right)+{\mathbf{mr}}\left(12\right)\right)\times {\mathbf{mr}}\left(14\right)$,
then
 if ${\mathbf{nmr}}=14$, ${\mathbf{nwa}}\ge K\times \left(K+2\right)$;
 if ${\mathbf{nmr}}=7$, ${\mathbf{nwa}}\ge 1$.
 10: $\mathbf{b}\left({\mathbf{nb}}\right)$ – Real (Kind=nag_wp) arrayOutput

On exit: the filtered output series. If the ARIMA model for the time series was known, and hence
${Q}_{y}^{\prime}$ backforecasts were supplied in
y, then
b contains
${Q}_{y}^{\prime}$ ‘filtered’ backforecasts followed by the filtered series. Otherwise, the filtered series begins at the start of
b just as the original series began at the start of
y. In either case, if the value of the series at time
$t$ is held in
${\mathbf{y}}\left(t\right)$, then the filtered value at time
$t$ is held in
${\mathbf{b}}\left(t\right)$.
 11: $\mathbf{nb}$ – IntegerInput

On entry: the dimension of the array
b as declared in the (sub)program from which
g13baf is called. In addition to holding the returned filtered series,
b is also used as an intermediate work array if the ARIMA model for the time series was known.
Constraints:
 if ${\mathbf{nmr}}=14$, ${\mathbf{nb}}\ge {\mathbf{ny}}+\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left({K}_{3},{K}_{1}+{K}_{2}\right)$;
 if ${\mathbf{nmr}}=7$, ${\mathbf{nb}}\ge {\mathbf{ny}}$.
Where
 ${K}_{1}={\mathbf{mr}}\left(1\right)+{\mathbf{mr}}\left(4\right)\times {\mathbf{mr}}\left(7\right)$;
 ${K}_{2}={\mathbf{mr}}\left(2\right)+{\mathbf{mr}}\left(5\right)\times {\mathbf{mr}}\left(7\right)$;
 ${K}_{3}={\mathbf{mr}}\left(3\right)+{\mathbf{mr}}\left(6\right)\times {\mathbf{mr}}\left(7\right)$.
 12: $\mathbf{ifail}$ – IntegerInput/Output

On entry:
ifail must be set to
$0$,
$1\text{ or}1$. If you are unfamiliar with this argument you should refer to
Section 3.4 in How to Use the NAG Library and its Documentation for details.
For environments where it might be inappropriate to halt program execution when an error is detected, the value
$1\text{ or}1$ is recommended. If the output of error messages is undesirable, then the value
$1$ is recommended. Otherwise, if you are not familiar with this argument, the recommended value is
$0$.
When the value $\mathbf{1}\text{ or}\mathbf{1}$ is used it is essential to test the value of ifail on exit.
On exit:
${\mathbf{ifail}}={\mathbf{0}}$ unless the routine detects an error or a warning has been flagged (see
Section 6).
6
Error Indicators and Warnings
If on entry
${\mathbf{ifail}}=0$ or
$1$, explanatory error messages are output on the current error message unit (as defined by
x04aaf).
Errors or warnings detected by the routine:
 ${\mathbf{ifail}}=1$

On entry,  ${\mathbf{nmr}}\ne 7$ and ${\mathbf{nmr}}\ne 14$. 
 ${\mathbf{ifail}}=2$

On entry,  the orders vector mr does not satisfy the constraints given in Section 5. 
 ${\mathbf{ifail}}=3$

 ${\mathbf{ifail}}=4$

On entry,  ny is too small to successfully carry out the requested filtering, (see Section 5). 
 ${\mathbf{ifail}}=5$

On entry,  the work array wa is too small. 
 ${\mathbf{ifail}}=6$

On entry,  the array b is too small. 
 ${\mathbf{ifail}}=7$

The orders vector for the filtering model is invalid.
 ${\mathbf{ifail}}=8$

The orders vector for the ARIMA model is invalid. (Only occurs if ${\mathbf{nmr}}=14$.)
 ${\mathbf{ifail}}=9$

The initial values of the filtered series are indeterminate for the given models.
 ${\mathbf{ifail}}=99$
An unexpected error has been triggered by this routine. Please
contact
NAG.
See
Section 3.9 in How to Use the NAG Library and its Documentation for further information.
 ${\mathbf{ifail}}=399$
Your licence key may have expired or may not have been installed correctly.
See
Section 3.8 in How to Use the NAG Library and its Documentation for further information.
 ${\mathbf{ifail}}=999$
Dynamic memory allocation failed.
See
Section 3.7 in How to Use the NAG Library and its Documentation for further information.
7
Accuracy
Accuracy and stability are high except when the MA parameters are close to the invertibility boundary.
8
Parallelism and Performance
g13baf is threaded by NAG for parallel execution in multithreaded implementations of the NAG Library.
g13baf makes calls to BLAS and/or LAPACK routines, which may be threaded within the vendor library used by this implementation. Consult the documentation for the vendor library for further information.
Please consult the
X06 Chapter Introduction for information on how to control and interrogate the OpenMP environment used within this routine. Please also consult the
Users' Note for your implementation for any additional implementationspecific information.
If an ARIMA model is supplied,
a local workspace array
of fixed
length is
allocated internally by
g13baf. The total size of
this array
amounts to
$K$ integer
elements, where
$K$ is the expression defined in the description of the argument
wa.
The time taken by
g13baf is approximately proportional to
with an appreciable fixed increase if an ARIMA model is supplied for the time series.
10
Example
This example reads a time series of length
$296$. It reads the univariate ARIMA
$\left(4,0,2,0,0,0,0\right)$ model and the ARIMA filtering
$\left(3,0,0,0,0,0,0\right)$ model for the series. Two initial backforecasts are required and these are calculated by a call to
g13ajf
.
The backforecasts are inserted at the start of the series and
g13baf is called to perform the calculations.
10.1
Program Text
Program Text (g13bafe.f90)
10.2
Program Data
Program Data (g13bafe.d)
10.3
Program Results
Program Results (g13bafe.r)