library(caracas)Consider this model: \[ x_i = a x_{i-1} + e_i, \quad i=1, \dots, 3 \] and \(x_0=e_0\). All terms \(e_0, \dots, e_3\) are independent and \(N(0,v^2)\) distributed. Let \(e=(e_0, \dots, e_3)\) and \(x=(x_0, \dots x_3)\). Hence \(e \sim N(0, v^2 I)\). Isolating error terms gives \[ e= \left[\begin{matrix}e_{0}\\e_{1}\\e_{2}\\e_{3}\end{matrix}\right] = \left[\begin{matrix}1 & 0 & 0 & 0\\- a & 1 & 0 & 0\\0 & - a & 1 & 0\\0 & 0 & - a & 1\end{matrix}\right] \left[\begin{matrix}x_{0}\\x_{1}\\x_{2}\\x_{3}\end{matrix}\right] = L_1 x \]
N <- 3
L1 <- diag(4)
L1[cbind(1 + (1:N), 1:N)] <- "-a"
L1 <- as_sym(L1)Since \(\mathbf{Var}(e)=v^2 I\) we have \(\mathbf{Var}(e)=v^2 I=L \mathbf{Var}(x) L'\) so the covariance matrix of \(x\) is \(V_1=\mathbf{Var}(x) = v^2 L^- (L^-)'\) while the concentration matrix (the inverse covariances matrix) is \(K=v^{-2}L' L\).
def_sym(v2)
L1inv <- inv(L1)
V1 <- v2 * L1inv %*% t(L1inv)
K1 <- (t(L1) %*% L1) / v2\[\begin{align} K_1 &= \left[\begin{matrix}\frac{a^{2} + 1}{v_{2}} & - \frac{a}{v_{2}} & 0 & 0\\- \frac{a}{v_{2}} & \frac{a^{2} + 1}{v_{2}} & - \frac{a}{v_{2}} & 0\\0 & - \frac{a}{v_{2}} & \frac{a^{2} + 1}{v_{2}} & - \frac{a}{v_{2}}\\0 & 0 & - \frac{a}{v_{2}} & \frac{1}{v_{2}}\end{matrix}\right] \\ V_1 &= \left[\begin{matrix}v_{2} & a v_{2} & a^{2} v_{2} & a^{3} v_{2}\\a v_{2} & v_{2} \left(a^{2} + 1\right) & v_{2} \left(a^{3} + a\right) & v_{2} \left(a^{4} + a^{2}\right)\\a^{2} v_{2} & v_{2} \left(a^{3} + a\right) & v_{2} \left(a^{4} + a^{2} + 1\right) & v_{2} \left(a^{5} + a^{3} + a\right)\\a^{3} v_{2} & v_{2} \left(a^{4} + a^{2}\right) & v_{2} \left(a^{5} + a^{3} + a\right) & v_{2} \left(a^{6} + a^{4} + a^{2} + 1\right)\end{matrix}\right] \end{align}\]
Augment the \(AR(1)\) process above with \(y_i=b x_i + u_i\) for \(i=1,2,3\). Suppose \(u_i\sim N(0, w^2)\) and all \(u_i\) are independent and inpendent of \(e\). Then \((e,u)\) can be expressed in terms of \((x,y)\) as \[ (e,u) = \left[\begin{matrix}e_{0}\\e_{1}\\e_{2}\\e_{3}\\u_{1}\\u_{2}\\u_{3}\end{matrix}\right] = \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0\\- a & 1 & 0 & 0 & 0 & 0 & 0\\0 & - a & 1 & 0 & 0 & 0 & 0\\0 & 0 & - a & 1 & 0 & 0 & 0\\0 & - b & 0 & 0 & 1 & 0 & 0\\0 & 0 & - b & 0 & 0 & 1 & 0\\0 & 0 & 0 & - b & 0 & 0 & 1\end{matrix}\right] \left[\begin{matrix}x_{0}\\x_{1}\\x_{2}\\x_{3}\\y_{1}\\y_{2}\\y_{3}\end{matrix}\right] = L_2 (x,y) \] where
N <- 3
L2 <- diag("1", 1 + 2*N)
L2[cbind(1 + (1:N), 1:N)] <- "-a"
L2[cbind(1 + N + (1:N), 1 + 1:N)] <- "-b"
L2 <- as_sym(L2)Veu <- diag(1, 7)
diag(Veu)[1:4] <- "v2"
diag(Veu)[5:7] <- "w2"
Veu
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7]
#> [1,] "v2" "0" "0" "0" "0" "0" "0"
#> [2,] "0" "v2" "0" "0" "0" "0" "0"
#> [3,] "0" "0" "v2" "0" "0" "0" "0"
#> [4,] "0" "0" "0" "v2" "0" "0" "0"
#> [5,] "0" "0" "0" "0" "w2" "0" "0"
#> [6,] "0" "0" "0" "0" "0" "w2" "0"
#> [7,] "0" "0" "0" "0" "0" "0" "w2"
Veu <- as_sym(Veu)
Veu
#> [caracas]: ⎡v₂ 0 0 0 0 0 0 ⎤
#> ⎢ ⎥
#> ⎢0 v₂ 0 0 0 0 0 ⎥
#> ⎢ ⎥
#> ⎢0 0 v₂ 0 0 0 0 ⎥
#> ⎢ ⎥
#> ⎢0 0 0 v₂ 0 0 0 ⎥
#> ⎢ ⎥
#> ⎢0 0 0 0 w₂ 0 0 ⎥
#> ⎢ ⎥
#> ⎢0 0 0 0 0 w₂ 0 ⎥
#> ⎢ ⎥
#> ⎣0 0 0 0 0 0 w₂⎦
L2inv <- inv(L2)
V2 <- L2inv %*% Veu %*% t(L2inv)
K2 <- t(L2) %*% inv(Veu) %*% L2\[\begin{align} K_2 &= \left[\begin{matrix}\frac{a^{2}}{v_{2}} + \frac{1}{v_{2}} & - \frac{a}{v_{2}} & 0 & 0 & 0 & 0 & 0\\- \frac{a}{v_{2}} & \frac{a^{2}}{v_{2}} + \frac{b^{2}}{w_{2}} + \frac{1}{v_{2}} & - \frac{a}{v_{2}} & 0 & - \frac{b}{w_{2}} & 0 & 0\\0 & - \frac{a}{v_{2}} & \frac{a^{2}}{v_{2}} + \frac{b^{2}}{w_{2}} + \frac{1}{v_{2}} & - \frac{a}{v_{2}} & 0 & - \frac{b}{w_{2}} & 0\\0 & 0 & - \frac{a}{v_{2}} & \frac{b^{2}}{w_{2}} + \frac{1}{v_{2}} & 0 & 0 & - \frac{b}{w_{2}}\\0 & - \frac{b}{w_{2}} & 0 & 0 & \frac{1}{w_{2}} & 0 & 0\\0 & 0 & - \frac{b}{w_{2}} & 0 & 0 & \frac{1}{w_{2}} & 0\\0 & 0 & 0 & - \frac{b}{w_{2}} & 0 & 0 & \frac{1}{w_{2}}\end{matrix}\right] \\ V_2 &= \left[\begin{matrix}v_{2} & a v_{2} & a^{2} v_{2} & a^{3} v_{2} & a b v_{2} & a^{2} b v_{2} & a^{3} b v_{2}\\a v_{2} & a^{2} v_{2} + v_{2} & a^{3} v_{2} + a v_{2} & a^{4} v_{2} + a^{2} v_{2} & a^{2} b v_{2} + b v_{2} & a^{3} b v_{2} + a b v_{2} & a^{4} b v_{2} + a^{2} b v_{2}\\a^{2} v_{2} & a^{3} v_{2} + a v_{2} & a^{4} v_{2} + a^{2} v_{2} + v_{2} & a^{5} v_{2} + a^{3} v_{2} + a v_{2} & a^{3} b v_{2} + a b v_{2} & a^{4} b v_{2} + a^{2} b v_{2} + b v_{2} & a^{5} b v_{2} + a^{3} b v_{2} + a b v_{2}\\a^{3} v_{2} & a^{4} v_{2} + a^{2} v_{2} & a^{5} v_{2} + a^{3} v_{2} + a v_{2} & a^{6} v_{2} + a^{4} v_{2} + a^{2} v_{2} + v_{2} & a^{4} b v_{2} + a^{2} b v_{2} & a^{5} b v_{2} + a^{3} b v_{2} + a b v_{2} & a^{6} b v_{2} + a^{4} b v_{2} + a^{2} b v_{2} + b v_{2}\\a b v_{2} & a^{2} b v_{2} + b v_{2} & a^{3} b v_{2} + a b v_{2} & a^{4} b v_{2} + a^{2} b v_{2} & a^{2} b^{2} v_{2} + b^{2} v_{2} + w_{2} & a^{3} b^{2} v_{2} + a b^{2} v_{2} & a^{4} b^{2} v_{2} + a^{2} b^{2} v_{2}\\a^{2} b v_{2} & a^{3} b v_{2} + a b v_{2} & a^{4} b v_{2} + a^{2} b v_{2} + b v_{2} & a^{5} b v_{2} + a^{3} b v_{2} + a b v_{2} & a^{3} b^{2} v_{2} + a b^{2} v_{2} & a^{4} b^{2} v_{2} + a^{2} b^{2} v_{2} + b^{2} v_{2} + w_{2} & a^{5} b^{2} v_{2} + a^{3} b^{2} v_{2} + a b^{2} v_{2}\\a^{3} b v_{2} & a^{4} b v_{2} + a^{2} b v_{2} & a^{5} b v_{2} + a^{3} b v_{2} + a b v_{2} & a^{6} b v_{2} + a^{4} b v_{2} + a^{2} b v_{2} + b v_{2} & a^{4} b^{2} v_{2} + a^{2} b^{2} v_{2} & a^{5} b^{2} v_{2} + a^{3} b^{2} v_{2} + a b^{2} v_{2} & a^{6} b^{2} v_{2} + a^{4} b^{2} v_{2} + a^{2} b^{2} v_{2} + b^{2} v_{2} + w_{2}\end{matrix}\right] \end{align}\]