Defining models#

Models outside the built-in components are defined in one of three ways.

kind

how

estimation

MappedModel

declare which entries each parameter sets

in Fortran; works with fit_many

MLEModel

subclass and write update(params)

calls Python at every evaluation

Fortran ssm_model_t

extend the type

in Fortran

Declared models#

A MappedModel is set up in three steps: the fixed matrices, as mod[name] = value; the initialization, with the initialize_* methods; and the parameters. map() lets a parameter set an entry, cov() a covariance block, and constrain() sets the transforms. An ARMA(1, 1) in DK’s form (DK §3.4):

>>> y = np.random.default_rng(8).standard_normal(200)
>>> arma = ss.MappedModel(y, k_states=2, k_posdef=1, k_params=3,
...                       param_names=["ar.L1", "ma.L1", "sigma2"],
...                       start_params=[0.0, 0.0, 1.0])
>>> arma["design"] = [[1.0, 0.0]]
>>> arma["transition"] = [[0.0, 1.0], [0.0, 0.0]]
>>> arma["selection"] = [[1.0], [0.0]]
>>> arma.initialize_stationary()
>>> _ = arma.map(0, "transition", 0, 0)          # T[0, 0] = phi
>>> _ = arma.map(1, "selection", 1, 0)           # R[1, 0] = theta
>>> _ = arma.map(2, "state_cov", 0, 0)           # Q = sigma2
>>> _ = arma.constrain(0, "stationary").constrain(1, "invertible").constrain(2, "positive")
>>> res = arma.fit()
>>> res.param_names
['ar.L1', 'ma.L1', 'sigma2']

It gives the same likelihood as ss.ARIMA(order=(1, 0, 1)). Declared models run entirely in Fortran, so they are as fast as the built-in ones and fit_many() can fit them in parallel.

Starting from components#

The built-in components can supply the fixed part of a declared model. representation() returns a structural model’s representation at given parameters: the components’ system matrices and initialization, with their states in the order the components were given. Passed in place of the data, it becomes a MappedModel’s starting matrices and initialization; mapped entries overwrite them at each evaluation. Mapping the structural model’s own parameters reproduces it:

>>> y = np.loadtxt("data/nile.csv", delimiter=",", skiprows=1)[:, 1]
>>> st = ss.StructuralModel(y, [ss.Irregular(), ss.Trend()])
>>> st.param_names
['sigma2.irregular', 'sigma2.level', 'sigma2.slope']
>>> rep = st.representation(st.start_params)
>>> mod = ss.MappedModel(rep, k_params=3, param_names=st.param_names,
...                      start_params=st.start_params)
>>> _ = mod.map(0, "obs_cov", 0, 0)       # H = sigma2.irregular
>>> _ = mod.map(1, "state_cov", 0, 0)     # Q[0, 0] = sigma2.level
>>> _ = mod.map(2, "state_cov", 1, 1)     # Q[1, 1] = sigma2.slope
>>> _ = mod.constrain([0, 1, 2], "positive")
>>> bool(np.isclose(mod.fit().llf, st.fit().llf))
True

From there, further map calls add parameters of your own; give the model more parameters and names to match. rep["transition"] and the other matrices show the layout to map into.

A map sets an entry to a multiple of one parameter, which covers the irregular, level, trend, seasonal and regression components and non-seasonal ARIMA. The cycle (\(\rho \cos\lambda\) in \(T\)) and seasonal ARIMA (products of coefficients) are not of that form. For those, use an MLEModel whose update copies the matrices of st.representation(params) and changes what it needs; building the representation at each evaluation costs some speed.

Python models#

MLEModel follows statsmodels: the subclass sets the fixed matrices in __init__ and overrides update, start_params and, optionally, the transforms and names. See the example in MLEModel.

The library calls update at every likelihood evaluation, and the score calls it 2k more times per gradient for its matrix derivatives. For the local level model this is still about twice as fast as statsmodels (see Performance). An exception raised in update propagates to the caller, and the model remains usable. Python models cannot be fitted with fit_many.

Fortran models#

In Fortran, extend ssm_model_t, set up rep in a constructor, and implement update and start_params; see statespace_model. New components for structural models extend component_t (statespace_components).