Rolling forecast with StateSpaceModels.jl

General Usage
,
1

Hello,

I have a question on using StateSpaceModels.jl to do a rolling window prediction of a time series.

This means:

  1. At time t using past data y_1...y_t, predict h steps ahead (\hat{y}_{t+1|t} to \hat{y}_{t+h|t}). This includes fitting the time series model (hyper parameters in the doc of StateSpaceModels.jl).
  2. Then time moves forward so that the actual value y_{t+1} gets available. Update the prediction filter with this new data point to forecast over the new horizon (\hat{y}_{t+2|t+1} to \hat{y}_{t+h+1 |t+1}).

And this second step (which is meant to be iterated) should be done:

  • without refitting the model hyperparameters over the updated timeseries y_1...y_{t+1}
  • without having the filter computation rerun over the entire updated timeseries

i.e. I would like to just move the Kalman filter one step further with the new observation y_{t+1}

Also as a variant of this question, how to use a fitted state-space model to do forecasting on a different time series, keeping all hyperparameters (probably including initial state estimates).

Data and code example

For the illustration, here is a code which creates a noisy sine wave and makes a forecast at some point of the time series (t=kf). It works, but my goal is now to move this forecast instant forward.

using StateSpaceModels
using Random

# noisy data gen
m = 20 # periodicity
nm = 10 # number of periods

t = range(0, m*nm, step=1)
n = length(t)

# noise
σ = 0.3
Random.seed!(0)
e = σ*randn(n)

y = sin.(2*pi*t/m) + e;

Results:

sine-series
sine-series1750×625 85 KB

Forecast setup

h = m*2 # Forecast horizon: 2 periods
kf = n-m # Instant at which to forecast: close to the end

# time series fragment up to forecast instant, for fitting:
tfit = t[1:kf]
yfit = y[1:kf]
# time vector of forecast
tf = t[kf] .+ range(1, h, step=1)
# true value, if it exists:
yf_true = [kf+hh<= n ? y[kf+hh] : NaN for hh=1:h]
tf

Forecast with Exponential Smoothing with Seasonality ETS(A,N,A)

mod = ExponentialSmoothing(yfit, seasonal=m)
fit!(mod)
forec = forecast(mod2, h); # plot forec.expected_value

which works as expected:

forecast_season
forecast_season1500×625 73 KB

Only minor isse: I get the warning

Warning: The optimization process converged but the Hessian matrix is not positive definite. This means that StateSpaceModels.jl cannot estimate the distribution of the hyperparameters…

But perhaps this relate to me using a deterministic trend, so that the estimation problem is ill-defined? In particular, in the result the smoothing level and smoothing_seasonal parameters are close to (0.0001), so I suppose that the initial state estimates (level and seasonal levels) have a strong weight in the forecast. and perhaps these initial seasonal levels corresponds to the “seasonal averages” over the time series, I didn’t check. Anyway, my question is not on this!