PyMCStateSpace#

class pymc_experimental.statespace.core.PyMCStateSpace(k_endog: int, k_states: int, k_posdef: int, filter_type: str = 'standard', verbose: bool = True, measurement_error: bool = False)[source]#

Base class for Linear Gaussian Statespace models in PyMC.

Holds a PytensorRepresentation and KalmanFilter, and provides a mapping between a PyMC model and the statespace model.

Parameters:

k_endog (int) – The number of endogenous variables (observed time series).
k_states (int) – The number of state variables.
k_posdef (int) – The number of shocks in the model
filter_type (str, optional) – The type of Kalman filter to use. Valid options are “standard”, “univariate”, “single”, “cholesky”, and “steady_state”. For more information, see the docs for each filter. Default is “standard”.
verbose (bool, optional) – If True, displays information about the initialized model. Defaults to True.
measurement_error (bool, optional) – If true, the model contains measurement error. Needed by post-estimation sampling methods to decide how to compute the observation errors. If False, these errors are deterministically zero; if True, they are sampled from a multivariate normal.

Notes

Based on the statsmodels statespace implementation statsmodels/statsmodels, described in [1].

All statespace models inherit from this base class, which is responsible for providing an interface between a PyMC model and a PytensorRepresentation of a linear statespace model. This is done via the update method, which takes as input a vector of PyMC random variables and assigns them to their correct positions inside the underlying PytensorRepresentation. Construction of the parameter vector, called theta, is done automatically, but depend on the names provided in the param_names property.

To implement a new statespace model, one needs to:

Overload the param_names property to return a list of parameter names.
Overload the update method to put these parameters into their respective statespace matrices

In addition, a number of additional properties can be overloaded to provide users with additional resources when writing their PyMC models. For details, see the attributes section of the docs for this class.

Finally, this class holds post-estimation methods common to all statespace models, which do not need to be overloaded when writing a custom statespace model.

Examples

The local level model is a simple statespace model. It is a Gaussian random walk with a drift term that itself also follows a Gaussian random walk, as described by the following two equations:

\[\begin{split}\begin{align} y_{t} &= y_{t-1} + x_t + \nu_t \tag{1} \\ x_{t} &= x_{t-1} + \eta_t \tag{2} \end{align}\end{split}\]

Where \(y_t\) is the observed data, and \(x_t\) is an unobserved trend term. The model has two unknown parameters, the variances on the two innovations, \(sigma_\nu\) and \(sigma_\eta\). Take the hidden state vector to be \(\begin{bmatrix} y_t & x_t \end{bmatrix}^T\) and the shock vector \(\varepsilon_t = \begin{bmatrix} \nu_t & \eta_t \end{bmatrix}^T\). Then this model can be cast into state-space form with the following matrices:

\[\begin{split}\begin{align} T &= \begin{bmatrix}1 & 1 \\ 0 & 1 \end{bmatrix} & R &= \begin{bmatrix}1 & 0 \\ 0 & 1 \end{bmatrix} & Q &= \begin{bmatrix} \sigma_\nu & 0 \\ 0 & \sigma_\eta \end{bmatrix} & Z &= \begin{bmatrix} 1 & 0 \end{bmatrix} \end{align}\end{split}\]

With the remaining statespace matrices as zero matrices of the appropriate sizes. The model has two states, two shocks, and one observed state. Knowing all this, a very simple local level model can be implemented as follows:

from pymc_experimental.statespace.core import PyMCStateSpace
import numpy as np

class LocalLevel(PyMCStateSpace):
    def __init__():
        # Initialize the superclass. This creates the PytensorRepresentation and the Kalman Filter
        super().__init__(k_endog=1, k_states=2, k_posdef=2)

        # Declare the non-zero, non-parameterized matrices
        self.ssm['transition', :, :] = np.array([[1.0, 1.0], [0.0, 1.0]])
        self.ssm['selection', :, :] = np.eye(2)
        self.ssm['design', :, :] = np.array([[1.0, 0.0]])

    @property
    def param_names(self):
        return ['x0', 'P0', 'sigma_nu', 'sigma_eta']

    def update(self, theta, mode=None):
        # Since the param_names are ['x0', 'P0', 'sigma_nu', 'sigma_eta'], theta will come in as
        # [x0.ravel(), P0.ravel(), sigma_nu, sigma_eta]
        # It will have length 2 + 4 + 1 + 1 = 8

        x0 = theta[:2]
        P0 = theta[2:6].reshape(2,2)
        sigma_nu = theta[6]
        sigma_eta = theta[7]

        # Assign parameters to their correct locations
        self.ssm['initial_state', :] = x0
        self.ssm['initial_state_cov', :, :] = P0
        self.ssm['state_cov', 0, 0] = sigma_nu
        self.ssm['state_cov', 1, 1] = sigma_eta

After defining priors over the named parameters P0, x0, sigma_eta, and sigma_nu, we can sample from this model:

import pymc as pm

ll = LocalLevel()

with pm.Model() as mod:
    x0 = pm.Normal('x0', shape=(2,))
    P0_diag = pm.Exponential('P0_diag', 1, shape=(2,))
    P0 = pm.Deterministic('P0', pt.diag(P0_diag))
    sigma_nu = pm.Exponential('sigma_nu', 1)
    sigma_eta = pm.Exponential('sigma_eta', 1)

    ll.build_statespace_graph(data = data)
    idata = pm.sample()

References

__init__(k_endog: int, k_states: int, k_posdef: int, filter_type: str = 'standard', verbose: bool = True, measurement_error: bool = False)[source]#

Methods

`__init__`(k_endog, k_states, k_posdef[, ...])
`add_default_priors`()	Add default priors to the active PyMC model context
`add_exogenous`(exog)	Add an exogenous process to the statespace model
`build_statespace_graph`(data[, ...])	Given a parameter vector theta, constructs the full computational graph describing the state space model and the associated log probability of the data.
`forecast`(idata, start[, periods, end, ...])	Generate forecasts of state space model trajectories into the future.
`impulse_response_function`(idata[, n_steps, ...])	Generate impulse response functions (IRF) from state space model dynamics.
`make_and_register_data`(name, shape[, dtype])	Helper function to create a pytensor symbolic variable and register it in the _name_to_data dictionary
`make_and_register_variable`(name, shape[, dtype])	Helper function to create a pytensor symbolic variable and register it in the _name_to_variable dictionary
`make_symbolic_graph`()	The purpose of the make_symbolic_graph function is to hide tedious parameter allocations from the user.
`sample_conditional_posterior`(idata[, ...])	Sample from the conditional posterior; that is, given parameter draws from the posterior distribution, compute Kalman filtered trajectories.
`sample_conditional_prior`(idata[, random_seed])	Sample from the conditional prior; that is, given parameter draws from the prior distribution, compute Kalman filtered trajectories.
`sample_unconditional_posterior`(idata[, ...])	Draw unconditional sample trajectories according to state space dynamics, using random samples from the posterior distribution over model parameters.
`sample_unconditional_prior`(idata[, steps, ...])	Draw unconditional sample trajectories according to state space dynamics, using random samples from the prior distribution over model parameters.
`unpack_statespace`()	Helper function to quickly obtain all statespace matrices in the standard order.

Attributes

`coords`	PyMC model coordinates
`data_info`	Information about Data variables that need to be declared in the PyMC model block.
`data_names`	Names of data variables expected by the model.
`default_priors`	Dictionary of parameter names and callable functions to construct default priors for the model
`observed_states`	A k_endog length list of strings, associated with the model's observed states
`param_dims`	Dictionary of named dimensions for each model parameter
`param_info`	Information about parameters needed to declare priors
`param_names`	Names of model parameters
`shock_names`	A k_posdef length list of strings, associated with the model's shock processes
`state_names`	A k_states length list of strings, associated with the model's hidden states