# How to include seasonal effects into the system matrices of a state space model

I am working on learning state space models and am leaning heavily on this very helpful documentation. However, I'm really confused about the best way to include both a seasonal effect and dynamic regressors into the model.

Let's further assuming that I want to use a simulation smoother to estimate the coefficients on the dynamic regressors.

My initial thought was to do something like the following (with the seasonal effect specified as dummy variables that sum to 1). That is to say, first regressing the dependent variable on the seasonal effect, then subtracting that effect and passing the remainder into the state space model.


y ~ normal(seasonal_effect, sigma)
y_prime = y - seasonal_effect
y_prime ~ state_space_model(...system matrices go here...)



The issue that I'm having with this specification is that:

1. It doesn't estimate any kind of seasonal effect even though my model clearly has one
2. It's estimating the dynamic linear model and the seasonal effect separately when they should seemingly be estimated together

From the available literature it seems that I should be able to code the seasonal effect directly into the system matrices, but I have no idea how.

1. If I put the seasonal dummies into the Z (or F, depending on the literature) design matrix, they are treated as a dynamic regressor, unless...
2. I restrict the state selection covariance matrix to only allow the non-seasonal dummies to evolve over time, but then the state estimates that I get from the seasonal dummies stay fixed at the values that I pass into the simulation smoother (as a1)

The libraries and models that I'm working with are very very clear once you have the system matrices specified, but I am lost on how to include seasonal effects into the system matrices of the dynamic linear model

Any help would be appreciated!

Here is an example from KFAS with two covariates and seasonal of length 12:

library(KFAS)
model_drivers <- SSModel(log(drivers) ~
SSMseasonal(period = 12,  Q = 0) +
log(PetrolPrice) + law, data = Seatbelts, H = NA)


First time point:

model_drivers$Z[,,1] (Intercept) log(PetrolPrice) law sea_dummy1 sea_dummy2 sea_dummy3 sea_dummy4 sea_dummy5 sea_dummy6 sea_dummy7 sea_dummy8 sea_dummy9 sea_dummy10 sea_dummy11 1.0000 -2.2733 0.0000 1.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000  And matrix T: model_drivers$T
, , 1

(Intercept) log(PetrolPrice) law sea_dummy1 sea_dummy2 sea_dummy3 sea_dummy4 sea_dummy5 sea_dummy6 sea_dummy7 sea_dummy8 sea_dummy9 sea_dummy10 sea_dummy11
(Intercept)                1                0   0          0          0          0          0          0          0          0          0          0           0           0
log(PetrolPrice)           0                1   0          0          0          0          0          0          0          0          0          0           0           0
law                        0                0   1          0          0          0          0          0          0          0          0          0           0           0
sea_dummy1                 0                0   0         -1         -1         -1         -1         -1         -1         -1         -1         -1          -1          -1
sea_dummy2                 0                0   0          1          0          0          0          0          0          0          0          0           0           0
sea_dummy3                 0                0   0          0          1          0          0          0          0          0          0          0           0           0
sea_dummy4                 0                0   0          0          0          1          0          0          0          0          0          0           0           0
sea_dummy5                 0                0   0          0          0          0          1          0          0          0          0          0           0           0
sea_dummy6                 0                0   0          0          0          0          0          1          0          0          0          0           0           0
sea_dummy7                 0                0   0          0          0          0          0          0          1          0          0          0           0           0
sea_dummy8                 0                0   0          0          0          0          0          0          0          1          0          0           0           0
sea_dummy9                 0                0   0          0          0          0          0          0          0          0          1          0           0           0
sea_dummy10                0                0   0          0          0          0          0          0          0          0          0          1           0           0
sea_dummy11                0                0   0          0          0          0          0          0          0          0          0          0           1           0


You can probably figure out the logic here, essentially the first row of seasonal part is full of values -1 and then there is values of 1 just below the diagonal of the seasonal part. And Z just contains one in the first entry and zeros on the rest.

Then you also have to remember to add some uncertainty for the seasonal states, i.e. prior state_1 ~ N(prior_mean, prior_variance).

With this model and wide priors for the variance parameter H, you should get the same results with KFAS as with some Bayesian solutions, so you can use KFAS to check the results even if you want to do fully Bayesian inference. And there is also an R package bssm which you can use for these kind of models directly for fully Bayesian inference (check functions bsm and run_mcmc).