DeepSORT(目标跟踪算法)中卡尔曼增益的理解
flyfish
先用最简单的例子来理解卡尔曼增益
公式 (1)
首先,通过多次测量一个物理量,并使用取平均值的方式来计算其真实值:
x ^ k = 1 k ( z 1 + z 2 + ⋯ + z k ) \hat{x}_k = \frac{1}{k} (z_1 + z_2 + \dots + z_k) x^k=k1(z1+z2+⋯+zk)
其中, x ^ k \hat{x}_k x^k 是对真实值的估计, z i z_i zi 是第 i i i 次测量的测量值。
公式 (2)
对公式 (1) 进行简单变换可以得到:
x ^ k = 1 k ( z 1 + z 2 + ⋯ + z k ) \hat{x}_k = \frac{1}{k} (z_1 + z_2 + \dots + z_k) x^k=k1(z1+z2+⋯+zk)
= 1 k ( ( z 1 + z 2 + ⋯ + z k − 1 ) + z k ) = \frac{1}{k} \left( (z_1 + z_2 + \dots + z_{k-1}) + z_k \right) =k1((z1+z2+⋯+zk−1)+zk)
= 1 k ( k − 1 ) x ^ k − 1 + 1 k z k = \frac{1}{k} (k-1) \hat{x}_{k-1} + \frac{1}{k} z_k =k1(k−1)x^k−1+k1zk
= x ^ k − 1 + 1 k ( z k − x ^ k − 1 ) = \hat{x}_{k-1} + \frac{1}{k} (z_k - \hat{x}_{k-1}) =x^k−1+k1(zk−x^k−1)
通过递归的方式计算当前的估计值 x ^ k \hat{x}_k x^k,其中包含了上一次的估计值 x ^ k − 1 \hat{x}_{k-1} x^k−1 和当前的测量值 z k z_k zk。
公式 (3)
可以把公式 (2) 中的后一项改写为:
x ^ k = x ^ k − 1 + K k ( z k − x ^ k − 1 ) \hat{x}_k = \hat{x}_{k-1} + K_k (z_k - \hat{x}_{k-1}) x^k=x^k−1+Kk(zk−x^k−1)
其中:
K k = 1 k K_k = \frac{1}{k} Kk=k1
这个 K k K_k Kk 就相当于卡尔曼滤波器中的卡尔曼增益(Kalman Gain)。
详细说明
- 初始估计值:
x ^ 1 = z 1 \hat{x}_1 = z_1 x^1=z1
初始时,对真实值的估计就是第一次测量的测量值。 - 递归更新:
递归公式: x ^ k = x ^ k − 1 + 1 k ( z k − x ^ k − 1 ) \hat{x}_k = \hat{x}_{k-1} + \frac{1}{k} (z_k - \hat{x}_{k-1}) x^k=x^k−1+k1(zk−x^k−1)当前的估计值 x ^ k \hat{x}_k x^k 是通过上一次的估计值 x ^ k − 1 \hat{x}_{k-1} x^k−1 加上调整量 1 k ( z k − x ^ k − 1 ) \frac{1}{k} (z_k - \hat{x}_{k-1}) k1(zk−x^k−1) 得到的。
- z k − x ^ k − 1 z_k - \hat{x}_{k-1} zk−x^k−1 是当前测量值与上一次估计值的差。
- 1 k \frac{1}{k} k1 是卡尔曼增益 K k K_k Kk,表示当前测量值对估计值的影响权重。
- 卡尔曼增益:
卡尔曼增益 K k K_k Kk 的作用是平衡先验估计和当前测量之间的权重。在这个简单的例子中, K k K_k Kk 随着测量次数 k k k 的增加逐渐减小,使得新的测量对估计的影响逐渐减弱。
通过这种递归方式,估计值会随着测量次数的增加逐渐收敛到真实值。随着 k k k 的增加,估计值 x ^ k \hat{x}_k x^k 会越来越稳定,变动也会越来越小。
卡尔曼滤波器在实际应用中,卡尔曼增益 K k K_k Kk 通常不是简单的 1 k \frac{1}{k} k1,而是根据测量误差和估计误差动态调整的。
多目标跟踪算法中的不同的数据融合方案
下图展示了样本1、样本2及其对应的正态分布曲线,以及根据两个不同公式计算的合并后正态分布曲线
代码如下
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation, PillowWriter
import numpy as np
# 设置中文字体
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = Falsefrom scipy.stats import normimport numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm# 定义两个样本数据
n1 = 50
n2 = 60
mu1 = 5
mu2 = 8
sigma1 = 1.5
sigma2 = 2.0# 生成两个样本
sample1 = np.random.normal(mu1, sigma1, n1)
sample2 = np.random.normal(mu2, sigma2, n2)# 计算样本的方差
sigma1_hat = np.var(sample1, ddof=1)
sigma2_hat = np.var(sample2, ddof=1)# 计算样本的均值
mu1_hat = np.mean(sample1)
mu2_hat = np.mean(sample2)# 第一个公式
sigma_new_hat_1 = (n1 * sigma1_hat + n2 * sigma2_hat + (n1 * n2 / (n1 + n2)) * (mu1_hat - mu2_hat)**2) / (n1 + n2)# 第二个公式
sigma_new_hat_2 = (n1 * sigma1_hat + n2 * sigma2_hat) / (n1 + n2)# 打印结果
print(f"样本1的方差: {sigma1_hat}")
print(f"样本2的方差: {sigma2_hat}")
print(f"使用第一个公式合并后的方差: {sigma_new_hat_1}")
print(f"使用第二个公式合并后的方差: {sigma_new_hat_2}")# 绘图
plt.figure(figsize=(10, 6))# 样本1的分布
plt.hist(sample1, bins=30, alpha=0.5, density=True, label='样本1')
# 样本2的分布
plt.hist(sample2, bins=30, alpha=0.5, density=True, label='样本2')# 添加正态分布曲线
x = np.linspace(min(np.min(sample1), np.min(sample2)), max(np.max(sample1), np.max(sample2)), 1000)
pdf1 = norm.pdf(x, mu1_hat, np.sqrt(sigma1_hat))
pdf2 = norm.pdf(x, mu2_hat, np.sqrt(sigma2_hat))
plt.plot(x, pdf1, 'blue', linewidth=2, label='样本1正态分布')
plt.plot(x, pdf2, 'orange', linewidth=2, label='样本2正态分布')# 合并后的均值
mu_combined = (n1 * mu1_hat + n2 * mu2_hat) / (n1 + n2)# 添加合并后的正态分布曲线
pdf_combined_1 = norm.pdf(x, mu_combined, np.sqrt(sigma_new_hat_1))
pdf_combined_2 = norm.pdf(x, mu_combined, np.sqrt(sigma_new_hat_2))
plt.plot(x, pdf_combined_1, 'green', linestyle='dashdot', linewidth=2, label='合并后正态分布(第一个公式)')
plt.plot(x, pdf_combined_2, 'red', linestyle='dashdot', linewidth=2, label='合并后正态分布(第二个公式)')# 添加均值和方差线
plt.axvline(mu1_hat, color='blue', linestyle='dashed', linewidth=2, label=f'样本1均值: {mu1_hat:.2f}')
plt.axvline(mu2_hat, color='orange', linestyle='dashed', linewidth=2, label=f'样本2均值: {mu2_hat:.2f}')
plt.axvline(mu_combined, color='purple', linestyle='dashed', linewidth=2, label=f'合并后均值: {mu_combined:.2f}')
plt.axvline(mu_combined + np.sqrt(sigma_new_hat_1), color='green', linestyle='dotted', linewidth=2, label=f'合并后方差(第一个公式): {sigma_new_hat_1:.2f}')
plt.axvline(mu_combined + np.sqrt(sigma_new_hat_2), color='red', linestyle='dotted', linewidth=2, label=f'合并后方差(第二个公式): {sigma_new_hat_2:.2f}')plt.legend()
plt.title('样本分布与合并后方差比较')
plt.xlabel('值')
plt.ylabel('概率密度')
plt.show()
融合两个样本的方差估计量(sample variance estimators),但它们考虑了不同的因素,因此得到了不同的结果。
让我们逐个分析这两个公式。
第一个公式
σ new 2 = n 1 σ 1 2 + n 2 σ 2 2 + ( n 1 n 2 n 1 + n 2 ) ( μ 1 − μ 2 ) 2 n 1 + n 2 \sigma_{\text{new}}^2 = \frac{n_1 \sigma_1^2 + n_2 \sigma_2^2 + \left( \frac{n_1 n_2}{n_1 + n_2} \right) (\mu_1 - \mu_2)^2}{n_1 + n_2} σnew2=n1+n2n1σ12+n2σ22+(n1+n2n1n2)(μ1−μ2)2
这个公式包含了一个额外的项 ( n 1 n 2 n 1 + n 2 ) ( μ 1 − μ 2 ) 2 \left( \frac{n_1 n_2}{n_1 + n_2} \right) (\mu_1 - \mu_2)^2 (n1+n2n1n2)(μ1−μ2)2,它考虑了两个样本均值 μ 1 \mu_1 μ1 和 μ 2 \mu_2 μ2 之间的差异。这个公式实际上是在计算两个独立样本的联合方差,适用于当两个样本来自不同总体且均值不同的情况。
额外的项 ( n 1 n 2 n 1 + n 2 ) ( μ 1 − μ 2 ) 2 \left( \frac{n_1 n_2}{n_1 + n_2} \right) (\mu_1 - \mu_2)^2 (n1+n2n1n2)(μ1−μ2)2 反映了两个样本均值之间的差异对整体方差的影响。这是因为当均值差异较大时,合并后的方差应该增大。
第一个公式中的额外项是:
( n 1 n 2 n 1 + n 2 ) ( μ 1 − μ 2 ) 2 \left( \frac{n_1 n_2}{n_1 + n_2} \right) (\mu_1 - \mu_2)^2 (n1+n2n1n2)(μ1−μ2)2
这个项反映了两个样本均值 μ 1 \mu_1 μ1 和 μ 2 \mu_2 μ2 之间的差异对合并后的方差的影响。具体来看:
- μ 1 \mu_1 μ1 和 μ 2 \mu_2 μ2 分别是两个样本的均值。
- n 1 n_1 n1 和 n 2 n_2 n2 分别是两个样本的样本量。
- σ 1 2 \sigma_1^2 σ12 和 σ 2 2 \sigma_2^2 σ22 分别是两个样本的方差估计。
公式中的这项 ( n 1 n 2 n 1 + n 2 ) ( μ 1 − μ 2 ) 2 \left( \frac{n_1 n_2}{n_1 + n_2} \right) (\mu_1 - \mu_2)^2 (n1+n2n1n2)(μ1−μ2)2 通过以下方式影响合并后的方差:
- 均值差异的平方 ( μ 1 − μ 2 ) 2 (\mu_1 - \mu_2)^2 (μ1−μ2)2:
- 如果两个样本的均值差异很大,这个项的值会很大,意味着两个样本之间的差异不能忽视,合并后的方差应该更大。
- 如果两个样本的均值很接近,这个项的值会很小,接近于零,这样合并后的方差主要由各自样本的方差决定。
- 权重因子 n 1 n 2 n 1 + n 2 \frac{n_1 n_2}{n_1 + n_2} n1+n2n1n2:
- 这个因子是一个加权平均系数,考虑了两个样本的大小。样本量越大,对方差的贡献越大。
- 当两个样本量相等时,因子等于样本量的一半。当一个样本量远大于另一个时,这个因子更接近于较小的样本量。
解释为什么会有这个额外项
这个额外项来源于两个样本均值差异导致的总体方差增加。通过加上这项,可以更准确地反映两个样本融合后的整体方差。它考虑了均值差异对方差的贡献,是一种更精确的方差合并方法。
在卡尔曼滤波中,虽然直接计算方差时没有显式使用这个公式,但卡尔曼增益的计算中也隐含了类似的思想。卡尔曼滤波通过状态预测和测量更新的步骤,将预测和测量的误差综合起来,调整状态估计和误差协方差。
第二个公式
σ new 2 = n 1 σ 1 2 + n 2 σ 2 2 n 1 + n 2 \sigma_{\text{new}}^2 = \frac{n_1 \sigma_1^2 + n_2 \sigma_2^2}{n_1 + n_2} σnew2=n1+n2n1σ12+n2σ22
这个公式则简单得多,仅仅是加权平均两个样本的方差估计。它假设两个样本来自相同总体,即它们的均值是相同的或者差异可以忽略不计。因此,在这种情况下,不需要考虑均值之间的差异对方差的影响。
- n 1 n_1 n1 和 n 2 n_2 n2 是两个样本的样本量。
- σ 1 2 \sigma_1^2 σ12 和 σ 2 2 \sigma_2^2 σ22 是两个样本的方差估计。
这个公式适用于两个样本均值相同或者均值差异很小的情况。
综述
- 第一个公式更通用,适用于两个样本来自不同总体且均值可能不同的情况,因为它包含了均值差异的项。
- 第二个公式更简化,适用于两个样本来自相同总体且均值相同或近似相同的情况,因为它忽略了均值差异的影响。
选择哪个公式取决于样本来源的假设。如果样本来自不同总体且均值可能不同,使用第一个公式。如果样本来自相同总体且均值相同或差异可以忽略,使用第二个公式。
卡尔曼增益的计算更接近于第二个公式的思想,即在预测和测量之间进行加权平均。不过,由于卡尔曼滤波中的状态更新和误差协方差更新步骤会考虑到预测误差和测量误差之间的差异,这一点与第一个公式中的额外项有所相似。总的来说,卡尔曼滤波是一种动态的、递归的估计过程,它综合了这两种融合方法的思想。
卡尔曼滤波是用于递归估计动态系统状态的经典算法。卡尔曼滤波中的卡尔曼增益(Kalman gain)在更新步骤中起关键作用,用于平衡预测与测量之间的误差。为了理解卡尔曼增益与上面的两个公式之间的关系,我们需要仔细分析卡尔曼滤波的更新公式。
在卡尔曼滤波中,估计误差协方差矩阵的更新公式如下:
- 预测步骤(Prediction Step): P pred = F P prev F T + Q \mathbf{P}_{\text{pred}} = \mathbf{F} \mathbf{P}_{\text{prev}} \mathbf{F}^T + \mathbf{Q} Ppred=FPprevFT+Q其中, P pred \mathbf{P}_{\text{pred}} Ppred 是预测误差协方差矩阵, F \mathbf{F} F 是状态转移矩阵, P prev \mathbf{P}_{\text{prev}} Pprev 是之前的误差协方差矩阵, Q \mathbf{Q} Q 是过程噪声协方差矩阵。
- 更新步骤(Update Step): K = P pred H T ( H P pred H T + R ) − 1 \mathbf{K} = \mathbf{P}_{\text{pred}} \mathbf{H}^T (\mathbf{H} \mathbf{P}_{\text{pred}} \mathbf{H}^T + \mathbf{R})^{-1} K=PpredHT(HPpredHT+R)−1其中, K \mathbf{K} K 是卡尔曼增益, H \mathbf{H} H 是测量矩阵, R \mathbf{R} R 是测量噪声协方差矩阵。
- 状态更新(State Update): x new = x pred + K ( z − H x pred ) \mathbf{x}_{\text{new}} = \mathbf{x}_{\text{pred}} + \mathbf{K} (\mathbf{z} - \mathbf{H} \mathbf{x}_{\text{pred}}) xnew=xpred+K(z−Hxpred)其中, x new \mathbf{x}_{\text{new}} xnew 是更新后的状态估计, x pred \mathbf{x}_{\text{pred}} xpred 是预测的状态估计, z \mathbf{z} z 是实际测量值。
- 误差协方差更新(Error Covariance Update): P new = ( I − K H ) P pred \mathbf{P}_{\text{new}} = (\mathbf{I} - \mathbf{K} \mathbf{H}) \mathbf{P}_{\text{pred}} Pnew=(I−KH)Ppred
在这些步骤中,卡尔曼增益 K \mathbf{K} K 的计算中并没有直接涉及两个方差融合公式。但是,我们可以从整体误差协方差矩阵的更新来间接理解卡尔曼滤波中的误差协方差矩阵与这两个公式之间的联系。
与两个方差融合公式的关系
- 与第一个公式的关系: σ new 2 = n 1 σ 1 2 + n 2 σ 2 2 + ( n 1 n 2 n 1 + n 2 ) ( μ 1 − μ 2 ) 2 n 1 + n 2 \sigma_{\text{new}}^2 = \frac{n_1 \sigma_1^2 + n_2 \sigma_2^2 + \left( \frac{n_1 n_2}{n_1 + n_2} \right) (\mu_1 - \mu_2)^2}{n_1 + n_2} σnew2=n1+n2n1σ12+n2σ22+(n1+n2n1n2)(μ1−μ2)2这个公式中的额外项反映了均值差异对方差的影响。在卡尔曼滤波中,预测步骤和更新步骤中的状态估计和测量值的差异(即 z − H x pred \mathbf{z} - \mathbf{H} \mathbf{x}_{\text{pred}} z−Hxpred)会通过卡尔曼增益 K \mathbf{K} K 调整到状态估计和误差协方差中,类似于第一个公式中的额外项。
- 与第二个公式的关系: σ new 2 = n 1 σ 1 2 + n 2 σ 2 2 n 1 + n 2 \sigma_{\text{new}}^2 = \frac{n_1 \sigma_1^2 + n_2 \sigma_2^2}{n_1 + n_2} σnew2=n1+n2n1σ12+n2σ22这个公式是简单的加权平均,假设两个样本的均值相同。在卡尔曼滤波中,卡尔曼增益 K \mathbf{K} K 的计算反映了测量噪声和过程噪声的权重,因此具有类似加权平均的效果。
从估计误差的最小化出发,来推导状态更新方程。
1. 问题描述
我们有一个预测的状态 x ^ k ∣ k − 1 \hat{x}_{k|k-1} x^k∣k−1 和一个新的观测值 z k z_k zk。我们的目标是更新状态估计,使得更新后的估计 x ^ k ∣ k \hat{x}_{k|k} x^k∣k 更加准确。
2. 估计误差
定义估计误差:
e k = x k − x ^ k ∣ k e_k = x_k - \hat{x}_{k|k} ek=xk−x^k∣k
我们的目标是最小化估计误差的协方差:
P k ∣ k = E [ ( e k ) ( e k ) T ] P_{k|k} = E[(e_k)(e_k)^T] Pk∣k=E[(ek)(ek)T]
3. 更新后的状态估计
假设更新后的状态估计是线性组合:
x ^ k ∣ k = x ^ k ∣ k − 1 + K k ( z k − H k x ^ k ∣ k − 1 ) \hat{x}_{k|k} = \hat{x}_{k|k-1} + K_k (z_k - H_k \hat{x}_{k|k-1}) x^k∣k=x^k∣k−1+Kk(zk−Hkx^k∣k−1)
其中, K k K_k Kk 是卡尔曼增益矩阵, z k − H k x ^ k ∣ k − 1 z_k - H_k \hat{x}_{k|k-1} zk−Hkx^k∣k−1 是观测残差(即预测与观测之间的差异)。
4. 误差代入
代入误差表达式:
e k = x k − ( x ^ k ∣ k − 1 + K k ( z k − H k x ^ k ∣ k − 1 ) ) e_k = x_k - (\hat{x}_{k|k-1} + K_k (z_k - H_k \hat{x}_{k|k-1})) ek=xk−(x^k∣k−1+Kk(zk−Hkx^k∣k−1))
改变下得到:
e k = ( x k − x ^ k ∣ k − 1 ) − K k ( z k − H k x ^ k ∣ k − 1 ) e_k = (x_k - \hat{x}_{k|k-1}) - K_k (z_k - H_k \hat{x}_{k|k-1}) ek=(xk−x^k∣k−1)−Kk(zk−Hkx^k∣k−1)
5. 观测模型
根据观测模型:
z k = H k x k + v k z_k = H_k x_k + v_k zk=Hkxk+vk
其中 v k v_k vk 是观测噪声,具有零均值和协方差 R k R_k Rk。
代入观测模型:
e k = ( x k − x ^ k ∣ k − 1 ) − K k ( H k x k + v k − H k x ^ k ∣ k − 1 ) e_k = (x_k - \hat{x}_{k|k-1}) - K_k (H_k x_k + v_k - H_k \hat{x}_{k|k-1}) ek=(xk−x^k∣k−1)−Kk(Hkxk+vk−Hkx^k∣k−1)
进一步简化:
e k = ( I − K k H k ) ( x k − x ^ k ∣ k − 1 ) − K k v k e_k = (I - K_k H_k)(x_k - \hat{x}_{k|k-1}) - K_k v_k ek=(I−KkHk)(xk−x^k∣k−1)−Kkvk
6. 误差协方差
计算误差协方差:
P k ∣ k = E [ e k e k T ] P_{k|k} = E[e_k e_k^T] Pk∣k=E[ekekT]
代入误差表达式:
P k ∣ k = E [ ( ( I − K k H k ) ( x k − x ^ k ∣ k − 1 ) − K k v k ) ( ( I − K k H k ) ( x k − x ^ k ∣ k − 1 ) − K k v k ) T ] P_{k|k} = E[((I - K_k H_k)(x_k - \hat{x}_{k|k-1}) - K_k v_k)((I - K_k H_k)(x_k - \hat{x}_{k|k-1}) - K_k v_k)^T] Pk∣k=E[((I−KkHk)(xk−x^k∣k−1)−Kkvk)((I−KkHk)(xk−x^k∣k−1)−Kkvk)T]
展开并利用统计特性(即 x k − x ^ k ∣ k − 1 x_k - \hat{x}_{k|k-1} xk−x^k∣k−1 和 v k v_k vk 互不相关,且 E [ v k ] = 0 E[v_k] = 0 E[vk]=0):
P k ∣ k = ( I − K k H k ) P k ∣ k − 1 ( I − K k H k ) T + K k R k K k T P_{k|k} = (I - K_k H_k) P_{k|k-1} (I - K_k H_k)^T + K_k R_k K_k^T Pk∣k=(I−KkHk)Pk∣k−1(I−KkHk)T+KkRkKkT
7. 最小化误差协方差
为了最小化误差协方差,我们选择 K k K_k Kk 使得 P k ∣ k P_{k|k} Pk∣k 尽可能小。通过数学推导和优化(这一步一般是通过求解最优 K k K_k Kk 来实现),可以得到卡尔曼增益矩阵的具体形式:
K k = P k ∣ k − 1 H k T ( H k P k ∣ k − 1 H k T + R k ) − 1 K_k = P_{k|k-1} H_k^T (H_k P_{k|k-1} H_k^T + R_k)^{-1} Kk=Pk∣k−1HkT(HkPk∣k−1HkT+Rk)−1
代入 K k K_k Kk 的具体形式,得到状态更新方程:
x ^ k ∣ k = x ^ k ∣ k − 1 + K k ( z k − H k x ^ k ∣ k − 1 ) \hat{x}_{k|k} = \hat{x}_{k|k-1} + K_k (z_k - H_k \hat{x}_{k|k-1}) x^k∣k=x^k∣k−1+Kk(zk−Hkx^k∣k−1)
8. 最后
最开始的,从平均数开始计算的式子
x ^ k = x ^ k − 1 + K k ( z k − x ^ k − 1 ) \hat{x}_k = \hat{x}_{k-1} + K_k (z_k - \hat{x}_{k-1}) x^k=x^k−1+Kk(zk−x^k−1) 与
x ^ k ∣ k = x ^ k ∣ k − 1 + K k ( z k − H k x ^ k ∣ k − 1 ) \hat{x}_{k|k} = \hat{x}_{k|k-1} + K_k (z_k - H_k \hat{x}_{k|k-1}) x^k∣k=x^k∣k−1+Kk(zk−Hkx^k∣k−1)是不是很像?