辨识实例:最小二乘法拟合圆的标准方程

辨识实例:最小二乘法拟合圆的标准方程

以下以拟合圆的方程来演示最小二乘二阶辨识的实例,已知圆的通式为
( x − a ) 2 + ( y − b ) 2 = r 2 \begin{equation} \begin{split} (x-a)^2 + (y-b)^2 = r^2 \end{split} \end{equation} (xa)2+(yb)2=r2
由于最小二乘法不能直接辨识二阶以上的方程,所以我们需要对该方程进行处理
( x − a ) 2 + ( y − b ) 2 = r 2 x 2 − 2 a x + a 2 + y 2 − 2 b y + b 2 = r 2 x 2 + y 2 − 2 a x − 2 b y + a 2 + b 2 − r 2 = 0 \begin{equation} \begin{split} (x-a)^2 + (y-b)^2 &= r^2 \\ x^2 - 2ax + a^2 + y^2 - 2by + b^2 &= r^2 \\ x^2 + y^2 - 2ax - 2by + a^2+ b^2 - r^2 &= 0 \end{split} \end{equation} (xa)2+(yb)2x22ax+a2+y22by+b2x2+y22ax2by+a2+b2r2=r2=r2=0
上式中 a 2 a^2 a2 b 2 b^2 b2以及 r 2 r^2 r2在应用时都是常数项,所以可以将其当作一个常数项 c c c来看待。待辨识方程即为
f ( a , b , c ) = x 2 + y 2 − 2 a x − 2 b y + c \begin{equation} \begin{split} f(a,b,c) = x^2 + y^2 - 2ax - 2by + c \end{split} \end{equation} f(a,b,c)=x2+y22ax2by+c
这里将上式化简为最小二乘的矩阵求解形式如下
f ( θ ) = Z m = [ − 2 x − 2 y E ] [ a b c ] = H m θ ^ \begin{equation} \begin{split} f(\theta) = Z_m &= \begin{bmatrix} -2x & -2y & E \end{bmatrix} \begin{bmatrix} a \\ b \\ c \end{bmatrix} \\ &=H_m \hat \theta \end{split} \end{equation} f(θ)=Zm=[2x2yE] abc =Hmθ^
式中 E E E为单位矩阵,其他字母均为矩阵,需要保证维度一致。
使用最小二乘的求解公式对其进行求解
θ ^ = [ a ^ b ^ c ^ ] = H m − 1 Z m \begin{equation} \begin{split} \hat \theta = \begin{bmatrix} \hat a \\ \hat b \\ \hat c \end{bmatrix}= H_m^{-1}Z_m \end{split} \end{equation} θ^= a^b^c^ =Hm1Zm
求解得到的 a ^ \hat a a^ b ^ \hat b b^即是辨识得到的圆的圆心坐标 ( x , y ) (x,y) (x,y),而半径 r r r则由于我们对其进行了合并处理,所以需要进一步处理才能得到,如下
a 2 + b 2 + r 2 = c r 2 = c − a 2 − b 2 r = c − a 2 − b 2 \begin{equation} \begin{split} a^2 + b^2 + r^2 &= c \\ r^2 &= c - a^2 - b^2 \\ r &= \sqrt{c - a^2 - b^2} \end{split} \end{equation} a2+b2+r2r2r=c=ca2b2=ca2b2
这里以圆心为 ( 0 , 0 ) (0,0) (0,0),半径 r = 3 r = 3 r=3m作为圆的真实参数。以下是该实例对应的matlab代码

% 生成没有噪音的圆的轨迹点
t = linspace(0, 2*pi, 100);
r = 3;
x = r*cos(t);
y = r*sin(t);
figure;
plot(x,y); hold on;
axis equal
% 采集带噪音的真实坐标数据
z_x = x + 0.2*randn(1,length(t));
z_y = y + 0.2*randn(1,length(t));
plot(z_x,z_y,'r.');
% 构造最小二乘形式
zm = -(z_x.^2 + z_y.^2)';
H = [-2*z_x', -2*z_y', ones(1,length(z_x))'];
% 最小二乘求解
hat = H\zm;
x_hat = hat(1);
y_hat = hat(2);
r_hat = sqrt(-hat(3)-hat(1)^2-hat(2)^2);
fprintf("圆心位置为:(%f, %f)\n",x_hat,y_hat);
fprintf("圆心半径为:%f m\n",r_hat);
%画图验证
x_fit = x_hat + r_hat*cos(t);
y_fit = y_hat + r_hat*sin(t);
plot(x_fit, y_fit,'r-');
legend('真实轨迹','采集数据','拟合轨迹','Location','best')

在这里插入图片描述在这里插入图片描述

从结果我们可以看出

  • 最小二乘并不是无偏的,它是将所有的误差平摊到每个采集的数据点上,它只能起到类似滤波器的效果,对于那些对精度要求较高的参数(小数点后两位以上),使用最小二乘较难满足精度要求
  • 最小二乘是最简单的辨识方法,它是基于假设噪音为高斯(正态)分布时的白噪音时才能获得较好的效果,如果不是白噪音,则可能需要使用到卡尔曼滤波、神经网络辨识、粒子群辨识、粒子滤波等高级辨识方法
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值