最小二乘拟合实验报告

最小二乘拟合实验报告南昌工程学院《计算方法》实验报告课程名称计算方法系院理学院专业信息与计算科学班级12级一班学生姓名魏志辉学号2012101316《最小二乘求解》引言在科学实验和生产实践中，经常要从一组实验数据出发，寻求函数y=f（x）的一个近似表达式y=φ(x)，称为经验公式，从几何上来看，这就是一个曲线拟合的问题。多项式的插值虽然在一定程度上解决了由函数表求函数近似表达式的问题，但用它来解决这里的问题，是有明显的缺陷的。首先，由实验提供的数据往往有测试误差。如果要求近似曲线y=φ(x)严格地通过所给的每个数据点，就会使曲线保留原来的测试误差，因输入及m，n生成中间矩阵C输入及m，n生成中间矩阵C生成法方程组的系数矩阵生成法方程组的右端向量解法方程组得输出i=1，2，…，mj=2，3，…，n+1 整体流程图 生成矩阵C流程图程序代码及注释%最小二乘拟合%a为线性拟合中的常数，b为一次项系数%t为均方误差，maxi为最大偏差function[a,b,t,maxi]=polyfit(x0,y0,n)m=length(x0);p=length(y0);%x0和y0长度不等时，报错ifm~=pfprintf('Error!Pleaseinputagain!\n');end%生成中间矩阵Cfori=1:mC(i,1)=1;forj=2:n+1C(i,j)=x0(i)*C(i,j-1);endend%生成系数矩阵AA=C'*C;%将题目中的y的每项求自然对数后得到y1fori=1:my1(i)=log(y0(i));end%生成法方程组的右端向量BB=C'*y1';%求解拟合系数X=A\B;%题中，a=exp（X(1)），b=X(2)a=exp(X(1));b=X(2);%先求偏差平方和，再求均方误差sum=0;fork=1:my2(k)=a*exp(b*x0(k));l(k)=y0(k)-y2(k);sum=sum+l(k).^2;endt=sqrt(sum);%最大偏差为偏差矩阵中绝对值最大的一项maxi=max(max(abs(l)));end算例分析1、测试示例>>x=1:18;>>y=[0.8982.383.071.842.021.942.222.774.024.765.466.5310.916.522.535.750.661.681.8];>>[abtmax]=polyfit(x,y,1)Error!Pleaseinputagain!a=0.7185b=0.2227t=31.6255max=22.06312、计算过程（1）首先输入已知点>>x=1:19;>>y=[0.8982.383.071.842.021.942.222.774.024.765.466.5310.916.522.535.750.661.681.8];（2）输出结果>>[abtmax]=polyfit(x,y,1)a=0.6814b=0.2306t=38.3255max=27.3047其中，a，b，t，max分别为常数项，一次项系数，均方误差，最大偏差。讨论与结论时间复杂度：>>tic;[abtmax]=polyfit(x,y,1);tocElapsedtimeis0.859861seconds.说明该算法具有一定的复杂性。直观展示：输入以下命令：>>x=1:19;>>y=[0.8982.383.071.842.021.942.222.774.024.765.466.5310.916.522.535.750.661.681.8];>>fori=1:19y0(i)=log(y(i));end>>x0=0:0.01:20;>>y1=0.6814*exp(0.2306*x0);>>plot(x,y,'+')>>holdon>>plot(x0,y1)>>plot(x,y0,'+')>>y2=log(0.6814)+0.2306*x0;>>holdon>>plot(