Introduction

We report all the formulae used in the main computations that take place in the package as a reference for the users and developers.

Basic model for the HP filter with jumps

The basic model is based on the state-space representation $\begin{aligned} y_{t} & = α_{t}^{(1)} + ε_{t}, & ε_{t} \sim N I D (0, σ_{ε}^{2}) \\ α_{t + 1}^{(1)} & = α_{t}^{(1)} + α_{t}^{(2)} + η_{t}, & η_{t} \sim N I D (0, σ_{η, t}^{2}) \\ α_{t + 1}^{(2)} & = α_{t}^{(2)} + ζ_{t}, & ζ_{t} \sim N I D (0, σ_{ζ, t}^{2}) \end{aligned}$ with initial conditions $[\begin{matrix} α_{1}^{(1)} \\ α_{1}^{(2)} \end{matrix}] \sim N ([\begin{matrix} a_{1}^{(1)} \\ a_{1}^{(2)} \end{matrix}], [\begin{matrix} p_{1}^{(11)} & p_{1}^{(12)} \\ p_{1}^{(12)} & p_{1}^{(22)} \end{matrix}]) .$

Scalar Kalman filtering recursions

The Kalman filtering recursions, written in scalar form (to gain computational speed and insights) are the following (for generality we let also the variance of the measurement error vary over time).

The initial innovation, its variance and the Kalman gains are: $\begin{aligned} i_{1} & = y_{1} - a_{1}^{(1)} \\ f_{1} & = p_{1}^{(11)} + σ_{ε, 1}^{2} \\ k_{1}^{(1)} & = (p_{1}^{(11)} + p_{1}^{(12)}) / f_{1} \\ k_{1}^{(2)} & = p_{1}^{(12)} / f_{1} \end{aligned}$ For t = 1, 2, …, n − 1 the recursions are $\begin{aligned} a_{t + 1}^{(1)} & = a_{t}^{(1)} + a_{t}^{(2)} + k_{t}^{(1)} i_{t} \\ a_{t + 1}^{(2)} & = a_{t}^{(2)} + k_{t}^{(2)} i_{t} \\ p_{t + 1}^{(11)} & = p_{t}^{(11)} + 2 p_{t}^{(12)} + p_{t}^{(22)} + σ_{η, t}^{2} - k_{t}^{(1)} k_{t}^{(1)} f_{t} \\ p_{t + 1}^{(12)} & = p_{t}^{(12)} + p_{t}^{(22)} - k_{t}^{(1)} k_{t}^{(2)} f_{t} \\ p_{t + 1}^{(22)} & = p_{t}^{(22)} + σ_{ζ, t}^{2} - k_{t}^{(2)} k_{t}^{(2)} f_{t} \\ i_{t + 1} & = y_{t + 1} - a_{t + 1}^{(1)} \\ f_{t + 1} & = p_{t + 1}^{11} + σ_{ε, t + 1}^{2} \\ k_{t + 1}^{(1)} & = (p_{t + 1}^{(11)} + p_{t + 1}^{(12)}) / f_{t + 1} \\ k_{t + 1}^{(2)} & = p_{t + 1}^{(12)} / f_{t + 1} \end{aligned}$

Modifications when missing observations are present

When one or more values of y_t are missing, then the only modifications to the above recursions are the following. $\begin{aligned} i_{t + 1} & = 0 \\ f_{t + 1} & = \infty \\ k_{t + 1}^{(1)} & = 0 \\ k_{t + 1}^{(1)} & = 0 \end{aligned} .$

Diffuse initial conditions

Since the two state variables are nonstationary, their initialization should be diffuse: $[\begin{matrix} α_{1}^{(1)} \\ α_{1}^{(2)} \end{matrix}] \sim N ([\begin{matrix} 0 \\ 0 \end{matrix}], [\begin{matrix} v & 0 \\ 0 & v \end{matrix}]),$ with v → ∞.

As it will be clear from the computations below, when v is infinite, the mean squared errors of a_t⁽¹⁾ and a_t⁽²⁾, and the variances of the innovations are infinite for t = 1, 2, while from t = 3 on they are finite.

Let us carry out the computations for t = 1, 2, 3 and then take the limit for v → ∞.

t = 1

$\begin{aligned} a_{1}^{(1)} & = 0 \\ a_{1}^{(2)} & = 0 \\ p_{1}^{(11)} & = v \\ p_{1}^{(12)} & = 0 \\ p_{1}^{(22)} & = v \\ i_{1} & = y_{1} \\ f_{1} & = v + σ_{ε, 1}^{2} \\ k_{1}^{(1)} & = v / (v + σ_{ε, 1}^{2}) \\ k_{1}^{(2)} & = 0 \end{aligned}$

t = 2

$\begin{aligned} a_{2}^{(1)} & = y_{1} \\ a_{2}^{(2)} & = 0 \\ p_{2}^{(11)} & = 2 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} \\ p_{2}^{(12)} & = v \\ p_{2}^{(22)} & = v + σ_{ζ, 1}^{2} \\ i_{2} & = y_{2} - \frac{v}{v + σ_{ε}^{2}} y_{1} \to y_{2} - y_{1} \\ f_{2} & = 2 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε}^{2}} + σ_{ε, 2}^{2} \\ k_{2}^{(1)} & = \frac{3 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}}}{2 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} + σ_{ε, 2}^{2}} \to 2 \\ k_{2}^{(2)} & = \frac{v}{2 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} + σ_{ε, 2}^{2}} \to 1 \end{aligned}$

t = 3

$\begin{aligned} a_{3}^{(1)} & = \frac{v^{2}}{v + σ_{ε, 1}^{2}} y_{1} + \frac{3 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}}}{2 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} + σ_{ε, 2}^{2}} (y_{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} y_{1}) \to 2 y_{2} - y_{1} \\ a_{3}^{(2)} & = \frac{3 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}}}{2 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} + σ_{ε, 2}^{2}} (y_{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} y_{1}) \to y_{2} - y_{1} \\ p_{3}^{(11)} & = 5 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} + σ_{ζ, 1}^{2} + σ_{η, 2}^{2} - \frac{{(3 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}})}^{2}}{2 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} + σ_{ε, 2}^{2}} \to σ_{η, 1}^{2} + σ_{η, 2}^{2} + σ_{ζ, 1}^{2} \\ p_{3}^{(12)} & = 2 v + σ_{ζ, 1}^{2} - \frac{v (3 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}})}{2 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} + σ_{ε, 2}^{2}} \to σ_{ζ, 1}^{2} \\ p_{3}^{(22)} & = v + σ_{ζ, 1}^{2} + σ_{ζ, 2}^{2} - \frac{v^{2}}{2 v + σ_{η, 1}^{2} - \frac{v^{2}}{v + σ_{ε, 1}^{2}} + σ_{ε, 2}^{2}} \to σ_{ζ, 1}^{2} + σ_{ζ, 2}^{2} \\ i_{3} & \to y_{3} - 2 y_{2} + y_{1} \\ f_{3} & \to σ_{η, 1}^{2} + σ_{η, 2}^{2} + σ_{ζ, 1}^{2} + σ_{ε, 3}^{2} \\ k_{3}^{(1)} & \to \frac{σ_{η, 1}^{2} + σ_{η, 2}^{2} + 2 σ_{ζ, 1}^{2}}{σ_{η, 1}^{2} + σ_{η, 2}^{2} + σ_{ζ, 1}^{2} + σ_{ε, 3}^{2}} \\ k_{3}^{(2)} & \to \frac{σ_{ζ, 1}^{2}}{σ_{η, 1}^{2} + σ_{η, 2}^{2} + σ_{ζ, 1}^{2} + σ_{ε, 3}^{2}} \end{aligned}$

Smoothing

The smoothing recursions start from t = n and work backwards down to t = 1. The following quantities are auxiliar to compute the smoothed values of α_t⁽¹⁾ and their MSE. $\begin{aligned} r_{n + 1}^{(1)} & = 0 \\ r_{n + 1}^{(2)} & = 0 \\ n_{n + 1}^{(11)} & = 0 \\ n_{n + 1}^{(12)} & = 0 \\ n_{n + 1}^{(22)} & = 0 \\ e_{n} & = i_{n} / f_{n} \\ d_{n} & = 1 / f_{n} \end{aligned}$ For t = n, n − 1, …, 1, compute $\begin{aligned} r_{t}^{(1)} & = i_{t} / f_{t} + (1 - k_{t}^{(1)}) r_{t + 1}^{(1)} - k_{t}^{(2)} r_{t + 1}^{(2)} \\ r_{t}^{(2)} & = r_{t + 1}^{(1)} + r_{t + 1}^{(2)} \\ n_{t}^{(11)} & = {(1 - k_{t}^{(1)})}^{2} n_{t + 1}^{(11)} - 2 (1 - k_{t}^{(1)}) k_{t}^{(2)} n_{t + 1}^{(12)} + k_{t}^{(2)} k_{t}^{(2)} n_{t + 1}^{(22)} + 1 / f_{t} \\ n_{t}^{(12)} & = (1 - k_{t}^{(1)}) (n_{t + 1}^{(11)} + n_{t + 1}^{(12)}) - k_{t}^{(2)} (n_{t + 1}^{(12)} + n_{t + 1}^{(22)}) \\ n_{t}^{(22)} & = n_{t + 1}^{(11)} + 2 n_{t + 1}^{(12)} + n_{t + 1}^{(22)} \\ e_{t - 1} & = i_{t - 1} / f_{t - 1} - k_{t - 1}^{(1)} r_{t}^{(1)} - k_{t - 1}^{(2)} r_{t}^{(2)} \\ d_{t - 1} & = 1 / f_{t - 1} + k_{t - 1}^{(1)} k_{t - 1}^{(1)} n_{t}^{(11)} + 2 k_{t - 1}^{(1)} k_{t - 1}^{(2)} n_{t}^{(12)} + k_{t - 1}^{(2)} k_{t - 1}^{(2)} n_{t}^{(22)} \end{aligned}$ The smoothed values of α_t⁽¹⁾, that is the Hodrick-Prescott filtered time series, and their mean squared errors are given by $\begin{aligned} a_{t | n}^{(1)} & = a_{t}^{(1)} + p_{t}^{(11)} r_{t}^{(1)} + p_{t}^{(12)} r_{t}^{(2)} \\ p_{t | n}^{(11)} & = p_{t}^{(11)} - p_{t}^{(11)} p_{t}^{(11)} n_{t}^{(11)} - 2 p_{t}^{(11)} p_{t}^{(12)} n_{t}^{(12)} - p_{t}^{(12)} p_{t}^{(12)} n_{t}^{(22)} \end{aligned}$

Weights for computing the effective degrees of freedom

Since the smoother is linear in the observations, the vector of smoothed α_t⁽¹⁾, say s, is just a linear transformation of the vector of observations, y: s = Wy. The number of effective degrees of freedom is the trace of the weighting matrix W (cf. Hastie, Tibshirani and Friedman, 2009, The Elements of Statistical Learning, Section 5.4.1). The formulae for computing such weights in a general state-space form can be found in Koopman and Harvey (2003) Journal of Economic Dynamics and Control, vol. 27. In our framework, the diagonal elements of the matrix W are given by $\begin{aligned} w_{t t} & = p_{t}^{(11)} (1 / f_{t} + k_{t}^{(1)} k_{t}^{(1)} n_{t}^{(11)} + 2 k_{t}^{(1)} k_{t}^{(2)} n_{t}^{(12)} + k_{t}^{(2)} k_{t}^{(2)} n_{t}^{(22)} - k_{t}^{(1)} n_{t}^{(11)} - k_{t}^{(2)} n_{t}^{(12)}) \\ - p_{t}^{(12)} (k_{t}^{(1)} (n_{t}^{(11)} + n_{t}^{(12)}) + k_{t}^{(2)} (n_{t}^{(12)} + n_{t}^{(22)})) \end{aligned}$

Analytical scores

The log-likelihood must be maximised with respect to a very large number of parameters (n + 3). Thus, providing the numerical optimiser with analytical scores is important for stability and speed. Since all of our parameters are related to quantities in the disturbance covariance matrices, we can adapt the results in Koopman and Shephard (1992, Biometrika vol. 79).

Recall that our (slightly re-parametrised) model is $\begin{aligned} y_{t} & = α_{t}^{(1)} + ε_{t}, & ε_{t} \sim N I D (0, σ_{ε}^{2}) \\ α_{t + 1}^{(1)} & = α_{t}^{(1)} + α_{t}^{(2)} + η_{t}, & η_{t} \sim N I D (0, σ_{t}^{2}) \\ α_{t + 1}^{(2)} & = α_{t}^{(2)} + ζ_{t}, & ζ_{t} \sim N I D (0, σ^{2} + γ^{2} σ_{t}^{2}) \end{aligned}$ where the parameters to estimate are σ_ε, σ, γ, and the sequence {σ_t}_{t = 1, …, n}, which are all non-negative. Notice that in this parametrisation λ = σ_ε²/σ².

λ free

If λ is not fixed and ℓ(θ) represents the log-likelihood function, with θ vector all of the parameters, then $\begin{aligned} \frac{\partial ℓ}{\partial σ_{ε}} & = σ_{ε} \sum_{t = 1}^{n} (e_{t} e_{t} - d_{t}) \\ \frac{\partial ℓ}{\partial σ} & = σ \sum_{t = 1}^{n} (r_{t}^{(2)} r_{t}^{(2)} - n_{t}^{(22)}) \\ \frac{\partial ℓ}{\partial γ} & = γ \sum_{t = 1}^{n} (r_{t}^{(2)} r_{t}^{(2)} - n_{t}^{(22)}) σ_{t}^{2} \\ \frac{\partial ℓ}{\partial σ_{t}} & = (r_{t}^{(1)} r_{t}^{(1)} - n_{t}^{(11)} + (r_{t}^{(2)} r_{t}^{(2)} - n_{t}^{(22)}) γ^{2}) σ_{t}^{2} \end{aligned}$

Generally, constrained optimisation problems also need the derivatives of the constraining function, which in our case is $g(\boldsymbol{\theta}) = \sum_{t=1}^n \sigma_t$. The solution to the regularised maximum likelihood problem must satisfy g(θ) ≤ M. The derivatives are trivial: $\frac{\partial g}{\partial σ_{ε}} = 0, \frac{\partial g}{\partial σ} = 0, \frac{\partial g}{\partial γ} = 0, \frac{\partial g}{\partial σ_{t}} = 1.$

λ fixed

If λ is fixed, σ_ε² = λσ² and, in the log-likelihood function ℓ(θ), the vector of parameters θ does not contain λ or σ_ε². The derivatives are now $\begin{aligned} \frac{\partial ℓ}{\partial σ} & = σ \sum_{t = 1}^{n} (r_{t}^{(2)} r_{t}^{(2)} - n_{t}^{(22)}) - σ λ \sum_{t = 1}^{n} (e_{t} e_{t} - d_{t}) \\ \frac{\partial ℓ}{\partial γ} & = γ \sum_{t = 1}^{n} (r_{t}^{(2)} r_{t}^{(2)} - n_{t}^{(22)}) σ_{t} \\ \frac{\partial ℓ}{\partial σ_{t}} & = (r_{t}^{(1)} r_{t}^{(1)} - n_{t}^{(11)} + (r_{t}^{(2)} r_{t}^{(2)} - n_{t}^{(22)}) γ^{2}) σ_{t}^{2} \end{aligned}$ The derivatives of the constraining function are $\frac{\partial g}{\partial σ} = 0, \frac{\partial g}{\partial γ} = 0, \frac{\partial g}{\partial σ_{t}} = 1.$

- Introduction
- Basic model for the HP filter with jumps

Formulae