A good example of deriving a likelihood function is the normal distribution:

The PDF of the normal distribution is:

$$ f(x;\mu, \sigma) = \frac{1}{\sigma\sqrt{2\pi}}exp[\frac{(x - \mu)^2}{2\sigma^2}] $$

For a series of observations, if we can consider them i.i.d. we see that for observation $(x_1, ..., x_p)$ we can write this as:

$L(\mu, \sigma) = \prod_{i=1}^n f(x_i; \mu, \sigma)$

Which is our likelihood function.

and logging this gives

$LL(\mu, \sigma) = \sum_{i=1}^n ln(f(x_i; \mu, \sigma))$

We can further simplify this but it is unnecessary. This is the log-likelihood function that we can maximize by minimizing its negative and solving with BFGS or some other non-linear solver. We can use automatic differentiation to get it's gradients for $\mu$ and $\sigma$.

Given that I would like to solve this using numerical methods - the book I am reading (Hamilton) uses prediction error decomposition to show the likelihood function for AR(p). It seems like I could proceed in a straightforward manner and derive the log likelihood from the density of the observations:

$Y_T \tilde{} N(\mu_p, \sigma^2V_p)$

Where $V_p$ is the variance-covariance matrix of $Y_T$, $\sigma^2$ is the first entry of the variance-covariance matrix, and $\mu_p$ is the vector of means given by:

$\mu_p = \frac{c}{1 - \phi_1 - \phi_2 - ... - \phi_p}$

Why does every paper/book/etc use prediction error decomposition instead of writing the entire formula out like my example of the normal distribution?


Here is my derivation of the MLE for AR(p):

Since $Y_t$ is not i.i.d. (it's autoregressive) we can't use the same simplifying assumptions as above. Instead however, we can use the fact that we can decompose any joint distribution into a product of some initial marginal distribution and a series of conditional distributions.

Consider a general stationary AR(p) process:

$Y_t = \phi_1y_{t-1} + ... \phi_py_{t-p} + \epsilon_t$

$\epsilon_t \tilde{} N(0, \sigma^2)$

$\mathbf{\theta} = (c, \phi_1, ... \phi_p, \sigma^2)$

As mentioned above, we can then write our MLE like so:

$f(Y_T; \mathbf{\theta}) = f(y_1,...y_p; \mathbf{\theta})\prod_{i=p+1}^Tf(y_t|y_{t-1}, ..., y_{t-p}; \mathbf{\theta})$

Since we're collecting a vector of $p$ $y$'s rather than a scalar $1$ for an AR(1) process we now are dealing with a multivariate normal distribution for

$f(y_1,...y_p; \mathbf{\theta})$

Who has a variance-covariance matrix $\sigma^2V_p$ and mean vector $\mathbf{\mu_p}$ with each element represented as

$\frac{c}{1 - \phi_1 - \phi_2 - ... - \phi_p}$. Knowing this, we can simply plug in to a standard multivariate normal distribution and simplify.

$f(\mathbf{y}; \mathbf{\theta}) = (2\pi)^{-p/2}(\sigma^2)^{-p/2}det(V_p)^{-1/2}exp[\frac{(\mathbf{y}-\mathbf{\mu_p})'V_p^{-1}(\mathbf{y}-\mathbf{\mu_p})}{2\sigma^2}]$

where $\mathbf{y}$ is simply a vector $(y_1, ...,y_p)$ of our $p$ $y$'s.

So now we've solved the first part (the marginal distribution)

The last thing to do is solve for the conditional distribution for the next $t-p$ values of $Y_t$. Conditioning $Y_t$ on $y_{t-1}, ..., y_{t-p}$ we see that

$\mu = c+\phi_1y_{t-1}+\phi_2y_{t-2},...,+\phi_py_{t-p}$ $\sigma^2 = \sigma^2$

This is represented as a univariate normal distribution, so

$f_{y_t|y_{t-1},...y_{t-p}}(y_t|y_{t-1}, ..., y_{t-p};\mathbf{\theta}) = \frac{1}{\sigma\sqrt{2\pi}}exp(\frac{-(y_t - c -\phi_1y_{t-1}-\phi_2y_{t-2},...,-\phi_py_{t-p})^2}{2\sigma^2})$

and now we can pull these together nicely:

$f(Y_T; \mathbf{\theta}) = f(y_1,...y_p; \mathbf{\theta})\prod_{i=p+1}^Tf(y_t|y_{t-1}, ..., y_{t-p}; \mathbf{\theta}) = (2\pi)^{-p/2}(\sigma^2)^{-p/2}det(V_p)^{-1/2}exp[\frac{(\mathbf{y}-\mathbf{\mu_p})'V_p^{-1}(\mathbf{y}-\mathbf{\mu_p})}{2\sigma^2}]\prod_{t=p+1}^T\frac{1}{\sigma\sqrt{2\pi}}exp(\frac{-(y_t - c -\phi_1y_{t-1}-\phi_2y_{t-2},...,-\phi_py_{t-p})^2}{2\sigma^2})$

To get the log likelihood we can find it by logging both sides:

$LL(Y_T; \mathbf{\theta}) = \ln{(2\pi)^{-p/2}(\sigma^2)^{-p/2}det(V_p)^{-1/2}exp[\frac{(\mathbf{y}-\mathbf{\mu_p})'V_p^{-1}(\mathbf{y}-\mathbf{\mu_p})}{2\sigma^2}]} + \sum_{i=p+1}^T\ln{\frac{1}{\sigma\sqrt{2\pi}}exp(\frac{-(y_t - c -\phi_1y_{t-1}-\phi_2y_{t-2},...,-\phi_py_{t-p})^2}{2\sigma^2})}$

Which can be simplified:

$$ = -\frac{p}{2}\ln{2\pi} - \frac{p}{2}\ln{\sigma^2} - \frac{1}{2}\ln{det(V_p)} + \frac{(\mathbf{y}-\mathbf{\mu_p})'V_p^{-1}(\mathbf{y}-\mathbf{\mu_p})}{2\sigma^2} + \sum_{t=p+1}^T\ln{\frac{1}{\sqrt{2\pi}\sigma} - \frac{(y_t - c -\phi_1y_{t-1}-\phi_2y_{t-2},...,-\phi_py_{t-p})^2}{2\sigma^2})}$$

and then further simplified using the algebra around summations:

$$ = -\frac{p}{2}\ln{2\pi} - \frac{p}{2}\ln{\sigma^2} - \frac{1}{2}\ln{det(V_p)} + \frac{(\mathbf{y}-\mathbf{\mu_p})'V_p^{-1}(\mathbf{y}-\mathbf{\mu_p})}{2\sigma^2} - \frac{T-p}{2}\ln{2\pi} - \frac{T-p}{2}\ln{\sigma^2} - \sum_{t=p+1}^T\frac{(y_t - c -\phi_1y_{t-1}-\phi_2y_{t-2},...,-\phi_py_{t-p})^2}{2\sigma^2})$$

Having derived this - I see HOW it is derived, but I do not see where Prediction Error Decomposition comes in - because the PED version of this seems much, much simpler.

Any advice would be helpful - thank you!

  • 1
    $\begingroup$ Thank you for the elaborate derivation. I will try to think about it. $\endgroup$ – Richard Hardy Dec 31 '16 at 13:15
  • $\begingroup$ @richardhardy no problem, looking forward to your ideas. $\endgroup$ – user124589 Jan 1 '17 at 0:11
  • 1
    $\begingroup$ @RichardHardy it seems exceptionally difficult to find good literature regarding the PED version for the exact likelihood. Most sources I am seeing are essentially saying that because inverting a TxT covariance matrix is computationally intractable for a long enough time series, MLE on the conditional likelihood reduces to OLS on the parameters, so people seem to do that. However - that doesn't answer the question. let me know if you figure anything out on the PED version, I'm still digging. $\endgroup$ – user124589 Jan 2 '17 at 21:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy