👨🎓 博主简介:博士研究生
🔬 超级学长:超级学长@实验室(提供各种程序开发、实验复现与论文指导)
📧 个人邮箱:easy_optics@126.com
💬 个人微信:easy_optics
🐧 个人企鹅:754357517
摘要
主要介绍了MATLAB中的曲线拟合方法,涵盖多项式拟合、加权最小方差拟合及非线性曲线拟合。在多项式拟合中,函数polyfit()可通过最小二乘法找到合适多项式系数,不同阶次拟合效果不同,阶次最高不超length(x)-1。加权最小方差拟合根据数据准确度赋予不同加权值,更符合拟合初衷,文中还给出其原理及求解公式,并通过实例展示拟合结果。对于非线性曲线拟合,已知输入输出向量及函数关系但未知系数向量时,可利用lsqcurvefit函数求解,同时介绍了该函数多种调用格式,最后通过具体实例阐述其应用及结果。
一、引言
在科学和工程领域, 曲线拟合的主要功能是寻求平滑的曲线来最好地表现带有噪声的测量数据, 从这些测量数据中寻求两个函数变量之间的关系或者变化趋势, 最后得到曲线拟合的函数表达式y =f (x ) 。
使用多项式进行数据拟合会出现数据振荡, 而使用 Spline 插值的方法可以得到很好的平滑效果, 但是关于该插值方法有太多的参数, 不适合曲线拟合的方法。
同时, 由于在进行曲线拟合的时候, 已经认为所有测量数据中包含噪声, 因此, 最后的拟合曲线并不要求通过每一个已知数据点, 衡量拟合数据的标准则是整体数据拟合的误差最小。
一般情况下, MATLAB 的曲线拟合方法用的是“最小方差”函数,其中方差的数值是拟合曲线和已知数据之间的垂直距离。
二、 多项式拟合
在 MATLAB 中, 函数polyfit() 采用最小二乘法对给定的数据进行多项式拟合, 得到该多项式的系数。 该函数的调用方式如下:
■ polyfit(x,y,n) ———找到次数为n 的多项式系数, 对于数据集合{ (xi, yi) } , 满足差的平方和最小。
■ [p,E]=polyfit(x,y,n) ———返回同上的多项式p 和矩阵E。 多项式系数在向量p 中, 矩阵 E 用在polyval 函数中来计算误差。
运行程序后, 得到的输出结果如图1所示。由图可以看出, 使用5 阶多项式拟合时, 得到的结果比较差。
当采用9阶多项式拟合时,得到的结果与原始数据符合得比较好。当使用函数polyfit()进行拟合时,多项式的阶次最大不超过length(x)-1。
三、加权最小方差拟合原理及实例
所谓加权最小方差(WLS) , 就是根据基础数据本身各自的准确度的不同, 在拟合的时候给每个数据以不同的加权数值。 这种方法比前面所介绍的单纯最小方差方法要更加符合拟合的初衷。
对应 N 阶多项式的拟合公式, 所需要求解的拟合系数需要求解线性方程组, 其中线性方程组的系数矩阵和需要求解的拟合系数矩阵分别为:
A = [ x 1 N ⋯ x 1 ⋯ 1 x 2 N ⋯ x 2 ⋯ 1 ⋮ ⋱ ⋮ x m N ⋯ x m ⋯ 1 ] , θ = [ θ n θ n − 1 ⋮ θ 1 ] \begin{aligned}\mathbf{A}=&\begin{bmatrix}x_1^N&\cdots&x_1\cdots1\\x_2^N&\cdots&x_2\cdots1\\\vdots&\ddots&\vdots\\x_m^N&\cdots&x_m\cdots1\end{bmatrix},\quad\boldsymbol{\theta}=&\begin{bmatrix}\theta_n\\\theta_{n-1}\\\vdots\\\theta_1\end{bmatrix}\end{aligned} A= x1Nx2N⋮xmN⋯⋯⋱⋯x1⋯1x2⋯1⋮xm⋯1 ,θ= θnθn−1⋮θ1
使用加权最小方差方法求解得到拟合系数为:
θ m n = [ θ m n n θ m n − 1 n ⋮ θ 1 n ] = [ A T M A ] − 1 A T M y \theta_m^n=\begin{bmatrix}\theta_{mn}^n\\\theta_{mn-1}^n\\\vdots\\\theta_1^n\end{bmatrix}=\begin{bmatrix}A^{\mathrm{T}}MA\end{bmatrix}^{-1}A^{\mathrm{T}}My θmn= θmnnθmn−1n⋮θ1n =[ATMA]−1ATMy
其对应的加权最小方差为表达式
J m = [ A θ − y ] T W [ A θ − y ] J_m=\begin{bmatrix}A&\theta-y\end{bmatrix}^\mathrm{T}W\begin{bmatrix}A&\theta-y\end{bmatrix} Jm=[Aθ−y]TW[Aθ−y]
根据 WLS 数据拟合方法, 自行编写使用 WLS 方法拟合数据的 M 函数,然后使用 WLS 方法进行数据拟合。
得到的拟合结果如图2所示。
四、非线性曲线拟合
非线性曲线拟合是已知输入向量 xdata、 输出向量ydata, 并知道输入与输出的函数关系为ydata=F(x, xdata) , 但不清楚系数向量x。 进行曲线拟合即求x 使得下式成立:
min x 1 2 ∥ F ( x , x d a t a ) − y d a t a ∥ 2 2 = 1 2 ∑ i ( F ( x , x d a t a i ) − y d a t a i ) 2 \min_x\frac{1}{2}\parallel F(x,xdata)-ydata\parallel_2^2=\frac{1}{2}\sum_i(F(x,xdata_i)-ydata_i)^2 xmin21∥F(x,xdata)−ydata∥22=21i∑(F(x,xdatai)−ydatai)2
在 MATLAB 中, 可以使用函数lsqcurvefit 解决此类问题。 其调用格式如下:
■ x = lsqcurvefit(fun,x0,xdata,ydata) ———x0 为初始解向量; xdata、ydata 为满足关系ydata=F(x, xdata) 的数据。
■ x = lsqcurvefit(fun,x0,xdata,ydata,lb,ub) ———lb、ub 为解向量的下界和上界lb ≤x ≤ub , 若没有指定界, 则lb=[ ],ub=[ ]。
■ x = lsqcurvefit(fun,x0,xdata,ydata,lb,ub,options) ———options 为指定的优化参数。
■ [x,resnorm] = lsqcurvefit(…) ———resnorm 是在x 处残差的平方和。
■ [x,resnorm,residual] = lsqcurvefit(…) ———residual 为在x 处的残差。
■ [x,resnorm,residual,exitflag] = lsqcurvefit ( … ) ———exitflag 为 终 止 迭 代 的条件。
■ [x,resnorm,residual,exitflag,output] = lsqcurvefit(… ) ———output 为输出的优化信息。
已知输入向量xdata 和输出向量ydata, 且长度都是n, 使用最小二乘非线性拟合函数:
y d a t a ( i ) = x ( 1 ) ∙ x d a t a ( i ) 2 + x ( 2 ) ∙ sin ( x d a t a ( i ) ) + x ( 3 ) ∙ x d a t a ( i ) 3 ydata(i)=x(1)\bullet xdata(i)^2+x(2)\bullet\sin(xdata(i))+x(3)\bullet xdata(i)^3 ydata(i)=x(1)∙xdata(i)2+x(2)∙sin(xdata(i))+x(3)∙xdata(i)3
目标函数为:
min x 1 2 ∑ i = 1 n ( F ( x , x d a t a i ) − y d a t a i ) 2 \min_x\frac{1}{2}\sum_{i=1}^n(F(x,xdata_i)-ydata_i)^2 xmin21i=1∑n(F(x,xdatai)−ydatai)2
F ( x , x d a t a ) = x ( 1 ) ∙ x d a t a 2 + x ( 2 ) ∙ sin ( x d a t a ) + x ( 3 ) ∙ x d a t a 3 F(x,xdata)=x(1)\bullet xdata^2+x(2)\bullet\sin(xdata)+x(3)\bullet xdata^3 F(x,xdata)=x(1)∙xdata2+x(2)∙sin(xdata)+x(3)∙xdata3
初始解向量定位x0= [0.3, 0.4, 0.1]。
首先建立拟合函数文件
function F = func_fit (x,xdata)
F = x(1) * xdata.^2 + x(2) * sin(xdata) + x(3) * xdata.^3;
编写函数拟合代码
xdata = [3.6 7.7 9.3 4.1 8.6 2.8 1.3 7.9 10.0 5.4];
ydata = [16.5 150.6 263.1 24.7 208.5 9.9 2.7 163.9 325.0 54.3];
x0 = [10, 10, 10];
[x, resnorm] = lsqcurvefit(@func_fit, x0, xdata, ydata)
结果为:
x =
0.2269 0.3385 0.3022
resnorm =
6.2950
即函数在x =0.2269、x =0.3385、x =0.3022 处残差的平方和均为6.295。
Matlab程序获取
上述matlab程序,可从以下链接处获取:
MATLAB 曲线拟合方法全解析:多项式、加权最小方差与非线性拟合
博主简介:擅长智能优化算法、信号处理、图像处理、机器视觉、深度学习、神经网络等领域Matlab仿真以及实验数据分析等,matlab代码问题、商业合作、课题选题与科研指导等均可私信交流。