MSE（均方误差）计算封装Matlab函数

背景

本人最近需要写多个仿真，需要大量用到MSE（均方误差）计算，于是干脆将MSE运算封装为函数，后续使用直接进行调用即可。

函数代码

%Project: 均方误差函数
%Author: Jace
%Data: 2021/11/01
%====================函数体====================
function [MSE]=MSE(Dim,Step,N,xkf,x)
    %====================分配空间========================
    MSE=zeros(Dim,N);
    MSEC=zeros(1,Dim);
    %====================迭代过程====================
    for n=1:Dim
        %--------维度循环--------
        for k=Step:N
            %--------时刻循环--------
            MSEC(n)=0;%误差平方和变量清零
            for i = k-(Step-1):k
                %--------加和循环--------
                MSEC(n)=MSEC(n)+(xkf(n,i)-x(n,i))^2;%误差平方和       
            end
            MSE(n,k)=MSEC(n)/Step;%取平均
        end
    end
end

调用方法

[MSE]=MSE(Dim,Step,N,xkf,x)

对应

[MSE矩阵]=MSE(状态维度Dim,MSE长度Step,总时长N,估计矩阵xkf,状态矩阵x)

注意：

1. 维度默认为状态维度，可以直接计算出每个状态值估计与真实的MSE，得到的MSE矩阵为Dim_nN维。如果不需要完全计算，可以自定义Dim数值。如Dim=1，则只计算第一个状态值的MSE，相应算得的MSE也只有1N维；
2. 输入的估计矩阵xkf和状态矩阵x都是估计算法迭代计算之后的结果矩阵，维度应该是Dim_n*N维；
3. 由于前Step长度不足计算，因此得到的数值为0。

调用测试函数

%Project: 基本二维Kalman测试函数
%Author: Jace
%Data: 2021/11/02
%====================准备====================
close all;
clear all;
clc;
%====================设定全局参数====================
%--------全局参数--------
N=100;%设定采样点数，即持续时长
%--------设定维度--------
Dim_n=2;%状态维度
Dim_m=2;%量测维度
%--------系统模型参数--------
A=[1.002,0;0,0.998];%状态转移矩阵
H=[1,0;0,1];%局部量测1量测矩阵
Gamma=1;

%--------噪声相关参数--------
P0=0.01;%初始状态噪声协方差

Q=0.01*eye(Dim_n);%设定系统噪声
R=0.1*eye(Dim_m);%设定观测噪声
w=sqrt(Q)*randn(Dim_n,N);
v=sqrt(R)*randn(Dim_m,N);

%====================分配空间========================
%--------系统参数--------
x=zeros(Dim_n,N);
z=zeros(Dim_m,N);%量测值

%--------Kalman过程参数--------
p=zeros(Dim_n,Dim_n,N);
xkf=zeros(Dim_n,N);%估计状态

%====================初始化====================
%--------系统参数初始化--------
x(:,1)=[10+P0*randn(1);20+P0*randn(1)];%物体初始真实状态值
z(:,1)=H*x(:,1)+v(:,1);%观测真实值初始值

%--------估计参数初始化--------
p(:,:,1)=0.1*eye(Dim_n);%误差协方差初始值
xkf(:,1)=x(:,1);%全局估计状态初始化

%--------特定矩阵初始化--------
In=eye(Dim_n);%2*2单位矩阵

%====================迭代过程====================
for k=2:N
    
    %系统模型
    x(:,k)=A*x(:,k-1)+w(:,k);
    %量测模型，标量
    z(:,k)=H*x(:,k)+v(:,k);
    %====================标准Kalman过程====================
    [p(:,:,k),xkf(:,k)]=Lkf(2,A,H,Gamma,Q,R,p(:,:,k-1),xkf(:,k-1),z(:,k));
end
%====================MSE计算====================
Step=10;
[MSE]=MSE(Dim_n,Step,N,xkf,x);
%====================绘图====================
%MSE
figure;
hold on,box on;
plot(MSE(1,:),'-r.');
plot(MSE(2,:),'-g.');
legend('MSE1','MSE2');
xlabel('采样时间');ylabel('数值');
title('MSE');