pymc-experimental
pymc-experimental copied to clipboard
Add ETS model
Adds BayesianETS as an available statespace model.
According to this book, this is great baseline model for forecasting. I don't have a lot of experience with it, but statsmodels implements it, so I did too.
Todo:
- [ ] Write an example notebook
- [ ] Add a class method utility to decompose the time series into level, trend, and seasonal components
- [ ] Test time series generation against statsmodels
- [ ] Check whether the first state is correctly interpreted in the Kalman filter
The last point is a but subtle. From the statsmodel implementation:
One consequence is that the "initial state" corresponds to the "filtered" state at time t=0, but this is different from the usual state space initialization used in Statsmodels, which initializes the model with the "predicted" state at time t=1.
I have no idea whether this is important or not for the present implementation, I need to look at it more closely.
Some comments: In general I just follow the statsmodel implementation, and I'm testing that all the statespace matrices match those of statsmodels. The one difference is in the selection matrix. Looking at the book equations, for example:
$$ \begin{align} y_t &= \ell_t + \epsilon_t \ \ell_t &= \ell_{t-1} + \alpha \epsilon_t \end{align} $$
The general statespace equation is:
$$ \begin{align} x_t &= T x_{t-1} + c + R \epsilon_t \ y_t &= Z x_t + d + H \end{align} $$
This implies states:
$$x_t = \begin{bmatrix} \epsilon_t \ \ell_t \end{bmatrix}$$
So:
$$T = \begin{bmatrix} 0 & 0 \ 0 & 1 \end{bmatrix}, R = \begin{bmatrix} 1 \ \alpha \end{bmatrix} $$
This all seems straightforward, but statsmodels sets the (0, 0) element of $R$ to be $1 - \alpha$. That really doesn't seem right, so I didn't do it.
So statsmodels was right of course. My equations above were slightly wrong. We want to have:
$$ \begin{align} y_t &= \ell_{t-1} + \epsilon_t \ \ell_t &= \ell_{t-1} + \alpha \epsilon_t \end{align} $$
This could be done by making $y_t$ a hidden state, then:
$$T = \begin{bmatrix} 0 & 1 \ 0 & 1 \end{bmatrix}, Z_t = \begin{bmatrix} 1 & 0 \end{bmatrix}$$
But we can also get back to the timing that I had above, with $y_t = \ell_t + \epsilon_t$, if we do some algebra. Specifically, $\ell_{t-1} = \ell_t - \alpha \epsilon_t$, so $y_t = \ell_t + \epsilon_t - \alpha \epsilon_t = \ell_t + (1 - \alpha) \epsilon_t$
So that's where the $R$ adjustments in statsmodels were coming from. The same thing happens when you bring in a seasonal component; it becomes $y_t = \ell_t + s_{t-m}+ (1 - \alpha - \gamma) \epsilon_t$
I also added an alternative parameterization that I thing might be nicer for priors. Under the $\alpha, \beta, \gamma$ parameterization, you need $\alpha \in (0, 1)$, $\beta \in (0, \alpha)$, and $\gamma \in (0, 1 - \alpha)$. Which is pretty hairy to set priors for. If instead you go with:
$$ \begin{align} y_t &= \ell_{t-1} + s_{t-m-1} + \epsilon_t \ \ell_t &= \ell_{t-1} + b_{t-1} + \alpha \epsilon_t \ b_t &= b_{t-1} + \alpha \beta^\star \epsilon_t \ s_{t-m} &= s_{t-m-1} + (1 - \alpha) \gamma^\star \epsilon_t \end{align} $$
Then all of $\alpha, \beta^\star, \gamma^\star$ are between 0 and 1, which is nice. So that's the default now, and you can ask for the other one with use_transformed_parameterization=True.
Check out this pull request on ![]()
See visual diffs & provide feedback on Jupyter Notebooks.
Powered by ReviewNB
Typo: 3 instances of "forcasting" (e missing) in the notebook
Oh wow, I have seen folks express ETS as a states pace model, but never in a way that so clear and straightforward - very nice!