最小二乘法

最小二乘法 本实验的目的是:简单介绍“最小二乘法”的思想～利用Mathmatica求最小二乘解。 在科学研究和实际工作中～常常会遇到这样的问题:给定两个变量x, y的m组实验数据～如何从中找出这两个变量间的函数关系的近似解析达(x,y),(x,y),？,(x,y)1122mm 式,也称为经验公式,～使得能对x与y之间的除了实验数据外的对应情况作出某种判断。 这样的问题一般可以分为两类:一类是对要对与之间所存在的对应规律一无所知～xy 这时要从实验数据中找出切合实际的近似解析表达式是相当困难的～俗称这类问题为黑箱问题,另一类是依据对问题所作的分析～通过数学建模或者通过整理归纳实验数据～能够判定 y,f(x,a)出x与y之间满足或大体上满足某种类型的函数关系式～其中 是个待定的参数～这些参数的值可以通过组实验数据来确定,一般nma,(a,a,？,a)12n 要求,～这类问题称为灰箱问题。解决灰箱问题的原则通常是使拟合函数在处的值xm,ni n2[f(x,a),y]与实验数值的偏差平方和最小～即取得最小值。这种在方差意义下对实,ii,1i y,f(x,a)验数据实现最佳拟合的方法称为“最小二乘法”～称为最小二乘解～a,a,？,a12n 称为拟合函数。下面对于“最小二乘法”通过两个例子来说明。 1(线性拟合问题 例1(在某化工厂生产过程中～为研究温度x,单位:摄氏度,对收率,产量,y (%)的影响～可测得一组数据如下表所示～试根据这些数据建立以x与y之间的拟合函数。 温度x 100 110 120 130 140 150 160 170 180 190 收率y(%) 45 51 54 61 66 70 74 78 85 89 解:,1,确定函数的类型 为此～应对给定的数据以x为横坐标～y为纵坐标作出散点图。输入如下语句: x Table100,10,i,i,0,9; 45,51,54,61,66,70,74,78,85,89;y Tablexi,yi,i,1,10xy ListPlotxy,PlotStyle PointSize0.015 运行后可得到数据表和图1。 从图1种可以看出这些点近似地落在一条直线周围～可以认为在y和x之间存在着线性关系～之所以不完全落在直线上～是因为观察数据本身存在的误差。下面我们就用最小二乘法求出与这些数据点最接近的直线方程。 90 80 100,45 , 110,51 ,70 120,54 , 130,61 , 60 140,66 , 150,70 , 160,74 , 170,78 , 120140160180 180,85 , 190,89 图1 ,2,利用Mathmatica求最小二乘解 y,ax，b设直线方程式为～其中a, b是待定系数～于是～可以把拟合函数在处的xi值与观察值的偏差平方和表示为关于待定系数a, b的二元函数: n2Q(a,b),[(ax，b),y]～ ,ii,1i 按照二元函数求极值的理论～其最小值应满足方程组: nn,,Q,,，,,，,,2[(axb)y]x0[(axb)y]x0,iii,iii,,,,,ai,1i,1～即 。 ,,nn,Q,,,，,,,，,,2[(axb)y](1)0[(axb)y]0,,iiii,,,bi,1i,1,, 于是～输入以下Mathematica语句便可求得系数a, b: qa_,b_: Sumaxi,b,yi^2,i,1,10 0,Dqa,b,b 0,a,bSolveDqa,b,a 797452运行后可得解: a , ,b1650165 也可以用“NSolve”命令得到a, b的小数解:。 a 0.48303,b ,2.73939 y,0.48303x,2.73979这样～我们得到方程～这就是前面所说的经验公式。 在此例题中～观察的数据点大致呈一直线状～如果在实际问题中数据点落在某一条曲线周围～则要根据曲线的形状来确定拟合函数的类型。比如考虑用x的n次多项式 2n来拟合～这成为多项式拟合。此时仍可用最小二乘法～P(x),a，ax，ax，？，axn012n n2n，1Q(a,a,a,？,a),[P(x),y]考虑元函数取最小的必要条件～令此函数对,012nnii,1i n，1各个参数的偏导等于0～解一个元的方程组便可求得这些参数的最小二乘解。 二、可化为线性拟合的曲线拟合问题 在许多场合下～拟合函数不具有线性形式～但是由实际经验或相关的学科理论～能够提供拟合函数的可取类型～而且可以通过适当的变量代换将拟合函数线性化～同样可以建立经验公式。 1 bxY,lny,X,xY,lna，bX,1,模型可以用变量替换将函数化为线性函数:。 y,ae y,a，blnxY,y,X,lnx可以用变量替换将函数化为线性函数。 ,2,模型 ,(y),a，bxz,,(y),3,模型可以用变量替换将其化为和之间的线性拟合。 zx例2(研究黏虫的生长过程～可测的一组数据下表所示。 温度t 11.8 14.7 15.4 16.5 17.1 18.3 19.8 20.3 历期N 30.4 15 13.8 12.7 10.7 7.5 6.8 5.7 k其中历期是指卵块孵化成幼虫的天数。昆虫学家认为在与之间有关系式:N,～NNtt,c其中k, c为常数。试求最小二乘解。 11cy,,x,t,a,,b,,解:作变换～由此把N, t的关系式化为了关于x, y的线性关Nkk t,c1c系式:。与例1相同～输入以下语句求解参数a, b: y,,x,,,ax，bkkk t 11.8,14.7,15.4,16.5,17.1,18.1,19.8,20.3; w 30.4,15.0,13.8,12.7,10.7,7.5,6.8,5.7; y 1w; x t; xy Tablexi,yi,i,1,8; Sumaxi,b,yi^2,i,1,8qa_,b_: SolveDqa,b,a 0,Dqa,b,b 0,a,b运行后的结果为:。最后为了求出参数k, c以及经 a 0.016481,b ,0.175433 验公式～还需要回代变量～为此输入命令: A a,b.%; 1A1,1k ,k,A1,2c 运行后得到参数k, c的值为:k=60.6758; :c=10.6445。最后～拟合曲线方程为: 60.6758N, 。 t,10.6445 为了比较得到的拟合函数和已知的数据点～我们再在同一坐标下绘出数据点的散点图及拟合函数的图形～输入语句为: 2 data Tableti,wi,i,1,8; t1 ListPlotdata, PlotStyle PointSize0.02,150 DisplayFunction Identity;125 fx_: kx,c;100 t2 Plotfx,x,11,25,75 AxesOrigin 11,0,50 DisplayFunction Identity;25 Showt1,t2,DisplayFunction 141618202224 $DisplayFunction 图2 运行结果如图2所示。 实验习题 1(为测定刀具的磨损速度～每隔一小时测量一次刀具的厚度～由此得到以下数据: 时间t 0 1 2 3 4 5 6 7 厚度y 27.0 26.8 26.5 26.3 26.1 25.7 25.3 24.8 试根据这组数据建立y与t之间的拟合函数。 2(一种合金在某种添加剂的不同浓度下进行实验～得到如下数据: 浓度x 10.0 15.0 20.0 25.0 30.0 抗压强度y 27.0 26.8 26.5 26.3 26.1 2已知函数y与x的关系适合模型:～试用最小二乘法确定系数a,b,c～并y,a，bx，cx 求出拟合曲线。 3(在研究化学反应速度时～得到下列数据: x3 6 9 12 15 18 21 24 i y57.6 41.9 31.0 22.7 16.6 12.2 8.9 6.5 i 其中表示实验中作记录的时间～表示在相应时刻反应混合物中物质的量～试根据这些xyii 数据建立经验公式。 3