# Estimating or specifying parameters in state space models¶

In this notebook we show how to fix specific values of certain parameters in statsmodels’ state space models while estimating others.

In general, state space models allow users to:

Estimate all parameters by maximum likelihood

Fix some parameters and estimate the rest

Fix all parameters (so that no parameters are estimated)

```
[1]:
```

```
%matplotlib inline
from importlib import reload
import numpy as np
import pandas as pd
import statsmodels.api as sm
import matplotlib.pyplot as plt
from pandas_datareader.data import DataReader
```

```
/tmp/ipykernel_3703/4283343036.py:5: DeprecationWarning:
Pyarrow will become a required dependency of pandas in the next major release of pandas (pandas 3.0),
(to allow more performant data types, such as the Arrow string type, and better interoperability with other libraries)
but was not found to be installed on your system.
If this would cause problems for you,
please provide us feedback at https://github.com/pandas-dev/pandas/issues/54466
import pandas as pd
```

To illustrate, we will use the Consumer Price Index for Apparel, which has a time-varying level and a strong seasonal component.

```
[2]:
```

```
endog = DataReader('CPIAPPNS', 'fred', start='1980').asfreq('MS')
endog.plot(figsize=(15, 3));
```

It is well known (e.g. Harvey and Jaeger [1993]) that the HP filter output can be generated by an unobserved components model given certain restrictions on the parameters.

The unobserved components model is:

For the trend to match the output of the HP filter, the parameters must be set as follows:

where \(\lambda\) is the parameter of the associated HP filter. For the monthly data that we use here, it is usually recommended that \(\lambda = 129600\).

```
[3]:
```

```
# Run the HP filter with lambda = 129600
hp_cycle, hp_trend = sm.tsa.filters.hpfilter(endog, lamb=129600)
# The unobserved components model above is the local linear trend, or "lltrend", specification
mod = sm.tsa.UnobservedComponents(endog, 'lltrend')
print(mod.param_names)
```

```
['sigma2.irregular', 'sigma2.level', 'sigma2.trend']
```

The parameters of the unobserved components model (UCM) are written as:

\(\sigma_\varepsilon^2 = \text{sigma2.irregular}\)

\(\sigma_\eta^2 = \text{sigma2.level}\)

\(\sigma_\zeta^2 = \text{sigma2.trend}\)

To satisfy the above restrictions, we will set \((\sigma_\varepsilon^2, \sigma_\eta^2, \sigma_\zeta^2) = (1, 0, 1 / 129600)\).

Since we are fixing all parameters here, we do not need to use the `fit`

method at all, since that method is used to perform maximum likelihood estimation. Instead, we can directly run the Kalman filter and smoother at our chosen parameters using the `smooth`

method.

```
[4]:
```

```
res = mod.smooth([1., 0, 1. / 129600])
print(res.summary())
```

```
Unobserved Components Results
==============================================================================
Dep. Variable: CPIAPPNS No. Observations: 529
Model: local linear trend Log Likelihood -2970.031
Date: Tue, 13 Feb 2024 AIC 5946.062
Time: 23:51:02 BIC 5958.864
Sample: 01-01-1980 HQIC 5951.074
- 01-01-2024
Covariance Type: opg
====================================================================================
coef std err z P>|z| [0.025 0.975]
------------------------------------------------------------------------------------
sigma2.irregular 1.0000 0.009 115.115 0.000 0.983 1.017
sigma2.level 0 0.000 0 1.000 -0.000 0.000
sigma2.trend 7.716e-06 2.01e-07 38.356 0.000 7.32e-06 8.11e-06
===================================================================================
Ljung-Box (L1) (Q): 251.07 Jarque-Bera (JB): 1.65
Prob(Q): 0.00 Prob(JB): 0.44
Heteroskedasticity (H): 2.33 Skew: 0.05
Prob(H) (two-sided): 0.00 Kurtosis: 2.75
===================================================================================
Warnings:
[1] Covariance matrix calculated using the outer product of gradients (complex-step).
```

The estimate that corresponds to the HP filter’s trend estimate is given by the smoothed estimate of the `level`

(which is \(\mu_t\) in the notation above):

```
[5]:
```

```
ucm_trend = pd.Series(res.level.smoothed, index=endog.index)
```

It is easy to see that the estimate of the smoothed level from the UCM is equal to the output of the HP filter:

```
[6]:
```

```
fig, ax = plt.subplots(figsize=(15, 3))
ax.plot(hp_trend, label='HP estimate')
ax.plot(ucm_trend, label='UCM estimate')
ax.legend();
```

## Adding a seasonal component¶

However, unobserved components models are more flexible than the HP filter. For example, the data shown above is clearly seasonal, but with time-varying seasonal effects (the seasonality is much weaker at the beginning than at the end). One of the benefits of the unobserved components framework is that we can add a stochastic seasonal component. In this case, we will estimate the variance of the seasonal component by maximum likelihood while still including the restriction on the parameters implied above so that the trend corresponds to the HP filter concept.

Adding the stochastic seasonal component adds one new parameter, `sigma2.seasonal`

.

```
[7]:
```

```
# Construct a local linear trend model with a stochastic seasonal component of period 1 year
mod = sm.tsa.UnobservedComponents(endog, 'lltrend', seasonal=12, stochastic_seasonal=True)
print(mod.param_names)
```

```
['sigma2.irregular', 'sigma2.level', 'sigma2.trend', 'sigma2.seasonal']
```

In this case, we will continue to restrict the first three parameters as described above, but we want to estimate the value of `sigma2.seasonal`

by maximum likelihood. Therefore, we will use the `fit`

method along with the `fix_params`

context manager.

The `fix_params`

method takes a dictionary of parameters names and associated values. Within the generated context, those parameters will be used in all cases. In the case of the `fit`

method, only the parameters that were not fixed will be estimated.

```
[8]:
```

```
# Here we restrict the first three parameters to specific values
with mod.fix_params({'sigma2.irregular': 1, 'sigma2.level': 0, 'sigma2.trend': 1. / 129600}):
# Now we fit any remaining parameters, which in this case
# is just `sigma2.seasonal`
res_restricted = mod.fit()
```

```
RUNNING THE L-BFGS-B CODE
* * *
Machine precision = 2.220D-16
N = 1 M = 10
At X0 0 variables are exactly at the bounds
At iterate 0 f= 3.89432D+00 |proj g|= 2.24624D-01
* * *
Tit = total number of iterations
Tnf = total number of function evaluations
Tnint = total number of segments explored during Cauchy searches
Skip = number of BFGS updates skipped
Nact = number of active bounds at final generalized Cauchy point
Projg = norm of the final projected gradient
F = final function value
* * *
N Tit Tnf Tnint Skip Nact Projg F
1 3 9 1 0 0 6.685D-06 3.272D+00
F = 3.2717006946573406
CONVERGENCE: NORM_OF_PROJECTED_GRADIENT_<=_PGTOL
```

```
This problem is unconstrained.
```

Alternatively, we could have simply used the `fit_constrained`

method, which also accepts a dictionary of constraints:

```
[9]:
```

```
res_restricted = mod.fit_constrained({'sigma2.irregular': 1, 'sigma2.level': 0, 'sigma2.trend': 1. / 129600})
```

```
RUNNING THE L-BFGS-B CODE
* * *
Machine precision = 2.220D-16
N = 1 M = 10
At X0 0 variables are exactly at the bounds
At iterate 0 f= 3.89432D+00 |proj g|= 2.24624D-01
* * *
Tit = total number of iterations
Tnf = total number of function evaluations
Tnint = total number of segments explored during Cauchy searches
Skip = number of BFGS updates skipped
Nact = number of active bounds at final generalized Cauchy point
Projg = norm of the final projected gradient
F = final function value
* * *
N Tit Tnf Tnint Skip Nact Projg F
1 3 9 1 0 0 6.685D-06 3.272D+00
F = 3.2717006946573406
CONVERGENCE: NORM_OF_PROJECTED_GRADIENT_<=_PGTOL
```

```
This problem is unconstrained.
```

The summary output includes all parameters, but indicates that the first three were fixed (and so were not estimated).

```
[10]:
```

```
print(res_restricted.summary())
```

```
Unobserved Components Results
=====================================================================================
Dep. Variable: CPIAPPNS No. Observations: 529
Model: local linear trend Log Likelihood -1730.730
+ stochastic seasonal(12) AIC 3463.459
Date: Tue, 13 Feb 2024 BIC 3467.705
Time: 23:51:03 HQIC 3465.123
Sample: 01-01-1980
- 01-01-2024
Covariance Type: opg
============================================================================================
coef std err z P>|z| [0.025 0.975]
--------------------------------------------------------------------------------------------
sigma2.irregular (fixed) 1.0000 nan nan nan nan nan
sigma2.level (fixed) 0 nan nan nan nan nan
sigma2.trend (fixed) 7.716e-06 nan nan nan nan nan
sigma2.seasonal 0.0920 0.007 12.660 0.000 0.078 0.106
===================================================================================
Ljung-Box (L1) (Q): 453.74 Jarque-Bera (JB): 39.93
Prob(Q): 0.00 Prob(JB): 0.00
Heteroskedasticity (H): 2.58 Skew: 0.32
Prob(H) (two-sided): 0.00 Kurtosis: 4.20
===================================================================================
Warnings:
[1] Covariance matrix calculated using the outer product of gradients (complex-step).
```

For comparison, we construct the unrestricted maximum likelihood estimates (MLE). In this case, the estimate of the level will no longer correspond to the HP filter concept.

```
[11]:
```

```
res_unrestricted = mod.fit()
```

```
RUNNING THE L-BFGS-B CODE
* * *
Machine precision = 2.220D-16
N = 4 M = 10
At X0 0 variables are exactly at the bounds
At iterate 0 f= 3.63589D+00 |proj g|= 1.07139D-01
```

```
This problem is unconstrained.
```

```
At iterate 5 f= 1.98480D+00 |proj g|= 5.13198D-01
ys=-9.333E-01 -gs= 3.188E-01 BFGS update SKIPPED
```

```
At iterate 10 f= 1.57685D+00 |proj g|= 4.30096D-01
At iterate 15 f= 1.48416D+00 |proj g|= 1.96427D-01
```

```
At iterate 20 f= 1.41376D+00 |proj g|= 3.11958D-01
```

```
At iterate 25 f= 1.38725D+00 |proj g|= 2.38467D-02
At iterate 30 f= 1.38723D+00 |proj g|= 7.81475D-06
* * *
Tit = total number of iterations
Tnf = total number of function evaluations
Tnint = total number of segments explored during Cauchy searches
Skip = number of BFGS updates skipped
Nact = number of active bounds at final generalized Cauchy point
Projg = norm of the final projected gradient
F = final function value
* * *
N Tit Tnf Tnint Skip Nact Projg F
4 30 58 1 1 0 7.815D-06 1.387D+00
F = 1.3872326906027199
CONVERGENCE: NORM_OF_PROJECTED_GRADIENT_<=_PGTOL
```

Finally, we can retrieve the smoothed estimates of the trend and seasonal components.

```
[12]:
```

```
# Construct the smoothed level estimates
unrestricted_trend = pd.Series(res_unrestricted.level.smoothed, index=endog.index)
restricted_trend = pd.Series(res_restricted.level.smoothed, index=endog.index)
# Construct the smoothed estimates of the seasonal pattern
unrestricted_seasonal = pd.Series(res_unrestricted.seasonal.smoothed, index=endog.index)
restricted_seasonal = pd.Series(res_restricted.seasonal.smoothed, index=endog.index)
```

Comparing the estimated level, it is clear that the seasonal UCM with fixed parameters still produces a trend that corresponds very closely (although no longer exactly) to the HP filter output.

Meanwhile, the estimated level from the model with no parameter restrictions (the MLE model) is much less smooth than these.

```
[13]:
```

```
fig, ax = plt.subplots(figsize=(15, 3))
ax.plot(unrestricted_trend, label='MLE, with seasonal')
ax.plot(restricted_trend, label='Fixed parameters, with seasonal')
ax.plot(hp_trend, label='HP filter, no seasonal')
ax.legend();
```

Finally, the UCM with the parameter restrictions is still able to pick up the time-varying seasonal component quite well.

```
[14]:
```

```
fig, ax = plt.subplots(figsize=(15, 3))
ax.plot(unrestricted_seasonal, label='MLE')
ax.plot(restricted_seasonal, label='Fixed parameters')
ax.legend();
```