g13 Chapter Contents
g13 Chapter Introduction
NAG C Library Manual

# NAG Library Function Documentnag_tsa_varma_diagnostic (g13dsc)

## 1  Purpose

nag_tsa_varma_diagnostic (g13dsc) is a diagnostic checking function suitable for use after fitting a vector ARMA model to a multivariate time series using nag_tsa_varma_estimate (g13ddc). The residual cross-correlation matrices are returned along with an estimate of their asymptotic standard errors and correlations. Also, nag_tsa_varma_diagnostic (g13dsc) calculates the modified Li–McLeod portmanteau statistic and its significance level for testing model adequacy.

## 2  Specification

 #include #include
 void nag_tsa_varma_diagnostic (Integer k, Integer n, const double v[], Integer kmax, Integer ip, Integer iq, Integer m, const double par[], const Nag_Boolean parhld[], double qq[], Integer ishow, const char *outfile, double r0[], double r[], double rcm[], Integer pdrcm, double *chi, Integer *idf, double *siglev, NagError *fail)

## 3  Description

Let ${W}_{\mathit{t}}={\left({w}_{1\mathit{t}},{w}_{2\mathit{t}},\dots ,{w}_{\mathit{k}\mathit{t}}\right)}^{\mathrm{T}}$, for $\mathit{t}=1,2,\dots ,n$, denote a vector of $k$ time series which is assumed to follow a multivariate ARMA model of the form
 $Wt-μ= ϕ1Wt-1-μ+ϕ2Wt-2-μ+⋯+ϕpWt-p-μ +εt-θ1εt-1-θ2εt-2-⋯-θqεt-q,$ (1)
where ${\epsilon }_{\mathit{t}}={\left({\epsilon }_{1\mathit{t}},{\epsilon }_{2\mathit{t}},\dots ,{\epsilon }_{k\mathit{t}}\right)}^{\mathrm{T}}$, for $\mathit{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 $\Sigma$. The components of ${\epsilon }_{t}$ are assumed to be uncorrelated at non-simultaneous lags. The ${\varphi }_{i}$ and ${\theta }_{j}$ are $k$ by $k$ matrices of parameters. $\left\{{\varphi }_{\mathit{i}}\right\}$, for $\mathit{i}=1,2,\dots ,p$, are called the autoregressive (AR) parameter matrices, and $\left\{{\theta }_{\mathit{i}}\right\}$, for $\mathit{i}=1,2,\dots ,q$, the moving average (MA) parameter matrices. The parameters in the model are thus the $p$ ($k$ by $k$) $\varphi$-matrices, the $q$ ($k$ by $k$) $\theta$-matrices, the mean vector $\mu$ and the residual error covariance matrix $\Sigma$. Let
 $Aϕ= ϕ1 I 0 . . . 0 ϕ2 0 I 0 . . 0 . . . . . . ϕp-1 0 . . . 0 I ϕp 0 . . . 0 0 pk×pk and Bθ= θ1 I 0 . . . 0 θ2 0 I 0 . . 0 . . . . . . θq-1 0 . . . I θq 0 . . . . 0 qk×qk$
where $I$ denotes the $k$ by $k$ identity matrix.
The ARMA model (1) is said to be stationary if the eigenvalues of $A\left(\varphi \right)$ lie inside the unit circle, and invertible if the eigenvalues of $B\left(\theta \right)$ lie inside the unit circle. The ARMA model is assumed to be both stationary and invertible. Note that some of the elements of the $\varphi$- and/or $\theta$-matrices may have been fixed at pre-specified values (for example by calling nag_tsa_varma_estimate (g13ddc)).
The estimated residual cross-correlation matrix at lag $l$ is defined to the $k$ by $k$ matrix ${\stackrel{^}{R}}_{l}$ whose $\left(i,j\right)$th element is computed as
 $r ^ i j l = ∑ t = l + 1 n ε ^ i t - l - ε - i ε ^ j t - ε - j ∑ t = 1 n ε ^ i t - ε - i 2 ∑ t = 1 n ε ^ j t - ε - j 2 , l =0,1,…,i​ and ​j=1,2,…,k ,$
where ${\stackrel{^}{\epsilon }}_{it}$ denotes an estimate of the $t$th residual for the $i$th series ${\epsilon }_{it}$ and ${\stackrel{-}{\epsilon }}_{i}=\sum _{t=1}^{n}{\stackrel{^}{\epsilon }}_{it}/n$. (Note that ${\stackrel{^}{R}}_{l}$ is an estimate of $E\left({\epsilon }_{t-l}{\epsilon }_{t}^{\mathrm{T}}\right)$, where $E$ is the expected value.)
A modified portmanteau statistic, ${Q}_{\left(m\right)}^{*}$, is calculated from the formula (see Li and McLeod (1981))
 $Qm * = k2 mm+1 2n + n ∑ l=1 m r^ lT R^ 0 -1 ⊗ R^ 0 -1 r^ l ,$
where $\otimes$ denotes Kronecker product, ${\stackrel{^}{R}}_{0}$ is the estimated residual cross-correlation matrix at lag zero and $\stackrel{^}{r}\left(l\right)=\mathrm{vec}\left({\stackrel{^}{R}}_{l}^{\mathrm{T}}\right)$, where $\text{vec}$ of a $k$ by $k$ matrix is a vector with the $\left(i,j\right)$th element in position $\left(i-1\right)k+j$. $m$ denotes the number of residual cross-correlation matrices computed. (Advice on the choice of $m$ is given in Section 8.2.) Let ${l}_{C}$ denote the total number of ‘free’ parameters in the ARMA model excluding the mean, $\mu$, and the residual error covariance matrix $\Sigma$. Then, under the hypothesis of model adequacy, ${Q}_{\left(m\right)}^{*}$, has an asymptotic ${\chi }^{2}$-distribution on $m{k}^{2}-{l}_{C}$ degrees of freedom.
Let $\underline{\stackrel{^}{r}}=\left(\mathrm{vec}\left({R}_{1}^{\mathrm{T}}\right),\mathrm{vec}\left({R}_{2}^{\mathrm{T}}\right),\dots ,\mathrm{vec}\left({R}_{m}^{\mathrm{T}}\right)\right)$ then the covariance matrix of $\underline{\stackrel{^}{r}}$ is given by
 $Varr̲^=Y-XXTGGTX-1XT/n,$
where $Y={I}_{m}\otimes \left(\Delta \otimes \Delta \right)$ and $G={I}_{m}\left(G{G}^{\mathrm{T}}\right)$. $\Delta$ is the dispersion matrix $\Sigma$ in correlation form and $G$ a nonsingular $k$ by $k$ matrix such that $G{G}^{\mathrm{T}}={\Delta }^{-1}$ and $G\Delta {G}^{\mathrm{T}}={I}_{k}$. The construction of the matrix $X$ is discussed in Li and McLeod (1981). (Note that the mean, $\mu$, plays no part in calculating $\mathrm{Var}\left(\stackrel{^}{\underline{r}}\right)$ and therefore is not required as input to nag_tsa_varma_diagnostic (g13dsc).)

## 4  References

Li W K and McLeod A I (1981) Distribution of the residual autocorrelations in multivariate ARMA time series models J. Roy. Statist. Soc. Ser. B 43 231–239

## 5  Arguments

The output quantities k, n, v, kmax, ip, iq, par, parhld and qq from nag_tsa_varma_estimate (g13ddc) are suitable for input to nag_tsa_varma_diagnostic (g13dsc).
1:     kIntegerInput
On entry: $k$, the number of residual time series.
Constraint: ${\mathbf{k}}\ge 1$.
2:     nIntegerInput
On entry: $n$, the number of observations in each residual series.
3:     v[${\mathbf{kmax}}×{\mathbf{n}}$]const doubleInput
On entry: ${\mathbf{v}}\left[{\mathbf{kmax}}×\left(\mathit{t}-1\right)+\mathit{i}-1\right]$ must contain an estimate of the $\mathit{i}$th component of ${\epsilon }_{\mathit{t}}$, for $\mathit{i}=1,2,\dots ,k$ and $\mathit{t}=1,2,\dots ,n$.
Constraints:
• no two rows of ${\mathbf{v}}$ may be identical;
• in each row there must be at least two distinct elements.
4:     kmaxIntegerInput
On entry: the first dimension of the arrays v, qq and r0
Constraint: ${\mathbf{kmax}}\ge {\mathbf{k}}$.
5:     ipIntegerInput
On entry: $p$, the number of AR parameter matrices.
Constraint: ${\mathbf{ip}}\ge 0$.
6:     iqIntegerInput
On entry: $q$, the number of MA parameter matrices.
Constraint: ${\mathbf{iq}}\ge 0$.
Note: ${\mathbf{ip}}={\mathbf{iq}}=0$ is not permitted.
7:     mIntegerInput
On entry: the value of $m$, the number of residual cross-correlation matrices to be computed. See Section 8.2 for advice on the choice of m.
Constraint: ${\mathbf{ip}}+{\mathbf{iq}}<{\mathbf{m}}<{\mathbf{n}}$.
8:     par[$\left({\mathbf{ip}}+{\mathbf{iq}}\right)×{\mathbf{k}}×{\mathbf{k}}$]const doubleInput
On entry: the parameter estimates read in row by row in the order ${\varphi }_{1},{\varphi }_{2},\dots ,{\varphi }_{p}$, ${\theta }_{1},{\theta }_{2},\dots ,{\theta }_{q}$.
Thus,
• if ${\mathbf{ip}}>0$, ${\mathbf{par}}\left[\left(\mathit{l}-1\right)×k×k+\left(\mathit{i}-1\right)×k+j-1\right]$ must be set equal to an estimate of the $\left(\mathit{i},j\right)$th element of ${\varphi }_{\mathit{l}}$, for $\mathit{l}=1,2,\dots ,p$ and $\mathit{i}=1,2,\dots ,k$;
• if ${\mathbf{iq}}\ge 0$, ${\mathbf{par}}\left[p×k×k+\left(\mathit{l}-1\right)×k×k+\left(\mathit{i}-1\right)×k+j-1\right]$ must be set equal to an estimate of the $\left(\mathit{i},j\right)$th element of ${\theta }_{\mathit{l}}$, for $\mathit{l}=1,2,\dots ,q$ and $\mathit{i}=1,2,\dots ,k$.
The first $p×k×k$ elements of par must satisfy the stationarity condition and the next $q×k×k$ elements of par must satisfy the invertibility condition.
9:     parhld[$\left({\mathbf{ip}}+{\mathbf{iq}}\right)×{\mathbf{k}}×{\mathbf{k}}$]const Nag_BooleanInput
On entry: ${\mathbf{parhld}}\left[\mathit{i}-1\right]$ must be set to Nag_TRUE if ${\mathbf{par}}\left[\mathit{i}-1\right]$ has been held constant at a pre-specified value and Nag_FALSE if ${\mathbf{par}}\left[\mathit{i}-1\right]$ is a free parameter, for $\mathit{i}=1,2,\dots ,\left(p+q\right)×k×k$.
10:   qq[${\mathbf{kmax}}×{\mathbf{k}}$]doubleInput/Output
On entry: ${\mathbf{qq}}\left[{\mathbf{kmax}}×\left(j-1\right)+i-1\right]$ is an efficient estimate of the $\left(i,j\right)$th element of $\Sigma$. The lower triangle only is needed.
Constraint: ${\mathbf{qq}}$ must be positive definite.
On exit: if NE_G13D_AR, NE_G13D_DIAG, NE_G13D_FACT, NE_G13D_ITER, NE_G13D_MA, NE_G13D_RES, NE_G13D_ZERO_VAR or NE_NOT_POS_DEF, then the upper triangle is set equal to the lower triangle.
11:   ishowIntegerInput
On entry: must be nonzero if the residual cross-correlation matrices $\left\{{\stackrel{^}{r}}_{ij}\left(l\right)\right\}$ and their standard errors $\left\{\mathrm{se}\left({\stackrel{^}{r}}_{ij}\left(l\right)\right)\right\}$, the modified portmanteau statistic with its significance and a summary table are to be printed. The summary table indicates which elements of the residual correlation matrices are significant at the $5%$ level in either a positive or negative direction; i.e., if ${\stackrel{^}{r}}_{ij}\left(l\right)>1.96×\mathrm{se}\left({\stackrel{^}{r}}_{ij}\left(l\right)\right)$ then a ‘$+$’ is printed, if ${\stackrel{^}{r}}_{ij}\left(l\right)<-1.96×\mathrm{se}\left({\stackrel{^}{r}}_{ij}\left(l\right)\right)$ then a ‘$-$’ is printed, otherwise a fullstop (.) is printed. The summary table is only printed if $k\le 6$ on entry.
The residual cross-correlation matrices, their standard errors and the modified portmanteau statistic with its significance are available also as output variables in r, rcm, chi, idf and siglev.
12:   outfileconst char *Input
On entry: the name of a file to which diagnostic output will be directed. If outfile is NULL the diagnostic output will be directed to standard output.
13:   r0[${\mathbf{kmax}}×{\mathbf{k}}$]doubleOutput
On exit: if $i\ne j$, then ${\mathbf{r0}}\left[{\mathbf{kmax}}×\left(j-1\right)+i-1\right]$ contains an estimate of the $\left(i,j\right)$th element of the residual cross-correlation matrix at lag zero, ${\stackrel{^}{R}}_{0}$. When $i=j$, ${\mathbf{r0}}\left[{\mathbf{kmax}}×\left(j-1\right)+i-1\right]$ contains the standard deviation of the $i$th residual series. If NE_G13D_RES or NE_G13D_ZERO_VAR on exit then the first k rows and columns of r0 are set to zero.
14:   r[${\mathbf{kmax}}×{\mathbf{kmax}}×{\mathbf{m}}$]doubleOutput
Note: the element ${\mathbf{R}}\left(l,i,j\right)$ is stored in the array element ${\mathbf{r}}\left[\left(j-1\right)×{\mathbf{kmax}}×{\mathbf{kmax}}+\left(i-1\right)×{\mathbf{kmax}}+l-1\right]$.
On exit: ${\mathbf{R}}\left(\mathit{l},\mathit{i},\mathit{j}\right)$ is an estimate of the $\left(\mathit{i},\mathit{j}\right)$th element of the residual cross-correlation matrix at lag $\mathit{l}$, for $\mathit{i}=1,2,\dots ,k$, $\mathit{j}=1,2,\dots ,k$ and $\mathit{l}=1,2,\dots ,m$. If NE_G13D_RES or NE_G13D_ZERO_VAR on exit then all elements of r are set to zero.
15:   rcm[${\mathbf{pdrcm}}×{\mathbf{m}}×{\mathbf{k}}×{\mathbf{k}}$]doubleOutput
On exit: the estimated standard errors and correlations of the elements in the array r. The correlation between ${\mathbf{R}}\left(l,i,j\right)$ and ${\mathbf{R}}\left({l}_{2},{i}_{2},{j}_{2}\right)$ is returned as ${\mathbf{rcm}}\left[{\mathbf{pdrcm}}×t+s\right]$ where $s=\left(l-1\right)×k×k+\left(j-1\right)×k+i$ and $t=\left({l}_{2}-1\right)×k×k+\left({j}_{2}-1\right)×k+{i}_{2}$ except that if $s=t$, then ${\mathbf{rcm}}\left[{\mathbf{pdrcm}}×t+s\right]$ contains the standard error of ${\mathbf{R}}\left(l,i,j\right)$. If on exit, NE_G13D_DIAG or NE_G13D_FACT, then all off-diagonal elements of RCM are set to zero and all diagonal elements are set to $1/\sqrt{n}$.
16:   pdrcmIntegerInput
On entry: the first dimension of the array rcm.
Constraint: ${\mathbf{pdrcm}}\ge {\mathbf{m}}×{\mathbf{k}}×{\mathbf{k}}$.
17:   chidouble *Output
On exit: the value of the modified portmanteau statistic, ${Q}_{\left(m\right)}^{*}$. If NE_G13D_RES or NE_G13D_ZERO_VAR on exit then chi is returned as zero.
18:   idfInteger *Output
On exit: the number of degrees of freedom of chi.
19:   siglevdouble *Output
On exit: the significance level of chi based on idf degrees of freedom. If NE_G13D_RES or NE_G13D_ZERO_VAR on exit, siglev is returned as one.
20:   failNagError *Input/Output
The NAG error argument (see Section 3.6 in the Essential Introduction).

## 6  Error Indicators and Warnings

NE_ALLOC_FAIL
Dynamic memory allocation failed.
On entry, argument $〈\mathit{\text{value}}〉$ had an illegal value.
NE_G13D_AR
On entry, the AR parameter estimates are outside the stationarity region.
NE_G13D_ARMA
On entry, ${\mathbf{ip}}=0$ and ${\mathbf{iq}}=0$.
NE_G13D_DIAG
The matrix rcm could not be computed because one of its diagonal elements was found to be non-positive.
NE_G13D_FACT
On entry, the AR operator has a factor in common with the MA operator.
NE_G13D_ITER
Excessive iterations needed to find zeros of determinental polynomials.
NE_G13D_MA
On entry, the MA parameter matrices are outside the invertibility region.
NE_G13D_RES
On entry, at least two of the residual series are identical.
NE_G13D_ZERO_VAR
On entry, at least one of the residual series in the array v has near-zero variance.
NE_INT
On entry, ${\mathbf{ip}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{ip}}\ge 0$.
On entry, ${\mathbf{iq}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{iq}}\ge 0$.
On entry, ${\mathbf{k}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{k}}\ge 1$.
NE_INT_2
On entry, ${\mathbf{kmax}}=〈\mathit{\text{value}}〉$ and ${\mathbf{k}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{kmax}}\ge {\mathbf{k}}$.
On entry, ${\mathbf{m}}=〈\mathit{\text{value}}〉$ and ${\mathbf{n}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{m}}<{\mathbf{n}}$.
NE_INT_3
On entry, ${\mathbf{m}}=〈\mathit{\text{value}}〉$, ${\mathbf{ip}}=〈\mathit{\text{value}}〉$ and ${\mathbf{iq}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{m}}>{\mathbf{ip}}+{\mathbf{iq}}$.
On entry, ${\mathbf{pdrcm}}=〈\mathit{\text{value}}〉$, ${\mathbf{m}}=〈\mathit{\text{value}}〉$ and ${\mathbf{k}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{pdrcm}}\ge {\mathbf{m}}×{\mathbf{k}}×{\mathbf{k}}$.
NE_INTERNAL_ERROR
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact NAG for assistance.
NE_NOT_CLOSE_FILE
Cannot close file $〈\mathit{\text{value}}〉$.
NE_NOT_POS_DEF
On entry, the covariance matrix qq is not positive definite.
NE_NOT_WRITE_FILE
Cannot open file $〈\mathit{\text{value}}〉$ for writing.

## 7  Accuracy

The computations are believed to be stable.

### 8.1  Timing

The time taken by nag_tsa_varma_diagnostic (g13dsc) depends upon the number of residual cross-correlation matrices to be computed, $m$, and the number of time series, $k$.

### 8.2  Choice of $m$

The number of residual cross-correlation matrices to be computed, $m$, should be chosen to ensure that when the ARMA model (1) is written as either an infinite order autoregressive process, i.e.,
 $Wt-μ=∑j=1∞πjWt-j-μ+εt$
or as an infinite order moving average process, i.e.,
 $Wt-μ=∑j= 1∞ψjεt-j+εt$
then the two sequences of $k$ by $k$ matrices $\left\{{\pi }_{1},{\pi }_{2},\dots \right\}$ and $\left\{{\psi }_{1},{\psi }_{2},\dots \right\}$ are such that ${\pi }_{j}$ and ${\psi }_{j}$ are approximately zero for $j>m$. An overestimate of $m$ is therefore preferable to an under-estimate of $m$. In many instances the choice $m=10$ will suffice. In practice, to be on the safe side, you should try setting $m=20$.

### 8.3  Checking a ‘White Noise’ Model

If you have fitted the ‘white noise’ model
 $Wt-μ=εt$
then nag_tsa_varma_diagnostic (g13dsc) should be entered with $p=1$, $q=0$, and the first ${k}^{2}$ elements of par and parhld set to zero and Nag_TRUE respectively.

### 8.4  Approximate Standard Errors

When NE_G13D_DIAG or NE_G13D_FACT all the standard errors in rcm are set to $1/\sqrt{n}$. This is the asymptotic standard error of ${\stackrel{^}{r}}_{ij}\left(l\right)$ when all the autoregressive and moving average parameters are assumed to be known rather than estimated.

### 8.5  Alternative Tests

${\stackrel{^}{R}}_{0}$ is useful in testing for instantaneous causality. If you wish to carry out a likelihood ratio test then the covariance matrix at lag zero $\left({\stackrel{^}{C}}_{0}\right)$ can be used. It can be recovered from ${\stackrel{^}{R}}_{0}$ by setting
 $C^0i,j =R^0i,j×R^0i,i×R^0j,j, for ​i≠j =R^0i,j×R^0i,j, for ​i=j$

## 9  Example

This example fits a bivariate AR(1) model to two series each of length $48$. $\mu$ has been estimated but ${\varphi }_{1}\left(2,1\right)$ has been constrained to be zero. Ten residual cross-correlation matrices are to be computed.

### 9.1  Program Text

Program Text (g13dsce.c)

### 9.2  Program Data

Program Data (g13dsce.d)

### 9.3  Program Results

Program Results (g13dsce.r)