昆工研究生数值分析报告-曲线拟合最小二乘法

数值分析上机实验课题报告《课题四曲线拟合的最小二乘法》专业：测绘工程姓名：钟凯学号：201722010102017年12月课题四曲线拟合的最小二乘法一、目的和意义在实际问题中测得的实验数据有时需寻求较简单函数逼近来解,曲线拟合的最小二乘法在解决这类问题的数据处理和误差分析中应用非常广泛,可提高数据处理的效率和精确度,已成为这类问题数据处理的重要的比较可靠的技术手段。在科学实验的统计方法研究中，往往要从一组实验数据,iixy0,1,2,,im中，寻找自变量x与因变量y之间的函数关系yFx。由于观测数据往往不准确，因此不要求yFx经过所有点,iixy，而只要求在给定点ix上误差而只要求所在所有给定点ix上的误差()iiiFxy0,1,2,,im按某种标准最小。若记012,,,,Tm，就是要求向量的范数最小。如果用最大范数，计算上困难较大，通常采用欧式范数2作为误差度量的标准。Fx的函数类型往往与实验的物理背景以及数据的实际分布有关，它一般含有某些待定参数。如果Fx是所有待定参数的线性函数，那么相应的问题称为线性最小二乘问题，否则称为非线性最小二乘问题。最小二乘法还是实验数据参数估计的重要工具。这是因为这种方法比其他方法更容易理解，即使在其他方法失效的情况下，用最小二乘法还能提供解答，而且从统计学的观点分析，用该方法求得各项估计具有最优统计特征，因此这一方法也是系统识别的重要基础。线性最小二乘问题可以借助多元微分学知识通过求解法方程组得到解答。用最小二乘法求拟合曲线时，首先要确定Sx的形式。这不单纯是数学问题，还与所研究问题的运动规律以及所得观测数据,iixy有关；通常要从问题的运动规律以及给定数据描图，确定Sx的形式，并通过实际计算选出较好的结果。为了使问题的提法更有一般性，通常把最小二乘法中的22都考虑为加权平方和2220miiiixSxfx这里0ix是,ab上的加权函数，它表示不同点,iixfx处的数据比重不同。二、问题与计算方法从随机的数据中找出其规律性，给出其近似表达式的问题，在生产实践和科学实验中大量存在，通常利用数据的最小二乘法求得拟合曲线。在某冶炼过程中，根据统计数据的含碳量与时间关系，试求含碳量y与时间t的拟合曲线。本题要求我们用33221tatatat对曲线进行拟合，这里2312311,,,,mttt故112110,12650iit,11312210,,544500iit,114220,24983750iit,11413310,,24983750iit,11532230,,1193362500iit,116330,58593218750iit1110,1365.55iiiyty,11220,54350.75iiiyty,11330,2379846.25iiiyty由于0,0,1,,nkjjjjadkn,可以利用此式算出拟合曲线的ia，即(t分）0510152025303540455055410y01.272.162.863.443.874.154.374.514.584.024.64111213111212223222331323333,,,,,,,,,,,,yaayay所以求得512.6610a，725.2910a，933.5210a，572932.66105.29103.5210tttt，误差为iiiyt0,1,11imax0.455ii，而均方误差为11200.51293892ii下表可见拟合出的每一点的误差以及均方误差。ty拟合值误差误差平方0000051.271.20215-0.067850.004603623102.162.16620.00623.844E-05152.862.918550.058550.003428103203.443.48560.04560.00207936253.873.893750.023750.000564063304.154.16940.01940.00037636354.374.33895-0.031050.000964102404.514.4288-0.08120.00659344454.584.46535-0.114650.013144622504.024.4750.4550.207025554.644.48415-0.155850.024289222均方误差0.51293892三、结构程序设计在本题使用VisualStudioC#.NET语言来编译，使用MicrosoftVisualStudio2010软件平台进行编译。1.软件使用：首先打开MicrosoftVisualStudio2010软件，点击新建项目，在弹出的窗体中选择VisualC#语言，并选择windows窗体应用程序，如图下图。这时我们已经新建了一个项目，并自动生成窗体Form1，然后，在右边的工具箱中添加各种编程所需的控件例如：按钮Button、文字lable、文本框Textbox等等，如图下图。这个时候我们通过添加按钮和输入文本框这些必要控件，已经建立好计算器的初步模型，如图下图。这个时候我们只需点击计算按钮，在弹出的窗体中输入代码即可。代码输入完毕后点击调试，便能够生成计算器在计算器程序中输入数值点击计算即可，如图下图。2.部分核心代码如下：（1）在窗体的Load事件里调用InitialDeal和Deal函数来处理数据。privatevoidResultReport_Load(objectsender,System.EventArgse){this.InitialDeal();this.Deal();}（2）初始化各个变量。privatevoidInitialDeal(){inti=0;while(tString.Length0){t[i++]=Convert.ToDouble(tString.Substring(0,tString.IndexOf(,)))*tUnit;tString=tString.Remove(0,tString.IndexOf(,)+1);}num=i;i=0;while(yString.Length0){y[i++]=Convert.ToDouble(yString.Substring(0,yString.IndexOf(,)))*yUnit;yString=yString.Remove(0,yString.IndexOf(,)+1);}i=0;while(inum){a11+=y[i]*y[i];i++;}this.labelA11.Text=a11.ToString(#.00E0;(#.00E0);0.00);i=0;while(inum){a12+=y[i]*y[i]*y[i];i++;}a21=a12;this.labelA12.Text=a12.ToString(#.00E0;(#.00E0);0.00);this.labelA21.Text=a21.ToString(#.00E0;(#.00E0);0.00);i=0;while(inum){a13+=y[i]*y[i]*y[i]*y[i];i++;}a31=a13;a22=a13;this.labelA31.Text=a31.ToString(#.00E0;(#.00E0);0.00);this.labelA13.Text=a13.ToString(#.00E0;(#.00E0);0.00);this.labelA22.Text=a22.ToString(#.00E0;(#.00E0);0.00);i=0;while(inum){a23+=y[i]*y[i]*y[i]*y[i]*y[i];i++;}a32=a23;this.labelA23.Text=a23.ToString(#.00E0;(#.00E0);0.00);this.labelA32.Text=a32.ToString(#.00E0;(#.00E0);0.00);i=0;while(inum){a33+=y[i]*y[i]*y[i]*y[i]*y[i]*y[i];i++;}this.labelA33.Text=a33.ToString(#.00E0;(#.00E0);0.00);i=0;while(inum){b1+=y[i]*t[i];i++;}this.labelB1.Text=b1.ToString(#.00E0;(#.00E0);0.00);i=0;while(inum){b2+=y[i]*y[i]*t[i];i++;}this.labelB2.Text=b2.ToString(#.00E0;(#.00E0);0.00);i=0;while(inum){b3+=y[i]*y[i]*y[i]*t[i];i++;}this.labelB3.Text=b3.ToString(#.00E0;(#.00E0);0.00);}（3）进行最小二乘法处理。privatevoidDeal(){doubleab11,ab12,ab13,ab21,ab22,ab23,ab31,ab32,ab33=0;doubleaDiterminal=0;aDiterminal=a11*a22*a33+a12*a23*a31+a13*a21*a32-a13*a22*a31-a23*a32*a11-a12*a21*a33;ab11=(a22*a33-a32*a23)/aDiterminal;ab12=(a13*a32-a12*a33)/aDiterminal;ab13=(a12*a23-a13*a22)/aDiterminal;ab21=(a23*a31-a21*a33)/aDiterminal;ab22=(a11*a33-a13*a31)/aDiterminal;ab23=(a13*a21-a11*a23)/aDiterminal;ab31=(a21*a32-a31*a22)/aDiterminal;ab32=(a12*a31-a11*a32)/aDiterminal;ab33=(a11*a22-a12*a21)/aDiterminal;（4）声明一个3行3列的两维数组，用于存放矩阵A并初始化矩阵Adouble[,]matrixA=newdouble[3,3]{{ab11,ab12,ab13},{ab21,ab22,ab23},{ab31,ab32,ab33}};（5）声明一个3行1列的两维数组，用于存放矩阵B并初始化矩阵Bdouble[,]matrixB=newdouble[3,1]{{b1},{b2},{b3}};（6）声明一个3行1列的两维数组，用于存放矩阵A和矩阵B的乘积double[,]matrixC=newdouble[3,1];//计算矩阵A的逆转//计算矩阵A和矩阵B的乘积for(inti=0;i3;i++){for(intj=0;j1;j++){//初始化矩阵CmatrixC[i,j]=0;//计算矩阵A和矩阵B的乘积，并把值存放在矩阵C中for(intk=0;k3;k++){matrixC[i,j]+=matrixA[i,k]*matrixB[k,j];}}}this.labela1Result.Text=matrixC[0,0].ToString(#.00E+0;#.00E-0;0.00);this.labela2Result.Tex