Estimating a TVP Model in Gretl

A state space model takes the form of the state (transition) equation \[\alpha_{t+1} = F_{t} \alpha_{t} + v_{t}\] and the observation (measurement) equation \[y_{t} = A_{t}^{\prime} x_{t} + H_{t}^{\prime} \alpha_{t} + w_{t}.\] \(y_{t}\) is the variable we observe and want to predict/explain. \(x_{t}\) is a vector of observed variables at time \(t\). \(\alpha_{t}\) is the value of the unobserved state vector at time \(t\). \(F\), \(A\), and \(H\) are matrices holding coefficients. \(F\), \(A\), and \(H\) can be either constant or time varying.

Gretl will do all the work of estimation for you once you’ve properly set up the problem, but due to the general nature of the state space model, it can be a challenge figuring out how to put everything into that form. This post walks you through an extension of one of the examples in the Gretl manual.

The Model

We want to estimate the following model, which is a modified version of the one in chapter 33 of the manual. The model is more general in order to add a variable to the observation equation. We have \[INFQ_{t} = \beta_0 + \beta_{1} URX_{t} + \beta_{2} INFT_{t} + \varepsilon_{t}\] where \(INFQ\) is the quarterly GDP deflator inflation rate, \(URX\) is the unemployment rate, and \(INFT\) is the growth of the consumption deflator.

If \(\beta_{0}\), \(\beta_{1}\), and \(\beta_{2}\) are constants, you can estimate the model using OLS, and do not need to concern yourself with the state space form or the Kalman filter. We want to allow for a time-varying Phillips curve relationship, so \(\beta_{1}\) follows the process \[\beta_{1, t+1} = \beta_{1, t} + \eta_{t}.\]

Load the Data

We’ll use a Euro-area macroeconomic dataset provided with Gretl:

open AWM.gdt --quiet
smpl 1974:1 1994:1

Setting Up the State Space Model

The state space is set up in the form of a bundle in Gretl. A bundle is similar to a C struct if you have programmed in C. The first thing we need to do is initialize the bundle. We have to pass the values of \(y\), \(H\), \(F\), and \(Q\) (the covariance matrix for \(v\) - with only one regressor in our model, it is the variance of \(v\)). We will use the ksetup function:

bundle B = ksetup(INFQ, 1, 1, 1)

The state vector includes only one variable, \(\beta_{1}\), so \(H\) has only one element, which we initialize with a 1. \(F\), which is the coefficient matrix for \(\alpha_{t}\), will therefore also have one element. \(Q\) has one element because there is only one state variable.

Let’s set up the variable with the time-varying coefficient. We need to add it to the bundle, and we’ll call it mX:

matrix B.mX = {URX}

The dependent variable when doing MLE won’t actually be \(INFQ\). It will be \(INFQ_{t} - \beta_{0}\) for the current guess of \(\beta_{0}\). Let’s create a matrix inside the bundle to hold \(INFQ\):

matrix B.depvar = {INFQ}

We also need to specify \(x\):

B.obsx = {INFT}

It is required to initialize the Kalman filter to some value. We don’t have much information, so we’ll specify to use a diffuse prior, which (hopefully) has no effect on the estimation.

B.diffuse = 1

Finally, because \(\beta_{1}\) is time varying, computing the error term of the observation equation requires the values \(\beta_{1t}\) and \(URX_{t}\) for every \(t\). The term \(H_{t}^{\prime} \alpha_{t}\) becomes \(URX_{t} \beta_{1t}\).

The state vector always changes for each \(t\), so Gretl will handle that, but \(H_{t}^{\prime}\) doesn’t need to change for each \(t\), so we handle that by writing a function to calculate the value of \(H_{t}^{\prime}\) and attaching that function to the bundle:

function void at_each_step(bundle *b)
  b.obsymat = transp(b.mX[b.t,])
end function

B.timevar_call = "at_each_step"

B.obsymat is the name for \(H\) inside the bundle. Now the function at_each_step will be called once for each observation during the Kalman filtering process. Note that b.t is the observation number at that time. Gretl automatically sets b.t to the appropriate value.

This completes the setup of the state space model. We have initialized the bundle, added all of the data, and added information about the time varying pieces.

MLE

In addition to the state vector, there are four parameters that need to be estimated, \(\beta_{0}\), \(\beta_{2}\), and the standard deviations of the observation and state equation error terms. We will call those parameters b0, b2, s_obs, and s_state inside the optimization routine. We have to specify initial values for each of them:

scalar b0 = mean(INFQ)
scalar b2 = 0
scalar s_obs = 0.1
scalar s_state = 0.1

Now implement the MLE routine. We update the values of the components of the state space model inside the bundle using the current choice of the four parameters, then apply the Kalman filter to the bundle to get the value of the likelihood function evaluated at that choice of the parameters. (After calling kfilter, the value of the likelihood function is stored in B.llt.)

mle LL = err ? NA : B.llt
    B.obsy = B.depvar - b0
    B.obsxmat = {b2}
    B.obsvar = s_obs^2
    B.statevar = s_state^2
    err = kfilter(&B)
    params b0 b2 s_obs s_state
end mle

Summarizing the Output

The primary goal of the estimation is to see the evolution of \(\widehat{\beta}_{1}\) over the sample. You need to run the Kalman smoother on the bundle to get the smoothed values of the state vector at each \(t\):

ksmooth(&B)

After ksmooth is done, there is a matrix B.state holding the values of the state vector and at each \(t\), and a matrix B.stvar holding the values of the var of the estimate of b1 at each \(t\). We can pull them out by extracting the first column of both matrices:

series tvar_b1hat = B.state[,1]
series tvar_b1se = sqrt(B.stvar[,1])

Make a nice-looking plot of \(\widehat{\beta}_{1}\) and a 95% confidence band:

gnuplot tvar_b1hat --time-series --with-lines --output=display \
  --band=tvar_b1hat,tvar_b1se,1.96 --band-style=fill

Complete Program

You can download the complete program here. Note that you’ll need a recent version of Gretl for the program to run. State space modeling was recently overhauled, so older versions of Gretl will not work.

Last Updated: October 19, 2017