曲线拟合——最小二乘法算法
- 格式:doc
- 大小:79.00 KB
- 文档页数:5
Lab04.曲线拟合的最小二乘法实验【实验目的和要求】1.让学生体验曲线拟合的最小二乘法,加深对曲线拟合的最小二乘法的理解;2.掌握函数ployfit和函数lsqcurvefit功能和使用方法,分别用这两个函数进行多项式拟合和非多项式拟合。
【实验内容】1.在Matlab命令窗口,用help命令查询函数polyfit和函数lsqcurvefit 功能和使用方法。
2.用多项式y=x3-6x2+5x-3,产生一组数据(xi,yi)(i=1,2,…,n),再在yi上添加随机干扰(可用rand产生(0,1)均匀分布随机数,或用randn产生N(0,1)均匀分布随机数),然后对xi和添加了随机干扰的yi用Matlab提供的函数ployfit用3次多项式拟合,将结果与原系数比较。
再作2或4次多项式拟合,分析所得结果。
3.用电压V=10伏的电池给电容器充电,电容器上t时刻的电压为,其中V0是电容器的初始电压,τ是充电常数。
对于下面的一组t,v数据,用Matlab提供的函数lsqcurvefit确定V0和τ。
t(秒) 0.5 1 2 3 4 5 7 9v(伏) 6.36 6.48 7.26 8.22 8.66 8.99 9.43 9.63 【实验仪器与软件】1.CPU主频在1GHz以上,内存在128Mb以上的PC;2.Matlab 6.0及以上版本。
实验讲评:实验成绩:评阅教师:200 年月日问题及算法分析:1、利用help命令,在MATLAB中查找polyfit和lsqcurvefit函数的用法。
2、在一组数据(xi,yi)(i=1,2,…,n)上,对yi上添加随机干扰,运用多项式拟合函数,对数据进行拟合(分别用2次,3次,4次拟合),分析拟合的效果。
3、根据t和V的关系画散点图,再根据给定的函数运用最小二乘拟合函数,确定其相应参数。
第一题:(1)>> help polyfitPOLYFIT Fit polynomial to data.P = POLYFIT(X,Y,N) finds the coefficients of a polynomial P(X) ofdegree N that fits the data Y best in a least-squares sense. P is arow vector of length N+1 containing the polynomial coefficients indescending powers, P(1)*X^N + P(2)*X^(N-1) +...+ P(N)*X + P(N+1).[P,S] = POLYFIT(X,Y,N) returns the polynomial coefficients P and astructure S for use with POLYVAL to obtain error estimates forpredictions. S contains fields for the triangular factor (R) from a QRdecomposition of the Vandermonde matrix of X, the degrees of freedom(df), and the norm of the residuals (normr). If the data Y are random,an estimate of the covariance matrix of P is(Rinv*Rinv')*normr^2/df,where Rinv is the inverse of R.[P,S,MU] = POLYFIT(X,Y,N) finds the coefficients of a polynomial inXHAT = (X-MU(1))/MU(2) where MU(1) = MEAN(X) and MU(2) = STD(X). Thiscentering and scaling transformation improves the numerical propertiesof both the polynomial and the fitting algorithm.Warning messages result if N is >= length(X), if X has repeated, ornearly repeated, points, or if X might need centering and scaling.Class support for inputs X,Y:float: double, singleSee also poly, polyval, roots.Reference page in Help browserdoc polyfit>>(2)>> help lsqcurvefitLSQCURVEFIT solves non-linear least squares problems.LSQCURVEFIT attempts to solve problems of the form:min sum {(FUN(X,XDATA)-YDATA).^2} where X, XDATA, YDATA and the valuesX returned by FUN can be vectors ormatrices.X=LSQCURVEFIT(FUN,X0,XDATA,YDATA) starts at X0 and finds coefficients Xto best fit the nonlinear functions in FUN to the data YDATA (in theleast-squares sense). FUN accepts inputs X and XDATA and returns avector (or matrix) of function values F, where F is the same size asYDATA, evaluated at X and XDATA. NOTE: FUN should returnFUN(X,XDATA)and not the sum-of-squares sum((FUN(X,XDATA)-YDATA).^2).((FUN(X,XDATA)-YDATA) is squared and summed implicitly in thealgorithm.)X=LSQCURVEFIT(FUN,X0,XDATA,YDATA,LB,UB) defines a set of lower andupper bounds on the design variables, X, so that the solution is in therange LB <= X <= UB. Use empty matrices for LB and UB if no boundsexist. Set LB(i) = -Inf if X(i) is unbounded below; set UB(i) = Inf ifX(i) is unbounded above.X=LSQCURVEFIT(FUN,X0,XDATA,YDATA,LB,UB,OPTIONS) minimizes with thedefault parameters replaced by values in the structure OPTIONS, anargument created with the OPTIMSET function. See OPTIMSET for details.Used options are Display, TolX, TolFun, DerivativeCheck, Diagnostics,FunValCheck, Jacobian, JacobMult, JacobPattern, LineSearchType,LevenbergMarquardt, MaxFunEvals, MaxIter, DiffMinChange andDiffMaxChange, LargeScale, MaxPCGIter, PrecondBandWidth, TolPCG,OutputFcn, and TypicalX. Use the Jacobian option to specify that FUNalso returns a second output argument J that is the Jacobian matrix atthe point X. If FUN returns a vector F of m components when X has length n, then J is an m-by-n matrix where J(i,j) is the partialderivative of F(i) with respect to x(j). (Note that the Jacobian J isthe transpose of the gradient of F.)[X,RESNORM]=LSQCURVEFIT(FUN,X0,XDATA,YDATA,...) returns the valueof thesquared 2-norm of the residual at X: sum {(FUN(X,XDATA)-YDATA).^2}.[X,RESNORM,RESIDUAL]=LSQCURVEFIT(FUN,X0,...) returns the value of residual,FUN(X,XDATA)-YDATA, at the solution X.[X,RESNORM,RESIDUAL,EXITFLAG]=LSQCURVEFIT(FUN,X0,XDATA,YDATA,...) returnsan EXITFLAG that describes the exit condition of LSQCURVEFIT. Possiblevalues of EXITFLAG and the corresponding exit conditions are1 LSQCURVEFIT converged to a solution X.2 Change in X smaller than the specified tolerance.3 Change in the residual smaller than the specified tolerance.4 Magnitude of search direction smaller than the specified tolerance.0 Maximum number of function evaluations or of iterations reached.-1 Algorithm terminated by the output function.-2 Bounds are inconsistent.-4 Line search cannot sufficiently decrease the residual alongthecurrent search direction.[X,RESNORM,RESIDUAL,EXITFLAG,OUTPUT]=LSQCURVEFIT(FUN,X0,XDATA,YDATA ,...)returns a structure OUTPUT with the number of iterations taken inOUTPUT.iterations, the number of function evaluations inOUTPUT.funcCount,the algorithm used in OUTPUT.algorithm, the number of CG iterations (ifused) in OUTPUT.cgiterations, the first-order optimality (if used)inOUTPUT.firstorderopt, and the exit message in OUTPUT.message.[X,RESNORM,RESIDUAL,EXITFLAG,OUTPUT,LAMBDA]=LSQCURVEFIT(FUN,X0,XDAT A,YDATA,...)returns the set of Lagrangian multipliers, LAMBDA, at the solution:LAMBDA.lower for LB and LAMBDA.upper for UB.[X,RESNORM,RESIDUAL,EXITFLAG,OUTPUT,LAMBDA,JACOBIAN]=LSQCURVEFIT(FU N,X0,XDATA,YDATA,...)returns the Jacobian of FUN at X.ExamplesFUN can be specified using @:xdata = [5;4;6]; % example xdataydata = 3*sin([5;4;6])+6; % example ydatax = lsqcurvefit(@myfun, [2 7], xdata, ydata)where myfun is a MATLAB function such as:function F = myfun(x,xdata)F = x(1)*sin(xdata)+x(2);FUN can also be an anonymous function:x = lsqcurvefit(@(x,xdata) x(1)*sin(xdata)+x(2),[2 7],xdata,ydata)If FUN is parameterized, you can use anonymous functions to capture theproblem-dependent parameters. Suppose you want to solve the curve-fittingproblem given in the function myfun, which is parameterized by its secondargument c. Here myfun is an M-file function such asfunction F = myfun(x,xdata,c)F = x(1)*exp(c*xdata)+x(2);To solve the curve-fitting problem for a specific value of c, first assignthe value to c. Then create a two-argument anonymous function that capturesthat value of c and calls myfun with three arguments. Finally, pass thisanonymous function to LSQCURVEFIT:xdata = [3; 1; 4]; % example xdataydata = 6*exp(-1.5*xdata)+3; % example ydatac = -1.5; % define parameterx = lsqcurvefit(@(x,xdata) myfun(x,xdata,c),[5;1],xdata,ydata) See also optimset, lsqnonlin, fsolve, @, inline.Reference page in Help browserdoc lsqcurvefit>>第二题:1 三次线性拟合clear allx=0:0.5:5;y=x.^3-6*x.^2+5*x-3;y1=y;for i=1:length(y)y1(i)=y1(i)+rand;enda=polyfit(x,y1,3);b=polyval(a,x);plot(x,y,'*',x,b),aa =1.0121 -6.1033 5.1933 -2.4782② 二次线性拟合clear allx=0:0.5:20;y=x.^3-6*x.^2+5*x-3;y1=y;for i=1:length(y)y1(i)=y1(i)+rand;enda=polyfit(x,y1,2);b=polyval(a,x);plot(x,y,'*',x,b),aa =23.9982 -232.0179 367.9756③ 四次线性拟合clear allx=0:0.5:20;y=x.^3-6*x.^2+5*x-3;y1=y;for j=1:length(y)y1(j)=y1(j)+rand;enda=polyfit(x,y1,4);b=polyval(a,x);plot(x,y,'*',x,b),aa =-0.0001 1.0038 -6.0561 5.2890 -2.8249 >>第三题:1 拟合曲线为:f(x)=定义函数:function f=fun(a,x)f=a(1)-(a(1)-a(2))*exp(-a(3)*x);主程序:clear allclcx=[0.5 1 2 3 4 5 7 9];y=[6.36 6.48 7.26 8.22 8.66 8.99 9.43 9.63];a0=[1 1 1];a=lsqcurvefit('fun',a0,x,y);y1=a(1)-(a(1)-a(2))*exp(-a(3)*x);plot(x,y,'r*',x,y1,'b')V1=a(2)tei=1/a(3)Optimization terminated: relative function value changing by less than OPTIONS.TolFun.。
§7 曲线拟合的最小二乘法7-1 一般的最小二乘逼近(曲线拟合的最小二乘法)最小二乘法的一般提法是:对给定的一组数据(,)(0,1,,)i i x y i m =L ,要求在函数类01{,,,}n ϕϕϕϕ=L 中找一个函数*()y S x =,使误差平方和22*222()001[()]min [()]m m m i i i i iS x i i i S x y S x y ϕδδ∈=====−=−∑∑∑其中 0011()()()()()n n S x a x a x a x n m ϕϕϕ=+++<L 带权的最小二乘法: 2220()[()()]mi i i i x S x f x δω==−∑其中()0x ω≥是[a,b ]上的权函数。
用最小二乘法求曲线拟合的问题,就是在()S x 中求一函数*()y S x =,使22δ取的最小。
它转化为求多元函数20100(,,,)()[()()]m n n i j j i i i j I a a a x a x f x ωϕ===−∑∑L 的极小点***01(,,,)n a a a L 问题。
由求多元函数极值的必要条件,有002()[()()]()0m n i j j i i k i i j k I x a x f x x a ωϕϕ==∂=−=∂∑∑ ),,1,0(n k L =若记 0(,)()()()m j k ij i k i i x x x ϕϕωϕϕ==∑ 0(,)()()()mk i i k i k i f x f x x d ϕωϕ==≡∑ ),,1,0(n k L = 则上式可改写为0(,)n kj j k j a d ϕϕ==∑ ),,1,0(n k L =这个方程称为法方程,矩阵形式.Ga d =其中 0101(,,,),(,,,)T T n n a a a a d d d d ==L L ,0001010111011(,)(,)(,)(,)(,)(,)(,)(,)(,)n n n n n n G ϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕ− = LL L L L L L 由于01,,,n ϕϕϕL 线性无关,故0G ≠,方程组存在唯一解*(0,1,,),k k a a k n ==L从而得到函数()f x 的最小二乘解为****0011()()()()n n S x a x a x a x ϕϕϕ=+++L 可证 *2200()[()()]()[()()]m mii i i i i i i x S x f x x S x f x ωω==−≤−∑∑ 故*()S x 使所求最小二乘解。
---------------------------------------------------------------最新资料推荐------------------------------------------------------ 数值分析论文--曲线拟合的最小二乘法曲线拟合的最小二乘法姓名:徐志超学号:2019730059 专业:材料工程学院:材料科学与工程学院科目:数值分析曲线拟合的最小二乘法一、目的和意义在物理实验中经常要观测两个有函数关系的物理量。
根据两个量的许多组观测数据来确定它们的函数曲线,这就是实验数据处理中的曲线拟合问题。
这类问题通常有两种情况:一种是两个观测量 x 与 y 之间的函数形式已知,但一些参数未知,需要确定未知参数的最佳估计值;另一种是 x 与 y 之间的函数形式还不知道,需要找出它们之间的经验公式。
后一种情况常假设 x 与 y 之间的关系是一个待定的多项式,多项式系数就是待定的未知参数,从而可采用类似于前一种情况的处理方法。
在两个观测量中,往往总有一个量精度比另一个高得多,为简单起见把精度较高的观测量看作没有误差,并把这个观测量选作x,而把所有的误差只认为是y 的误差。
设 x 和 y 的函数关系由理论公式 y=f(x; c1, c2, cm)1 / 13(0-0-1)给出,其中 c1, c2, cm 是 m 个要通过实验确定的参数。
对于每组观测数据(xi, yi) i=1, 2,, N。
都对应于 xy 平面上一个点。
若不存在测量误差,则这些数据点都准确落在理论曲线上。
只要选取m 组测量值代入式(0-0-1),便得到方程组yi=f (x;c1,c2,cm)(0-0-2)式中 i=1,2,, m.求 m 个方程的联立解即得 m 个参数的数值。
显然Nm 时,参数不能确定。
在 Nm 的情况下,式(0-0-2)成为矛盾方程组,不能直接用解方程的方法求得 m 个参数值,只能用曲线拟合的方法来处理。
Excel拟合曲线用的最小二乘法1. 介绍Excel作为一款常用的办公软件,被广泛应用于数据分析和处理,而拟合曲线是数据分析中常用的方法之一。
拟合曲线用的最小二乘法是一种常见的拟合方法,通过最小化数据点与拟合曲线之间的距离来找到最佳拟合曲线,从而对数据进行预测和分析。
在本文中,我将从深度和广度的角度来探讨Excel拟合曲线用的最小二乘法,带你深入探索这一主题。
2. 最小二乘法的原理在Excel中进行曲线拟合时,最小二乘法是一种常用的拟合方法。
其原理是通过最小化残差平方和来找到最佳拟合曲线。
残差是指每个数据点到拟合曲线的垂直距离,最小二乘法通过调整拟合曲线的参数,使得残差平方和最小化,从而得到最佳拟合曲线。
在Excel中,可以利用内置函数或插件来实现最小二乘法的曲线拟合,对于不同类型的曲线拟合,可以选择不同的拟合函数进行拟合。
3. Excel中的拟合曲线在Excel中进行拟合曲线时,首先需要将数据导入Excel,然后利用内置的数据分析工具或者插件来进行曲线拟合。
通过选择拟合函数、调整参数等操作,可以得到拟合曲线的相关信息,如拟合优度、参数估计值等。
可以根据拟合曲线的结果来对数据进行预测和分析,从而得到对应的结论和见解。
4. 个人观点与理解对于Excel拟合曲线用的最小二乘法,我认为这是一种简单而有效的数据分析方法。
它能够快速对数据进行拟合,并得到拟合曲线的相关信息,对于数据的预测和分析具有一定的帮助。
然而,也需要注意到拟合曲线并不一定能够准确描述数据的真实情况,需要结合实际背景和专业知识进行分析和判断。
在使用最小二乘法进行曲线拟合时,需要注意数据的可靠性和拟合结果的可信度,以避免出现不准确的结论和偏差的情况。
5. 总结通过本文的探讨,我们对Excel拟合曲线用的最小二乘法有了更深入的了解。
最小二乘法的原理、Excel中的实际操作以及个人观点与理解都得到了充分的展示和探讨。
在实际应用中,需要结合具体情况和专业知识来灵活运用最小二乘法进行曲线拟合,从而得到准确的分析和预测结果。
最小二乘法拟合曲线算法1、概述给定数据点P(x i ,y i ),其中i=1,2,…,n 。
求近似曲线y= φ(x)。
并且使得近似曲线与y=f(x)的偏差最小。
近似曲线在点pi 处的偏差δi = φ(x i )-y ,i=1,2,...,n 。
按偏差平方和最小的原则选取拟合曲线,这种方法为最小二乘法,偏差平方和公式为:min σ2 = (φ(x i −y i ))2ni =1n i =12、推导过程1)设拟合多项式为:y = a 0+ a 1x +⋯+a k x k2)各点到这条曲线的距离之和,即偏差平方如下:R 2= [y i −(a 0+ a 1x +⋯+a k x k )]2ni =13)多项式系数为学习对象,为了求得符合条件的系数值,对上面等式的a i 分别求导,得:−2 y i − a 0+ a 1x +⋯+a k x k =0ni=1−2 y i − a 0+ a 1x +⋯+a k x k x =0ni=1……−2 y i − a 0+ a 1x +⋯+a k x k x k =0ni=14)将等式移项化简,得:a 0+ a 1x +⋯+a k x k = y i ni =1n i =1a 0+ a 1x +⋯+a k x k x = y i x ni =1n i =1……a 0+ a 1x +⋯+a k x k x k = y i x k ni =1n i =15)依上式得矩阵为:x i0 ni=1⋯x i kni=1⋮⋱⋮x i k ni=1⋯x i2kni=1a0⋮a k=y i x i0ni=1⋮y i x i kni=1上边等式左边为1+K阶对称矩阵,解此矩阵方程即可得到曲线系数a k6)对于AX=B,A为对称矩阵,对称矩阵可以分解为一个下三角矩阵、一个上三角矩阵(下三角矩阵的转置)和一个对角线矩阵相乘。
即A=LDL T所以AX=LDL T X=B,令DL T X=Y -> LY=B,其中L为下三角矩阵,且已知,可求出Y。
python曲线拟合的最小二乘法最小二乘法是一种常用的曲线拟合方法。
它的核心思想是寻找一条曲线,使得曲线上的点到实际数据点的距离之和最小。
在Python中,可以使用scipy库中的optimize.curve_fit函数来实现最小二乘法的曲线拟合。
该函数接受两个参数:要拟合的函数和实际数据点。
首先,需要定义一个函数,表示要拟合的曲线形式。
然后,通过curve_fit函数传入该函数和实际数据点,即可得到最小二乘法的曲线拟合结果。
以下是一个简单的示例代码:```pythonimport numpy as npfrom scipy import optimizeimport matplotlib.pyplot as plt# 定义要拟合的曲线形式def func(x, a, b, c):return a * np.exp(-b * x) + c# 实际数据点xdata = np.array([0, 1, 2, 3, 4, 5])ydata = np.array([1, 2, 3, 4, 5, 6])# 曲线拟合params, params_covariance = optimize.curve_fit(func, xdata, ydata)# 拟合结果a, b, c = paramsprint("拟合参数:a =", a, ", b =", b, ", c =", c)# 画出拟合曲线和实际数据点plt.plot(xdata, ydata, 'bo', label='实际数据点')plt.plot(xdata, func(xdata, a, b, c), 'r-', label='拟合曲线') plt.xlabel('x')plt.ylabel('y')plt.legend()plt.show()```在这个示例代码中,我们定义了一个指数形式的曲线函数func,并使用实际数据点进行拟合。
曲线拟合--最小二乘法1:已知平面上四个点:(0,1)、(1,2.1)、(2,2.9)和(3,3.2),求出一条直线拟合这四个点,使得偏差平方和变为极小。
解:设直线方程为:0 1 0 01 2.1 1 2.12 2.9 4 5.83 3.2 9 9.6Sum=6 Sum=9.2 Sum=14 Sum=17.5 代入正规方程:,编程求解上方程组:>> eq1='14*A+6*B=17.5';>>eq2='6*A+4*B=9.2';>> [A,B]=solve(eq1,eq2,'A,B');>> disp(A)0.74>> disp(B)1.19所以直线方程为:2:已知数据如下表所示1 2 4 610 5 2 1试求(1)用抛物线拟合这些数据使得偏差平方和最小;(2)用型如的函数来拟合这些数据使得偏差平方和最小。
(3)比较这两种拟合结果。
解:(1)设抛物线方程为:1 10 1 1 1 10 102 5 4 8 16 10 20 4 2 16 64 256 8 326 1 36 216 1296 6 36 Sum=13 Sum=18 Sum=57 Sum=289 Sum=1569 Sum=34 Sum=98代入正规方程:得到系数A,B,C的方程组:编程求解上方程组:>>eq1='1569*A+289*B+57*C=98';>>eq2='289*A+57*B+13*C=34';>>eq3='57*A+13*B+4*C=18';>> [A,B,C]=solve(eq1,eq2,eq3,'A,B,C');>> disp(A); disp(B); disp(C)102/199-1048/1992848/199>> A=102/199; disp(A) 0.5126>> B=-1048/199; disp(B) -5.2663>> C=2848/199; disp(C) 14.3116所以得到抛物线的方程为:(2)设函数1 10 1 1 102 5 1/2 1/4 5/24 2 1/4 1/16 1/26 1 1/6 1/36 1/6Sum=13 Sum=18 Sum=23/12 Sum=193/144 Sum=79/6 得到系数A,B的方程组:编程求解上方程组:>> eq1='4*A+23*B/12=18';>>eq2='23*A/12+193*B/144=79/6';>> [A,B]=solve(eq1,eq2,'A,B');>> disp(A); disp(B)-160/243872/81>> A=-160/243; disp(A)-0.6584>> B=827/81; disp(B)10.2099所以得到的函数为:(3)比较(1)和(2)两种方法拟合的方程:编程画出抛物线的图像为:>> x=-2:0.1:12;>> y=0.5126*x.^2-5.2663*x+14.3116;plot(x,y);grid on(a)再编程画出的图像为:>> x=-2:0.1:12;>> y=-0.6584+10.2099*(x.^(-1));>> plot(x,y);grid on>> x=-1:0.01:1;>> y=-0.6584+10.2099*(x.^(-1));plot(x,y);grid on(b)比较两图像可知,图像(b)在点(0,0)处不连续。
最小二乘法拟合曲线公式
最小二乘法是一种常用的数学方法,可以用来拟合一条曲线,使得曲线上的点与实际观测值的误差最小化。
最小二乘法拟合曲线的公式为:
y = a + bx
其中,y 是因变量,x 是自变量,a 和 b 是拟合曲线的系数。
最小二乘法通过最小化误差平方和来确定 a 和 b 的值,即:
b = (n∑xy - ∑x∑y) / (n∑x^2 - (∑x)^2)
a = (∑y - b∑x) / n
其中,n 是数据点的个数,∑表示求和符号,x 和 y 分别表示自变量和因变量的值。
拟合曲线的误差可以通过计算残差平方和来评估,即:
SSR = ∑(y - )^2
其中,y 是实际观测值,是拟合曲线的预测值。
最小二乘法拟合曲线的优点在于可以用简单的数学公式表示,易于理解和应用。
- 1 -。
使用最小二乘法拟合直线数据
最小二乘法是一种拟合曲线或直线的常用方法。
它通过最小化误差的平方和来找到最佳拟合参数。
假设我们有一组点(x, y),我们希望找到一条曲线y = f(x) 来拟合这些点。
最小二乘法的基本思想是:找到一条曲线,使得所有点到这条曲线的垂直距离的平方和最小。
最小二乘法通常用于线性回归和非线性回归。
在非线性回归中,我们通常需要使用一些优化算法来找到最佳参数。
下面是一个简单的Python 代码示例,演示如何使用最小二乘法拟合一条直线:
import numpy as np
import matplotlib.pyplot as plt
# 生成一组随机数据
x = np.random.rand(50)
y = 2 * x + 1 + np.random.randn(50) * 0.1
# 使用最小二乘法拟合直线
p, _, _, _ = np.polyfit(x, y, 1, method='leastsq')
# 绘制原始数据和拟合直线
plt.scatter(x, y, color='blue')
plt.plot(x, p[0] * x + p[1], color='red')
plt.show()
在这个例子中,我们首先生成了一组随机数据(x, y)。
然后,我们使用np.polyfit函数来拟合一条直线。
最后,我们将原始数据和拟合直线绘制在同一张图上。
曲线拟合——最小二乘法算法
一、目的和要求
1)了解最小二乘法的基本原理,熟悉最小二乘算法;
2)掌握最小二乘进行曲线拟合的编程,通过程序解决实际问题。
二、实习内容
1)最小二乘进行多项式拟合的编程实现。
2)用完成的程序解决实际问题。
三、算法
1)输入数据节点数n ,拟合的多项式次数m ,循环输入各节点的数据x j , y j (j=0,1,…,n-1)
2)由x j 求S ;由x j ,y j 求T :
S k =
∑-=10n j k j x ( k=0,1,2, … 2*m ) T k = ∑-=1
0n j k j j x y ( k=0,1,2,… m )
3)由S 形成系数矩阵数组c i,j :c[i][j]=S[i+j] (i=0,1,2,…m, j=0,1,2,…,m);由T 形成系数矩阵增广部分c i,m+1:c[i][m+1]=T[i] (i=0,1,2,…m)
4)对线性方程组CA=T[或A C ],用列主元高斯消去法求解系数矩阵A=(a 0,a 1,…,a m )T
四、实验步骤
1)完成最小二乘法进行曲线拟合的程序设计及录入、编辑;
2)完成程序的编译和链接,并进行修改;
3)用书上P105例2的例子对程序进行验证,并进行修改;
4)用完成的程序求解下面的实际问题。
5)完成实验报告。
五、实验结果
1. 经编译、链接及例子验证结果正确的源程序:
#include<stdio.h>
#include<math.h>
#define Q 100
float CF(int,float);
main()
{
int i,j,n1,n,p,k,q;
float x[Q],y[Q],s[Q]={0},t[Q]={0},a[Q][Q]={0},l,sum=0;
/*以下是最小二乘的程序*/
printf("input 数据组数n");
scanf("%d",&n);
printf("input 拟合次数n1");
scanf("%d",&n1);
for(i=0;i<n;i++)
{
printf("x[%d]=",i);
scanf("%f",&x[i]);
printf("y[%d]=",i);
scanf("%f",&y[i]);
}
for(i=0;i<=2*n1;i++)
for(j=0;j<n;j++)
{
s[i]=s[i]+CF(i,x[j]);
if(i<=n1)
t[i]=t[i]+y[j]*CF(i,x[j]);
}
for(i=0;i<n1+1;i++)
for(j=0;j<n1+2;j++)
{
a[i][j]=s[i+j];
if(j==n1+1)
a[i][j]=t[i];
}
for(i=0;i<n1+1;i++)
for(j=0;j<n1+2;j++)
printf("a[%d][%d]=%f",i,j,a[i][j]); /*以下的是削去法的程序*/
for(j=0;j<=n1-1;j++)
{p=j;
for(i=j+1;i<=n1;i++)
{
if(fabs(a[j][j])<fabs(a[i][j]))
p=i;
}
if(p!=j)
for(i=j;i<=n1+1;i++)
{l=a[p][i];
a[p][i]=a[j][i];
a[j][i]=l;
}
for(k=j+1;k<=n1;k++)
{l=a[k][j]/a[j][j];
for(q=j;q<=n1+1;q++)
a[k][q]=a[k][q]-l*a[j][q];
}
}
for(i=0;i<n1+1;i++)
{for(j=0;j<n1+2;j++)
printf("a[%d][%d]=%f\n",i,j,a[i][j]);
printf("\n");}
x[n1]=a[n1][n1+1]/a[n1][n1];
for(i=n1-1;i>=0;i--)
{for(j=i+1;j<=n1;j++)
sum=a[i][j]*x[j]+sum;
x[i]=(a[i][n1+1]-sum)/a[i][i];
sum=0;
}
for(i=0;i<=n1;i++)
printf("x[%d]=%f\n",i,x[i]);
}
float CF(int i,float v)
{float a=1.0;
while(i--)a*=v;
return a;
}
2. 实例验证结果:
1)输入初始参数:
n=9,m=2
X:1 3 4 5 6 7 8 9 10
Y:10 5 4 2 1 1 2 3 4
2)结果输出:
1.实际应用
问题:
作物体运动的观测实验,得出以下实验测量数据,用最小二乘拟合求物体的运动方程。
时间t(秒) 0 0.9 1.9 3.0 3.9 5.0 距离s(cm) 0 10 30 50 80 110
解题步骤:
1)画草图
2)确定拟合方程次数为1:
用完成的程序输入数据,求取拟合方程中的未知数,得出方程:
S=-7.855.58+22.253763T
计算误差:
3)确定拟合方程次数为2:
用完成的程序输入数据,求取拟合方程中的未知数,得出方程:
S=-0.583370+11.081389T+2.248811T*T
计算误差:
六、分析和讨论
结合实际问题,进行拟合次数的分析和讨论:
七、心得(*可选)
①调试过程中遇到的问题和解决对策;②经验体会等。
以上都是我业余所写,程序肯定有很多不足之处希望大家多多指教。