Matlab函数逼近与拟合法.ppt

,单击此处编辑母版标题样式,*,*,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,科学计算与,MATLAB,主讲：唐建国,中南大学材料科学与工程学院,2013.09,第五讲 函数逼近与拟合法,内容提要,引言,函数逼近,傅里叶逼近,最小二乘法拟合,最小二乘法,多元线性拟合,非线性拟合,MATLAB,的拟合函数,小结,2025/12/13 周六,3,例：考察某种纤维的强度与其拉伸倍数的关系,下表是实际测定的,24,个纤维样品的强度与相应的拉伸倍数的记录,:,纤维强度随拉伸倍数增加而增加。,1,、引言,2025/12/13 周六,4,24,个点大致分布在一条直线附近。,故可认为强度,y,与拉伸倍数,x,的主要关系应为线性关系：,必须找到一种度量标准来衡量什么曲线最接近所有数据点。,2025/12/13 周六,5,2025/12/13 周六,6,在一个包含有很多数据点的区间内构造插值函数，必然使用高次多项式。而高次插值多项式是不稳定的。,由于数据本身存在误差，利用插值方法得到的插值多项式必然保留了所有的测量误差，导致插值函数与物理规律差异较大。,实验数据的拟合可以克服插值方法在处理这类问题中存在的缺点。,对这样的数据采用上一讲介绍的插值方法近似求描述物理规律的解析函数，必然存在下列缺点：,2025/12/13 周六,7,实验数据拟合的基本思想：,使近似函数尽量靠近数据点，而不要求近似函数一定通过所有数据点。,实验数据拟合可以在一定精度内找出反映物理量间客观函数关系的解析式。如果实验数据存在误差，这种做法可以部分抵消原来数据中的测量误差，从而使所得到的拟合函数更好地反映物理规律。,2025/12/13 周六,8,利用拟合可以解决两类物理问题：,物理规律已知，但描述物理规律的解析式中某些系数未知，可以利用实验方法获得了物理量之间的关系，通过拟合的方法，求出这些系数的近似值。,物理规律未知，利用实验方法获得了物理量之间的关系，通过拟合的方法，得到一个近似的解析式，用于描述物理规律。,拟合函数尽量靠近数据点如何实现？,2025/12/13 周六,9,2,、,函数逼近,在区间,a,b,上已知一连续函数,f(x,),，如果该函数表达式太过于复杂不利于进行计算机运算，就会利用一个简单函数去近似,f(x,),，这就是函数逼近问题。,如果,f(x,),的表达式未知，只知道描述,f(x,),的一条曲线，这就是曲线拟合问题。,和插值问题不同，逼近和拟合并不要求逼近函数在已知点上的值一定等于原函数的函数值，而是按照某种标准使得二者的差值达到最小。,2025/12/13 周六,10,逼近方法：,Chebyshev,(,切比雪夫,),逼近,：连续函数，多项式。,F=Chebyshev(y,k,x0),Legendre(,勒让德,),逼近：多项式。,F=Legendre(y,k,x0),Pade,(,帕德,),逼近：有理分式。,F=Pade(y,k,x0),傅里叶逼近：,周期函数，三角多项式。,连续周期函数，,A0,A,B=,FZZ(func,T,n,),离散周期函数，,c=,DFF(f,N,),2025/12/13 周六,11,Chebyshev,(,切比雪夫,),逼近,当一个连续函数定义在区间,-1,，,1,上时，可以展开成为切比雪夫级数。,2025/12/13 周六,12,function f=Chebyshev(y,k,x0),%,用切比雪夫多项式逼近已知函数,%,已知函数：,y,%,逼近已知函数所需项数：,k,%,逼近点的,x,坐标：,x0,%,求得的切比雪夫逼近多项式或在,x0,处,的逼近值：,f,syms,t;,T(1:k+1)=t;,T(1)=1;,T(2)=t;,c(1:k+1)=0.0;,c(1)=,int(subs(y,findsym(sym(y),sym(t,)*T(1)/sqrt(1-t2),t,-1,1)/pi;,c(2)=2*,int(subs(y,findsym(sym(y),sym(t,)*T(2)/sqrt(1-t2),t,-1,1)/pi;,f=c(1)+c(2)*t;,for i=3:k+1,T(i,)=2*t*T(i-1)-T(i-2);,c(i,)=2*,int(subs(y,findsym(sym(y),sym(t,)*T(i)/sqrt(1-t2),t,-1,1)/2;,f=f+,c(i,)*,T(i,);,f=vpa(f,6);,if(i,=k+1),if(nargin,=3),f=subs(f,t,x0);,else,f=vpa(f,6);,end,end,end,2025/12/13 周六,13,离散周期函数的傅里叶逼近,function c=,DFF(f,N,),%,用傅里叶级数逼近已知的离散周期函数,%,离散数据点：,f,%,展开项数：,N,%,离散傅里叶逼近系数：,c,c(1:N)=0;,for(m,=1:N),for(n,=1:N),c(m,)=,c(m)+f(n,)*exp(-i*m*n*2*pi/N);,end,c(m,)=,c(m,)/N;,end,2025/12/13 周六,14,例,N,1,2,3,4,5,6,Y,0.8415,0.9093,0.1411,-0.7568,-0.9589,-0.2794,y=0.8415 0.9093 0.1411-0.7568-0.9589-0.2794;,c=DFF(y,6),c=,Columns 1 through 4,-0.0926-0.5003i -0.0260-0.0194i -0.0251+0.0000i -0.0260+0.0194i,Columns 5 through 6,-0.0926+0.5003i -0.0172-0.0000i,2025/12/13 周六,15,例：基于神经网络的高炉铁水硅含量预报模型,根据,RBF,神经网络具有收敛速度快和全局优化的特点,建立了,RBF,网络模型,并将其应用对高炉铁水硅含量预报。监于铁水硅含量与炉缸温度之间的密切相关性,通过铁水硅含量来间接地反映炉内温度变化。采用,MATLAB,中的,Newrbe,函数进行函数逼近,对高炉一段连续时期内正常生产的数据的归一化处理后进行训练和仿真,提高了铁水硅含量预报的命中率。,径向基函数,RBF,神经网络是一个只有一个隐藏层的三层前馈神经网络结构。它不仅具有很强的非线性映射能力,并且具有收敛速度快,全局优化的特点,2025/12/13 周六,16,2025/12/13 周六,17,2025/12/13 周六,18,2025/12/13 周六,19,3.1,最小二乘法,首先，从一个简单的例子来讨论一元线性拟合与最小二乘法问题。,为了具有一般性，把上式改写为：,通过实验测量，求金属铜电阻温度系数,，金属电阻与温度关系如下：,3,、最小二乘法拟合,2025/12/13 周六,20,通过实验测得金属铜温度,x,与电阻,y,数据如下：,x,i,(),Y,i,(,),x,i,(),Y,i,(,),x,i,(),Y,i,(,),0,4.38,70,5.58,140,6.74,10,4.56,80,5.74,150,6.94,20,4.70,90,5.96,160,7.12,30,4.86,100,6.06,170,7.28,40,5.08,110,6.26,180,7.42,50,5.24,120,6.44,190,7.60,60,5.40,130,6.58,200,7.78,2025/12/13 周六,21,设一元线性拟合函数为：,将实验数据代入,拟合函数，得到方程组,21,个线性方程,矛盾方程组,2025/12/13 周六,22,由于以上矛盾方程组不能确定一组唯一的,A,0,和,A,1,，,也就是说，由方程组可求得,A,0,和,A,1,的多组解，那么究竟哪一组解最接近客观真实值呢,?,按照拟合的思想，应当使在每一个测量点拟合函数的函数值尽量接近测量值，这样的拟合函数才是满足要求的，即：,定义偏差：,2025/12/13 周六,23,按照拟合的思想，必须在每一个测量点的偏差都很小，如何达到这一要求？,但是由于偏差有正有负，求和时可能互相抵消，这并不能保证在每一个测量点的偏差都很小。,方法一：偏差之和最小,尽管这种方法可以保证在每一个测量点的偏差都很小，但这种方法数学处理比较困难。,方法二：偏差绝对值之和最小,2025/12/13 周六,24,这种方法既可以保证在每一个测量点的偏差都很小，又方便数学处理，所以这种方法是可行的。,方法三：偏差的平方和最小,-,最小二乘法,2025/12/13 周六,25,残差向量的各分量平方和记为：,最小二乘法：,以残差平方和最小问题的解来确定拟合函数的方法。,令,在回归分析中称为残差,(,i,=1,2,m,),残差向量：,2025/12/13 周六,26,由多元函数求极值的必要条件，有,可得,即,2025/12/13 周六,27,上式为由,n,+1,个方程组成的方程组，称正规方程组。,由,得,即,2025/12/13 周六,28,引入记号,则由内积的概念可知,显然内积满足交换律,正规方程组便可化为,2025/12/13 周六,29,将其表示成矩阵形式：,其系数矩阵为对称阵。,所以正规方程组的系数矩阵非奇异,即,根据,Crame,法则,正规方程组有唯一解，称其为最小二乘解。,2025/12/13 周六,30,2025/12/13 周六,31,实现流程图,function a,b=LZXEC(x,y),%,离散试验数据点的线性,最小二乘拟合,%,试验数据点的,x,坐标向量：,X,%,试验数据点的,y,坐标向量：,Y,%,拟合的一次项系数：,a,%,拟合的常数项：,b,if(length(x)=length(y),n=length(x);,else,disp(x,和,y,的维数不相等！,);,return;,end%,维数检查,A=zeros(2,2);,A(2,2)=n;,B=zeros(2,1);,for i=1:n,A(1,1)=A(1,1)+x(i)*x(i);,A(1,2)=A(1,2)+x(i);,B(1,1)=B(1,1)+x(i)*y(i);,B(2,1)=B(2,1)+y(i);,end,A(2,1)=A(1,2);,s=AB;,a=s(1);,b=s(2);,2025/12/13 周六,32,例,X,1,2,3,4,5,y,1.5,1.8,4,3.4,5.7,x=1:5;,y=1.5 1.8 4 3.4 5.7;,a,b=LZXEC(x,y),a=1.0000,b=0.2800,2025/12/13 周六,33,例,.,回到引言的实例,，从散点图可以看出，,纤维强度和拉伸倍数之间近似线性关系，故可选取线性函数,为拟合函数建立正规方程组，其基函数为,根据内积公式，可得,正规方程组为,2025/12/13 周六,34,解得,残差平方和：,拟合曲线与散点,的关系如右图,:,即为所求的最小二乘解。,故,2025/12/13 周六,35,例：金属铜温度,x,与电阻,y,，,线性拟合,matlab,程序,2025/12/13 周六,36,2025/12/13 周六,37,线性拟合在物理实验中应用十分广泛，例如弹性介质杨氏模量测量中应变与应力的关系，电阻电路中电流与电压的关系等。,有些物理量之间在一定范围内是线性关系，也可使用线性拟合的方法，只是要注意其适用范围。,还有一种情况是量物理量之间并不存在线性关系，但经过适当变换后可转化为线性关系。,2025/12/13 周六,38,常用的线性变换,函数,变换后的函数,2025/12/13 周六,39,3.2,多元线性拟合,影响一个物理量的因素，很多情况下不止一个，为了得到描述物理规律的近似函数，就必须采取多元线性拟合。,设物理量,y,随,x,1,，,x,2,，,，,x,k,等,k,个物理量而变化，即：,为了寻求描述物理规律的近似函数，通过实验测得,n,组数据,(,一般,nk),，利用拟合的方法近似函数。,2025/12/13 周六,40,测得,n,组实验数据如下：,i,1,2,3,n,x,1i,x,11,x,12,x,13,x,1n,x,2i,x,21,x,22,x,23,x,2n,x,3i,x,31,x,32,x,33,x,3n,.,.,.,.,.,x,ki,.,.,.,.,.,x,k1,.,.,.,.,.,x,k2,.,.,.,.,.,x,k3,.,.,.,.,.,X,kn,y,i,y,1,y,2,y,3,y,n,2025/12/13 周六,41,设近似,函数,为,与一元线性拟合的思路相同，由n组实验数据代入,上,式，得到n个方程式构成的k元线性方程组，用最小二乘原理确定函数系数A,0,A,1,A,k,使y,i,与Y,i,的偏差的平方和最小。,2025/12/13 周六,42,偏差为,即,(i=1,2,n),令,2025/12/13 周六,43,求极值,2025/12/13 周六,44,得到,2025/12/13 周六,45,化简整理后可得,将第一个方程式中的A,0,提出，代入其它各式后，关于A,0,A,1,A,2,A,k,的正规方程组则为,2025/12/13 周六,46,2025/12/13 周六,47,上式中,(r,s=1,2,k),2025/12/13 周六,48,物理学及各科学技术领域都普遍存在非线性函数关系，对此不能直接使用线性拟合的方法，对这类问题可以采取不同的方法解决。,方法一：变换为线性拟合,有些非线性函数可以经过变量替换，转化成线性函数关系，然后对替换变量进行线性拟合，最后再还原为原始的物理量。但不是所有的函数关系都可经过变量替换而转化成线性函数关系的。,方法二：多项式拟合,任何一个连续函数，在一个比较小的邻域内均可用多项式任意逼近。所以在物理学的许多问题中，不论物理量直接存在何种函数关系，都可用多项式进行数据拟合。,2025/12/13 周六,49,直线：,多项式：,有理函数：,指数：,3.3,非线性拟合,2025/12/13 周六,50,Log-linear:,Gaussians:,2025/12/13 周六,51,设有n组实验数据x,i,，y,i,，（i=1,2,n），用k次多项式拟合，设拟合方程为：,即：,偏差为,2025/12/13 周六,52,偏差的平方和,为：,取极小值,即,2025/12/13 周六,53,整理得,：,即,即得到了正则方程组,2025/12/13 周六,54,MATLAB,程序实现,function A=multifit(X,Y,m),%,离散试验数据点的多项式曲线拟合,%,试验数据点的,x,坐标向量：,X,%,试验数据点的,y,坐标向量：,Y,%,拟合多项式的次数：,m,%,拟合多项式的系数向量：,A,N=length(X);,M=length(Y);,if(N=M),disp(,数据点坐标不匹配！,);,return;,end,c(1:(2*m+1)=0;,b(1:(m+1)=0;,。,for j=1:(2*m+1)%,求出,c,和,b,for k=1:N,c(j)=c(j)+X(k)(j-1);,if(j x=1:4;,y=4 10 18 26;,A=multifit(x,y,2),A=0.0489 0.1612 0.5672,m=1:0.01:4;,n=-1.5+0.49.*m+0.5.*m.2;,l=0.0489+0.1612.*m+,0.5672.*m.2;,plot(x,y,-*,m,n,-y,m,l,-r),P=polyfit(x,y,2),P=0.5000 4.9000 -1.5000,2025/12/13 周六,60,4,、,matlab,拟合函数,线性拟合函数,格式：,p=linefit(x,y),说明：,x,，,y,输入同维数据向量,p,输出拟合多项式的系数向量,2025/12/13 周六,61,例：实验测量获得如下数据，求此物理规律的近似表达式,2025/12/13 周六,62,2025/12/13 周六,63,多项式拟合函数,格式：,p=polyfit(x,y,m),说明：,x,，,y,输入同维数据向量,m,拟合多项式次数,p,输出拟合多项式的系数 向量,2025/12/13 周六,64,例：多项式拟合,2025/12/13 周六,65,2025/12/13 周六,66,2025/12/13 周六,67,2025/12/13 周六,68,2025/12/13 周六,69,非线性函数拟合,格式：,a=lsqcurvefit(fun,a0,x,y),说明：,x,，,y,输入同维数据向量,fun,拟合函数文件,a,输出拟合函数的系数 向量,a0 a,的初值,2025/12/13 周六,70,例：已知,利用下列数据，求,a,b,c,d,2025/12/13 周六,71,2025/12/13 周六,72,2025/12/13 周六,73,2025/12/13 周六,74,2025/12/13 周六,75,例：布氏硬度图像自动测量,2025/12/13 周六,76,布氏硬度测量的关键是对压痕图像进行自动图像处理,以得到压痕的直径,代入公式从而得到布氏硬度值。最好的办法是抓取压痕的边界,对边界用最小二乘法拟合出一个最贴近的圆,求出该圆的像素直径,从而求出真实直径。,2025/12/13 周六,77,2025/12/13 周六,78,2025/12/13 周六,79,2025/12/13 周六,80,小 结,函数逼近和拟合与插值法的区别,线性最小二乘法的应用和计算,最小二乘法的程序设计,拟合函数的理解与应用,2025/12/13 周六,81,谢 谢！,