估值分析作业二

转自：https://blog.csdn.net/gangdanerya/article/details/105105611

问题：

要求：使用扩展卡尔曼滤波（EKF）的方法解决

 

卡尔曼滤波器只适用于线性系统模型，然而实际中的系统往往都是非线性模型。所以必须对卡尔曼滤波器进行修改。

EKF基于线性化系统的思想，将系统函数的非线性函数作一阶Taylor展开，得到线性化系统方程。扩展的卡尔曼滤波器算法就是适用于非线性系统的卡尔曼滤波器。它与经典的线性卡尔曼滤波器很相似，算法步骤和结构都相同。不同在于系统模型和矩阵A和H。.

 

 

 

代码：

N = 300;  %连续计算N个时刻

Q = 10;   %过程方差

R = 1;          %测量方差

v=sqrt(Q)*randn(1,N); %过程标准差

n=sqrt(R)*randn(1,N);     %测量标准差

P =eye(1);

 

x=zeros(1,N);%真实值

Xresult=zeros(1,N); %最优估计值

 

x(1,1)=0;

Xresult (1,1)=x(1,1);

z=zeros(1,N);           %状态测量值

z(1)=x(1,1)^2/20+v(1);

z1=zeros(1,N);        %测量值的估计

z1 (1,1)=z(1);

 

for k = 2 : N

    %模拟系统

    x(:,k) = 0.5 * x(:,k-1) + (25 * x(:,k-1) / (1 + x(:,k-1).^2)) + 8 * cos(1.2*k) + v(k-1); 

    z(k) = x(:,k).^2 / 20 + n(k);

    %EKF开始

    Xmid= 0.5* Xresult (:,k-1)+ 25* Xresult (:,k-1)/(1+ Xresult (:,k-1).^2) + 8 * cos(1.2*k);

    %预测值的估计

    z1= Xmid.^2/20; 

    %雅克比矩阵

    F = 0.5 + 25 * (1- Xresult.^2)/((1+ Xresult.^2).^2);

    H = Xmid /10;

    PP=F*P*F'+Q;   %过程方差预测

    Kk=PP*H'*inv(H*PP*H'+R);  %inv返回逆矩阵

    Xresult (k)= Xmid +Kk*(z(k)- z1); %修正值

    P=PP-Kk*H*PP;       %EKF方差

end

%画图

   t = 1 : N;

   figure;

   plot(t,x(1,t),'c',t, Xresult (1,t),'b');

   legend('真实值','EKF估计值');

 

 

运行结果：

当N=30时：

 

 

 

当N=300时：

 

 

 

观测最后算法运行得到的结果，发现当时间规模较小的时候，EKF的估计值和真实值相差较小，而当时间规模较大时，通过EKF拟合而来的曲线则与真实值偏差较大。