实验数据分析方法_误差理论与最小二乘法.ppt

1、实验数据分析方法_Chap.5,1,第二部分 实验数据的统计分析 第五章 误差理论与最小二乘法 第六章 回归分析 第七章 多变量分析 第八章 功率谱与周期分析,实验数据分析方法,教材：天文数据处理方法：丁月蓉 编著 主要参考书：实验的数学处理：李惕陪 著,教学方法：基本理论 + 具体实例 + 上机实习(课后),实验数据分析方法_Chap.5,2,第五章 误差理论与最小二乘法,天文学的诸多理论是以天文观测为基础的，如地球自转理论、人造卫星运动理论等都离不开天文观测。人们通过对某一天文量(静态的或动态的)的直接或间接观测，获得大量的数据。而任何观测都不可避免的含有误差。因此，当我们在利用观测结果时

2、，必须分析这些数据的可靠程度：只有当它们的误差在我们允许的范围之内时，我们才能放心大胆的去使用它，否则则不能使用。,误差的研究无论是对生产实践还是基础理论研究都有着重要意义！,实验数据分析方法_Chap.5,3,例1：由于牛顿在其最初计算中使用了具有较大误差的地球半径值，使得他测得的月球加速度的值和理论计算值相差约10，因而推迟了20年发表他的引力理论！例2：爱因斯坦广义相对论的观测证明：1916年爱因斯坦在德国物理学纪事上发表了具有划时代意义的重要文献广义相对论基础。文章指出，当光线行经太阳附近时，光线产生弯曲，其弯曲曲率预计为 = 1.”75，而1911年他用经典方法得到=0.”9，相差两

3、倍。如果观测能测得在1.”75附近，这将证明他的广义相对论是正确的，如果测得的值是在经典值附近，则将否定其理论。幸好1919年英国天文学家爱丁顿爵士在西非几内亚湾的普林西比岛的日全食观测中测得1.”610.”30；与此同时有人在巴西东北海岸外索伯雷尔的日食观测中测得1.”980.”12。这两个结果与广义相对论的预言值相近，远大于经典理论值，强有力的证明了广义相对论的正确性! 如果他们当时的观测误差很大，置信度很低，以致于和理论值相差甚远，那么也就很难由此来验证这个理论了。由此可见，观测和误差分析对基础理论的研究起了一个不可估量的作用!,实验数据分析方法_Chap.5,4,最小二乘法是用来处理具

4、有误差的观测数据的一种有效的方法，也是最早用于天文观测资料处理的一种数学工具。早在l794年，高斯为了利用小行星坐标的多次观测准确地推算小行星的轨道，第一次应用了最小二乘法。1805年勒让德应用测量平差方法确定了彗星的轨道和地球子午线弧长。1809年高斯又推证了误差的概率定律，从而使最小二乘法高度完善化，成为数据处理中应用最广的一个分支。随着概率统计学和矩阵理论的发展以及电子计算机的广泛应用，最小二乘法进入了近代数据处理方法的行列。,实验数据分析方法_Chap.5,5,误差是实验科学术语，指测量结果偏离真值的程度。对任何一个物理量进行的测量都不可能得出一个绝对准确的数值，即使采用测量技术所能达

5、到的最完善的方法，测出的数值也和真实值存在差异，这种测量值和真实值的差异称为误差。(from Wiki) 误差按其表达形式分：绝对误差、相对误差 误差按其性质及产生原因分：系统误差、随机误差、过失（人为）误差 误差不仅存在于测量值中，计算时采用近似的理论模型，计算中一些理论常数的不准确以及数值计算中取位的多少等也会在计算结果中产生误差。,5.1 误差的定义与分类,实验数据分析方法_Chap.5,6,5.1.1 绝对误差和相对误差,一个量值的给出值的绝对误差定义为该量值的给出值与其真值之差，或用公式表示为：绝对误差给出值-真值公式中的给出值如果是被测量的观测结果，则相应的误差为观测误差；如果给出

6、值是某量的计算近似值，则相应的误差为计算近似值的误差。式中的真值是被测量本身的真实大小，它是一个理想的概念：一般说来，真值是未知的，通常用约定值来代替。例如某一系统的天文常数也可看作相应量值的真值。从绝对误差的定义式不难看出，绝对误差和被测量具有相同的量纲。因此，若说一颗星其位置误差为0.“1，测时的记录误差为0.“0001，都是指的绝对误差。,实验数据分析方法_Chap.5,7,我们把误差的反号值定义为修正值，则可得：真值给出值 - 误差给出值 + 修正值这表明，带有误差的给出值加上修正值后可消除或减小误差的影响。在有些情况下用绝对误差来表示测量的精度是不恰当的：如目前卫星激光测距的准确度(

7、测量值与被测量真值之间的偏离程度)已达cm级，卫星的距离一般为103km量级；但如果我们测定的是恒星的距离(这里指离太阳在20pc以内的恒星)，用三角视差法一般可准确到0.”02，相当于2pc的测距误差，显然它和卫星的测距误差是无法直接比较的！但如果我们引入相对误差的概念，它们的测距误差就有了可比性。,实验数据分析方法_Chap.5,8,被测量的绝对误差与其真值a之比定义为这个量的相对误差，并用下式表示：当误差较小时，相对误差式中真值a可用给定值代替。对于上面的例子，它们测距的相对误差分别为110和110-1。 即三角视差测量的相对误差反而要比卫星激光测距的相对误差小！,实验数据分析方法_Ch

8、ap.5,9,由观测的环境因素差异、仪器性能、不同的观测者等因素造成的按某一确定的规律变化的误差称为系统误差。系统误差的大小和符号在多次重复观测中几乎相同，通常使观测值往一个方向偏离。另外，这种误差可以归结为某一因素或某几个因素的函数，而这种函数通常可以用解析公式表达出来。人们总是设法找出代表系统误差的解析表达式，然后在观测结果中扣除。由某些难以控制的随机因素造成的，绝对值和符号的变化时大时小、时正时负，以不可预测的方式变化的误差称为随机误差。虽然就其个体而言，随机误差没有规律、不可预料，但就其总体而言，随着观测次数的增加，它又服从某种统计规律。下面我们将从概率论的角度出发讨论随机误差所满足的

9、统计规律。,5.1.2 系统误差、随机误差和过失误差,实验数据分析方法_Chap.5,10,古典误差理论认为，随机误差服从正态分布，因此我们可以用正态分布密度曲线来表征随机误差，随机误差的分布密度曲线可表为：其被称为高斯误差方程，其相应图形也常被称为高斯误差曲线。式中 称为精密度指数，=x-a, 为随机误差的均方差。 高斯误差方程的一般表达式：,实验数据分析方法_Chap.5,11,随机误差有下列统计特征，当观测样本足够大时： (1) 绝对值相等、符号相反的正负误差近于相等。因此，随机误差的算术平均值随着观测次数的增加愈来愈小，以零为极限。 (2) 误差的概率与误差的大小有关，绝对值小的误差出

10、现的概率比绝对值大的误差出现的概率大，绝对值很大的误差出现的概率很小。,根据随机误差的这些特征，当不存在系统误差的影响时，多次测量结果的平均值将更接近于真值。随机误差产生的原因很多，观测时环境因素的微小变化，设备中的热噪声等都是产生随机误差的重要原因。,实验数据分析方法_Chap.5,12,实际上，系统误差和随机误差之间并没有明显的界限 有时，我们把一些具有复杂规律但暂末掌握的系统误差都当作随机误差处理。而随着人们对误差及其规律的认识的加深，就有可能把这些以往认识不到因而归之于随机误差的这类误差确认为系统误差。反之，在一个较短时期内可能呈现出某种规律，故而归为系统误差，但经过一段较长时间的观测

11、，发现这种变化规律破坏了，并呈现出随机性，这就是说，随着时间的推移，两种不同性质的误差有可能互相转化。 过失（人为）误差是指测量结果与事实明显不符的一种误差。如观测时对错星或观测过程中望远镜/记录仪器的小故障等过失原因造成的结果异常。这种误差一般比较容易发现，而且只要观测人员认真细致，基本上是可以避免的。,实验数据分析方法_Chap.5,13,数据处理中一个很重要的方面是评定一列观测值的可靠程度。它是指观测结果与真值的一致程度，是观测结果中系统误差和随机误差大小的综合度量，常用准确度这个词来表征。在消除了系统误差之后，观测的可靠程度由随机误差的大小来衡量。一列观测值精度高低必须从全列观测值的误

12、差来衡量，而不能只根据个别值的误差来判断。另外，观测的目的是要从一列观测值中确定(直接地或间接地)被测量的真值，但由于观测手段和观测次数的限制，真值实际上是测不到的，只能得到它的一个近似值或估计值。在天文学中通常把最接近于被测量的真值的一个近似值称为它们的最或然值，因此，数据处理的又一个重要的问题是给出被测量的最或然值及其精度。最或然值的精度是衡量观测结果的精度和处理方法有效性的综合指标。,5.2 观测精度,实验数据分析方法_Chap.5,14,标准偏差(又称均方误差)是用来衡量一列观测值精度高低的一个较好指标。 设 为被测量的一组观测值，a为被测量的真值，且xi中只包含随机误差，则 称为xi

13、的真误差，我们定义真误差的平方的算术平均值的平方根为这列观测值的标准偏差或标准误差，天文上又常称之为中误差，并用表示，即：这里定义的标准误差和统计学中从方差的正平方根定义的标准差是一致的，因为从概率论的角度来说，xi的真值可用其数学期望表示。,5.2.1 精度标准,实验数据分析方法_Chap.5,15,下面我们来说明标准偏差的大小为什么可以用来衡量一列观测值的精度高低：由正态分布的性质可知，观测值xi在（a，a+）区间上的概率，或说i出现在（，+）范围内的概率为68.3%, 已知1 2. 则区间（a1，a+1）小于（a2 , a+2），也就是说 = 1的观测数据在a 周围的分布较密集，而 =

14、的观测值在a 周围的分布较分散，即标准偏差 的大小可以衡量一列观测值在真值周围分布的密度程度，而这种密集程度是具有概率含义的，即误差在（，+）内的置信水平是68.3%。,实验数据分析方法_Chap.5,16,下表列出了一些常用的置信水平误差限：,可见，误差落在3中的概率为99.7,亦即绝对值大于3的误差仅有0.3%,这显然是一个小概率事件。所以在有限次观测中，误差值大于3的观测值可能含有过失误差，应考虑舍去该观测值; 当然, 也有可能这个值并不含有过失误差, 如舍去它会犯“弃真”错误，但这种误差的最大概率也只有0.3%。这种取舍观测值的原则称为拉依达准则或简称为3 准则。,实验数据分析方法_C

15、hap.5,17,高斯函数的性质,实验数据分析方法_Chap.5,18,在比较两个观测结果时，应在相同的置信水平上比较它们的误差限，误差限较小的观测较精确，为了说明观测的精度，通常把观测结果报导为(置信水平)。凡是没有注明置信水平的，一般均指 68.3，相应的误差限即为标准误差。在上述各式中，真值a(或x)通常是未知的，因此真误差也是未知的，通常用被测量的最或然值或真值的估计值代替真值，观测值与其最或然值之差称为观测值的残差或离差。标准误差不取决于观测中个别误差的符号，对观测值中较大误差和较小误差比较灵敏，是表示精度的较好方法。实际应用中，有时也常用平均误差离差绝对值的算术平均值来表示精度；也

16、有时采用概率误差：即绝对值比它大的误差和绝对值比它小的误差出现的可能性一样大，将误差绝对值按大小顺序排列，序列的中位数即为概率误差。平均误差和概率误差只有当N较大时才较可靠。天体物理中还经常采用半峰宽度来表示观测的精度，所谓半峰宽度，即观测值分布曲线在极大值半高度处的全宽（Full Width at Half Maximum)。,实验数据分析方法_Chap.5,19,在很多实际问题中, 待求量往往不能直接观测得到，但它们可通过对其它量的观测，再利用它们之间的函数关系换算求得：这种情况就称为间接观测。间接观测在天文观测中是普遍存在的，例如：在人造卫星的定轨预报中要测的是卫星在某一历元的轨道根数，

17、但它们不能直接测得而只能通过测定卫星的赤经、赤纬换算而得到。对于间接观测的情况，应首先由直接观测量求出间接观测量的最或然值，然后由直接观测量的精度估计出间接观测量的精度。通常用下面的式子表示间接观测量y与m个直接观测量xk (k = 1 m)的关系：,5.2.2 误差传递公式,实验数据分析方法_Chap.5,20,为了求得间接观测时误差传递的关系，需要对上式进行线性化处理 如果直接观测量的误差相对于它们的观测值来说是较小的量，则非线性函数可以在各个观测值的邻近点上展开成泰勒级数，然后取误差的一阶项而略去一切高阶误差项：式中 为观测量xk的离差，我们把它记为k。若对xk(k=1m)各进行了N次观

18、测，设间接观测量任一次观测的离差为yyy0，y0f(x10，x20，xm0), 将 yy+y0, k=xkxk0 代入上式，可得：直接观测量xk的误差以的形式出现在间接观测量y的误差中，或说间接观测量y的误差是m个直接观测量的误差加权和，权重因子 称为y的误差传递系数。,实验数据分析方法_Chap.5,21,设m个直接观测量的标准偏差为 ，根据标准偏差的定义及随机变量方差的运算法则，可得间接观测量y的标准偏差为：式中kj为第k个观测量与第j个观测量的相关系数。当各个直接观测量相互独立时，有kj=0, 则有：上式通常称为独立观测量的误差合成定理。 若间接观测量与直接观测量的关系为线性关系时，即：

19、 。则有：此式即为线性情况下的标准偏差传递公式。,实验数据分析方法_Chap.5,22,例：利用IRAF进行测光时，其会根据误差传递以如下的公式给出测光误差：,根据信噪比的定义：S/N = Flux/Err，故1/MerrS/N，即IRAF里给出的测光误差的倒数即为信噪比。除了信噪比会引起测光误差外，还有很多其他的因素也会带来误差，如减本底、除平场、减暗流等过程都会带来附加的误差：一般平场的精度可以达到千分之五左右。,目标源的测光误差可以按如下形式给出：,Eref为多颗比较星测光误差的平均值，Eobj为目标源测光误差，Eothers为其他误差，根据不同的情况确定，比如误差小于千分之五的时候“其

20、他误差”就可能需要包括平场误差，再比如比较星的定标误差等。,实验数据分析方法_Chap.5,23,观测精度的高低是由观测条件决定的，它包括观测的手段、仪器的精度、观测的次数、观测者技术熟练的程度等, 因此我们按观测时的条件把观测分成两大类: 如果某一列观测是在完全相同的条件下进行的，则为等精度观测，所得到的序列称为等精度观测列；如果某一列观测是在不同的条件下进行的，称为非等精度观测，相应的观测序列为非等精度观测列。等精度观测列的标准偏差对于等精度观测列，可以用全列观测值的标准偏差来衡量这列观测值的精度。但是，由于观测值的真误差一般是未知的，为此通常用观测值的残差代替真误差。而对于一列等精度观测

21、值来说，被测量的最或然值就是这列观测值的算术平均值 ，则有残差 ，而真误差为:,5.2.3 等精度观测和非等精度观测,实验数据分析方法_Chap.5,24,为算术平均值的真误差，对上式两边求平方和，得：并有： 由线性情况下的标准偏差传递公式，并将算术平均值的标准偏差 代入上式则得：整理后得到一等精度观测列用残差表示的标准偏差公式(这里用高斯符号 表示求和)：,实验数据分析方法_Chap.5,25,权与非等精度观测列处理非等精度观测序列的情况在天文学中是很普遍的：例如利用观测星表编制基本星表就是一个典型的例子。各种星表中的星位都具有误差；即使是在同一星表中，它所包含的星位也不都具有相同的标准偏差

22、。 它们大多数和观测次数的多少有关，故而大多数星表中有一栏同时列出了各恒星观测的次数，相应的精度随所用的观测数目的增加而增加。因此，在编制基本星表时，需根据它们精度的高低区别对待。在数据处理中，通常用数值pi表示对某一观测结果xi的重视程度，并称之为权。观测值精度的高低是和其误差大小密切相关的：误差越大, 观测值精度就越低，对它的重视程度也应相应减小。在观测值只包含随机误差的情况下，通常定义权与标准偏差的平方成反比。,实验数据分析方法_Chap.5,26,设非等精度观测列的标准偏差分别为1,2,N ,通常把和最大的标准偏差对应的观测值的权定为1，设1maxi (i = 1N )，则标准偏差为i

23、的观测值xi的权为：不难看出p1=1，故x1被称为单位权观测值。对非等精度观测序列被测量的最或然值需要加权平均, 即:标准偏差公式为:权只是从相对意义上表示一个量的精确程度：我们同样可以取和最小的i对应的观测值为单位权观测值；这时虽然各个观测值权的数值和原来不同了，但这些观测值权的比值并未改变。有时为了使所有观测值的权均为整数，可以根据要求选取单位权观测值。,实验数据分析方法_Chap.5,27,由于被测量的真值在有限次观测中是无法得到的，数据处理的任务是通过对被测量的有限次观测求出被测量的最接近于真值的量，即被测量的最或然值。,5.3 直接观测量的最或然值及其精度,5.3.1 最小二乘准则,

24、最小二乘法是求解被测量最或然值的基本方法。按照最或然值的定义，它是最接近于真值的值。设一组观测值为x1, x2, , xN ，待求的最或然值为x*，则它们的残差为i = xi - x* (i =1 N), 最小二乘准则就是选择x*，使得残差平方和为最小。即x*必须满足：,实验数据分析方法_Chap.5,28,对于一列等精度观测列，设由最小二乘准则求出的最或然值为x*，由N个观测值可得N个残差方程：i = xi - x* (i = 1 N)根据最小二乘准则，最或然值x*应满足：由极值原理，有：于是得：设观测值的标准偏差为，则由上式并利用标准偏差的传递公式得：,5.3.2 等精度观测列的最或然值及

25、精度,多次观测取平均可以减小观测结果的随机误差！,实验数据分析方法_Chap.5,29,设x1, x2, , xN为一非等精度观测列，x*为被测量的最或然值，由于各个xi的精度不同，不能像处理等精度观测列那样直接应用 来求解x*，而必须先将它转化为等精度观测列，再利用等精度观测列的最小二乘准则来求最或然值及其精度。 设观测值xi的权为pi，可以证明，只要将每个观测值乘以相应的权的平方根，就可以把原来的非等精度观测列转化为一等精度观测列 ，与之对应的残差序列为 。由最小二乘准则有：,5.3.3 非等精度观测列的最或然值及精度,则非等精度观测列的加权平均值为,实验数据分析方法_Chap.5,30,

26、非等精度观测列的最或然值的标准偏差为：由于非等精度观测列中每个观测值的标准偏差可表示为 ，则上式又可写为：其中为单位权标准偏差，它可按等精度观测列的标准偏差公式计算，但它对应的残差是 ，最后得：,实例,实验数据分析方法_Chap.5,31,间接观测中一种较普遍的情况是观测量为待求量的线性函数。设对直接观测量进行了N次观测，待求的未知量为xk (k = 1 m)，则可得N个观测方程：如果li没有误差且各方程是独立的，则由其中m(mN)个方程可以解出m个未知量的真值。 但实际上观测值总会有误差。如果我们用未知量的最或然值代入上式, 则观测量li与待求量的最或然值的关系可表示成如下的方程组：,5.4

27、 间接观测量的最或然值及其精度,5.4.1 误差方程,式中1,2,N 分别 为l1,l2,lN 的残差。,实验数据分析方法_Chap.5,32,通常称以上方程组为误差方程或条件方程，在这个方程组中有N个方程，m+N个未知量, 即使不考虑vi 的影响，也不能找出严格满足所有方程的解，更何况残差i必须要考虑，但它又是未知的。因此，要求出未知量必须要有附加条件，而使用最小二乘准则能得到这个方程圆满的解。根据最小二乘准则，在等精度观测列的情况下，未知量的最或然值是使残差平方和最小的那些值，即由极值原理，xk (k 1 m) 应满足：,5.4.2 正态方程,实验数据分析方法_Chap.5,33,即：经过

28、简单整理并引用高斯符号，则由此可得到线性方程组常称以上方程组为正态方程或法方程。,实验数据分析方法_Chap.5,34,间接观测另一种常见的情况是观测值是待求量的非线性函数。例如，人造卫星的轨道改正中，观测量是某一历元卫星的球面坐标，待求量是相应历元的六个轨道根数，它们之间的关系是很复杂的非线性关系；利用甚长基线(VLBI)观测测定地球自转参数，观测量是来自射电源同一波前到达VLBI两个测站的钟面时之差即几何延迟，待求量是地球自转参数，它们之间的关系也是很复杂的非线性关系；又如，利用食双星的光变曲线确定其轨道要素是目前测定食双星轨道要素的惟一方法，而食双星的光变曲线不仅和轨道根数有关，还依赖于

29、其它一些因素：包括两颗子星的大小、光度、形状等因此，利用光变曲线得到食双星的轨道要素(称为食双星的测光轨道解)是一个典型的复杂非线性间接观测问题。,实验数据分析方法_Chap.5,35,观测量yi与待求量xk (k = 1 m)之间的非线性关系可写为 设x0k为xk的近似值(或初值)，并用xk表示xk与其近似值之差。则由上式可以算出已知待求量近似值的函数y0k，并记yi yi y0i，对上式在x0k (k 1 m)上进行泰勒展开，并略去xk的二次及二次以上的项，这样可得：其中 (k1m)当x0k给定时为己知系数，下面我们用bik(k1m)表示。因为观测值yi有误差，因此必须考虑yi中的误差，故

30、而得到误差方程：,实验数据分析方法_Chap.5,36,利用最小二乘准则，可得到法方程：解此方程得到xk(k1m)，分别加上近似值x0k(k1m)，就可得待求量的最或然值。当|xk|较大时，可将得到的xk代替原来的近似值x0k重新算出系数bik 和yi并解法方程得到新的xk。这种过程可以反复迭代，直到最后的|xk|值小于给定的误差限为止，这时最后得到的xk即为所求。这种算法常被称为高斯牛顿法或泰勒展开法，此法在求解过程中需反复迭代和修正，逐次迭代的结果将使最后的xk更接近真解。当初值选得较好时，随着迭代次数的增加，修正值|k|将越来越小，即为迭代“收敛”；否则称迭代“发散”：迭代得到的新值可能

31、比原来的值更远离真解，而这种情况在实际应用中时有发生，所以初值的选取是至关重要的。,实验数据分析方法_Chap.5,37,为了求最或然值的标准偏差，必须要知道它们与观测值li的标准偏差之间的关系以及li的标准偏差；要求li的标准偏差，首先要求出li的残差，而这只要将从法方程解得的未知量的最或然值代入误差方程便可得到。(由残差求标准偏差的公式推导请详见书中叙述) 观测值的标准偏差为：其中，Nm称为自由度，意思是指求解m个未知量只需在m个不同条件下测得m个观测值；但现有Nm个测得值，故而多测了Nm个值。从上面的推导可知，用最小二乘法求解未知量时，为了得到较小的标准偏差；通常要求N - m越大越好。

32、,5.4.3 最或然值的标准偏差,=,实验数据分析方法_Chap.5,38,有了观测值的标准偏差后, 就可以求m个最或然值的标准偏差。设m个最或然值的标准偏差为 对应的权分别为 ，则由非等精度观测列的标准偏差公式 可以得到：式中 按观测值标准偏差公式计算。pxk的计算可借助法方程求得，即只要将法方程右端项b1l, b2l,bml 改为1,0,0,0，解此法方程得到的x1即为 ；若把法方程右端项分别改为0,1,0则由 可解得px2。依次类推 ,实验数据分析方法_Chap.5,39,5.5 最小二乘曲线拟合,天文工作中常遇到达样两个问题，其一是：y和x是可被观测的天文量, 且y是x的函数，它们的函

33、数关系由公式(曲线)：y f (x，ck) (k1 m) 给出，但式中含有m个未知参数ck (k 1 m )。我们的任务是根据y和x的N组观测值寻求参数ck的最佳估计k，进而得到以上公式(曲线)具体形式的最佳估计；另一问题是：y和x之间的函数形式未知，而需要利用对y和x的观测求出y和x之间关系的一个经验公式(或经验曲线)。,实验数据分析方法_Chap.5,40,由于观测值总含有误差，通常只能用曲线拟合的方法由y和x的观测值(yi, xi)i=1N，求得理论曲线或经验曲线中参数的估计值。曲线拟合的特点在于，被确定的曲线原则上并不特别要求真正通过给定的所有观测点，而只要尽可能在绝大多数观测点附近通

34、过。这对于含有误差的观测来说较之过所有点的曲线拟合更合理，并有利于减小对未知数据进行预测时的偏差*。 确定表达式中的参数是曲线拟合中的基本问题。另外，经验公式的确定又是参数估计的基础，但它与客观实际联系紧密，必须结合专业知识并依据经验才能得到较好的解决。,实验数据分析方法_Chap.5,41,实验数据分析方法_Chap.5,42,最小二乘法（又称最小平方法）是一种数学优化技术，它通过最小化误差的平方和找到一组数据的最佳函数匹配；其是用最简的方法求得一些绝对不可知的真值，而令误差平方之和为最小；最小二乘法通常用于曲线拟合。很多其他的优化问题也可通过最小化能量或最大化熵用最小二乘形式表达。 180

35、1意大利天文学家朱赛普皮亚齐发现了第一颗小行星谷神星，在40天的跟踪观测后，谷神星运行至太阳背后。皮亚齐失去了谷神星的位置。随后全世界的科学家通过皮亚齐的观测数据开始了寻找谷神星的行动。但是大多数的计算都没有结果，只有当时年仅24岁的高斯成功计算出了谷神星的轨道，奥地利天文学家海因里希奥尔伯斯在高斯计算出的轨道上重新发现了谷神星，从此高斯闻名世界。他的这个最小二乘的方法发表在1809年的著作天体运动论中。法国科学家勒让德也于1806年独立发明最小二乘法。1829年，高斯提供了这个方法较其它方法为优的证明：最小二乘法在很大方面上优化效果强于其它方法，被称为高斯-莫卡夫定理。,实验数据分析方法_C

36、hap.5,43,理论曲线(或经验公式)中参数的估计问题可用如下的数学语言描述：若y是关于自变量x和待定参数ck(k=1m)的形式已知的函数：yf(x,c)。 今给出(x, y)的N对观测值(xi, yi) (i=1N)，要确定参数ck(k=1m)，使某个目标函数 取极值(极大值或极小值)。因此曲线拟合就是对目标函数进行最优化计算，寻求使目标函数d取极值的一组参数值。目标函数的具体形式可根据具体问题的要求来选取, 可以在非最小二乘意义下确定c使得:,5.5.1 目标函数和最优化,实验数据分析方法_Chap.5,44,达到极小。也可以在最小二乘意义下求解c，即使目标函数： 达到极小。 我们称这种

37、选取各观测点的残差平方和作为目标函数的拟合为最小二乘曲线拟合最小二乘曲线拟合用拟合的2量：作为目标函数。寻求使2最小的参数c作为参数的估计值。其中pi为观测值yi的权重因子：,5.5.2 最小二乘曲线拟合,实验数据分析方法_Chap.5,45,满足最小二乘准则的参数值可由下列方程组解出，即由：解此参数的最小二乘估计 k (k=1m)。 线性情况： 理论曲线是未知参数的线性情况时，它的一般形式可表示为对于N组观测值(xi, yi)，把线性函数代入上述方程组，则可得到未知参数c的线性方程组：,实验数据分析方法_Chap.5,46,例如在m=2且为等精度的情况下（2个未知参数），方程组化为：,已知：

38、y=y0(x)+c1f1(x)+c2f2(x) c1f1(xi)f1(xi)+c2 f2(xi)f1(xi)=yi-y0(xi)f1(xi) c1f1(xi)f2(xi)+c2 f2(xi)f2(xi)=yi-y0(xi)f2(xi)解之便可以得到c1和c2 的最佳估计值,实验数据分析方法_Chap.5,47,把参数估计值代入理论关系式，可以得到对应各个自变量xi的y的估计值：线性情况最典型的例子是：这是标准的线性模型，形式简单。但是有些看来较复杂的模型，常常可以通过变量代换的方法简化成这样的形式。下面我们给出几个例子：,例1： 是一个多项式模型，尽管观测值y对自变量而言是非线性的，但它对参数

39、是线性的，因此仍属线性问题。只要作变量代换：,则多项式即可化为标准的线性形式,实验数据分析方法_Chap.5,48,例2： 观测量y对自变量x及参数均为非线性，但通过变量代换仍可 化为线性问题来处理。即对两边取导数，得令 ，得例3：,这是标准的直线模型，解出C0，Cl后，用逆变换求c0，c1：,式中Aj，j (j1，2)分别为周期函数的振幅和初相位，它们都是拟合过程中待估计的参数pj为已知的周期。,实验数据分析方法_Chap.5,49,这个函数形式是非线性的，但我们亦可以通过变量变换将其转化为线性的：,这是以c1, c2, c3, c4参数的标准化模型。由线性情况的最小二乘拟合的参数估计公式解

40、得参数c1,c2,c3 ,c4 后可得周期函数的拟合参数,将它们代入周期函数公式中即得周期函数拟合曲线 变量变换的方法可以把看来较复杂的模型化简，且变换既适用于待定参数也适用于观测量和自变量。这种能通过变量代换的方法化为线性模型的理论或经验公式称为广义线性模型。,实验数据分析方法_Chap.5,50,思考：若把观测量y和x进行调换，最终由上式得到的最小二乘拟合结果是否不变？, 对dy 对dx 对(dx2+dy2)1/2,实例：测定星系中心大质量黑洞的质量,Using RBLR the central mass is:V is the BLR clouds velocity (either fr

41、om FWHM or sLINE)f is a dimensionless factor that depends on the geometry and kinematics of the BLR.如何测定 RBLR ?,Finding the central (black hole) mass is one of the “holy grails” of reverberation mapping in the past decade. (but the sample might be biased.),实验数据分析方法_Chap.5,52,Continuum luminosity var

42、y.BLR respond to the variations (via photoionization).,测定RBLR：Reverberation Mapping,The entire BLR does not respond at the same time. A cloud at a distance R from the central source and angle q to the line of sight will appear to respond after a time:,实验数据分析方法_Chap.5,53,Time,Light curves,Line Flux,C

43、ontinuum Flux,Hb,Kaspi et al. 2000,实验数据分析方法_Chap.5,54,BLR size (RBLR) vs. Luminosity Both are fundamental measured quantities.Peterson et al. (2004) compiled all studies to date. 35 objects with Balmer (mainly Hb) lines time lag. Characteristic BLR size = Time Lag * speed of light. Luminosities in t

44、he Optical, UV, and X-rays. BLR size from averaging all Balmer lines time lags per object.,BLR Size Luminosity Relation,测定 RBLR 进而计算黑洞质量的更普适方法,实验数据分析方法_Chap.5,55,Linear Regression,Uncertainties in both quantities and Intrinsic scatter in the relation Two regression methods:,1. FITEXY from Press et a

45、l. (1992) implemented by Tremaine et al. (2002).,2. BCES (Bivariate Correlated Errors and intrinsic Scatter) by Akritas & Bershady (1996).,and also outlier points,实验数据分析方法_Chap.5,56,Hb RBLR Optical luminosity (5100 A),RBLR lLl(5100 ) (0.690.05),Kaspi et al. 2005,实验数据分析方法_Chap.5,57,Hb RBLR UV luminos

46、ity (1450 A),RBLR lLl(1450 ) (0.560.05),Kaspi et al. 2005,实验数据分析方法_Chap.5,58, 非线性情况： 例：温度为T，面积为A的黑体，其辐射波长为的能量可用下式表示：不难看出，这个理论公式对参数而言也是非线性的，而且也不能通过变量变换转化为线性形式。 对于理论公式是待定参数c的非线性函数的情况，最小二乘解亦应满足准则；但由于相应的2量也是c的非线性函数，因此无法得到关于的解析解。一般情况下，对非线性问题可使用泰勒展开使理论公式线性化，用逐次迭代法求解。,实验数据分析方法_Chap.5,59,把函数yf (x,c)在参数初值 c(

47、0)(c1(0), c2(0), , cm(0) 附近作泰勒展开，略去二次以上的高阶项，得到关于参数c的改正值 的线性函数：,利用线性拟合公式求出参数c的一级改正值：将上述公式中的c(0)换成c(1)，又可以得到c的二级近似值c(2)c(1)+(2)。如此反复迭代计算，由r级近似值c(r)求r +1级近似值c(r+1)的公式为,若从迭代的第r步到第r + 1步参数c的改正值变化很小，即,|i(r) | ；或最小2量减小得很少，即,实验数据分析方法_Chap.5,60,则可停止迭代。这时把参数的第r + 1级近似值作为最佳估计值。这里为所要求的精度，并要求1。 最小2量可以用来检验理论曲线的函数

48、形式是否合理：在观测值服从正态分布的情况下，由2分布定理知 服从自由度为N的2分布，若将其中的参数c用它们的最小二乘估计 代入，则可得最小2量：可以证明：N-m为拟合的自由度。若E（2min） N-m，表明理论曲线与观测值有显著矛盾。,实验数据分析方法_Chap.5,61,在前一小节，我们介绍了一般情况的最小二乘拟合法，即只要目标函数有明显的表达式，就可以用微分法或迭代法求得理论或经验公式中的待定参数。在实际问题中，有的理论公式很复杂，甚至写不出解析表达式，故而目标函数的导数难于求解，这时就不能用迭代法求解参数，而应用最优化方法直接在计算机上求解目标函数的极值，才能得到比较满意的参数估计结果。直接应用优化技术求极值问题的具体方法很多。下面简单介绍较常用且有效的网格搜索法和随机搜索法。,5.5.3 最优化求解,