-
Notifications
You must be signed in to change notification settings - Fork 15
Description
When using SSMseasonal with trigonometric seasonality and a reduced number of harmonics, the current setup of the package is a bit counterintuitive if you want to set a1 and P1. As an example, I've provided code here for the sunspots.month data which is present in base R. This time series has an extended seasonal cycle, occurring about every 11 years. Hence if we want to model the seasonality, we will most likely need the Fourier representation, and likely don't need/want the full harmonics. So in the example below, I've written the KFAS code to limit to the first five harmonics.
library(KFAS)
plot(sunspot.month)
y <- as.vector(sunspot.month)
KFAS_mod <- SSModel(y ~ SSMtrend(degree=1, Q=list(diag(1)), a1=c(0), P1= c(diag(3))) +
SSMseasonal(period=12*11, Q=as.matrix(0), sea.type="trigonometric", harmonics=1:5))
This works well. However, if we want to set a1 and P1, we might think that we need to give 10 values, since there are 10 harmonics, i.e. we would write SSMseasonal(period=12*11, Q=as.matrix(0), sea.type="trigonometric", harmonics=1:5, a1=rep(0,10), P1= c(diag(3,10))
But this returns an error:
a1 must be a (m x 1) matrix where m is the number of states
To avoid the error, we have to provide an a1 and P1 that would cover the full set of harmonics, even though we are only using a reduced set. So, this works: SSMseasonal(period=12*11, Q=as.matrix(0), sea.type="trigonometric", harmonics=1:5, a1=rep(0,(12*11)-1), P1= c(diag(3,(12*11)-1))
As I said, I think this is a bit counterintuitive and in my mind makes more sense if it worked like the first example. This is definitely a corner case, but thought it worth raising nonetheless.