# Conditional mean and co-variance in $VAR(p)$ conditional on one lag only

Suppose I have a $$p$$'th order vector auto regression

$$\vec Z_t = F_1\vec Z_{t-1}+F_2\vec Z_{t-2} + \cdots +F_p \vec Z_{t - p} + \vec \epsilon_t,\qquad \vec\epsilon_t\sim N_q(\vec0,Q)$$

where $$Z_t\in\mathbb{R}^q$$. Then we can put this into the state-space form

\begin{aligned} \vec X_t &= (\vec Z_t,\vec Z_{t-1},\dots, \vec Z_{t - p + 1}) \\ \vec X_t &= F\vec X_{t-1} + R\vec\epsilon_t\\ F &= \begin{pmatrix} F_1 & \cdots & \cdots & F_{p-1} & F_p \\ I_q & 0 & \cdots & 0 & 0 \\ 0 & I_q & \ddots & \vdots & \vdots \\ \vdots & \ddots & \ddots &0 & 0 \\ 0 & \cdots & 0 & I_q & 0 \end{pmatrix} & R &= \begin{pmatrix} I_q \\ 0 \\ \vdots \\ 0 \end{pmatrix} \end{aligned}

It follows that the conditional mean is

\begin{aligned} E(\vec X_{t + p} \mid \vec X_t = \vec x)&= E(\vec E(\vec X_{t + p} \mid \vec X_{t + p - 1}) \mid \vec X_t = \vec x) \\ &= E(F\vec X_{t + p -1} \mid \vec X_t = \vec x) = \cdots = F^p\vec x \end{aligned}

and for the conditional covariance, we use that

\begin{aligned} \vec Z_{t+h} &\propto \sum_{i = 1}^h G(h-i)\vec \epsilon_{t + i} \\ G(k) &=\left\{\begin{matrix} I_q & k = 0 \\ \sum_{i = 1}^{\min (k, p)}F_iG(k-i) & k > 0 \end{matrix}\right. \end{aligned}

Thus,

\begin{aligned} \text{Var}(\vec Z_{t+h}\mid \vec X_t) &= \sum_{i=1}^h G(h-i)QG(h-i)^\top \\ \text{cov}(\vec Z_{t+h}, \vec Z_{t+l}\mid \vec X_t) &= \text{cov}\left( \sum_{i = 1}^h G(h-i)\vec \epsilon_{t + i}, \sum_{i = 1}^l G(l-i)\vec \epsilon_{t + i} \mid \vec X_t\right) \\ &= \text{cov}\left( \sum_{i = 1}^l G(h-i)\vec \epsilon_{t + i}, \sum_{i = 1}^l G(l-i)\vec \epsilon_{t + i} \mid \vec X_t\right) \\ &= \sum_{i = 1}^l G(h-i)QG(l-i)^\top \end{aligned}

where I assume that $$h > l$$. I hope the above is correct though this is not the question. My question is what are similar expressions if we only condition on $$Z_t$$ and not $$X_t= (\vec Z_t,\vec Z_{t-1}, \vec Z_{t - p + 1})$$? I.e., what are

$$E(\vec X_{t+p}\mid \vec Z_t),\qquad \text{Var}(\vec Z_{t+h}\mid \vec Z_t),\qquad \text{Cov}(\vec Z_{t+h},\vec Z_{t+l}\mid \vec Z_t),\qquad h > l > 0$$

# Update

Say

$$\vec X_0 \sim N(\vec \mu_0, Q_0)$$

which we either set to be the stationary distribution or a distribution selected for convenience. Then

$$\vec Z_k = \sum_{i = 1}^k G(k-i)\vec \epsilon_i + R^\top F^k\vec X_0$$

and thus,

\begin{align*} E(\vec Z_k) &= R^\top F^k\vec\mu_0 \\ \text{Var}(\vec Z_k) &= \sum_{i = 1}^k G(k-i)QG(k-i)^\top + R^\top F^kQ_0F^{k\top} R \\ \text{Cov}(\vec Z_k, \vec Z_l) &= \sum_{i = 1}^l G(k-i)QG(l-i)^\top + R^\top F^kQ_0F^{l\top} R, & k &> l > 0 \end{align*}

Using the above, we can compute the joint mean and covariance matrix and find that

\begin{aligned} \begin{pmatrix} \vec X_{t+p} \\ \vec Z_t \end{pmatrix} = \begin{pmatrix} \vec Z_{t + p} \\ \vec Z_{t + p - 1} \\ \vdots \\ \vec Z_t \end{pmatrix} & \sim N\left(\vec\mu, \Sigma \right) \\ \vec\mu &= (E(\vec X_{t+p})^\top, E(\vec Z_t)^\top)^\top\\ \Sigma &= \begin{pmatrix} \text{Var}(\vec X_{t+p}) & \text{Cov}(\vec X_{t+p}, \vec Z_t) \\ \text{Cov}(\vec Z_t,\vec X_{t+p}) & \text{Var}(\vec Z_t) \end{pmatrix} \end{aligned}

From which it is follows that

\begin{aligned} \vec X_{t+p} \mid \vec Z_t = \vec z &\sim N(k_\vec z, K_\vec z) \\ k_\vec z &= E(\vec X_{t+p}) + \text{Cov}(\vec X_{t+p}, \vec Z_t) \text{Var}(\vec Z_t)^{-1}(\vec z - E(\vec Z_t)) \\ K_\vec z &= \text{Var}(\vec X_{t+p}) - \text{Cov}(\vec X_{t+p}, \vec Z_t) \text{Var}(\vec Z_t)^{-1} \text{Cov}(\vec Z_t, \vec X_{t+p}) \end{aligned}

Is the above correct? Can I simplify the final expression further to something that is faster to compute or neater? I gather the latter is simple when we have stationary model since the unconditional means and covariances are independent of time. What about in the non-stationary case where we assume some distribution for $$\vec X_0$$?

Suppose that we are in the univariate case with $$E(z_t) = 0$$ and the process is stationary. Then

\begin{aligned} \vec X_{t+p} \mid Z_t = z & \sim N(\vec\beta z, K_{Z}) \\ \vec\beta &=\text{Cov}(\vec X_{t+p}, Z_t) \text{Var}(Z_t)^{-1} \\ K_\vec z &= \text{Var}(\vec X_{t+p}) - \text{Cov}(\vec X_{t+p},Z_t) \text{Var}(Z_t)^{-1} \text{Cov}(Z_t, \vec X_{t+p}) \end{aligned}

and this is easy to check. E.g., consider a univariate $$AR(3)$$ model

$$z_t = \phi_1z_{t-1} + \phi_2z_{t-2}+\phi_3z_{t-3} + \epsilon_t, \qquad\epsilon_t\sim N(0,\sigma^2)$$

then here is an example where we can confirm the result

# assign loadings
p <- 3L
F. <- matrix(0., p, p)
F.[1, ] <- c(0.1, 0.5, 0.25)                      # phi_1, phi_2, phi_3
F.[2:p, 1:(p - 1L)] <- diag(p - 1L)

# get time-invariant covariance matrix
sig <- .2                                         # sigma
get_Q_0 <- function(F., sig){
p <- nrow(F.)
Q <- matrix(c(sig^2, numeric(p * p - 1L)), p, p)
eg <- eigen(F.)
las <- eg$$values U <- eg$$vectors
T. <- solve(U, t(solve(U, Q)))
Z <- T./(1 - tcrossprod(las))
Re(tcrossprod(U %*% Z, U))
}
Q_0 <- get_Q_0(F., sig)

F_aug <- matrix(0., p + 1L, p + 1L)
F_aug[1:p, 1:p] <- F.
F_aug[p + 1L, p] <- 1
F_aug
#R      [,1] [,2] [,3] [,4]
#R [1,]  0.1  0.5 0.25    0
#R [2,]  1.0  0.0 0.00    0
#R [3,]  0.0  1.0 0.00    0
#R [4,]  0.0  0.0 1.00    0

# get Sigma
Sigma <- get_Q_0(F_aug, sig)

# coefficients
beta <- Sigma[1:p, p + 1L] * Sigma[p + 1L, p + 1L]^(-1)
beta
#R [1] 0.592 0.691 0.545

# compute conditional covariance matrix
K <- Sigma[1:p, 1:p] -
tcrossprod(Sigma[1:p, p + 1L],
Sigma[p + 1L, p + 1L]^(-1) *  Sigma[p + 1L, 1:p])
K
#R        [,1]   [,2]   [,3]
#R [1,] 0.0575 0.0121 0.0326
#R [2,] 0.0121 0.0462 0.0149
#R [3,] 0.0326 0.0149 0.0622

#####
# simulate Z and X and fit linear regression model
set.seed(73408077)
n_sim <- 100000L
tmp <- crossprod(chol(Q_0), matrix(rnorm(p * n_sim), p))
tmp <- rbind(F.[1,] %*% tmp + rnorm(n_sim, sd = sig), tmp)
X <- t(tmp)[, 1:p]
Z <- tmp[p + 1L, ]

fit <- lm(X ~ Z - 1)
coef(fit)
#R    [,1] [,2]  [,3]
#R Z 0.593 0.69 0.546
beta # should be
#R [1] 0.592 0.691 0.545

# estimated conditional covariance matrix
crossprod(residuals(fit)) / (n_sim - p)
#R        [,1]   [,2]   [,3]
#R [1,] 0.0573 0.0121 0.0326
#R [2,] 0.0121 0.0464 0.0149
#R [3,] 0.0326 0.0149 0.0620
K # should be
#R        [,1]   [,2]   [,3]
#R [1,] 0.0575 0.0121 0.0326
#R [2,] 0.0121 0.0462 0.0149
#R [3,] 0.0326 0.0149 0.0622

• Thanks for coming back and answering your own question for the benefit of all other users who might be facing the same problem. – Richard Hardy Mar 23 at 20:48