Уравнение наблюдения
$$y_t=\mu_t+\sum_{j=1}^ {k}\beta_jx_{jt}+\varepsilon_t$$
Уравнение уровня
$$\mu_{t +1}=\mu_t+\delta_t\nu_t+\eta_t$$
Уравнение тренда
$$\nu_{ t+1}=\nu_t+\zeta_t$$
с $\delta_t$ как разницей во времени между наблюдениями или приращениями времени, поскольку наблюдения не равноудалены во времени.
У нас есть общая форма пространства состояний.
$$y_t=Z_t\alpha_t+\varepsilon_t$$ с $\varepsilon_t\sim\mathbb{N}(0 ,H_t)$
$$\alpha_{t+1}=T_t\alpha_t+\eta_t)$$ с $\eta_t\sim\mathbb{N}(0,Q_t)$
Если оставить в стороне экзогенные переменные, ниже приведены матрицы соответствующей модели пространства состояний.
Матрица проектирования $Z$
$$\begin{bmatrix}1&0\end{bmatrix}$$
Уравнение состояний $\alpha$
$$\begin{bmatrix}\mu_{t+1}\\nu_{t+1}\end{bmatrix}=\begin{bmatrix}1&\delta_t\0&1\ end{bmatrix}\begin{bmatrix}\mu_t\\nu_t\end{bmatrix}+\begin{bmatrix}\eta_t\\zeta_t\end{bmatrix}$$
Матрица перехода $T$
$$\begin{bmatrix}1&\delta_t\\0&1\end{bmatrix}$$
Матрица наблюдения ковариации $H$
$$\begin{bmatrix}\sigma^2_{\varepsilon}\end{bmatrix}$$
Матрица ковариационных состояний $Q$
$$\begin{bmatrix}\sigma^2_{\eta}&0\ \0&\sigma^2_{\zeta}\end{bmatrix}$$
Как вставить экзогенные переменные в матричное определение пространства состояний?
Я создаю экзогенный фрейм данных следующим образом
Код: Выделить всё
f_year = DeterministicProcess(X.index, period=period, fourier=harmonics).in_sample()
X_exog = pd.concat([X, f_year], axis=1)
Код: Выделить всё
kwargs = {'delta_t': delta_time[start : end]}
mod = models.DiscreteTimeLocalLinearTrend(Y[start : end], X[start : end], **kwargs)
mod_fit = mod.fit(maxiter=500, method='powell', disp=0)
kwargs = {'delta_t': delta_time[end : end + steps_predictions]}
y_pred = mod_fit.forecast(steps=steps_predictions, exog=X[end : end + steps_predictions], **kwargs)
Код: Выделить всё
class DiscreteTimeLocalLinearTrend(sm.tsa.statespace.MLEModel):
def __init__(self, endog, exog, **kwargs):
# Model order
n_exog = exog.shape[1]
k_states = 2 + n_exog
k_posdef = k_states
exog = exog.values
delta_t = kwargs.get('delta_t').values
# Initialize the statespace
super(DiscreteTimeLocalLinearTrend, self).__init__(
endog=endog, exog=exog,
k_states=k_states, k_posdef=k_posdef,
initialization='approximate_diffuse',
loglikelihood_burn=k_states,
**kwargs
)
# Design matrix
self.ssm['design'] = np.zeros((1, k_states, len(endog)))
self.ssm['design'][0, :2, :] = 1
self.ssm['design'][0, 2:, :] = exog.T
# Transition matrix
self.ssm['transition'] = np.zeros((k_states, k_states, len(endog)))
for t in range(len(endog)):
self.ssm['transition'][:2, :2, t] = np.array([[1, delta_t[t]], [0, 1]])
self.ssm['transition'][0, 2:, t] = exog[t]
# Selection matrix
self.ssm['selection'] = np.eye(k_states)
# State covariance matrix
self.ssm['state_cov'] = np.eye(k_posdef)
@property
def param_names(self):
return ['sigma2.measurement', 'sigma2.level', 'sigma2.trend']
@property
def start_params(self):
return [np.std(self.endog)]*3
def transform_params(self, unconstrained):
return unconstrained**2
def untransform_params(self, constrained):
return constrained**0.5
def clone(self, endog, exog=None, **kwargs):
return self._clone_from_init_kwds(endog, exog=exog, **kwargs)
def update(self, params, *args, **kwargs):
delta_t = kwargs.get('delta_t').values
params = super(DiscreteTimeLocalLinearTrend, self).update(params, *args, **kwargs)
# Observation covariance matrix
self.ssm['obs_cov',0,0] = params[0]
# State covariance matrix
self.ssm['state_cov'] = np.zeros((self.ssm.k_states, self.ssm.k_states))
np.fill_diagonal(self.ssm['state_cov'], np.real(params[1:self.ssm.k_states+1]))
# Design matrix
self.ssm['design'][0, 2:, :] = self.exog.T
# Transition matrix
for t in range(len(self.endog)):
self.ssm['transition'][:2, :2, t] = np.array([[1, delta_t[t]], [0, 1]])
Код: Выделить всё
# Create extended model
if extend_kwargs is None:
extend_kwargs = {}
Код: Выделить всё
# Create extended model
if extend_kwargs is None:
extend_kwargs = kwargs
После решения этой проблемы я получите ошибку в методе update, поскольку delta_t не найден в kwargs. Если я прокомментирую матрицу перехода кода в методе обновления, код запускается, но на графике прогноза я замечаю, что влияние экзогенных переменных не учитывается.
Оставляем код обновления матрицы перехода закомментирован и используется инициализация матрицы перехода как обычный локальный линейный тренд из
Код: Выделить всё
for t in range(len(endog)):
self.ssm['transition'][:2, :2, t] = np.array([[1, delta_t[t]], [0, 1]])
self.ssm['transition'][0, 2:, t] = exog[t]
Код: Выделить всё
self.ssm['transition'] = np.eye(k_states)
self.ssm['transition'][:2, :2] = np.array([[1, 1],
[0, 1]])
Где я ошибаюсь?
Подробнее здесь: https://stackoverflow.com/questions/791 ... ion-matrix
Мобильная версия