线性最小均方误差算法（LMSE），最小二乘法（LS）

目录

• 背景
• 正交投影引理
• LMSE算法
• LS算法
• 直线拟合

背景

??对于一个系统，在给予一定的输入，那么通常都会产生相对应的输出。在实际的系统中，这样的输出必然伴随着噪声，这样被噪声污染的输出通常是传感器的输出信号，也叫观测信号

??同时，如果系统的模型是清晰的，我们可以通过严格的理论计算来得到真实值，通过作差的方式变把噪声去除了。

??然而在实际的系统中，我们对整个系统的物理模型通常是未知的或者有一些参数未知，也可能是模型不准确，某些参数值有一定的偏差。因此，我们为了得到真实的信号，就需要利用观测信号估计系统的模型，尽可能地去除噪声的影响（其实是信号处理的思想）

数学描述

??$X=(x_1..x_n{)}^T$是观测信号，$f(u，\theta )$是系统的模型，$u,\theta$ 分别是系统的输入与未知参数，$n$是噪声。

??整个观测的过程可以描述为：具有的未知参数的系统在输入$u$的激励下，产生了伴随噪声$n$的输出（观测信号）$X$。

??即：$X=f(u,\theta )+n$

??而LMSE，LS算法要做的工作是利用$X$,去尽可能地得到一个$\hat{\theta }$，使$\hat{\theta }$接近真实的$\theta$，而算法的条件是整个系统模型$f(u,\theta )$得是线性的。

??即：$X=H\theta +n$,也称为线性观测模型。

以上是数学描述。

??既然是要通过已知的观测信号$X$来得到真实参数的一个进可能近似$\hat{\theta }$，那么在LMSE与LS算法中，我们是把$\hat{\theta }$看作是观测$X$的线性组合。

即：$\hat{\theta }=a+BX$,在空间里面表现为真实值在$X$上的投影。

正交投影引理

??定义：若$\theta ,X$是空间中随机矢量，那么$\theta$在$X$上的投影定义为$\theta ,X$的内积，记作$OP<\theta |X>$

引理Ⅰ：若$\theta ,X$是空间中随机矢量，$\theta$在$X$上的投影唯一。

引理Ⅱ：正交投影满足线性性

$OP[A_1{\theta }_1+A_2{\theta }_2|X]=A_{1}OP[{\theta }_1|X]+A_{2}OP[{\theta }_2|X]$

以上两个引理比较简单，第三个引理在LMSE递推算法中至关重要。

引理Ⅲ：

记$x(k)=[x(k?1)x_k{]}^T$，这里面$x(k)，x(k?1)，x_k$都是随机矢量。那么，对于随机矢量$s$，有：

$OP[s|x(k)]=OP[s|x(k?1)]+OP[\tilde{s}|{\tilde{x}}_k]$

其中，$\tilde{s}=s?OP[s|s(k?1)],{\tilde{x}}_k=x_k?OP[x_k|x(k?1)]$

因为要与随机变量的统计特征联系起来，所以可以推倒出：

$OP[\tilde{s}|{\tilde{x}}_k]=E(\tilde{s}{{\tilde{x}}_k}^T)[E({\tilde{x}}_k{{\tilde{x}}_k}^T){]}^{?1}{\tilde{x}}_k$；具体的推导过程

可以察看参考文献。

LMSE算法

先验条件;

$\theta$：均值${\mu }_{\theta }$,协方差矩阵$C_{\theta }$已知

$X$:均值$X_{\theta }$,协方差矩阵$C_X$已知

互协方差矩阵$C_{\theta X}$已知。

解析法

要想使$\hat{\theta }$与$\theta$尽可能接近，只需求解函数$E[(\theta ?\hat{\theta }{)}^T(\theta ?\hat{\theta })]$的最小值，这里面$\theta$是真实值，是一
个常数。函数可以理解为内积后取平均。

解析法思想很简单，求导数取极值即可。

令$\hat{\theta }=a+BX$

$E[(\theta ?\hat{\theta }{)}^T(\theta ?\hat{\theta })]$=$E[(\theta ?a?BX{)}^T(\theta ?a?BX)]$

分别对$a,B$求偏导=0

得到：
$a_L={\mu }_{\theta }?C_{\theta X}{C}_X^{?1}{\mu }_X$

$B_L=C_{\theta X}{C}_X^{?1}$

代入原函数中得到;

$\hat{\theta }={\mu }_{\theta }+C_{\theta X}{C}_X^{?1}(X?{\mu }_X)$

若观测与噪声独立，且噪声的统计特征已知，解析公式可以进一步简化为：

$\hat{\theta }={\mu }_{\theta }+C_{\theta }{H}^T(H{C}_{\theta }{H}^T+C_n{)}^{?1}(X?H{\mu }_{\theta }?{\mu }_n)$

迭代法

基本思想：需要迭代$k$次，且$\hat{\theta }(k)=OP[\theta (k)|x(k)]$，且$X(k)=[x(k?1)x_k{]}^T$，在这里，$x_k$是新的一个观测值。

那么，根据正交投影引理Ⅲ，有：

$OP[\theta (k)|x(k)]=OP[\theta (k?1)|x(k?1)]+OP[\tilde{\theta }(k)|\tilde{x}(k)]$

同理，根据正交投影引理Ⅲ：

$\tilde{\theta }(k)=\theta ?OP[\theta |X(k?1)],{\tilde{x}}_k=x_k?OP[x_k|X(k?1)]$ ,这一步是更新的增量，相当于把第$k$次观测中与前$k?1$次观测相关的信息去掉，留下新的信息。

进一步推倒可得到递推算法;

${\hat{\theta }}_0={\mu }_{\theta }$

$M_0=C_{\theta }$

$for1:K$

$K_k=M_{k?1}{H_k}^T(H_k{M}_{k?1}{H_k}^T+C_{n_k}{)}^{?1}$

${\hat{\theta }}_k={\hat{\theta }}_{k?1}+K_k(X_k?H_k{\hat{\theta }}_{k?1}?{\mu }_{n_k})$

$M_k=(I?K_k{H}_k)M_{k?1}$

$end$

LS算法

LS算法的特点是不需要各类参数以及噪声的先验信息，相当于直接把真实值看成参数，并且直接在$X$上面进行投影。然后误差均方最小。

解析法

解析法与LSME算法一样，只是目标函数变成了$(X?H\hat{\theta }{)}^T(X?H\hat{\theta })$

求导=0得:

$\hat{\theta }=(H^{T}H{)}^{?1}{H}^{T}X$

迭代法

令系统输出的真实值为$X_r,$，观测值为$X$

递推原理如下:

$OP[X_r|X(k)]=OP[X_r|X(k?1)]+OP[M|N]$

$M=X_r?OP[X_r|X(k?1)]$

$N=x_k?OP[x_k|X(k?1)]$

直线拟合

%%%本程序基于LMSE算法，利用20个点的数据拟合直线。%%产生数据
clear;clc;
N=50;
x=1:1:N; 
k=1;
b=0.5;n=normrnd(0.6,1.5,[N,1]);  %N(0.6,0.3) 20*1的噪声
y=k*x'+b+n;%%%把数据转换为观测方程
X=y;
AT=1:N;
%H=[X,ones(N,1)];
%问题一 H矩阵是模型矩阵,里面不可能有观测量X
H=[AT',ones(N,1)];%%给定初始值
u_theta=[0.5;0.8];     %theta均值矩阵
C_theta=[0.1,0;0,0.4]';%theta协方差矩阵u_n=ones(N,1);           %n的均值矩阵
C_n=diag(0.3*ones(N,1)); %n的协方差矩阵%%LMSE计算解析解：
%假设观测与噪声独立的%theta=u_theta+C_theta*H'/(H*C_theta*H'+C_n)*(X-H*u_theta-u_n);
%问题二 多减了u_n
theta=u_theta+C_theta*H'/(H*C_theta*H'+C_n)*(X-H*u_theta);%%%下面展示利用递推算法;
%%给定初始条件M=C_theta;
theta_1=u_theta;
%以上相当于迭代第一次for i=2:N%K=M*H(1:i,:)'/(H(1:i,:)*M*H(1:i,:)'+C_n(1:i,1:i));%计算K矩阵%问题三 递推过程中，H矩阵一直都是1*2的，元素值是随采样时刻变化的%C_n是2*1的H_1=[i;1]';%1行2列的递推H矩阵K=M*H_1'/(H_1*M*H_1'+C_n(i,i));%2行1列的矩阵theta_1=theta_1+K*(X(i,1)-H_1*theta_1);%更新theta矩阵%更新M矩阵[k,l]=size(M);M=(eye(k)-K*H_1)*M;end%%下面利用最小二乘法
theta_2=(H'*H)\H'*X;%%作图对比
figure(1);
plot(x',y,'*');  %观测点hold on;
x=1:0.1:N;
H_=[x',ones(length(x),1)];
y_=H_*theta;plot(x',y_,'r');   %LMSE解析法hold on;y_1=H_*theta_1;
plot(x',y_1,'g');  %LMSE迭代法 hold on;
y_2=H_*theta_2;
plot(x',y_2,'b');  %LS算法legend('测量曲线','LMSE解析法','LMSE迭代法','LS');
hold off

效果：

三种算法拟合效果几乎一样。

参考文献：

赵树杰, 赵建勋. 信号检测与估计理论[M]. 电子工业出版社, 2013.