非随机参数估计的克拉美罗界（Cramer-Rao Bounds）

最新推荐文章于 2025-03-20 21:55:23 发布

GHelpU

最新推荐文章于 2025-03-20 21:55:23 发布

阅读量8.1k

点赞数 16

分类专栏：信号处理文章标签： 1024程序员节 CRB 克拉美罗界

本文链接：https://blog.csdn.net/weixin_42647783/article/details/109238707

版权

信号处理专栏收录该内容

10 篇文章

订阅专栏

更新日志：

2020.10.24：发布

1 简介

克拉美罗界用以研究参数估计问题的一套理论。参数估计是现代信号处理领域的一个重头戏，人为观测到的信号往往是物理模型产生的，为了研究信号的产生过程，我们会设计各种各样的数学模型以近似物理模型。这些数学模型往往是由模型的结构和参数组成，若给定了模型的结构，而模型的参数还未知的话，就可以利用观测数据通过各种算法估计这些参数（这就叫后验）。这些估计参数的方法就很多了，可以分成两大类，一类是有偏估计（估计值的期望与待估参数的真值之间有系统误差），一类是无偏估计（估计值的期望与待估参数的真值之间无系统误差）。克拉美罗界只适用于无偏估计，它认为，任何无偏估计量的方差不会低于克拉美罗界，也就是说，设计的参数估计算法再强，给出的无偏估计量的方差不可能低于下界，因此也称其为Cramer-Rao Lower Bounds（CRLB），意为克拉美罗下界。

本篇博文主要研究非随机参数估计的克拉美罗下界，非随机参数（nonrandom parameter）是指参数的值是固定的，不随观测的变化而变化。

一个模型只要给定，估计参数的下界就固定了，所以我们在设计算法时候，设计的结果往往会通过多次试验计算估计的方差，看与下界的距离有多少，就知道算法好不好了。当然，如果是有偏估计的话，那方差完全有可能是低于下界的，但是有偏估计一般就没什么意义了，比如待估计的参数 $a = 5$ ，然后我每次估计时候都认为是8，那方差就是0了，但是这样的估计就失去意义了，所以这里只考虑无偏估计。

举个简单的例子，现在有一个电压真值为5V的接头，我们希望通过电压表测量这个电压（也就是待估计的非随机参数），但是电压表不准确，测量值老是跳，结果等于这个真值5V加上零均值的高斯白噪声。为了估计电压，我们采用这样的方法：每测得一个值，就认为这是接头的电压了。因此，比如总共测了100次，可以得到100个估计值，这100个估计值可能都接近于5V（如果测量次数接近于无穷的话，那么这些估计值的平均是等于5V的，因为这是个无偏估计方法），但是估计值都各不相等，那么这100个值的方差必然大于克拉美罗界。

那么到底怎么定义克拉美罗界呢？如下。（注意，本博文中加粗的符号表示向量或者矩阵，不加粗的表示标量）

2 克拉美罗界的定义

用向量 $\bf{x}$ 表示待估计的参数， $\bf{z}$ 表示对参数的测量， $\hat{\bf{x}}(\bf{z})$ 表示根据测量给出的参数估计结果，叫估计子。

无偏估计子 $\hat{\bf{x}}(\bf{z})$ 给出的结果往往是围绕着参数真值 $\bf{x}$ 波动，为了描述波动的程度，定义误差的方差，即平均平方误差（mean-square error matrix）为
${\bf \Lambda}_{\epsilon}={\rm MSE}(\hat{\bf{x}}({\bf{z}})) \triangleq {\bf E}\{[\hat{\bf x}({\bf z})-{\bf x}][\hat{\bf x}({\bf z})-{\bf x}] ^T\}$

用 $p({\bf z}|{\bf x})$ 表示测量的概率密度函数，定义尺寸为 $\times K$ 的Fisher信息矩阵（Fisher’s information matrix ，FIM）的矩阵 ${\bf J}_F ({\bf x})$ 如下：
$\begin{aligned} {\bf J}_F ({\bf x}) &\triangleq {\bf E}\{[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})][\triangledown_{\bf x}\ln p({\bf z}|{\bf x})] ^T\} \\ &= -{\bf E}\{\triangledown_{\bf x}[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})] ^T\} \end{aligned}$

其中矩阵 ${\bf J}_F ({\bf x})$ 内部的每一个元素可以写成
$\begin{aligned} {J}_{F_{ij}} ({\bf x}) &\triangleq {\bf E}\{\frac{\partial \ln p({\bf z}|{\bf x})}{\partial x_i}\frac{\partial \ln p({\bf z}|{\bf x})}{\partial x_j}\} \\ &= -{\bf E}\{\frac{\partial^2 \ln p({\bf z}|{\bf x})}{\partial x_i \partial x_j}\} \end{aligned}$

那么克拉美罗界具有如下性质：
${\bf \Lambda}_{\epsilon}-{\bf J}_F^{-1} ({\bf x})\geqslant {\bf 0}$

这里 $\geqslant$ 表示左边的矩阵式子是非负定的。

另外，这一不等式也表明了，左边的矩阵式子所有的子矩阵都是非负定的，其中对角线上的第 $i$ 元素满足
$\sigma_{\epsilon_i}^2 \triangleq {\rm Var}[\hat{x}_i({\bf z})-x_i]\geqslant J_F^{ii}({\bf x})$

当且仅当
$\hat{{\bf x}}({\bf z})-{\bf x}={\bf k}({\bf x})[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]$

其中 ${\bf k}({\bf x})={\bf J}_F^{-1}({\bf x})$

上面两个不等式子的等号成立。

3 示例

3.1 理论案例

仍然以简介中电压测量为例，接头电压固定为 $x$ ，即为待估计的量，而电压测量时，测量受到均值为 $\mu=0$ 、方差为 $\sigma^2$ 的加性高斯白噪声 $n$ 的影响，即
$z = x + n$

因此，测量的概率密度函数表示为
$p(z|x)=\frac{1}{\sqrt{2\pi}\sigma} \exp\left[-\frac{(z-x)^2}{2\sigma^2}\right]$

那么该模型的Fisher信息矩阵（在本例中其实是标量不算是矩阵，因此下面没有用粗体）为
$J_F(x)=-\frac{\partial^2 \ln [p(z|x)]}{\partial x^2}=\frac{\partial^2}{\partial x^2}\frac{(z-x)^2}{2\sigma^2}=\frac{1}{\sigma^2}$

那么克拉美罗下界为
${\Lambda}_{\epsilon} = {\rm MSE}{(\hat{x}({z}))}\geqslant \sigma^2$

也就是说，设计的算法只要是无偏的，估计子的方差就一定比 $\sigma^2$ 大。

3.2 实验验证

若我们的估计子取为观测，即 $\hat{x}({z})=z$ 。

那么接下来看看该估计子的均值：
$E[\hat{x}(z)]=\int_{-\infty}^{\infty}zp(z|x)dz=x$

说明该估计子是无偏的，然后计算该估计子的方差
${\rm Var}[\hat{x}(z)]=\int_{-\infty}^{\infty}(z-x)^2p(z|x)dz=\sigma^2$

估计子的方差恰好等于CRB，说明这种估计方法是最好的。

做个简单的实验验证一下，接头电压为5V，然后我们选取噪声的方差分别为 0.5 ~ 5 ，每一个噪声参数，都做10000次重复实验，然后计算这10000次实验给出的估计子的方差，再与克拉美罗下界对比。仿真代码如下，

// A Matlab Code
clear; close all; clc;
x = 5;              % 待估计的电压 5V
Times = 10000;      % 每个噪声模型下,都做10000次仿真
miu = 0;            % 噪声的均值为0
index = 1;          % 索引
for sigma = 0.5:0.5:5
    x_hat = zeros(Times, 1);
    for t = 1:Times
        n = normrnd(miu, sigma, 1, 1);  % 仿真高斯白噪声
        z = x + n;                      % 观测结果
        x_hat(t, 1) = z;                % 认为观测的值即为电压值
    end
    MSE(index, 1) = var(x_hat);         % 计算每个方差为sigma的噪声对应的估计子的方差,即为MSE
    CRB(index, 1) = sigma^2;            % 计算每个方差为sigma的噪声对应的克拉美罗下界
    index = index + 1;
end

figure;hold on;
plot(0.5:0.5:5, MSE, 'LineWidth',2, 'Marker', 'o');
plot(0.5:0.5:5, CRB, 'LineWidth',2, 'color', [119, 172, 48]/255);
legend('MSE', 'CRB');
title('估计子的方差与克拉美罗下界的对比');

仿真结果如下图，

实验结果非常直观了，MSE非常接近CRB，与我们的理论是一致的，当然一些偏差是合理的，因为只是进行了10000次实验，难免会有一点点小偏差，如果是无穷次实验，那就绝对一致了。

这是只观测一次的情况下，最优的办法了，所以平时我们基本电压表测出怎样的结果，就认为实际的电压是啥样子的。如果读者有个疑问，我们不是会连续测很多次，然后取平均作为最终的估计值，那这个时候做很多很多次试验，只要每次试验的测量次数足够多，那这些试验的估计子的方差不就非常接近于零了，那怎么会大于CRB呢？其实不然，上面我们算的CRB都是标量的，就是只得到一次结果便给出估计，如果测量次数非常多的话，那CRB就是不是上面这么3.1节这么算的了，读者可以自行尝试一下，也欢迎读者在评论区给出答案哈~

下面的证明如果读者不感冒就可以跳过了，如果想进一步深入研究克拉美罗界的理论的话，建议好好品尝一下。

4 证明

4.1 不等式的证明

因为 $\hat{x}_i({\bf z})$ 是无偏的，所以有
$\int_{-\infty}^{\infty}[\hat{x}_i({\bf z})-x_i]p({\bf z}|{\bf x})d{\bf z}=0$

对 $x_j$ 求偏导，有
$-\delta_{ij}\int_{-\infty}^{\infty}p({\bf z}|{\bf x})d{\bf z}+\int_{-\infty}^{\infty}[\hat{x}_i({\bf z})-x_i]\frac{\partial p({\bf z}|{\bf x})}{\partial x_j}d{\bf z}=0$

其中 $\delta_{ij}$ 当 $i = j$ 时取值1，其他情况取值0。显然上式的第一个积分等于1，而第二个积分式子为
$\frac{\partial p({\bf z}|{\bf x})}{\partial x_j}=\frac{\partial \ln p({\bf z}|{\bf x})}{\partial x_j}p({\bf z}|{\bf x})$

因此可以重新写为
$\int_{-\infty}^{\infty}\left\{[\hat{x}_i({\bf z})-x_i]\sqrt{p({\bf z}|{\bf x})}\right\}\left\{ \sqrt{p({\bf z}|{\bf x})}\frac{\partial \ln p({\bf z}|{\bf x})}{\partial x_j}\right\}d{\bf z}=\delta_{ij}$

或者写成矩阵的形式为，其中 ${\bf I}$ 表示单位阵
$\int_{-\infty}^{\infty}\left\{[\hat{\bf x}({\bf z})-{\bf x}]\sqrt{p({\bf z}|{\bf x})}\right\}\left\{ \sqrt{p({\bf z}|{\bf x})}[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]^T\right\}d{\bf z}={\bf I}$

为了用上Schwarz不等式，定义两个任意的 $\times 1$ 大小的向量 ${\bf b}_1$ 和 ${\bf b}_2$ ，且这两个向量与测量 ${\bf z}$ 无关，现在，对上面的等式左乘 ${\bf b}_1^T$ ，右乘 ${\bf b}_2$ ，得到
$\int_{-\infty}^{\infty}\left\{{\bf b}_1^T[\hat{\bf x}({\bf z})-{\bf x}]\sqrt{p({\bf z}|{\bf x})}\right\}\left\{ \sqrt{p({\bf z}|{\bf x})}[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]^T{\bf b}_2\right\}d{\bf z}={\bf b}_1^T{\bf b}_2$

这样以来，上面式子花括号中的量都变成标量了，利用Schwarz不等式，可以得到
$({\bf b}_1^T{\bf b}_2)^2 \leqslant \int_{-\infty}^{\infty}{\bf b}_1^T[\hat{\bf x}({\bf z})-{\bf x}][\hat{\bf x}({\bf z})-{\bf x}]^T {\bf b}_1 {p({\bf z}|{\bf x})} d{\bf z}\\ \cdot \int_{-\infty}^{\infty} {\bf b}_2^T [\triangledown_{\bf x}\ln p({\bf z}|{\bf x})][\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]^T{\bf b}_2 {p({\bf z}|{\bf x})}d{\bf z}$

等号成立的条件是当且仅当
${\bf b}_1^T[\hat{\bf x}({\bf z})-{\bf x}]=C({\bf x})[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]^T{\bf b}_2$

其中 $C({\bf x})$ 表示一个只与 ${\bf x}$ 有关、与 ${\bf z}$ 无关的量。如果把这里的两个积分公式分别换成 ${\bf \Lambda}_{\epsilon}$ 和 ${\bf J}_F({\bf x})$ ，那么这个式子就可以简化为
$({\bf b}_1^T{\bf b}_2)^2 \leqslant [{\bf b}_1^T{\bf \Lambda}_{\epsilon} {\bf b}_1] \cdot [{\bf b}_2^T {\bf J}_F({\bf x}){\bf b}_2]$

考虑到以上式子中 ${\bf b}_1$ 和 ${\bf b}_2$ 都是任意选取的，那么可以令
${\bf b}_2={\bf J}_F^{-1}({\bf x}){\bf b}_1$

这样上一个不等式就可以写成
$[{\bf b}_1^T{\bf \Lambda}_{\epsilon} {\bf b}_1] \cdot [{\bf b}_1^T {\bf J}_F^{-1}({\bf A}){\bf b}_1] \geqslant [{\bf b}_1^T {\bf J}_F^{-1}({\bf x}){\bf b}_1] ^2$

于是就可以简化为
${\bf b}_1^T{\bf \Lambda}_{\epsilon} {\bf b}_1 \geqslant {\bf b}_1^T {\bf J}_F^{-1}({\bf x}){\bf b}_1$

或者写成
${\bf b}_1^T[{\bf \Lambda}_{\epsilon}-{\bf J}_F^{-1}({\bf x})] {\bf b}_1 \geqslant 0$

由于前面已经说明了 ${\bf b}_1$ 的任意性，在任意选取 ${\bf b}_1$ 时上面这一不等式都必然成立，因此
${\bf \Lambda}_{\epsilon}-{\bf J}_F^{-1}({\bf x}) \geqslant {\bf 0}$

这里的 $\geqslant$ 仍然表示左边的矩阵式子是非负定的。

4.2 等号成立条件

再次利用 4.1 小节中的两个式子（直接摘自上小节），第二个式子为等号成立的条件
${\bf b}_2={\bf J}_F^{-1}({\bf x}){\bf b}_1$

${\bf b}_1^T[\hat{\bf x}({\bf z})-{\bf x}]=C({\bf x})[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]^T{\bf b}_2$

把第一个式子代入到第二个式子，则等号成立的条件变为
${\bf b}_1^T[\hat{\bf x}({\bf z})-{\bf x}]=C({\bf x})[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]^T{\bf J}_F^{-1}({\bf x}){\bf b}_1$

考虑到 ${\bf b}_1$ 的任意性，则可以进一步简化为
$[\hat{\bf x}({\bf z})-{\bf x}]^T=C({\bf x})[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]^T{\bf J}_F^{-1}({\bf x})$

重写该式子为
$[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]^T=[\hat{\bf x}({\bf z})-{\bf x}]^T\frac{{\bf J}_F({\bf x})}{C({\bf x})}$

式子两边同时对 ${\bf x}$ 求偏导，得到
$\triangledown_{\bf x}[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]^T=-{\bf I}\frac{{\bf J}_F({\bf x})}{C({\bf x})}+[\hat{\bf x}({\bf z})-{\bf x}]^T\triangledown_{\bf x}\left[ \frac{{\bf J}_F({\bf x})}{C({\bf x})} \right]$

这样，对等式两边同时乘以 -1，然后同时取期望，由于估计量是无偏的，则等号右边的第二项的期望为 0 ，那么可以得到
${\bf J}_F({\bf x})=\frac{1}{C({\bf x})}{\bf J}_F({\bf x})$

因此等号成立的条件对应于 $C({\bf x})=1$ ，于是等号成立的条件等价于
$\hat{\bf x}({\bf z})-{\bf x}={\bf J}_F^{-1}({\bf x})[\triangledown_{\bf x}\ln p({\bf z}|{\bf x})]$

或者写成
$\triangledown_{\bf x}\ln p({\bf z}|{\bf x})={\bf J}_F({\bf x})[\hat{\bf x}({\bf z})-{\bf x}]$

这就是克拉美罗界不等式中等号成立的条件了，当无偏估计子满足这一条件时，该估计子的方差达到克拉美罗界，便也是最好的了。

4.3 补充性质

有些时候，我们需要估计的是以上所述的 $\times 1$ 大小的参数向量 ${\bf x}$ 的函数，而并非 ${\bf x}$ 本身。假设要估计的量是 ${\bf d}=\bf \gamma(\bf x)$ ，这里 ${\bf d}=[d_1,d_2,...,d_M]^T$ 。 $M$ 和 $K$ 并没有关系，方程可能是非线性的。估计误差表示为
$\hat{\bf d}({\bf z})-\gamma(\bf x) \triangleq {\bf d}_{\epsilon}$

如果假设估计是无偏的，估计误差的协方差矩阵用 ${\bf \Lambda}_{\epsilon}(\bf d)$ 表示，利用上述性质，可以证明出如下几条性质：

性质1：估计参数d的克拉美罗界
${\bf \Lambda}_{\epsilon}(\bf d) -\left \{\triangledown_{\bf x}[\gamma^T({\bf x})]\right\}^T{\bf J}_F^{-1}({\bf x})\left \{\triangledown_{\bf x}[\gamma^T({\bf x})]\right\}\geqslant {\bf 0}$
性质2：对上面矩阵形式克拉美罗界中单独元素的拆解
${\rm Var}(d_{{\epsilon}_i}) \geqslant \sum_{k=1}^K \sum_{j=1}^K \frac{\partial \gamma_i({\bf x})}{\partial x_k} J_F^{kj}({\bf x})\frac{\partial \gamma_i({\bf x})}{\partial x_j}$
性质3：上述性质在线性模型下可以进一步简化，若有以下线性关系
$\gamma({\bf x})\triangleq {\bf \Gamma}{\bf x}$

其中 ${\bf \Gamma}$ 是一个 $\times K$ 的矩阵。若估计量是无偏的，则有
${\bf \Lambda}_{\epsilon}(\bf d) -{\bf \Gamma}^T{\bf J}_F^{-1}({\bf x}){\bf \Gamma}\geqslant {\bf 0}$