最小二乘法圆拟合及matlab程序
- 格式:ppt
- 大小:51.50 KB
- 文档页数:10
曲线拟合(curve-fitting ):工程实践中,用测量到的一些离散的数据},...2,1,0),,{(m i y x i i =求一个近似的函数)(x ϕ来拟合这组数据,要求所得的拟合曲线能最好的反映数据的基本趋势(即使)(x ϕ最好地逼近()x f ,而不必满足插值原则。
因此没必要取)(i x ϕ=i y ,只要使i i i y x -=)(ϕδ尽可能地小)。
原理:给定数据点},...2,1,0),,{(m i y x i i =。
求近似曲线)(x ϕ。
并且使得近似曲线与()x f 的偏差最小。
近似曲线在该点处的偏差i i i y x -=)(ϕδ,i=1,2,...,m 。
常见的曲线拟合方法:1.使偏差绝对值之和最小2.使偏差绝对值最大的最小3.使偏差平方和最小最小二乘法:按偏差平方和最小的原则选取拟合曲线,并且采取二项式方程为拟合曲线的方法,称为最小二乘法。
推导过程:1. 设拟合多项式为:kk x a x a a x +++=...)(10ϕ2. 各点到这条曲线的距离之和,即偏差平方和如下:3. 问题转化为求待定系数0a ...k a 对等式右边求i a 偏导数,因而我们得到了:.......4、 把这些等式化简并表示成矩阵的形式,就可以得到下面的矩阵:5. 将这个范德蒙得矩阵化简后可得到:6. 也就是说X*A=Y ,那么A = (X'*X)-1*X'*Y ,便得到了系数矩阵A ,同时,我们也就得到了拟合曲线。
MATLAB实现:MATLAB提供了polyfit()函数命令进行最小二乘曲线拟合。
调用格式:p=polyfit(x,y,n)[p,s]= polyfit(x,y,n)[p,s,mu]=polyfit(x,y,n)x,y为数据点,n为多项式阶数,返回p为幂次从高到低的多项式系数向量p。
x 必须是单调的。
矩阵s包括R(对x进行QR分解的三角元素)、df(自由度)、normr(残差)用于生成预测值的误差估计。
方便大家使用的最小二乘法曲线拟合的Matlab程序非常方便用户使用,直接按提示操作即可;这里我演示一个例子:(红色部分为用户输入部分,其余为程序运行的结果,结果图为Untitled.fig,Untitled2.fig) 请以向量的形式输入x,y.x=[1,2,3,4]y=[3,4,5,6]通过下面的交互式图形,你可以事先估计一下你要拟合的多项式的阶数,方便下面的计算.polytool()是交互式函数,在图形上方[Degree]框中输入阶数,右击左下角的[Export]输出图形回车打开polytool交互式界面回车继续进行拟合输入多项式拟合的阶数m = 4Warning: Polynomial is not unique; degree >= number of data points. > In polyfit at 72In zxecf at 64输出多项式的各项系数a = 0.0200000000000001a = -0.2000000000000008a = 0.7000000000000022a = 0.0000000000000000a = 2.4799999999999973输出多项式的有关信息 SR: [4x5 double]df: 0normr: 2.3915e-015Warning: Zero degrees of freedom implies infinite error bounds.> In polyval at 104In polyconf at 92In zxecf at 69观测数据拟合数据x y yh1.0000 3.0000 3.00002.0000 4.0000 4.00003 5 54.0000 6.0000 6.0000剩余平方和 Q = 0.000000标准误差 Sigma = 0.000000相关指数 RR = 1.000000请输入你所需要拟合的数据点,若没有请按回车键结束程序.输入插值点x0 = 3输出插值点拟合函数值 y0 = 5.0000>>结果:untitled.figuntitled2.fig一些matlab优化算法代码的分享代码的目录如下:欢迎讨论1.约束优化问题:minRosen(Rosen梯度法求解约束多维函数的极值)(算法还有bug) minPF(外点罚函数法解线性等式约束)minGeneralPF(外点罚函数法解一般等式约束)minNF(内点罚函数法)minMixFun(混合罚函数法)minJSMixFun(混合罚函数加速法)minFactor(乘子法)minconPS(坐标轮换法)(算法还有bug)minconSimpSearch(复合形法)2.非线性最小二乘优化问题minMGN(修正G-N法)3.线性规划:CmpSimpleMthd(完整单纯形法)4.整数规划(含0-1规划)DividePlane(割平面法)ZeroOneprog(枚举法)5.二次规划QuadLagR(拉格朗日法)ActivedeSet(起作用集法)6.辅助函数(在一些函数中会调用)minNT(牛顿法求多元函数的极值)Funval(求目标函数的值)minMNT(修正的牛顿法求多元函数极值)minHJ(黄金分割法求一维函数的极值)7.高级优化算法1)粒子群优化算法(求解无约束优化问题)1>PSO(基本粒子群算法)2>YSPSO(待压缩因子的粒子群算法)3>LinWPSO(线性递减权重粒子群优化算法)4>SAPSO(自适应权重粒子群优化算法)5>RandWSPO(随机权重粒子群优化算法)6>LnCPSO(同步变化的学习因子)7>AsyLnCPSO(异步变化的学习因子)(算法还有bug)8>SecPSO(用二阶粒子群优化算法求解无约束优化问题)9>SecVibratPSO(用二阶振荡粒子群优化算法求解五约束优化问题)10>CLSPSO(用混沌群粒子优化算法求解无约束优化问题)11>SelPSO(基于选择的粒子群优化算法)12>BreedPSO(基于交叉遗传的粒子群优化算法)13>SimuAPSO(基于模拟退火的粒子群优化算法)2)遗传算法1>myGA(基本遗传算法解决一维约束规划问题)2>SBOGA(顺序选择遗传算法求解一维无约束优化问题)3>NormFitGA(动态线性标定适应值的遗传算法求解一维无约束优化问题)4>GMGA(大变异遗传算法求解一维无约束优化问题)5>AdapGA(自适应遗传算法求解一维无约束优化问题)6>DblGEGA(双切点遗传算法求解一维无约束优化问题)7>MMAdapGA(多变异位自适应遗传算法求解一维无约束优化问题)自己编写的马尔科夫链程序A 代表一组数据序列一维数组本程序的操作对象也是如此t=length(A); % 计算序列“A”的总状态数B=unique(A); % 序列“A”的独立状态数顺序,“E”E=sort(B,'ascend');a=0;b=0;c=0;d=0;for j=1:1:ttLocalization=find(A==E(j)); % 序列“A”中找到其独立状态“E”的位置for i=1:1:length(Localization)if Localization(i)+1>tbreak; % 范围限定elseif A(Localization(i)+1)== E(1)a=a+1;elseif A(Localization(i)+1)== E(2)b=b+1;elseif A(Localization(i)+1)== E(3)c=c+1;% 依此类推,取决于独立状态“E”的个数elsed=d+1;endendT(j,1:tt)=[a,b,c,d]; % “T”为占位矩阵endTT=T;for u=2:1:ttTT(u,:)= T(u,:)- T(u-1,:);endTT; % 至此,得到转移频数矩阵Y=sum(TT,2);for uu=1:1:ttTR(uu,:)= TT(uu,:)./Y(uu,1);endTR % 最终得到马尔科夫转移频率/概率矩阵% 观测序列马尔科夫性质的检验:N=numel(TT);uuu=1;Col=sum(TT,2); % 对列求和Row=sum(TT,1); % 对行求和Total=sum(Row); % 频数总和for i=1:1:ttfor j=1:1:ttxx(uuu,1)=sum((TT(i,j)-(Row(i)*Col(j))./Total).^2./( (Row(i)*Col(j)). /Total));uuu=uuu+1; % 计算统计量x2endendxx=sum(xx)。
2曲线拟合的线性最小二乘法及其MATLAB?序例2给出一组数据点(X i, y i)列入表2中,试用线性最小二乘法求拟合曲线, 估计其误差,作出拟合曲线•解 (1)在MATLAB工作窗口输入程序>> x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.126.50 68.04];plot(x,y , 'r*'),lege nd( '实验数据(xi,yi)' )xlabel( 'x' ), ylabel( 'y'),title( '数据点(xi,yi) 的散点图’)运行后屏幕显示数据的散点图(略)(3)编写下列MATLAB程序计算f(x)在(X j,yj处的函数值,即输入程序>> syms al a2 a3 a4x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];fi=a1.*x.A3+ a2.*x.A2+ a3.*x+ a4运行后屏幕显示关于a1,a2, a3和a4的线性方程组fi =[ -125/8*a1+25/4*a2-5/2*a3+a4,-4913/1000*a1+289/100*a2-17/10*a3+a4,-1331/1000*a1+121/100*a2-11/10*a3+a4, -64/125*a1+16/25*a2-4/5*a3+a4,a4, 1/1000*a1+1/100*a2+1/10*a3+a4,27/8*a1+9/4*a2+3/2*a3+a4,19683/1000*a1+729/100*a2+27/10*a3+a4, 5832/125*a1+324/25*a2+18/5*a3+a4]编写构造误差平方和的MATLAB程序>> y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04];fi=[-125/8*a1+25/4*a2-5/2*a3+a4, -4913/1000*a1+289/100*a2-17/10*a3+a4, -1331/1000*a1+121/100*a2-11/10*a3+a4,-64/125*a1+16/25*a2-4/5*a3+a4, a4,1/1000*a1+1/100*a2+1/10*a3+a4,27/8*a1+9/4*a2+3/2*a3+a4,19683/1000*a1+729/100*a2+27/10*a3+a4,5832/125*a1+324/25*a2+18/5*a3+a4];fy=fi-y; fy2=fy.A2; J=sum(fy.A2)运行后屏幕显示误差平方和如下J=(-125/8*a1+25/4*a2-5/2*a3+a4+1929/10)A2+(-4913/1000*a1+289/100*a2-17/10*a3+a4+171/2)A2+(-1331/1000*a1+121/100*a2-11/10*a3+a4+723/20F2+(-64/125*a1+16/25*a2-4/5*a3+a4+663/25F2+(a4+91/10F2+(1/1000*a1+1/100*a2+1/10*a3+a4+843/100)A2+(27/8*a1+9/4*a2+3/2*a3+a4+328/25)A2+(19683/1000 *a1+729/100*a2+27/10*a3+a4-13/2F2+(5832/125*a1+324/25*a2 +18/5*a3+a4-1701/25F2为求31,32,33,34使J达到最小,只需利用极值的必要条件-丄 0 a k (k 1,2,3,4),得到关于31,32,33,34的线性方程组,这可以由下面的MATLAB程序完成,即输入程序>> syms a1 a2 a3 a4 J=(-125/8*a1+25/4*32-5/2*a3+34+1929/10)A2+(-4913/1000*a1+289/100*32-17/10*33+34...+171/2)A2+(-1331/1000*a1+1 21/100*a2-11/10*a3+34+723/20)A2+(-64/125*31+16/25*32-4/5*a3+34+663/25)A2+(34+91/10)A2+(1/1000*31+1/100*32+1/10*a3+a4+843/100)A2+(27/8*31+9/4*32+3/2*a3+34+328/25)A2+(19683/ 1000*a1+729/100*32+27/10*a3+34-13/2)A2+(5832/125*31+324/2 5*a2+18/5*a3+a4-1701/25)A2;Ja1=diff(J,a1); Ja2=diff(J,a2); Ja3=diff(J,a3);Ja4=diff(J,a4);Ja11=simple(Ja1), Ja21=simple(Ja2), Ja31=simple(Ja3),Ja41=simple(Ja4),运行后屏幕显示J分别对31, 32 ,33 ,34的偏导数如下Ja1仁56918107/10000*31+32097579/25000*32+1377283/2500*33+23667/250*34-8442429/625J321 =32097579/25000*31+1377283/2500*32+23667/250*33 +67*34+767319/625 J331 = 1377283/2500*31+23667/250*32+67*33+18/5*34-232638/125J341 = 23667/250*31+67*32+18/5*33+18*34+14859/25解线性方程组J311 =0,J321 =0,J331 =0,J341 =0,输入下列程序>>A=[56918107/10000, 32097579/25000, 1377283/2500,23667/250; 32097579/25000, 1377283/2500, 23667/250, 67; 1377283/2500, 23667/250, 67, 18/5; 23667/250, 67, 18/5, 18];B=[8442429/625, -767319/625, 232638/125, -14859/25];C=B/A, f=poly2sym(C)运行后屏幕显示拟合函数f及其系数C如下C = 5.0911 -14.1905 6.4102 -8.2574f=716503695845759/140737488355328*xA3-7988544102557579/562949953421312*xA2+1804307491277693/281474976710656*x -4648521160813215/562949953421312 故所求的拟合曲线为f (x) 5.0911 x314.1905 x2 6.4102 x 8.2574 .(4)编写下面的MATLAB 程序估计其误差,并作出拟合曲线和数据的图形.输入程序>> xi=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6]; y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04]; n=length(xi);f=5.0911.*xi.A3-14.1905.*xi.A2+6.4102.*xi -8.2574;x=-2.5:0.01: 3.6;F=5.0911.*x.A3-14.1905.*x.A2+6.4102.*x -8.2574;fy=abs(f-y); fy2=fy.A2; Ew=max(fy),E1=sum(fy)/n, E2=sqrt((sum(fy2))/n) plot(xi,y,'r*'), hold on, plot(x,F, 'b-'),hold off legend('数据点(xi,yi)','拟合曲线y=f(x)'),xlabel('x'), ylabel('y'),title(例2的数据点(xi,yi)和拟合曲线y=f(x)的图形')运行后屏幕显示数据(X i,yj与拟合函数f的最大误差平均误差E i和均方根误差E2及其数据点(X j,yj和拟合曲线y=f(x)的图形(略).Ew = E1 = E2 =3.105 4 0.903 4 1.240 96函数逼近及其MATLAB?序最佳均方逼近的MATLAB^程序function [yy1,a,WE]=zjjfbj(f,X,Y,xx) m=size(f);n=length(X);m=m(1);b=zeros(m,m); c=zeros(m,1);if n~=length(Y) error( 'X和丫的维数应该相同') end for j=1:m for k=1:mb(j,k)=0;for i=1:nb(j,k)=b(j,k)+feval(f(j,:),X(i))*feval(f(k,:),X(i));endendc(j)=0;for i=1:nc(j)=c(j)+feval(f(j,:),X(i))*Y(i);endenda=b\c;WE=0;for i=1:nff=0;for j=1:m ff=ff+a(j)*feval(f(j,:),X(i)); endWE=WE+(Y(i)-ff)*(Y(i)-ff);endif nargin==3return ;endyy=[];for i=1:ml=[];for j=1:length(xx) l=[l,feval(f(i,:),xx(j))]; end yy=[yy l'];endyy=yy*a; yy1=yy'; a=a';WE;例6. 1对数据X和Y,用函数y 1,y x, y x2进行逼近,用所得到的逼近函数计算在x 6.5处的函数值,并估计误差.其中X=(1 3 4 5 6 7 8 9); Y=(-11 -13 -11 -7 -1 7 17 29). 解在MATLA工作窗口输入程序>> X=[ 1 3 4 5 6 7 8 9]; Y=[-11 -13 -11 -7 -17 17 29];f=['fun0';'fun1';'fun2'];[yy,a,WE]=zjjfbj(f,X,Y,6.5) 运行后屏幕显示如下yy =2.75000000000003a =-7.00000000000010 -4.999999999999951.00000000000000WE =7.172323350269439e-027例 6.2 对数据X 和丫,用函数 y 1, y x, y x2,y cosx,y e x,y sinx进行逼近,其中X=(0 0.50 1.00 1.50 2.00 2.50 3.00 ),丫=(0 0.47940.8415 0.9815 0.9126 0.5985 0.1645 ) .解在MATLA工作窗口输入程序>> X=[ 0 0.50 1.00 1.50 2.00 2.50 3.00];丫=[0 0.4794 0.8415 0.9815 0.9126 0.1645];f=['fun0';'fun1';'fun2';'fun3';'fun4';'fun5'];xx=0:0.2:3;[yy,a,WE]=zjjfbj(f,X,Y,plot(X,Y,'ro',xx,yy,'b-')运行后屏幕显示如下(图略)yy = Columns 1 through 7-0.0005 0.2037 0.3939 0.5656 0.8348 0.9236Columns 8 through 140.9771 0.9926 0.9691 0.9069 0.6766 0.5191Columns 15 through 160.3444 0.1642 0.5985xx), 0.7141 0.8080a = 0.3828 0.4070 -0.3901 0.0765 -0.4598 0.5653 WE = 1.5769e-004 即,最佳逼近函数为y=0.3828+0.4070*x-0.3901*xA2+0.0765*exp(x) +0.5653*sin(x) .8随机数据点上的二元拟合及其MATLA 程序例 8 设节点 ( X,Y,Z ) 中的 X 和 Y 分别是在区间 [ 3,3] 和 [ 2.5,3.5]上的 5022个 随 机 数 , Z 是 函 数 Z=7-3 x 3e -x2 -y2 在 (X,Y ) 的 值 , 拟 合 点 ( X I ,Y I ) 中的 X I =-3:0.2:3, Y I =-2.5:0.2:3.5. 分别用二元拟合方法中最近邻内插法、三 角基线性内插法、三角基三次内插法和 MATLAB4 网格化坐标方法计算在 ( X I ,Y I ) 处的值,作出它们的图形,并与被拟和曲面进行比较 .解 (1 )最近邻内插法 .输入程序>> x=rand(50,1); y=rand(50,1); %生成50个一元均匀分布随机数x 和y , x , y .X=-3+(3-(-3))*x; %利用x 生成的随机变量.title( ' 用最近邻内插法拟合函数 z =7-3 xA3 exp(-xA2 - yA2) 的曲面和节点的图形 ' )%legend( ' 拟合曲面 ',' 节点 (xi,yi,zi)' )hold on%在当前图形上添加新图形面及其插值乙(略).(2)三角基线性内插法 输入程序>> x=rand(50,1); y=rand(50,1); %生成50个一元均匀分布随机数x 和y , x ,y .X=-3+(3-(-3))*x; %利用x 生成 上的随机变量.Y=-2.5+(3.5-(-2.5))*y;%利用 y 生成 上的随机变量 .Z=7-3* X.A3 .* exp(-X.A2 - Y.A2);%在每个随机点( X,Y )处计算Z 的值.-0.4598*cos(x)Y=-2.5+(3.5-(-2.5))*y;Z=7-3* X.A3 .* exp(-X.A2 - Y.A2);处计算Z 的值.X1=-3:0.2:3; Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); ZI=griddata(X,Y,Z,XI,YI, (XI,YI )处的插值 ZI.mesh(XI,YI, ZI)xlabel( 'x' ), ylabel( 'y'%利用 y 生成的随机变量 .%在每个随机点( X,Y )%将坐标( XI,YI )网格化 .'nearest' ) %计算在每个插值点 %作二元拟合图形 . ), zlabel( 'z' ),plot3(X,Y,Z, 'bo' ) (X,Y,Z). hold of 运行后屏幕显示用最近邻内插法拟合函数%用兰色小圆圈画出每个节点 %结束在当前图形上添加新图形 .22Z=7-3 x 3e-x2 -y2在两组不同节点处的曲X1=-3:0.2:3;title('用三角基线性内插法拟合函数z =7-3 x A3 exp(-x A2 -y A2)的曲面和节点的图形’) %legend( ' 拟合曲面',' hold on plot3(X,Y,Z, 'bo' ) (X,Y,Z). hold of22运行后屏幕显示用三角基线性内插法拟合函数Z=7-3 x 3e-x -y在两组不同节点处的曲面和节点的图形及其插值 乙(略).(3)三角基三次内插法 . 输入程序>> x=rand(50,1); y=rand(50,1); %生成50个一元均匀分布随机数x 和y , x ,y .title( ' 用三角基三次内插法拟合函数 z =7-3 xA3 exp(-xA2 -22运行后屏幕显示用三角基三次内插法拟合函数 Z=7-3 x 3e -x -y 在两组不同节点处的曲面和节点的图形及其插值 Z I (略).( 4 ) MATLAB 4网格化坐标方法 . 输入程序>> x=rand(50,1); y=rand(50,1); %生成50个一元均匀分布随机数x 和y , x ,y .X=-3+(3-(-3))*x;%利用x 生成 上的随机变量.Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); ZI=griddata(X,Y,Z,XI,YI,XI,YI )处的插值 ZI.mesh(XI,YI, ZI)xlabel( 'x' ), ylabel( 'y%将坐标( XI,YI )网格化 .'linear' ) %计算在每个插值点 %作二元拟合图形 . ), zlabel( 'z' ),节点 (xi,yi,zi)' )%在当前图形上添加新图形 .%用兰色小圆圈画出每个节点 %结束在当前图形上添加新图形 . X=-3+(3-(-3))*x; %利用x 生成上的随机变量.Y=-2.5+(3.5-(-2.5))*y;Z=7-3* X.A3 .* exp(-X.A2 - Y.A2);处计算Z 的值.X1=-3:0.2:3; Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); ZI=griddata(X,Y,Z,XI,YI,( XI,YI )处的插值 ZI.mesh(XI,YI, ZI) xlabel( 'x' ), ylabel( 'y'%利用y 生成上的随机变量•%在每个随机点( X,Y )%将坐标( XI,YI )网格化 .'cubic' ) %计算在每个插值点 %作二元拟合图形 . ), zlabel( 'z' ),yA2) 的曲面和节点的图形 ' )%legend( ' 拟合曲面','hold on 节点 (xi,yi,zi)' )%在当前图形上添加新图形plot3(X,Y,Z,'bo' ) (X,Y,Z).hold of%用兰色小圆圈画出每个节点%结束在当前图形上添加新图形 .22运行后屏幕显示用MATLAB 网格化坐标方法拟合函数Z=7-3 x 3e -x-y 在两组不同 节点处的曲面和节点的图形及其插值 ZI (略).22(5) 作被拟合曲面Z=7-3x 3e -x-y 和节点的图形. 输入程序>> x=ra nd(50,1); y=rand(50,1); %生成50个一元均匀分布随机数x 和y , x ,y .X=-3+(3-(-3))*x; %利用x 生成随机变量.Y=-2.5+(3.5-(-2.5))*y;%利用y 生成随机变量.Z=7-3* X.A3 .* exp(-X.A2 - Y.A2);%在每个随机点( X,Y )处计算Z 的值.X1=-3.:0.1:3.;Y1=-2.5:0.1:3.5;[XI,YI] = meshgrid(X1,Y1);%将坐标( XI,YI )网格化 .ZI=7-3* XI.A3 .* exp(-XI.A2 - YI.A2); mesh(XI,YI, ZI) %作二元拟合图形 .xlabel( 'x'), ylabel('y' ), zlabel( 'z' ),title( ' 被拟合函数 z =7-3 xA3 exp(-xA2 - yA2) 的曲面和节点的图形 ' )%legend('被拟合函数曲面','节点(xi,yi,zi)' )hold on%在当前图形上添加新图形 .plot3(X,Y,Z, 'bo' )%用兰色小圆圈画出每个节点 (X,Y,Z).hold of%结束在当前图形上添加新图形 .22运行后屏幕显示被拟合函数 Z=7-3 x 3e -x-y 的曲面和节点的图形及其函数值 ZI(略) .Y=-2.5+(3.5-(-2.5))*y;Z=7-3* X.A3 .* exp(-X.A2 - Y.A2);处计算Z 的值.X1=-3:0.2:3; Y1=-2.5:0.2:3.5; [XI,YI] = meshgrid(X1,Y1); ZI=griddata(X,Y,Z,XI,YI, (XI,YI )处的插值 ZI.mesh(XI,YI, ZI) xlabel( 'x' ), ylabel( 'y' %利用y 生成上的随机变量.%在每个随机点( X,Y )'v4'%将坐标( XI,YI )网格化 . ) %计算在每个插值点%作二元拟合图形 . ), zlabel( 'z' ),' 用 MATLAB 4 网格化坐标方法 拟合函数 z =7-3 xA3 的曲面和节点的图形 ' )%legend( ' 拟合曲面 ',' 节点 (xi,yi,zi)' ) hold on %在当前图形上添加新图形 . plot3(X,Y,Z, 'bo' ) %用兰色小圆圈画出每个节点(X,Y,Z).hold oftitle(exp(-x A2 - y A2)'bo' ) %结束在当前图形上添加新图形 .。
最小二乘法matlab程序最小二乘法是一种统计模型,它可以被用来拟合一元函数数据,或者拟合非线性曲线。
它的基本思想是找到一组参数,使得拟合的曲线与实际数据的差距最小。
本文将介绍如何使用Matlab实现一个最小二乘法的程序,并与现有的一些现成的最小二乘法的matlab程序进行比较,找出其优缺点。
首先,要使用最小二乘法拟合曲线,需要准备一组输入数据,一般可以将其表示为两个向量,分别是自变量x和因变量y。
这些数据可以是由测量和实验得到的,也可以是由人工输入的,但无论如何都要确保它们的准确性。
接下来,就可以使用Matlab输入数据进行处理,用最小二乘法计算出最拟合的曲线及其参数。
具体步骤主要分为三步:第一步是计算输入数据的均值和方差,包括自变量x和因变量y的均值和方差;第二步是计算自变量x和因变量y的关系,即最小二乘拟合曲线的系数;第三步是验证拟合的曲线的准确性,如果不满意,可以重新调整参数,以获得较好的拟合效果。
此外,Matlab除了提供自带的最小二乘法函数外,还支持第三方开发者开发现成的matlab程序,用于解决最小二乘法的问题。
这些程序中有一些是开源的,另一些则是出售的。
其中开源的有LEAST,CURVEFIT,CURVEFITTOOL等,而出售的有MATLAB Curve Fitting Toolbox,Optimization Toolbox和Statistics Toolbox等。
它们的突出特点是速度快,代码简洁,容易上手,适用于多种拟合类型。
然而,各种matlab程序也有自身的缺点,最明显的就是当输入数据非常庞大时,它们的计算能力就无法跟上,速度就会变慢。
此外,使用出售的matlab程序可能相对昂贵,而且有时需要安装某些复杂的库文件,这也是一种麻烦。
因此,使用最小二乘法拟合曲线时,可以参考现有的matlab程序,也可以自己编写matlab代码,同时要考虑到程序的可靠性、效率和可行性。
本文介绍的matlab程序的最大优势是它不需要依赖第三方的软件,而且能够满足大多数用户的需求,使得最小二乘法可以在短时间内被成功运用。
最⼩⼆乘法曲线拟合的Matlab程序⽅便⼤家使⽤的最⼩⼆乘法曲线拟合的Matlab程序⾮常⽅便⽤户使⽤,直接按提⽰操作即可;这⾥我演⽰⼀个例⼦:(红⾊部分为⽤户输⼊部分,其余为程序运⾏的结果,结果图为Untitled.fig,Untitled2.fig) 请以向量的形式输⼊x,y.x=[1,2,3,4]y=[3,4,5,6]通过下⾯的交互式图形,你可以事先估计⼀下你要拟合的多项式的阶数,⽅便下⾯的计算.polytool()是交互式函数,在图形上⽅[Degree]框中输⼊阶数,右击左下⾓的[Export]输出图形回车打开polytool交互式界⾯回车继续进⾏拟合输⼊多项式拟合的阶数m = 4Warning: Polynomial is not unique; degree >= number of data points. > In polyfit at 72In zxecf at 64输出多项式的各项系数a = 0.0200000000000001a = -0.2000000000000008a = 0.7000000000000022a = 0.0000000000000000a = 2.4799999999999973输出多项式的有关信息 SR: [4x5 double]df: 0normr: 2.3915e-015Warning: Zero degrees of freedom implies infinite error bounds.> In polyval at 104In polyconf at 92In zxecf at 69观测数据拟合数据x y yh1.0000 3.0000 3.00002.0000 4.0000 4.00003 5 54.0000 6.0000 6.0000剩余平⽅和 Q = 0.000000相关指数 RR = 1.000000请输⼊你所需要拟合的数据点,若没有请按回车键结束程序.输⼊插值点x0 = 3输出插值点拟合函数值 y0 = 5.0000>>结果:untitled.figuntitled2.fig⼀些matlab优化算法代码的分享代码的⽬录如下:欢迎讨论1.约束优化问题:minRosen(Rosen梯度法求解约束多维函数的极值)(算法还有bug) minPF(外点罚函数法解线性等式约束) minGeneralPF(外点罚函数法解⼀般等式约束)minNF(内点罚函数法)minMixFun(混合罚函数法)minJSMixFun(混合罚函数加速法)minFactor(乘⼦法)minconPS(坐标轮换法)(算法还有bug)minconSimpSearch(复合形法)2.⾮线性最⼩⼆乘优化问题minMGN(修正G-N法)3.线性规划:CmpSimpleMthd(完整单纯形法)4.整数规划(含0-1规划)DividePlane(割平⾯法)ZeroOneprog(枚举法)5.⼆次规划QuadLagR(拉格朗⽇法)ActivedeSet(起作⽤集法)6.辅助函数(在⼀些函数中会调⽤)minNT(⽜顿法求多元函数的极值)minMNT(修正的⽜顿法求多元函数极值)minHJ(黄⾦分割法求⼀维函数的极值)7.⾼级优化算法1)粒⼦群优化算法(求解⽆约束优化问题)1>PSO(基本粒⼦群算法)2>YSPSO(待压缩因⼦的粒⼦群算法)3>LinWPSO(线性递减权重粒⼦群优化算法)4>SAPSO(⾃适应权重粒⼦群优化算法)5>RandWSPO(随机权重粒⼦群优化算法)6>LnCPSO(同步变化的学习因⼦)7>AsyLnCPSO(异步变化的学习因⼦)(算法还有bug)8>SecPSO(⽤⼆阶粒⼦群优化算法求解⽆约束优化问题)9>SecVibratPSO(⽤⼆阶振荡粒⼦群优化算法求解五约束优化问题)10>CLSPSO(⽤混沌群粒⼦优化算法求解⽆约束优化问题)11>SelPSO(基于选择的粒⼦群优化算法)12>BreedPSO(基于交叉遗传的粒⼦群优化算法)13>SimuAPSO(基于模拟退⽕的粒⼦群优化算法)2)遗传算法1>myGA(基本遗传算法解决⼀维约束规划问题)2>SBOGA(顺序选择遗传算法求解⼀维⽆约束优化问题)3>NormFitGA(动态线性标定适应值的遗传算法求解⼀维⽆约束优化问题)4>GMGA(⼤变异遗传算法求解⼀维⽆约束优化问题)5>AdapGA(⾃适应遗传算法求解⼀维⽆约束优化问题)6>DblGEGA(双切点遗传算法求解⼀维⽆约束优化问题)7>MMAdapGA(多变异位⾃适应遗传算法求解⼀维⽆约束优化问题)⾃⼰编写的马尔科夫链程序A 代表⼀组数据序列⼀维数组本程序的操作对象也是如此t=length(A); % 计算序列“A”的总状态数B=unique(A); % 序列“A”的独⽴状态数顺序,“E”E=sort(B,'ascend');a=0;b=0;c=0;d=0;Localization=find(A==E(j)); % 序列“A”中找到其独⽴状态“E”的位置for i=1:1:length(Localization)if Localization(i)+1>tbreak; % 范围限定elseif A(Localization(i)+1)== E(1)a=a+1;elseif A(Localization(i)+1)== E(2)b=b+1;elseif A(Localization(i)+1)== E(3)c=c+1;% 依此类推,取决于独⽴状态“E”的个数elsed=d+1;endendT(j,1:tt)=[a,b,c,d]; % “T”为占位矩阵endTT=T;for u=2:1:ttTT(u,:)= T(u,:)- T(u-1,:);endTT; % ⾄此,得到转移频数矩阵Y=sum(TT,2);for uu=1:1:ttTR(uu,:)= TT(uu,:)./Y(uu,1);endTR % 最终得到马尔科夫转移频率/概率矩阵% 观测序列马尔科夫性质的检验:N=numel(TT);uuu=1;Col=sum(TT,2); % 对列求和Row=sum(TT,1); % 对⾏求和Total=sum(Row); % 频数总和for i=1:1:ttxx(uuu,1)=sum((TT(i,j)-(Row(i)*Col(j))./Total).^2./( (Row(i)*Col(j)). /Total)); uuu=uuu+1; % 计算统计量x2endendxx=sum(xx)。
最小二乘法曲线拟合的Matlab程序最小二乘法是一种常用的数学优化技术,它通过最小化误差的平方和来找到最佳函数匹配。
在曲线拟合中,最小二乘法被广泛使用来找到最佳拟合曲线。
下面的Matlab程序演示了如何使用最小二乘法进行曲线拟合。
% 输入数据x = [1, 2, 3, 4, 5];y = [2.2, 2.8, 3.6, 4.5, 5.1];% 构建矩阵A = [x(:), ones(size(x))]; % 使用x向量和单位矩阵构建矩阵A% 使用最小二乘法求解theta = (A' * A) \ (A' * y); % 利用最小二乘法的公式求解% 显示拟合曲线plot(x, theta(1) * x + theta(2), '-', 'LineWidth', 2); % 画出拟合曲线hold on; % 保持当前图像plot(x, y, 'o', 'MarkerSize', 10, 'MarkerFaceColor','b'); % 在图像上画出原始数据点xlabel('x'); % 设置x轴标签ylabel('y'); % 设置y轴标签legend('拟合曲线', '原始数据点'); % 设置图例这个程序首先定义了一组输入数据x和y。
然后,它构建了一个矩阵A,这个矩阵由输入数据x和单位矩阵构成。
然后,程序使用最小二乘法的公式来求解最佳拟合曲线的参数。
最后,程序画出拟合曲线和原始数据点。
这个程序使用的是线性最小二乘法,适用于一次曲线拟合。
如果你的数据更适合非线性模型,例如二次曲线或指数曲线,那么你需要使用非线性最小二乘法。
Matlab提供了lsqcurvefit函数,可以用于非线性曲线拟合。
例如:% 非线性模型 y = a * x^2 + b * x + cfun = @(theta, x) theta(1) * x.^2 + theta(2) * x +theta(3);guess = [1, 1, 1]; % 初始猜测值% 使用lsqcurvefit函数求解theta = lsqcurvefit(fun, guess, x, y);% 显示拟合曲线plot(x, fun(theta, x), '-', 'LineWidth', 2); % 画出拟合曲线hold on; % 保持当前图像plot(x, y, 'o', 'MarkerSize', 10, 'MarkerFaceColor','b'); % 在图像上画出原始数据点xlabel('x'); % 设置x轴标签ylabel('y'); % 设置y轴标签legend('拟合曲线', '原始数据点'); % 设置图例这个程序定义了一个非线性函数fun,然后使用lsqcurvefit函数来求解最佳拟合曲线的参数。
最小二乘法拟合和优化常用的 matlab 命令及其适用范围
最小二乘法拟合和优化常用的Matlab 命令及其适用范围如下:
1. polyfit:用于多项式拟合的函数,可以通过最小二乘法拟合
一组数据点的多项式曲线。
适用范围:适合进行多项式拟合的情况。
2. lsqcurvefit:用于非线性最小二乘法拟合的函数,可以通过
最小二乘法拟合一组数据点的非线性函数曲线。
适用范围:适合进行非线性函数拟合的情况。
3. fminsearch:用于寻找函数的最小值的函数,可以通过优化
算法寻找最适合数据点的参数值。
适用范围:适合进行简单的参数优化的情况。
4. fmincon:用于带约束条件的优化问题的函数,可以通过优
化算法寻找最适合数据点的参数值,并满足约束条件。
适用范围:适合进行带约束条件的参数优化的情况。
5. lsqlin:用于带等式约束的线性最小二乘法拟合的函数,可
以通过最小二乘法拟合一组数据点的线性函数曲线,并满足等式约束。
适用范围:适合进行带等式约束的线性函数拟合的情况。
通过这些 Matlab 命令,可以进行最小二乘法拟合和优化,并
得到最适合数据的拟合曲线或参数值。
具体使用哪个命令取决于数据的特点和问题的需求。
最小二乘法(附MATLAB代码)今天我主要是从如何使用MATLAB实现最小二乘法,首先给出今天重点使用的两个函数。
比如我想拟合下面这组数据x=[9,13,15,17,18.6,20,23,29,31.7,35];y=[-8,-6.45,-5.1,-4,-3,-1.95,-1.5,-0.4,0.2,-0.75];我先用matlab将这组离散点画出来,plot(x,y,'o')嗯,大概这个样子,这时我们想使用一次函数拟合上述曲线,可使用以下代码clearclcx=[9,13,15,17,18.6,20,23,29,31.7,35];y=[-8,-6.45,-5.1,-4,-3,-1.95,-1.5,-0.4,0.2,-0.75];coeff icient=polyfit(x,y,1); %用一次函数拟合曲线,想用几次函数拟合,就把n设成那个数y1=polyval(coefficient,x);%plot(x,y,'-',x,y1,'o'),这个地方原来'-'和'o'写反了,现已更正,可以得到正确的图形。
plot(x,y,'o',x,y1,'-')得到的结果是coefficient=[0.2989,-9.4107]所以得到的一次函数为y=0.2989*x-9.4107同理如果用二次函数拟合该曲线,得到的各项系数为coefficient=[-0.0157 1.0037 -16.2817]所以得到的二次函数为y=-0.0157*x^2+1.0037*x-16.2817其他阶数依此类推。
但是使用polyfit(x,y,n)函数有一个注意事项:举个例子,比如说我们想用9阶多项式拟合上述曲线时,我们发现拟合的曲线是正常的,得到的各项系数也是正常的但是当我们用10阶多项式拟合曲线时,此时各项系数如下,得到的曲线如下很明显出现了问题,所以使用polyfit(x,y,n)函数时要严格遵守上述事项。
最小二乘法及matlab 程序最小二乘法简介:最小二乘法(又称最小平方法)是一种数学优化技术。
它通过最小化误差的平方和寻找数据的最佳函数匹配。
利用最小二乘法可以简便地求得未知的数据,并使得这些求得的数据与实际数据之间误差的平方和为最小。
最小二乘法还可用于曲线拟合。
其他一些优化问题也可通过最小化能量或最大化熵用最小二乘法来表达。
最小二乘法的矩阵形式:Ax b = 若未知量的个数大于方程的个数,则方程无解,但在数值计算领域,我们通常是计算min ||||Ax b - ,解出其中的x 。
比较直观的做法是求解T T A Ax A b = ,但通常比较低效。
其中一种常见的解法是对 A 进行QR 分解(A=QR ),其中Q 是n k ⨯正交矩阵(Orthonormal Matrix ),R 是n k ⨯上三角矩阵(Upper Triangular Matrix ),则有:1min min min Ax b QRx b Rx Q b --=-=-。
Matlab 命令:一次函数线性拟合使用polyfit (x,y,1);多项式函数线性拟合使用 polyfit (x,y,n ),n 为次数;例如:x=[0.5,1.0,1.5,2.0,2.5,3.0]; y=[1.75,2.45,3.81,4.80,7.00,8.60]; p=polyfit(x,y,2) x1=0.5:0.5:3.0; y1=polyval(p,x1); plot(x,y,'*r',x1,y1,'-b') 计算结果为:p =0.5614 0.8287 1.1560即所得多项式为y=0.5614x^2+0.8287x+1.15560然后可以使用polyval在t处预测:y_hat=polyval(p,t)非线性函数使用:lsqcurvefit、nlinfit格式:lsqcurvefit(f,a,x,y)、nlinfit(x,y,f,a)f:符号函数句柄,如果是以m文件的形式调用的时候,别忘记加@.这里需要注意,f函数的返回值是和y匹对的,即拟合参数的标准是(f-y)^2取最小值,具体看下面的例子a:最开始预估的值(预拟合的未知参数的估计值)。
matlab最小二乘法拟合求参数
最小二乘法是一种数据拟合的常用方法,可以求得一组参数使得拟合函数与给定数据的残差平方和最小。
在Matlab中,可以通过以下步骤求解最小二乘法拟合的参数:
1. 输入数据:首先,将需要拟合的数据输入到Matlab中,例如,可以创建两个向量x和y来表示一组二维数据。
2. 选择拟合函数:根据数据的特点选择一个合适的拟合函数形式,例如,线性、二次、指数等。
假设选择线性拟合y = a*x + b。
3. 构建拟合方程:根据选择的拟合函数形式,构建拟合方程,即根据给定的数据和参数a、b,计算预测的y值。
4. 残差计算:计算预测值与实际值之间的差异,即残差。
可以使用Matlab的内置函数或者编写自定义函数来计算残差。
5. 残差平方和最小化:根据最小二乘法的原理,目标是使得残差平方和最小化。
可以使用Matlab的内置函数或者编写自定义函数来求解最小二乘法的参数。
6. 求解参数:使用最小化残差平方和的方法,求解拟合方程的参数。
在Matlab中,可以使用lsqcurvefit函数或者lsqnonlin函数等进行求解。
7. 结果评估:根据求解得到的参数,计算拟合方程在给定数据上的拟合度,可以计算相关系数等来评估拟合效果。
以上就是使用Matlab进行最小二乘法拟合求解参数的一般步骤。
具体的实现方法可以根据数据和拟合函数的不同进行调整和优化。
---------------------------------------------------------------最新资料推荐------------------------------------------------------ 用MatLab画图(最小二乘法做曲线拟合) 用 MatLab 画图(最小二乘法做曲线拟合) 帮朋友利用实验数据画图时,发现 MatLab 的确是画图的好工具,用它画的图比Excel光滑、精确。
利用一组数据要计算出这组数据对应的函数表达式从而得到相应图像,MatLab 的程序如下:x=[1 5 10 20 30 40 60 80] y=[15. 4 33. 9 42. 2 50. 556 62. 7 72 81. 1] plot(x, y, ‘ r*’ ) ; legend(‘ 实验数据(xi, yi) ‘ ) xlabel(‘ x’ ) , ylabel(‘ y’ ) , title(‘ 数据点(xi, yi) 的散点图’ ) syms a1 a2 a3 x=[15 10 20 30 40 60 80]; fi=a1. *x. +a2. *x+a3 y=[15. 4 33.9 42. 2 50. 5 56 62. 7 72 81. 1] fi =[a1+a2+a3,25*a1+5*a2+a2+(400*a1+20*a2+a3-101/2) +(900*a1+30*a2+a3-56) +(1600*a1+40*a2+a3-627/10) +(3600*a1+60*a2+a3-72)+(6400*a1+80*a2+a3-811/10) ; Ja1=diff(J, a1) ;Ja2=diff(J, a2) ; Ja3=diff(J, a3) ; Ja11=simple(Ja1) ,Ja21=simple(Ja2) , Ja31=simple(Ja3) A=[114921252, 1656252, 26052; 1656252, 26052, 492; 26052, 492, 16]; B=[9542429/5, 166129/5, 4138/5]; C=B/A, f=poly2sym(C) xi=[1 5 10 20 3040 60 80] ; y=[15. 4 33. 9 42. 2 50. 5 56 62. 7 72 81. 1]; n=length(xi) ; f=-0. 0086. *xi. +1. 3876. *xi+23. 1078;1 / 6x=1: 1/10: 80; F=-0. 0086. *x. +1. 3876. *x+23. 1078; fy=abs(f-y) ; fy2=fy. ; Ew=max(fy) , E1=sum(fy) /n,E2=sqrt((sum(fy2) ) /n) plot(xi, y, ‘ r*’ ) , hold on, plot(x, F, ‘ b-’ ) , hold off legend(‘ 数据点(xi, yi) ‘ , ‘ 拟合曲线f(x) = -0. 0086x +1. 3876x+23. 1078’ ) , xlabel(‘ x’ ) , ylabel(‘ y’ ) , title(‘ 实验数据点(xi, yi) 及拟合曲线f(x) ‘ ) 下图是程序运行后得到的:Su7Tw8VxaW ybXAcZB d#Cf!Eg%FhGj*Ik(Jl-Kn+M o0Np2Or3Ps4R t6Sv7Tw8VxaWzbXAcZBe#Cf! Eg%Fi Gj*Ik)Jl-Kn+Mo1Np2Or3Qs4Rt6Sv7 Uw8Vx aWzbYAc ZBe#Df!Eg%FiHj*Ik) Jm-Kn +Mo1Nq2Or3Qs 5Rt6Sv7Uw9VxaWzbYAdZBe#D f$Eg%F iHj(I k) Jm-Ln+Mo1Nq2Pr3Qs5Ru6S v7Uw9V yaWzbY AdZCe#Df$Eh%FiHj(Il) Jm-Ln0Mo1Nq2Pr4 Qs5Ru6Tv8Uw9VyaXzbYAdZCe !Df$Eh %GiHj (Il) Km-Ln0Mp1Nq2Pr4Qt5Ru 6Tv8U x9VyaXz cYAdZCe! Dg$Eh%Gi*Hj(Il) Km+Ln0M p1Oq2P r4Qt5Su6Tv8Ux9WyaXzcYBdZ Ce!Dg$Fh%Gi* Hk(Il) Km+Lo0Mp1Oq3Pr4Qt5 Su7Tv8Ux9Wyb XzcYBd#Ce!Dg$FhGi*Hk(Jl ) Km+L o0Np1Oq 3Ps4Rt 5Su7Tw8Ux9WybXAcY Bd#Cf!Dg$FhGj*Hk(Jl-Km+Lo0Np2Oq3Ps4 Rt6Su7Tw8Vx9 WybXAcZBd#Cf!Eg$FhGj*Ik (Jl-Kn +Lo0Np2Or3Ps4Rt6Sv7Tw8VxaWybXA cZBe#Cf!Eg%F hGj*Ik) Jl-K n+Mo0Np2Or3Q s4Rt6Sv 7Uw8V xaWzbXAcZBe# D f! Eg%FiGj* Ik) Jm- Kn+Mo1 Nq2Or3Qs5Rt6Sv7Uw9VxaWzb YAcZBe#Df$Eg %FiHj*Ik)---------------------------------------------------------------最新资料推荐------------------------------------------------------ Jm-Ln+Mo1Nq2Pr 3Qs5Ru 6Sv7Uw 9VyaWzbYAdZBe#Df$Eh%FiH j(Ik) J m-Ln0M o1Nq2Pr4Qs5Ru6Tv7Uw9VyaX zbYAdZC e#Df$ Eh%GiHj(Il) Jm-Ln0Mp1Nq2 Pr4Qt5Ru6Tv8 Uw9VyaXzcYAdZCe!Df$Eh%Gi *Hj(Il) Km+Ln 0Mp1Oq2Pr4Qt5Su6Tv8Ux9Vy aXzcYB dZCe!D g$Eh%Gi*Hk(Il) Km+Lo0Mp1O q3Pr4Qt5Su7Tv8Ux9WyaXzc Y Bd#Ce!Dg$Fh %Gi*Hk( Jl) Km +Lo0Np1Oq3Ps 4 Qt5Su7Tw8Ux 9WybXzcYBd#C f!Dg$FhGi*H k (Jl-Km+Lo0N p2Oq3Ps4Rt5S u7Tw8Vx9WybX AcYBd#Cf! Eg$ FhGj*Ik (Jl- Kn+Lo0Np2Or3 P s4Rt6Su7Tw8 VxaWybXA cZBd #Cf!Eg%FhGj * Ik) Jl-Kn+Mo 0Np2Or3Qs4Rt 6Sv7Tw8VxaWz bXAcZBe#Cf!E g%FiGj*Ik) J m-Kn+Mo1Np2O r 3Qs5Rt6Sv7U w8VxaWzbYAcZ Be#Df! Eg%Fi H j*Ik) Jm-Ln+ Mo1Nq2O r3Qs5 Ru6Sv7Uw9Vxa W zbYAdZBe#Df $Eh%Fi Hj(Ik ) Jm-Ln0Mo1Nq 2Pr3Qs5Ru6Tv 7Uw9Vya WzbYA dZCe#Df$Eh%G iHj(Il) Jm-L n0Mp1Nq2Pr4Q s5Ru6Tv8Uw9V ya XzbYAdZCe! Df$Eh%Gi*Hj( Il) Km-Ln0Mp1 Oq 2Pr4 Qt5Ru6Tv8Ux9Vy aXz cYAdZCe!Dg$E h%G i*Hk(Il) K m+Ln0Mp1O q3P r4Qt5Su6Tv8U x9WyaXzcYBd# Ce!Dg$Fh% Gi* Hk(Jl) Km+Lo0 Mp1Oq3Ps4Qt5 Su7Tv8Ux9Wyb XzcYBd#Cf! Dg $F hGi*Hk(Jl -Km+Lo0N p1Oq 3Ps4Rt5Su7Tw 8U x9WybXAcYB d#Cf!Eg$F hG j*Hk(Jl-Kn+Lo0Np2O q3Ps4Rt 6Su7Tw8Vx9WybXAcZBd#Cf!E g%FhGj*Ik(J l-Kn+Mo0Np2O r3Ps4Rt6Sv7T w8Vxa WzbXAcZ Be#Cf! Eg%Fi Gj*Ik) Jl-Kn+ Mo1Np2Or3Qs4 Rt6Sv7Uw8VxaWzbYAcZBe#Df !Eg%FiHj*Ik )3 / 6Jm-Kn+Mo1Nq2Or3Qs5Rt6Sv 7Uw9Vx aWzbYA dZBe#Df$Eg%FiHj(Ik) Jm-L n+Mo1Nq2Pr3Q s5Ru6Sv7Uw9VyaWzbYAdZCe# Df$Eh %FiHj( Il) Jm-Ln0Mo1Nq2Pr4Qs5Ru6 Tv8Uw9VyaXzb YAdZCe!Df$Eh%GiHj(Il) Km -Ln0Mp 1Nq2Pr 4Qt5Ru6Tv8Ux9VyaXzcYAdZC e!Dg$E h%Gi*H j(Il) Km+Ln0Mp1Oq2Pr4Qt5S u6Tv8U x9WyaX zcYBdZCe! Dg$Fh%Gi*Hk(Il) Km+Lo0Mp1Oq3 Pr4Qt5Su7Tv8Ux9WybXzcYBd #Ce!D g$FhGi *Hk(Jl ) Km+Lo0Np1Oq3Ps4R t5Su7Tw8Ux9T v7Uw9VyaXzbYAdZCe#Df$Eh% GiHj( Il) Jm- Ln0Mp1Nq2Pr4Qs5Ru6Tv8Uw9 VyaXzcY AdZCe !Df$Eh%Gi*Hj(Il) Km-Ln0Mp 1Oq2Pr 4Qt5Ru 6Tv8Ux9VyaXzcYBdZCe! Dg$E h%Gi*Hk(Il) K m+Ln0Mp1Oq3Pr4Qt5Su6Tv8U x9WyaX zcYBd# Ce!Dg$Fh%Gi*Hk(Jl) Km+Lo0 Mp1Oq3P s4Qt5 Su7Tw8Ux9Wyb X zcYBd#Cf!Dg $FhGi*Hk(Jl -Km+Lo0Np1Oq3Ps4Rt5Su7Tw 8Vx9Wy bXAcYB d#Cf! Eg$FhGj*Hk(Jl-Kn+L o0Np2O q3Ps4R t6Su7Tw8VxaWybXAcZBd#Cf! Eg%Fh Gj*Ik( Jl-Kn+Mo0Np2Or3Ps4Rt6Sv7 Tw8Vxa WzbXAc ZBe#Cf!Eg%FiGj*Ik) Jl-Kn +Mo1Np2Or3Qs 5Rt6Sv7Uw8Vx a WzbYAcZBe#D f! Eg%FiHj*Ik) Jm-Kn+Mo1 Nq 2Or3Qs5Ru6 Sv7Uw9Vx aWzb YAdZBe#Df$Eg %F iHj(Ik) Jm -Ln+Mo1N q2Pr 3Qs5Ru6Tv7Uw 9VyaWzbYAdZC e#Df$Eh%FiH j(Il) Jm-Ln0M o1Nq2Pr4Qs5R u6Tv8Uw9VyaX zbYAdZCe! Df$ E h%GiHj(Il) Km-Ln0Mp 1Oq2 Pr4Qt5Ru6Tv8 U x9VyaXzcYAd ZCe!Dg$E h%Gi *Hj(Il) Km+Ln 0Mp1Oq3Pr4Qt 5Su6Tv8Ux9Wy aXzcYBdZCe!D g $Fh%Gi*Hk(I l)---------------------------------------------------------------最新资料推荐------------------------------------------------------ Km+Lo0Mp1O q3Ps4Qt5Su7T v 8Ux9WybXzcY Bd#Ce! D g$Fh Gi*Hk(Jl) Km+ L o0Np1Oq3Ps4 Rt5Su7T w8Ux9 WybXAcYBd#Cf !Dg$FhGj*Hk (Jl-Km+ Lo0Np 2Oq3Ps4Rt6Su 7Tw8Vx9WybXA cZBd#Cf ! Eg$F hGj*Ik(Jl-K n+Mo0Np2Or3P s4Rt6Sv 7Tw8V xaWybXAcZBe# C f!Eg% FhGj*Ik) Jl-K n+Mo 1Np2Or3Qs4Rt 6Sv 7Uw8VxaWz bXAcZBe#D f!E g%FiGj*Ik) J m- Kn+Mo1Nq2O r3Qs5Rt6Sv7U w9VxaWzbYAcZ Be#Df$Eg%Fi Hj*Ik) Jm-Ln+ Mo1Nq2Pr3Qs5 Ru6Sv7Uw9Vya WzbYAdZBe#Df $Eh%FiHj(Il ) Jm -Ln0Mo1Nq 2Pr4Qs5Ru 6Tv 7Uw9VyaXzbYA dZCe#Df$Bd#C f! Eg%FhGj*I k(Jl-Kn+Mo0Np2O r3Qs4Rt 6Sv7Tw8VxaWzbXAcZBe#Cf!E g%FiGj*Ik) J l-Kn+Mo1Np2O r3Qs5Rt6Sv7U w8Vxa WzbYAcZ Be#Df! Eg%Fi Hj*Ik) Jm-Kn+ Mo1Nq2Or3Qs5 Ru6Sv7Uw9VxaWzbYAdZBe#Df $Eg%FiHj(Ik ) Jm-Ln0Mo1Nq2Pr3Qs5Ru6Tv 7Uw9Vy aWzbYA dZCe#Df$Eh%FiHj(Il) Jm-L n0Mp1Nq2Pr4Q s5Ru6Tv8Uw9VyaXzbYAdZCe! Df$Eh %GiHj( Il) Km-Ln0Mp1Oq2Pr4Qt5Ru6 Tv8Ux9VyaXzc YAdZCe!Dg$Eh%Gi*Hj(Il) Km +Ln0Mp 1Oq3Pr 4Qt5Su6Tv8Ux9WyaXzcYBdZC e!Dg$F h%Gi*H k(Jl) Km+Lo0Mp1Oq3Ps4Qt5S u7Ts5R u6Sv7U w9VyaWzbYAdZBe#Df$Eh%Fi Hj(Ik)Jm-Ln0 Mo1Nq2Pr4Qs5Ru6Tv7Uw9Vya XzbYA dZCe#Df $Eh%Gi Hj(Il) Jm-Ln0Mp1N q2Pr4Qt5Ru6T v8Uw9VyaXzcYAdZCe!Df$Eh% Gi*Hj( Il) Km- Ln0Mp1Oq2Pr4Qt5Su6Tv8Ux9 VyaXzcY BdZCe !Dg$Eh%Gi*Hk(Il) Km+Ln0Mp 1Oq3Pr4Qt5Su5 / 67Tv8Ux9WyaXzcYBd#Ce! Dg$F h%Gi*Hk (Jl) K m+Lo0Np1Oq3Ps4Qt5Su7Tw8U x9WybX zcYBd# Cf!Dg$FhGi*Hk(Jl-Km+Lo0 Np2Oq3Ps4Rt5 Su7Tw8Vx9WybXAcYBd#Cf!Eg $FhGj*Hk(Jl -Kn+Lo0Np2Or3Ps4Rt6Su7Tw 8VxaWy bXAcZB d#Cf! Eg%FhG j*Ik(Jl-Kn+M o0Np2O r3Qs4R t6Sr4Qt5Su6Tv8Ux9WyaXzcY BdZCe!Dg$Fh% Gi*Hk(Il) Km+Lo0Mp1Oq3Pr4 Qt5Su7Tv8Ux9 WybXzcYBd#Ce! Dg$FhGi*Hk (Jl) Km +Lo0Np 1Oq3Ps4Rt5Su7Tw8Ux9WybXA cYBd#Cf!Dg$FhGj*Hk(Jl- K m+Lo0Np2Oq3 Ps4Rt6Su7Tw8 Vx9WybXAcZBd # Cf!Eg$FhGj *Ik(Jl- Kn+Lo 0Np2Or3Ps4Rt 6Sv7Tw8VxaWy bXAcZBe#Cf!E g%FhGj*Ik) J l-Kn+Mo0Np2O r3Qs4Rt6Sv7U w8VxaWzbXAcZ B e#Df!Eg%Fi Gj*Ik) J m-Kn+ Mo1Nq2Or3Qs5 R t6Sv7Uw9Vxa WzbYAcZB e#Df $Eg%Ff! Dg$Fh Gi*Hk(Jl-Km +Lo0Np1Oq3Ps 4Rt5Su7Tw8Vx 9W ybXAcYBd#C f!Eg$Fh Gj*H k(Jl-Kn+Lo0N p2Oq3Ps4Rt6S u7Tw8Vxa WybX AcZBd#Cf!Eg% F hGj*Ik(Jl- Kn+Mo0N p2Or3 Ps4Rt6Sv7Tw8 V xaWzbXAcZBe #Cf!Eg% FiGj *Ik) Jl-Kn+Mo 1Np2Or3Qs5Rt 6Sv7Uw8V xaWz bYAcZBe#Df!E g%FiHj*Ik) J m-Kn+Mo1Nq2O r3Qs5Ru6Sv7U w9VxaW zbYAdZBe#Df$Eg %Fi Hj(Ik) Jm-Ln +M o1Nq2Pr3Qs 5Ru6Tv7U w9Vy aWzbYAdZCe#D f$Eh%FiHj(I l) Jm-Ln0Mo1N q2Pr4Qs5Or3P s4Rt6Su7Tw8V xaWybXAcZ Be# Cf! Eg%FhGj* I。
曲线拟合(curve-fitting ):工程实践中,用测量到的一些离散的数据},...2,1,0),,{(m i y x i i =求一个近似的函数)(x ϕ来拟合这组数据,要求所得的拟合曲线能最好的反映数据的基本趋势(即使)(x ϕ最好地逼近()x f ,而不必满足插值原则。
因此没必要取)(i x ϕ=i y ,只要使i i i y x -=)(ϕδ尽可能地小)。
原理:给定数据点},...2,1,0),,{(m i y x i i =。
求近似曲线)(x ϕ。
并且使得近似曲线与()x f 的偏差最小。
近似曲线在该点处的偏差i i i y x -=)(ϕδ,i=1,2,...,m 。
常见的曲线拟合方法:1.使偏差绝对值之和最小2.使偏差绝对值最大的最小3.使偏差平方和最小最小二乘法:按偏差平方和最小的原则选取拟合曲线,并且采取二项式方程为拟合曲线的方法,称为最小二乘法。
推导过程:1. 设拟合多项式为:2. 各点到这条曲线的距离之和,即偏差平方和如下:3. 问题转化为求待定系数0a ...k a 对等式右边求i a 偏导数,因而我们得到了: .......4、 把这些等式化简并表示成矩阵的形式,就可以得到下面的矩阵:5. 将这个范德蒙得矩阵化简后可得到:6. 也就是说X*A=Y ,那么A = (X'*X)-1*X'*Y ,便得到了系数矩阵A ,同时,我们也就得到了拟合曲线。
MATLAB 实现:MATLAB 提供了polyfit ()函数命令进行最小二乘曲线拟合。
调用格式:p=polyfit(x,y,n)[p,s]= polyfit(x,y,n)[p,s,mu]=polyfit(x,y,n)x,y 为数据点,n 为多项式阶数,返回p 为幂次从高到低的多项式系数向量p 。
x 必须是单调的。
矩阵s 包括R (对x 进行QR 分解的三角元素)、df(自由度)、normr(残差)用于生成预测值的误差估计。
2 曲线拟合的线性最小二乘法及其MATLAB 程序例2 给出一组数据点),(i i y x 列入表2中,试用线性最小二乘法求拟合曲线,估计其误差,作出拟合曲线.解 (1)在MATLAB 工作窗口输入程序>> x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.126.50 68.04];plot(x,y,'r*'),legend('实验数据(xi,yi)')xlabel('x'), ylabel('y'),title('数据点(xi,yi)的散点图')运行后屏幕显示数据的散点图(略).(3)编写下列MATLAB 程序计算)(x f 在),(i i y x 处的函数值,即输入程序>> syms a1 a2 a3 a4x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];fi=a1.*x.^3+ a2.*x.^2+ a3.*x+ a4运行后屏幕显示关于a 1,a 2, a 3和a 4的线性方程组fi =[ -125/8*a1+25/4*a2-5/2*a3+a4,-4913/1000*a1+289/100*a2-17/10*a3+a4,-1331/1000*a1+121/100*a2-11/10*a3+a4,-64/125*a1+16/25*a2-4/5*a3+a4,a4, 1/1000*a1+1/100*a2+1/10*a3+a4,27/8*a1+9/4*a2+3/2*a3+a4,19683/1000*a1+729/100*a2+27/10*a3+a4,5832/125*a1+324/25*a2+18/5*a3+a4]编写构造误差平方和的MATLAB 程序>> y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.126.50 68.04];fi=[-125/8*a1+25/4*a2-5/2*a3+a4,-4913/1000*a1+289/100*a2-17/10*a3+a4,-1331/1000*a1+121/100*a2-11/10*a3+a4,-64/125*a1+16/25*a2-4/5*a3+a4, a4, 1/1000*a1+1/100*a2+1/10*a3+a4,27/8*a1+9/4*a2+3/2*a3+a4,19683/1000*a1+729/100*a2+27/10*a3+a4,5832/125*a1+324/25*a2+18/5*a3+a4];fy=fi-y; fy2=fy.^2; J=sum(fy.^2)运行后屏幕显示误差平方和如下J=(-125/8*a1+25/4*a2-5/2*a3+a4+1929/10)^2+(-4913/1000*a1+289/100*a2-17/10*a3+a4+171/2)^2+(-1331/1000*a1+121/100*a2-11/10*a3+a4+723/20)^2+(-64/125*a1+16/25*a2-4/5*a3+a 4+663/25)^2+(a4+91/10)^2+(1/1000*a1+1/100*a2+1/10*a3+a4+843/100)^2+(27/8*a1+9/4*a2+3/2*a3+a4+328/25)^2+(19683/1000*a1+729/100*a2+27/10*a3+a4-13/2)^2+(5832/125*a1+324/25*a2+18/5*a3+a4-1701/25)^2为求4321,,,a a a a 使J 达到最小,只需利用极值的必要条件0=∂∂ka J )4,3,2,1(=k ,得到关于4321,,,a a a a 的线性方程组,这可以由下面的MATLAB 程序完成,即输入程序>> syms a1 a2 a3 a4J=(-125/8*a1+25/4*a2-5/2*a3+a4+1929/10)^2+(-4913/1000*a1+289/100*a2-17/10*a3+a4...+171/2)^2+(-1331/1000*a1+121/100*a2-11/10*a3+a4+723/20)^2+(-64/125*a1+16/25*a2-4/5*a3+a4+663/25)^2+(a4+91/10)^2+(1/1000*a1+1/100*a2+1/10*a3+a4+843/100)^2+(27/8*a1+9/4*a2+3/2*a3+a4+328/25)^2+(19683/1000*a1+729/100*a2+27/10*a3+a4-13/2)^2+(5832/125*a1+324/25*a2+18/5*a3+a4-1701/25)^2;Ja1=diff(J,a1); Ja2=diff(J,a2); Ja3=diff(J,a3); Ja4=diff(J,a4);Ja11=simple(Ja1), Ja21=simple(Ja2), Ja31=simple(Ja3), Ja41=simple(Ja4),运行后屏幕显示J 分别对a 1, a 2 ,a 3 ,a 4的偏导数如下Ja11=56918107/10000*a1+32097579/25000*a2+1377283/2500*a3+23667/250*a4-8442429/625Ja21 =32097579/25000*a1+1377283/2500*a2+23667/250*a3+67*a4+767319/625Ja31 =1377283/2500*a1+23667/250*a2+67*a3+18/5*a4-232638/125Ja41 =23667/250*a1+67*a2+18/5*a3+18*a4+14859/25解线性方程组Ja 11 =0,Ja 21 =0,Ja 31 =0,Ja 41 =0,输入下列程序>>A=[56918107/10000, 32097579/25000, 1377283/2500, 23667/250; 32097579/25000, 1377283/2500, 23667/250, 67; 1377283/2500, 23667/250, 67, 18/5; 23667/250, 67, 18/5, 18];B=[8442429/625, -767319/625, 232638/125, -14859/25]; C=B/A, f=poly2sym(C)运行后屏幕显示拟合函数f 及其系数C 如下C = 5.0911 -14.1905 6.4102 -8.2574f=716503695845759/140737488355328*x^3-7988544102557579/562949953421312*x^2+####77693/28####56*x-4648521160813215/562949953421312故所求的拟合曲线为8.25746.410214.19055.0911)(23-+-=x x x x f .(4)编写下面的MATLAB 程序估计其误差,并作出拟合曲线和数据的图形.输入程序>> xi=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6];y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04]; n=length(xi);f=5.0911.*xi.^3-14.1905.*xi.^2+6.4102.*xi -8.2574;x=-2.5:0.01: 3.6;F=5.0911.*x.^3-14.1905.*x.^2+6.4102.*x -8.2574;fy=abs(f-y); fy2=fy.^2; Ew=max(fy),E1=sum(fy)/n, E2=sqrt((sum(fy2))/n)plot(xi,y,'r*'), hold on, plot(x,F,'b-'), hold offlegend('数据点(xi,yi)','拟合曲线y=f(x)'),xlabel('x'), ylabel('y'),title('例2的数据点(xi,yi)和拟合曲线y=f(x)的图形')运行后屏幕显示数据),(i i y x 与拟合函数f 的最大误差E w ,平均误差E 1和均方根误差E 2及其数据点),(i i y x 和拟合曲线y =f (x )的图形(略).Ew = E1 = E2 =3.105 4 0.903 4 1.240 96 函数逼近及其MATLAB 程序最佳均方逼近的MATLAB 主程序function [yy1,a,WE]=zjjfbj(f,X,Y,xx)m=size(f);n=length(X);m=m(1);b=zeros(m,m);c=zeros(m,1);if n~=length(Y)error('X 和Y 的维数应该相同')endfor j=1:mfor k=1:mb(j,k)=0;for i=1:nb(j,k)=b(j,k)+feval(f(j,:),X(i))*feval(f(k,:),X(i));endendc(j)=0;for i=1:nc(j)=c(j)+feval(f(j,:),X(i))*Y(i);endenda=b\c;WE=0;for i=1:nff=0;for j=1:mff=ff+a(j)*feval(f(j,:),X(i));endWE=WE+(Y(i)-ff)*(Y(i)-ff);endif nargin==3return ;endyy=[];for i=1:ml=[];for j=1:length(xx)l=[l,feval(f(i,:),xx(j))];endyy=[yy l'];endyy=yy*a; yy1=yy'; a=a';WE;例6.1对数据X 和Y , 用函数2,,1x y x y y ===进行逼近,用所得到的逼近函数计算在 6.5=x 处的函数值,并估计误差.其中X =(1 3 4 5 6 7 8 9); Y =(-11 -13 -11 -7 -1 7 17 29). 解 在MATLAB 工作窗口输入程序>> X=[ 1 3 4 5 6 7 8 9]; Y=[-11 -13 -11 -7 -1 7 17 29];f=['fun0';'fun1';'fun2'];[yy,a,WE]=zjjfbj(f,X,Y,6.5)运行后屏幕显示如下yy =2.75000000000003a =-7.00000000000010 -4.99999999999995 1.00000000000000WE =7.172323350269439e-027例 6.2 对数据X 和Y ,用函数2,,1x y x y y ===,x y cos =,=y e x ,x y sin =进行逼近,其中X =(0 0.50 1.00 1.50 2.00 2.50 3.00),Y =(0 0.4794 0.8415 0.9815 0.9126 0.5985 0.1645).解 在MATLAB 工作窗口输入程序>> X=[ 0 0.50 1.00 1.50 2.00 2.50 3.00];Y=[0 0.4794 0.8415 0.9815 0.9126 0.5985 0.1645];f=['fun0';'fun1';'fun2';'fun3';'fun4';'fun5'];xx=0:0.2:3;[yy,a,WE]=zjjfbj(f,X,Y, xx), plot(X,Y,'ro',xx,yy,'b-')运行后屏幕显示如下(图略)yy = Columns 1 through 7-0.0005 0.2037 0.3939 0.5656 0.7141 0.8348 0.9236Columns 8 through 140.9771 0.9926 0.9691 0.9069 0.8080 0.6766 0.5191Columns 15 through 160.3444 0.1642a = 0.3828 0.4070 -0.3901 0.0765 -0.4598 0.5653WE = 1.5769e-004即,最佳逼近函数为y=0.3828+0.4070*x-0.3901*x^2+0.0765*exp(x) -0.4598*cos(x) +0.5653*sin(x).8 随机数据点上的二元拟合及其MATLAB 程序例8 设节点(X,Y,Z )中的X 和Y 分别是在区间]3,3[-和]5.3,5.2[-上的50个随机数,Z 是函数Z =7-3x 3e 22y - -x 在(X,Y )的值,拟合点(X I ,Y I )中的X I =-3:0.2:3, Y I =-2.5:0.2:3.5.分别用二元拟合方法中最近邻内插法、三角基线性内插法、三角基三次内插法和MATLAB 4网格化坐标方法计算在(X I ,Y I )处的值,作出它们的图形,并与被拟和曲面进行比较.解 (1)最近邻内插法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y , x ,y . X=-3+(3-(-3))*x;%利用x 生成的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点(X,Y )处计算Z 的值.X1=-3:0.2:3;Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标(XI,YI )网格化. ZI=griddata(X,Y,Z,XI,YI, 'nearest') %计算在每个插值点(XI,YI )处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用最近邻内插法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z).hold of %结束在当前图形上添加新图形. 运行后屏幕显示用最近邻内插法拟合函数Z =7-3x 3e 22y - -x 在两组不同节点处的曲面及其插值Z I (略).(2)三角基线性内插法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y , x ,y . X=-3+(3-(-3))*x;%利用x 生成 上的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成 上的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点(X,Y )处计算Z 的值.X1=-3:0.2:3;Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标(XI,YI )网格化. ZI=griddata(X,Y,Z,XI,YI, 'linear') %计算在每个插值点(XI,YI )处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用三角基线性内插法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z).hold of %结束在当前图形上添加新图形. 运行后屏幕显示用三角基线性内插法拟合函数Z =7-3x 3e 22y - -x 在两组不同节点处的曲面和节点的图形及其插值Z I (略).(3)三角基三次内插法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y , x ,y . X=-3+(3-(-3))*x;%利用x 生成 上的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成 上的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点(X,Y )处计算Z 的值.X1=-3:0.2:3;Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标(XI,YI )网格化. ZI=griddata(X,Y,Z,XI,YI, 'cubic') %计算在每个插值点(XI,YI )处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用三角基三次内插法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z).hold of %结束在当前图形上添加新图形. 运行后屏幕显示用三角基三次内插法拟合函数Z =7-3x 3e 22y - -x 在两组不同节点处的曲面和节点的图形及其插值Z I (略).(4)MATLAB 4网格化坐标方法.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y , x ,y .X=-3+(3-(-3))*x;%利用x 生成 上的随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成 上的随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点(X,Y )处计算Z 的值.X1=-3:0.2:3; Y1=-2.5:0.2:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标(XI,YI )网格化. ZI=griddata(X,Y,Z,XI,YI, 'v4') %计算在每个插值点(XI,YI )处的插值ZI.mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('用MATLAB 4网格化坐标方法拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('拟合曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z).hold of %结束在当前图形上添加新图形. 运行后屏幕显示用MATLAB 4网格化坐标方法拟合函数Z =7-3x 3e 22y - -x 在两组不同节点处的曲面和节点的图形及其插值ZI (略).(5)作被拟合曲面Z =7-3x 3e 22y - -x 和节点的图形.输入程序>> x=rand(50,1);y=rand(50,1); %生成50个一元均匀分布随机数x 和y , x ,y . X=-3+(3-(-3))*x;%利用x 生成随机变量.Y=-2.5+(3.5-(-2.5))*y; %利用y 生成随机变量.Z=7-3* X.^3 .* exp(-X.^2 - Y.^2); %在每个随机点(X,Y )处计算Z 的值.X1=-3.:0.1:3.;Y1=-2.5:0.1:3.5;[XI,YI] = meshgrid(X1,Y1); %将坐标(XI,YI )网格化. ZI=7-3* XI.^3 .* exp(-XI.^2 - YI.^2);mesh(XI,YI, ZI) %作二元拟合图形.xlabel('x'), ylabel('y'), zlabel('z'),title('被拟合函数z =7-3 x^3 exp(-x^2 - y^2) 的曲面和节点的图形')%legend('被拟合函数曲面','节点(xi,yi,zi)')hold on %在当前图形上添加新图形.plot3(X,Y,Z, 'bo') %用兰色小圆圈画出每个节点(X,Y,Z).hold of %结束在当前图形上添加新图形. 运行后屏幕显示被拟合函数Z =7-3x 3e22y - -x 的曲面和节点的图形及其函数值ZI(略).。
方便大家使用的最小二乘法曲线拟合的Matlab程序非常方便用户使用,直接按提示操作即可;这里我演示一个例子:(红色部分为用户输入部分,其余为程序运行的结果,结果图为Untitled.fig,Untitled2.fig) 请以向量的形式输入x,y.x=[1,2,3,4]y=[3,4,5,6]通过下面的交互式图形,你可以事先估计一下你要拟合的多项式的阶数,方便下面的计算.polytool()是交互式函数,在图形上方[Degree]框中输入阶数,右击左下角的[Export]输出图形回车打开polytool交互式界面回车继续进行拟合输入多项式拟合的阶数m = 4Warning: Polynomial is not unique; degree >= number of data points. > In polyfit at 72In zxecf at 64输出多项式的各项系数a = 0.0200000000000001a = -0.2000000000000008a = 0.7000000000000022a = 0.0000000000000000a = 2.4799999999999973输出多项式的有关信息 SR: [4x5 double]df: 0normr: 2.3915e-015Warning: Zero degrees of freedom implies infinite error bounds.> In polyval at 104In polyconf at 92In zxecf at 69观测数据拟合数据x y yh1.0000 3.0000 3.00002.0000 4.0000 4.00003 5 54.0000 6.0000 6.0000剩余平方和 Q = 0.000000标准误差 Sigma = 0.000000相关指数 RR = 1.000000请输入你所需要拟合的数据点,若没有请按回车键结束程序.输入插值点x0 = 3输出插值点拟合函数值 y0 = 5.0000>>结果:untitled.figuntitled2.fig一些matlab优化算法代码的分享代码的目录如下:欢迎讨论1.约束优化问题:minRosen(Rosen梯度法求解约束多维函数的极值)(算法还有bug) minPF(外点罚函数法解线性等式约束)minGeneralPF(外点罚函数法解一般等式约束)minNF(内点罚函数法)minMixFun(混合罚函数法)minJSMixFun(混合罚函数加速法)minFactor(乘子法)minconPS(坐标轮换法)(算法还有bug)minconSimpSearch(复合形法)2.非线性最小二乘优化问题minMGN(修正G-N法)3.线性规划:CmpSimpleMthd(完整单纯形法)4.整数规划(含0-1规划)DividePlane(割平面法)ZeroOneprog(枚举法)5.二次规划QuadLagR(拉格朗日法)ActivedeSet(起作用集法)6.辅助函数(在一些函数中会调用)minNT(牛顿法求多元函数的极值)Funval(求目标函数的值)minMNT(修正的牛顿法求多元函数极值)minHJ(黄金分割法求一维函数的极值)7.高级优化算法1)粒子群优化算法(求解无约束优化问题)1>PSO(基本粒子群算法)2>YSPSO(待压缩因子的粒子群算法)3>LinWPSO(线性递减权重粒子群优化算法)4>SAPSO(自适应权重粒子群优化算法)5>RandWSPO(随机权重粒子群优化算法)6>LnCPSO(同步变化的学习因子)7>AsyLnCPSO(异步变化的学习因子)(算法还有bug)8>SecPSO(用二阶粒子群优化算法求解无约束优化问题)9>SecVibratPSO(用二阶振荡粒子群优化算法求解五约束优化问题)10>CLSPSO(用混沌群粒子优化算法求解无约束优化问题)11>SelPSO(基于选择的粒子群优化算法)12>BreedPSO(基于交叉遗传的粒子群优化算法)13>SimuAPSO(基于模拟退火的粒子群优化算法)2)遗传算法1>myGA(基本遗传算法解决一维约束规划问题)2>SBOGA(顺序选择遗传算法求解一维无约束优化问题)3>NormFitGA(动态线性标定适应值的遗传算法求解一维无约束优化问题)4>GMGA(大变异遗传算法求解一维无约束优化问题)5>AdapGA(自适应遗传算法求解一维无约束优化问题)6>DblGEGA(双切点遗传算法求解一维无约束优化问题)7>MMAdapGA(多变异位自适应遗传算法求解一维无约束优化问题)自己编写的马尔科夫链程序A 代表一组数据序列一维数组本程序的操作对象也是如此t=length(A); % 计算序列“A”的总状态数B=unique(A); % 序列“A”的独立状态数顺序,“E”E=sort(B,'ascend');a=0;b=0;c=0;d=0;for j=1:1:ttLocalization=find(A==E(j)); % 序列“A”中找到其独立状态“E”的位置for i=1:1:length(Localization)if Localization(i)+1>tbreak; % 范围限定elseif A(Localization(i)+1)== E(1)a=a+1;elseif A(Localization(i)+1)== E(2)b=b+1;elseif A(Localization(i)+1)== E(3)c=c+1;% 依此类推,取决于独立状态“E”的个数elsed=d+1;endendT(j,1:tt)=[a,b,c,d]; % “T”为占位矩阵endTT=T;for u=2:1:ttTT(u,:)= T(u,:)- T(u-1,:);endTT; % 至此,得到转移频数矩阵Y=sum(TT,2);for uu=1:1:ttTR(uu,:)= TT(uu,:)./Y(uu,1);endTR % 最终得到马尔科夫转移频率/概率矩阵% 观测序列马尔科夫性质的检验:N=numel(TT);uuu=1;Col=sum(TT,2); % 对列求和Row=sum(TT,1); % 对行求和Total=sum(Row); % 频数总和for i=1:1:ttfor j=1:1:ttxx(uuu,1)=sum((TT(i,j)-(Row(i)*Col(j))./Total).^2./( (Row(i)*Col(j)). /Total));uuu=uuu+1; % 计算统计量x2endendxx=sum(xx)。