# Joint distribution of AR(1) model

Let

$$y_1 \sim \mathcal{N}(\phi_0,\,\sigma^2),$$ $$y_t|y_{t-1},\, \phi_0,\,\phi_1,\,\sigma^2 \sim \mathcal{N}(\phi_0+\phi_1(y_{t-1} - \phi_0),\,\sigma^2),$$ for $t=2,3,\cdots,T$.

I want to find the joint distribution of $y_1,\,y_2,\,\cdots,\,y_T$. I think that it is multivariate Normal with mean a vector of $\phi_0$ and precision matrix $Q= \frac{1}{\sigma^2} \begin{pmatrix} (1+\phi_1^2) & -\phi_1 & \cdots & 0 \\ -\phi_1 & (1+\phi_1^2) & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & 1 \end{pmatrix}.$

I would like to know if this is true and if there is a proof to read it.

Let us write the joint density as $$p(y_1,\ldots,y_T) = p(y_1)\,p(y_2\mid y_1) \, p(y_3 \mid y_2,y_1) \ldots \, p(y_T \mid y_{T-1},\ldots y_1).$$ Furthermore, since the process is AR(1), the past values influence future values only via the latest value, i.e., we have the Markov property $p(y_t \mid y_1,\ldots,y_{t-1}) = p(y_t \mid y_{t-1})$. Substituting this in the factorization, we get $$p(y_1,\ldots,y_T) = p(y_1)\, \prod_{i=2}^T p(y_i \mid y_{i-1}).$$ The marginal density of $y_1$ and the required conditional densities were given as assumptions. From now on, we shall ignore multiplicative constants (that are independent of $y$), since they are in the end be determined by the requirement that the joint density integrates to 1. $$\propto e^{-\frac{1}{2\sigma^2}(y_1 - \phi_0)^2} \times \prod_{i=2}^T e^{-\frac{1}{2\sigma^2}(y_i - \phi_0 - \phi_1\,(y_{i-1} - \phi_0))^2} = e^{-\frac{1}{2}\,E}$$ where $$E = \frac{1}{\sigma^2}(y_1 - \phi_0)^2 + \sum_{i=2}^T \frac{1}{\sigma^2}(y_i - \phi_0 - \phi_1(y_{i-1} - \phi_0))^2$$ $$= \frac{1}{\sigma^2}(y_1 - \phi_0)^2 + \sum_{i=2}^T \frac{1}{\sigma^2}\left((y_i - \phi_0)^2 - 2\,(y_i - \phi_0)\,\phi_1\,(y_{i-1} - \phi_0) + \phi_1^2 (y_{i-1} - \phi_0)^2 \right)$$ $$= \sum_{i=1}^{T-1}\frac{1}{\sigma^2}\,(1 + \phi_1^2)(y_i - \phi_0)^2 + \frac{1}{2\sigma^2} (y_T - \phi_0)^2 + \sum_{i=1}^{T-1}\frac{1}{\sigma^2}\,2\,(-\phi)\,(y_{i+1}-\phi_0)\,(y_i-\phi_0).$$ So the joint density is proportional to \begin{align} \mathrm{exp}\bigg(-\frac{1}{2}\,\sum_{i=1}^{T-1}\frac{1}{\sigma^2}\,(1 + \phi_1^2)(y_i - \phi_0)^2 + \frac{1}{2\sigma^2} (y_T - \phi_0)^2 \\+ \sum_{i=1}^{T-1}2\,\frac{1}{\sigma^2}\,(-\phi)\,(y_{i+1}-\phi_0)\,(y_i-\phi_0)\bigg), \end{align} Observe that this exponent is a quadratic form of a vector consisting of variables $(y_i - \phi_0)$. Thus we conclude that the joint density is a multivariate normal with means $E(y_i) = \phi_0$ and the precision matrix can be read from the previous expression, since we have $E = (y - \phi_0\mathbf{1})\,\Sigma^{-1}\,(y-\phi_0\,\mathbf{1})$. Namely,

• If $i=j$ and $i<T$, $\Sigma^{-1}_{ij} = (1 + \phi_1^2) / \sigma^2$
• If $i=j=T$, $\Sigma^{-1}_{ij} = 1 / \sigma^2$
• If $|i-j|=1$, $\Sigma^{-1}_{ij} = -\phi_1 / \sigma^2$
• If $|i-j|>1$, $\Sigma^{-1}_{ij} = 0$,

which is indeed the form of the precision matrix that was claimed in the question.

• So, what is the implied AR(1) model here? – Alecos Papadopoulos Jul 23 '15 at 11:40
• @AlecosPapadopoulos I'm not sure if I understand the comment, but $y_t= \phi_0 - \phi_1\,\phi_0 + \phi_1\,y_{t-1} + e_t,~e_t \sim N(0,\sigma^2)$ is the AR(1) model implied by the OP? – Juho Kokkala Jul 23 '15 at 11:44
• If this is the model, then, if also $e_t$ is uncorrelated with the past, then $$Var(y_2) = \phi_1^2Var(y_1) + \sigma^2 = (1+\phi_1^2)\sigma^2$$ and $$Var(y_3) = \phi_1^2Var(y_2) + \sigma^2 = [\phi_1^2(1+\phi_1^2)+1]\sigma^2$$. Etc. Does this conforms with the result you arrived at ? I haven't worked it through – Alecos Papadopoulos Jul 23 '15 at 11:59
• I agree those should be the variances, or there is an error. I haven't checked whether the obtained precision matrix corresponds with that (since my solution never works out the variance explicitly, but instead arrives directly at the precision). – Juho Kokkala Jul 23 '15 at 12:02
• @AlecosPapadopoulos At least numerically - picking random $\phi_1$ and $\sigma^2$ and inverting the precision matrix with $T=3$ I get diagonal elements of the covariance matrix that agree with the expressions in your comment. I suppose it is not too difficult to construct the entire covariance matrix and verify it is the inverse of the precision matrix stated here, but I haven't done that yet. – Juho Kokkala Jul 23 '15 at 12:09

The title of the question points towards a functional specification of the form

$$y_t= \phi_0 - \phi_1\,\phi_0 + \phi_1\,y_{t-1} + e_t,~e_t \sim N(0,\sigma^2), t>1, |\phi_1|<1$$

with the error term $e_t$ being i.i.d.

Given the assumption on the initial available observation (which does not necessarily represent the beginning of the process, just the first observation of the sample), we can determine that the process from then on is heteroskedastic. Specifically,

$${\rm Var}(y_2) = \phi_1^2{\rm Var}(y_1) + \sigma^2 = (1+\phi_1^2)\sigma^2$$

$${\rm Var}(y_3) = \phi_1^2{\rm Var}(y_2) + \sigma^2 = [\phi_1^2(1+\phi_1^2)+1]\sigma^2$$

$${\rm Var}(y_4) = \phi_1^2{\rm Var}(y_3) + \sigma^2 = [(\phi_1^2)^3+(\phi_1^2)^2+(\phi_1^2)+1]\sigma^2$$

The pattern is clear and, asymptotically, it leads to the familiar $\lim_{t \rightarrow \infty} {\rm Var}(y_t) = \sigma^2/(1-\phi^2)$. But only asymptotically. The joint distribution of the sample therefore will be a joint distribution of random variables with different and monotonically increasing (but bounded) variances of the marginal distributions. The expected value is common for all observations, and equal to $\phi_0$.

We now obtain the covariances for a sample of three observations $\{y_1,y_2, y_3\}$. We have

$${\rm Cov}(y_2,y_1) = E(y_2y_1) - \phi_0^2 = E\Big(\phi_0y_1 - \phi_1\phi_0y_1 + \phi_1y_1^2+ e_2y_1\Big) - \phi_0^2$$

$$=\phi_0^2 - \phi_1\phi_0^2 + \phi_1\big({\rm Var}(y_1) + \phi_0^2\big) - \phi_0^2 = - \phi_1\phi_0^2 + \phi_1\sigma^2 + \phi_1\phi_0^2$$

and so

$${\rm Cov}(y_2,y_1) = \phi_1\sigma^2,\;\;\; E(y_2y_1) = \phi_1\sigma^2 + \phi_0^2$$

Continuing,

$${\rm Cov}(y_3,y_1) = E(y_3y_1) - \phi_0^2 = E\Big(\phi_0y_1 - \phi_1\phi_0y_1 + \phi_1y_2y_1+ e_3y_1\Big) - \phi_0^2$$

$$= \phi_0^2 - \phi_1\phi_0^2 + \phi_1^2\sigma^2 + \phi_1\phi_0^2 - \phi_0^2$$

$$\implies {\rm Cov}(y_3,y_1) = \phi_1^2\sigma^2,\;\; E(y_3y_1) = \phi_1^2\sigma^2 + \phi_0^2$$

Finally,

$${\rm Cov}(y_3,y_2) = E(y_3y_2) - \phi_0^2 = E\Big(\phi_0y_3 - \phi_1\phi_0y_3 + \phi_1y_1y_3+ e_2y_3\Big) - \phi_0^2$$

$$=\phi_0^2 - \phi_1\phi_0^2 + \phi_1^3\sigma^2 + \phi_1\phi_0^2 +\phi_1\sigma^2-\phi_0^2$$

$$\implies {\rm Cov}(y_3,y_2) = \phi_1(1+\phi_1^2)\sigma^2$$

We observe that

$${\rm Cov}(y_3,y_2) \neq {\rm Cov}(y_2,y_1)$$

namely that the first-order autocovariance depends also on $t$. The covariance matrix of a sample of three observations is therefore

$${\rm Cov}(y_1,y_2,y_3)= \sigma^2 \begin{pmatrix} 1 & \phi_1 & \phi_1^2 \\ \phi_1 & (1+\phi_1^2) & \phi_1(1+\phi_1^2) \\ \phi_1^2 & \phi_1(1+\phi_1^2) & (1+\phi_1^2+\phi_1^4) \\ \end{pmatrix}$$

The inverse of this matrix (i.e. the precision matrix) (calculated online on this site) is given as

$$Q = {\rm Cov}^{-1}(y_1,y_2,y_3)= \frac {1}{\sigma^2} \begin{pmatrix} 1+\phi_1^2 & -\phi_1 & 0 \\ -\phi_1 & 1+\phi_1^2 & -\phi_1 \\ 0 & -\phi_1 & 1 \\ \end{pmatrix}$$

which is, the general result at @JuhoKokkala answer, for $T=3$.

• Thank you for your help, but I can't change the distribution of $y_t$. I am sure that the precision matrix I wrote is true, I don't know what is going on with the mean. – F.F. Jul 22 '15 at 17:24
• I am sorry, I did not understand that you mean the marginal distribution. I need to find it without any additional assumpions. – F.F. Jul 22 '15 at 17:51
• This answer assumes that the marginal variance of $y_t$ is $\sigma^2$ and then proceeds to give a covariance matrix with diagonal elements $\sigma^2 \, (1 - \phi_1^2)^{-1}$, which is a contradiction unless $\phi_1=0$. – Juho Kokkala Jul 22 '15 at 18:23
• Furthermore, I don't see how the joint distribution is not fully defined by giving marginal of $y_1$, conditionals $y_{t+1} \mid y_t$ and the Markov property. – Juho Kokkala Jul 22 '15 at 18:25
• There is a theorem stating if (y) is normal, and (x|y) is also normal then (x, y) are bivariate normal. So it is fine to assume this, and it is also true. – probabilityislogic Jul 23 '15 at 10:50