- html - 出于某种原因,IE8 对我的 Sass 文件中继承的 html5 CSS 不友好?
- JMeter 在响应断言中使用 span 标签的问题
- html - 在 :hover and :active? 上具有不同效果的 CSS 动画
- html - 相对于居中的 html 内容固定的 CSS 重复背景?
我正在使用 pykalman 模块中的 KalmanFilter,我想知道它如何处理缺失的观察结果。根据文档:
In real world systems, it is common to have sensors occasionally fail. The Kalman Filter, Kalman Smoother, and EM algorithm are all equipped to handle this scenario. To make use of it, one only need apply a NumPy mask to the measurement at the missing time step:
from numpy import ma X = ma.array([1,2,3]) X1 = ma.masked # hide measurement at time step 1 kf.em(X).smooth(X)
我们可以平滑输入时间序列。由于这是一个“附加”功能,我认为它不是自动完成的;那么在变量中包含 NaN 时的默认方法是什么?
这里解释了可能发生的理论方法;这也是 pykalman 所做的吗(在我看来这真的很棒):
Cross Validated - How to handle incomplete data in Kalman Filter?
最佳答案
我们看一下源码:
在 filter_update 函数中,pykalman 检查当前观察是否被屏蔽。
def filter_update(...)
# Make a masked observation if necessary
if observation is None:
n_dim_obs = observation_covariance.shape[0]
observation = np.ma.array(np.zeros(n_dim_obs))
observation.mask = True
else:
observation = np.ma.asarray(observation)
它不影响预测步骤。但校正步骤有两种选择。它发生在 _filter_correct 函数中。
def _filter_correct(...)
if not np.any(np.ma.getmask(observation)):
# the normal Kalman Filter math
else:
n_dim_state = predicted_state_covariance.shape[0]
n_dim_obs = observation_matrix.shape[0]
kalman_gain = np.zeros((n_dim_state, n_dim_obs))
# !!!! the corrected state takes the result of the prediction !!!!
corrected_state_mean = predicted_state_mean
corrected_state_covariance = predicted_state_covariance
如您所见,这正是理论上的方法。
这是一个简短的示例和可供使用的工作数据。
假设您有一个 GPS 接收器,并且您想在行走时跟踪自己。接收机具有良好的精度。为简化起见,假设您只向前走。
没有什么有趣的事情发生。由于 GPS 信号良好,过滤器可以很好地估计您的位置。有一段时间没有信号怎么办?
过滤器只能根据现有状态和系统动力学知识进行预测(见矩阵Q)。随着每个预测步骤的不确定性增加。估计位置周围的 1-Sigma 范围变大。一旦再次出现新的观察,状态就会得到纠正。
这里是代码和 data :
from pykalman import KalmanFilter
import numpy as np
import matplotlib.pyplot as plt
from numpy import ma
# enable or disable missing observations
use_mask = 1
# reading data (quick and dirty)
Time=[]
X=[]
for line in open('data/dataset_01.csv'):
f1, f2 = line.split(';')
Time.append(float(f1))
X.append(float(f2))
if (use_mask):
X = ma.asarray(X)
X[300:500] = ma.masked
# Filter Configuration
# time step
dt = Time[2] - Time[1]
# transition_matrix
F = [[1, dt, 0.5*dt*dt],
[0, 1, dt],
[0, 0, 1]]
# observation_matrix
H = [1, 0, 0]
# transition_covariance
Q = [[ 1, 0, 0],
[ 0, 1e-4, 0],
[ 0, 0, 1e-6]]
# observation_covariance
R = [0.04] # max error = 0.6m
# initial_state_mean
X0 = [0,
0,
0]
# initial_state_covariance
P0 = [[ 10, 0, 0],
[ 0, 1, 0],
[ 0, 0, 1]]
n_timesteps = len(Time)
n_dim_state = 3
filtered_state_means = np.zeros((n_timesteps, n_dim_state))
filtered_state_covariances = np.zeros((n_timesteps, n_dim_state, n_dim_state))
# Kalman-Filter initialization
kf = KalmanFilter(transition_matrices = F,
observation_matrices = H,
transition_covariance = Q,
observation_covariance = R,
initial_state_mean = X0,
initial_state_covariance = P0)
# iterative estimation for each new measurement
for t in range(n_timesteps):
if t == 0:
filtered_state_means[t] = X0
filtered_state_covariances[t] = P0
else:
filtered_state_means[t], filtered_state_covariances[t] = (
kf.filter_update(
filtered_state_means[t-1],
filtered_state_covariances[t-1],
observation = X[t])
)
position_sigma = np.sqrt(filtered_state_covariances[:, 0, 0]);
# plot of the resulted trajectory
plt.plot(Time, filtered_state_means[:, 0], "g-", label="Filtered position", markersize=1)
plt.plot(Time, filtered_state_means[:, 0] + position_sigma, "r--", label="+ sigma", markersize=1)
plt.plot(Time, filtered_state_means[:, 0] - position_sigma, "r--", label="- sigma", markersize=1)
plt.grid()
plt.legend(loc="upper left")
plt.xlabel("Time (s)")
plt.ylabel("Position (m)")
plt.show()
更新
如果你掩盖更长的时期(300:700),它看起来会更有趣。
正如你所看到的,位置又回来了。这是因为转移矩阵 F 绑定(bind)了位置、速度和加速度的预测。
如果你看一下速度状态,它会解释下降的位置。
在 300 秒的时间点,加速度卡住。速度以恒定斜率下降并穿过 0 值。在这个时间点之后,位置必须返回。
要绘制速度,请使用以下代码:
velocity_sigma = np.sqrt(filtered_state_covariances[:, 1, 1]);
# plot of the estimated velocity
plt.plot(Time, filtered_state_means[:, 1], "g-", label="Filtered velocity", markersize=1)
plt.plot(Time, filtered_state_means[:, 1] + velocity_sigma, "r--", label="+ sigma", markersize=1)
plt.plot(Time, filtered_state_means[:, 1] - velocity_sigma, "r--", label="- sigma", markersize=1)
plt.grid()
plt.legend(loc="upper left")
plt.xlabel("Time (s)")
plt.ylabel("Velocity (m/s)")
plt.show()
关于numpy - pykalman:(默认)处理缺失值,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/49662567/
在Pykalman的文档中,它说它只接受observation_matrices和transition_matrices参数的方阵。 有办法解决这个问题吗?我必须估计具有非方形观测矩阵的状态空间系统。
我正在使用 pykalman 模块中的 KalmanFilter,我想知道它如何处理缺失的观察结果。根据文档: In real world systems, it is common to have
有了一维测量数据,我想使用卡尔曼滤波器了解每个点的状态标准偏差。我的流程如下: from pykalman import KalmanFilter import numpy as np measure
这是我在 Stackoverflow 上的第一个问题,所以如果我措辞不当,我深表歉意。我正在编写代码以从 IMU 获取原始加速度数据,然后对其进行整合以更新对象的位置。目前这段代码每毫秒获取一个新的加
我可以运行 pykalman documentation 中给出的简单 pykalman 卡尔曼滤波器示例: import pykalman import numpy as np kf = pykal
我是一名优秀的程序员,十分优秀!