包括来自 statsmodels 的状态 space 模型中的参数

Question

根据之前 , and the helpful answer, I've subclassed the MLEModel to encapsulate the model. I'd like to allow for two parameters q1 and q2 so that the state noise covariance matrix is generalized as in Sarkka (2013) 的示例 4.3 构建模型（术语 re-arranged 我的约定）：

我想我会用下面的 update 方法来完成这个，但是我运行遇到 fit 方法的问题，因为它 returns UFuncTypeError: Cannot cast ufunc 'multiply' output from dtype('complex128') to dtype('float64') with casting rule 'same_kind'。我在这里错过了什么？

import numpy as np
import scipy.linalg as linalg
import statsmodels.api as sm

class Tracker2D(sm.tsa.statespace.MLEModel):
    """Position tracker in two dimensions with four states

    """
    start_params = [1.0, 1.0]
    param_names = ["q1", "q2"]

    def __init__(self, endog):
        super(Tracker2D, self).__init__(endog, k_states=4)
        self.endog = endog
        self._state_names = ["x1", "dx1/dt",
                             "x3", "dx3/dt"]
        # dt: sampling rate; s = standard deviation of the process noise
        # common to both dimensions
        dt, s = 0.1, 0.5
        # dynamic model matrices A and Q
        A2d = [[1, dt],
               [0, 1]]
        A = linalg.block_diag(A2d, A2d)
        Q2d = [[dt ** 3 / 3, dt ** 2 / 2],
               [dt ** 2 / 2, dt]]
        Q = linalg.block_diag(Q2d, Q2d)
        # measurement model matrices H and R
        H = np.array([[1, 0, 0, 0],
                      [0, 0, 1, 0]])
        R = s ** 2 * np.eye(2)
        self["design"] = H
        self["obs_cov"] = R
        self["transition"] = A
        self["selection"] = np.eye(4)
        self["state_cov"] = Q

    def update(self, params, **kwargs):
        self["state_cov", :2, :2] *= params[0]
        self["state_cov", 2:, 2:] *= params[1]


# Initialization
m0 = np.array([[0, 1, 0, -1]]).T  # state vector column vector
P0 = np.eye(4)                   # process covariance matrix

# With object Y below being the simulated measurements in downloadable
# data file from previous post
with open("measurements_2d.npy", "rb") as f:
    Y = np.load(f)

tracker2D = Tracker2D(pd.DataFrame(Y.T))
tracker2D.initialize_known((tracker2D["transition"] @ m0.flatten()),
                           (tracker2D["transition"] @ P0 @
                            tracker2D["transition"].T +
                            tracker2D["state_cov"]))
# Below throws the error
tracker2D.fit()

Answer 1

您收到的错误消息是关于尝试在 dtype=float 矩阵中设置复数值。你会得到同样的错误：

A = np.eye(2)
A *= 1.0j

错误出现在：

def update(self, params, **kwargs):
    self["state_cov", :2, :2] *= params[0]
    self["state_cov", 2:, 2:] *= params[1]

因为您正在原地修改“state_cov”。当params是一个复数向量，而现有的“state_cov”矩阵的dtype是float时，就会出现这个错误。 Statsmodels在计算参数的标准误差时会将参数向量设置为复数，因为它使用复步微分。

你可以使用像

这样的东西

def update(self, params, **kwargs):
    self["state_cov", :2, :2] = params[0] * self["state_cov", :2, :2]
    self["state_cov", 2:, 2:] = params[1] * self["state_cov", 2:, 2:]

虽然我应该指出，这不会给你我认为你想要的，因为它会根据以前的内容修改“state_cov”。相反，我认为你想要这样的东西：

class Tracker2D(sm.tsa.statespace.MLEModel):
    """Position tracker in two dimensions with four states

    """
    start_params = [1.0, 1.0]
    param_names = ["q1", "q2"]

    def __init__(self, endog):
        super(Tracker2D, self).__init__(endog, k_states=4)
        self.endog = endog
        self._state_names = ["x1", "dx1/dt",
                             "x3", "dx3/dt"]
        # dt: sampling rate; s = standard deviation of the process noise
        # common to both dimensions
        dt, s = 0.1, 0.5
        # dynamic model matrices A and Q
        A2d = [[1, dt],
               [0, 1]]
        A = linalg.block_diag(A2d, A2d)
        Q2d = [[dt ** 3 / 3, dt ** 2 / 2],
               [dt ** 2 / 2, dt]]

        # First we save the base Q matrix
        self.Q = linalg.block_diag(Q2d, Q2d)

        # measurement model matrices H and R
        H = np.array([[1, 0, 0, 0],
                      [0, 0, 1, 0]])
        R = s ** 2 * np.eye(2)
        self["design"] = H
        self["obs_cov"] = R
        self["transition"] = A
        self["selection"] = np.eye(4)
        self["state_cov"] = self.Q.copy()

    def update(self, params, **kwargs):
        # Now update the state cov based on the original Q
        # matrix, and set entire blocks of the matrix, rather
        # than modifying it in-place.
        self["state_cov", :2, :2] = params[0] * self.Q[:2, :2]
        self["state_cov", 2:, 2:] = params[1] * self.Q[2:, 2:]

包括来自 statsmodels 的状态 space 模型中的参数

Including parameters in state space model from statsmodels

python

statsmodels

state-space