Maximum Likelihood in a time series multi-population model

Question

Currently I am trying to implement a mortality rate model (Lee-Li), but I am having some trouble with the maximum likelihood procedure of a time series model given by:

\begin{align} K_{t+1} &= K_t + \theta + \epsilon_{t+1}\\ \kappa_{t+1} &= a \kappa_{t} + \delta_{t+1} \end{align}

Where $E_t := (\delta_t,\epsilon_t)$ are assumed to be independent and come from a bi variate normal distribution with mean $(0,0)$ and co-variance matrix $\mathbf C$.

I have data available for both $K_t$ and $\kappa_t$, and the idea is to maximize for $a$, $\theta$ and $\mathbf C$ using maximum likelihood. However I have no clue as how to approach this problem. I tried obtaining the MLE for $\mathbf C$ using the MLE for the covariance matrix of a bivariate normal distribution, but this seems to be wrong. Furthermore I have no idea on how to approach $\theta$ and a for this matter either, since I don't know how to set up a (log)-likelihood for this expression.

Any help would be appreciated.

Edit: The co-variance matrix $\mathbf C = \begin{bmatrix} \sigma^{2}_{\epsilon} & \rho \sigma_{\epsilon} \sigma_{\delta}& \\ \rho \sigma_{\epsilon} \sigma_{\delta} & \sigma^{2}_{\delta} \end{bmatrix}$ Where the correlation coefficient $\rho$ is non-zero. I have tried to do the following steps to derive the maximum likelihood for $\theta$ :

1: Condition the distribution of $K_{t+1}$ on past data $t$ and the parameters
$\mathbf{\Theta}=[a,\theta,\sigma^{2}_{\epsilon},\sigma^{2}_{\delta},\rho \sigma_{\epsilon} \sigma_{\delta}]$, which has a bi variate normal distribution

2: compute the conditional log-likelihood of this equation: \begin{align} \sum_{t}\log(L \left(\mathbf{\Theta} \right)) &= \log\left(\frac{1}{2 \pi \sigma_{\epsilon} \sigma_{\delta}\sqrt{ 1 - \rho^{2}}}\right) - \frac{z}{2(1 - \rho^{2}} \end{align} where \begin{align} z = \frac{(K_{t+1} - (K_{t} + \theta))^{2}}{\sigma_\epsilon^2} + \frac{(\kappa_{t+1} - a\kappa_{t}) ^{2}}{\sigma_\delta^2} - \frac{2\rho(K_{t+1} - (K_{t} + \theta))(\kappa_{t+1} - a\kappa_{t})}{\sigma_\epsilon \sigma_\delta} \end{align}

3: take the derivative with respect to $\theta$. Since $z$ is the only expression involving $\theta$ we can ignore the rest: \begin{align} \frac{\partial \log \left( L \right)}{\partial \theta} &= 0 \Leftrightarrow \\ \sum_{t}\frac{-2(K_{t+1} - (K_t+\theta))}{\sigma_{\epsilon}^2} + \sum_{t}\frac{2\rho(\kappa_{t+1} - a\kappa_t)}{\sigma_{\epsilon}\sigma_{\delta}}&=0 \Leftrightarrow \\ \sum_{t}\frac{K_{t+1} - (K_t+\theta)}{\sigma_{\epsilon}^2} &= \sum_{t}\frac{\rho(\kappa_{t+1} - a\kappa_t)}{\sigma_{\epsilon}\sigma_{\delta}} \Leftrightarrow \\ \sum_{t}K_{t+1} - K_t - n\theta &= \sigma_{\epsilon}^2\sum_{t}\frac{\rho(\kappa_{t+1} - a\kappa_t)}{\sigma_{\epsilon}\sigma_{\delta}} \Leftrightarrow \\ \hat{\theta_{MLE}} =\sum_{t}\frac{(K_{t+1} - K_t)}{n} &-\sigma_{\epsilon}^2\sum_{t}\frac{\rho(\kappa_{t+1} - a\kappa_t)}{n\sigma_{\epsilon}\sigma_{\delta}} \end{align}

Similar calculations for $a,\sigma_{\epsilon},\sigma_{\delta}$ yield \begin{align} \hat{a_{MLE}} &= \frac{\sum_{t} \kappa_{t+1} \kappa_{t} - \frac{\sigma_{\delta}}{\sigma_{\epsilon}} \rho \sum_{t} \kappa_{t}(K_{t+1} - K_{t} - \theta)}{ \sum_{t} \kappa_{t}^{2}} \\ \hat{\sigma_{\epsilon}^{MLE}} &= \frac{-\frac{\rho*(K_{t+1} - K_{t} - \theta)(\kappa_{t+1} - a*\kappa_{t})}{ \sigma_{\delta}} + \sqrt D_{\epsilon} )}{ 2n(1-\rho^2)} \\ \hat{\sigma_{\delta}^{MLE}} &= \frac{-\frac{\rho*(K_{t+1} - K_{t} - \theta)(\kappa_{t+1} - a*\kappa_{t})}{ \sigma_{\epsilon}} + \sqrt D_{\delta} )}{ 2n(1-\rho^2)} \end{align}

Where \begin{align} D_{\epsilon} &= \left(\frac{\rho \sum_{t} (K_{t+1} - K_{t} - \theta)(\kappa_{t+1} - a*\kappa_{t})}{\sigma_{\epsilon}}\right)^{2} + 4 n(1-\rho^2) \sum_{t} (\kappa_{t+1} - a\kappa_{t})^2 \\ D_{\delta} &=\left(\frac{\rho \sum_{t} (K_{t+1} - K_{t} - \theta)(\kappa_{t+1} - a*\kappa_{t})}{\sigma_{\epsilon}}\right)^{2} + 4 n(1-\rho^2) \sum_{t} (K_{t+1} - K_{t} - \theta)^2 \\ \end{align}

Finally, since the expression is too complicated analytically, one can use numerical techniques to solve the ML equations. I have used the quasi Newton-Raphson algorithm proposed by Goodman to determine the parameters.

Answer 1

If you define $\mathbf{K}\equiv [K,\kappa]$ and $\boldsymbol{\mu}\equiv [K_t+\theta,a\kappa_t]$, then your model is

$$ \mathbf{K}_{t+1}\vert(\mathbf{K}_t,a,\theta,\mathbf{C})\sim \mathrm{N}_{\boldsymbol{\mu},\mathbf{C}}$$

i.e. the distribution of your data at $t+1$ conditioned on both the data at $t$, and on your parameters, is a bivariate normal with the specified non-zero conditional mean. The equation above is the likelihood.

If your covariance matrix were isotropic ($\mathbf{C}=\sigma^2\mathbf{I}$) this would be simple linear regression (i.e. with $\mathbf{X}=\mathbf{K}_t$ and $\mathbf{Y}=\mathbf{K}_{t+1}$). So the only change is your have a general covariance matrix, so three components (e.g. two variances and a correlation coefficient).

For MLE you want to determine the parameter values that make the gradient of the log-likelihood zero. Your parameters are $\mathbf{\Theta}=[a,\theta,C_{K,K},C_{\kappa,\kappa},C_{K,\kappa}]$, where I have expanded out the diagonal and off-diagonal components of the symmetric covariance matrix. For each parameter you get an equation by setting the gradient to zero, i.e. $\frac{\partial L}{\partial \Theta_k}=0$ for $k=1,\ldots,5$, where $L$ is the negative log likelihood.

Note that for the first two parameters $a,\theta$ you can get the gradient from the chain rule, i.e. $$\frac{\partial L}{\partial a}=\frac{\partial L}{\partial \mu_{\kappa}}\frac{\partial \mu_{\kappa}}{\partial a}$$ and $$\frac{\partial L}{\partial \theta}=\frac{\partial L}{\partial \mu_K}\frac{\partial \mu_K}{\partial \theta}$$ where $\frac{\partial \mu_{\kappa}}{\partial a}=\kappa_t$ and $\frac{\partial \mu_K}{\partial \theta}=1$.