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

第二部分 实验数据的统计分析 n 第五章 误差理论与最小二乘法 n 第六章 回归分析 n 第七章 多变量分析 n 第八章 功率谱与周期分析 实验数据分析方法 教材 ： 天文数据处理方法 ：丁月蓉 编著 主要参考书 ： 实验的数学处理 ： 李惕陪 著 教学方法： 基本理论 + 具体实例 + 上机实习 (课后 ) 1实验数据分析方法实验数据分析方法 _Chap.5 第五章 误差理论与最小二乘法 n 天文学的诸多理论是以天文观测为基础的 ，如地球自转理论、人造卫星运动理论等都 离不开天文观测。人们通过对某一天文量 ( 静态的或动态的 )的直接或间接观测，获得 大量的数据。而任何观测都不可避免的含有 误差 。因此，当我们在利用观测结果时，必 须分析这些数据的 可靠程度 ：只有当它们的 误差在我们允许的范围之内时，我们才能放 心大胆的去使用它，否则则不能使用。 误差的研究无论是对生产实践还是基础理论误差的研究无论是对生产实践还是基础理论 研究都有着研究都有着 重要意义！重要意义！ 2实验数据分析方法实验数据分析方法 _Chap.5 n 例 1： 由于牛顿在其最初计算中使用了具有较大误差的 地球 半径值 ，使得他测得的月球加速度的值和理论计算值相差约 10，因而推迟了 20年发表他的引力理论！ n 例 2： 爱因斯坦广义相对论的观测证明： 1916年爱因斯坦在 德国 物理学纪事 上发表了具有划时代意义的重要文献 广义相对论基础 。文章指出，当光线行经太阳附近时，光 线产生弯曲，其弯曲曲率预计为 = 1.”75 ，而 1911年他用 经典方法得到 =0.”9 ，相差两倍。如果观测能测得 在 1.”75 附近，这将证明他的广义相对论是正确的，如果测得 的值是在经典值附近，则将否定其理论。幸好 1919年英国天 文学家爱丁顿爵士在西非几内亚湾的普林西比岛的日全食观 测中测得 1.”61 0.”30 ；与此同时有人在巴西东北海岸 外索伯雷尔的日食观测中测得 1.”98 0.”12 。这两个结 果与广义相对论的预言值相近，远大于经典理论值，强有力 的证明了广义相对论的正确性 ! 如果他们当时的观测误差很 大，置信度很低，以致于和理论值相差甚远，那么也就很难 由此来验证这个理论了。由此可见， 观测和误差分析对基础 理论的研究起了一个不可估量的作用 ! 3实验数据分析方法实验数据分析方法 _Chap.5 n 最小二乘法 是用来处理具有误差的观测数据 的一种有效的方法，也是最早用于天文观测资 料处理的一种数学工具。早在 l794年，高斯为 了利用小行星坐标的多次观测准确地推算小行 星的轨道，第一次应用了最小二乘法。 1805年 勒让德应用测量平差方法确定了彗星的轨道和 地球子午线弧长。 1809年高斯又推证了误差的 概率定律，从而使最小二乘法高度完善化，成 为数据处理中应用最广的一个分支。随着概率 统计学和矩阵理论的发展以及电子计算机的广 泛应用，最小二乘法进入了近代数据处理方法 的行列。 4实验数据分析方法实验数据分析方法 _Chap.5 n 误差 是实验科学术语，指 测量结 果偏离真值的程 度。对任何一个物理量进行的测量都不可能得出一 个绝对准确的数值，即使采用测量技术所能达到的 最完善的方法，测出的数值也和真实值存在差异， 这种测量值和真实值的差异称为误差。 (from Wiki) n 误差按其表达形式分： 绝对误差、相对误差 n 误差按其性质及产生原因分： 系统误差 、 随机误差 、 过失（人为）误差 n 误差不仅存在于测量值中，计算时采用近似的理论 模型，计算中一些理论常数的不准确以及数值计算 中取位的多少等也会在计算结果中产生误差。 5.1 误差的定义与分类 5实验数据分析方法实验数据分析方法 _Chap.5 5.1.1 绝对误差和相对误差 n 一个量值的给出值的 绝对误差 定义为该量值的给出 值与其真值之差，或用公式表示为： 绝对误差给出值 -真值 n 公式中的给出值如果是被测量的观测结果，则相 应的误差为观测误差；如果给出值是某量的计算近 似值，则相应的误差为计算近似值的误差。式中的 真值是被测量本身的真实大小，它是一个理想的概 念：一般说来，真值是未知的，通常用约定值来代 替。例如某一系统的天文常数也可看作相应量值的 真值。从绝对误差的定义式不难看出，绝对误差和 被测量具有相同的量纲。因此，若说一颗星其位置 误差为 0.“1，测时的记录误差为 0.“0001，都是指 的绝对误差。 6实验数据分析方法实验数据分析方法 _Chap.5 n 我们把误差的反号值定义为修正值，则可得： 真值给出值 - 误差给出值 + 修正值 这表明，带有误差的给出值加上修正值后可消除 或减小误差的影响。 n 在有些情况下用绝对误差来表示测量的精度是不 恰当的：如目前卫星激光测距的准确度 (测量值与 被测量真值之间的偏离程度 )已达 cm级，卫星的距 离一般为 103km量级；但如果我们测定的是恒星的 距离 (这里指离太阳在 20pc以内的恒星 )，用三角视 差法一般可准确到 0.”02 ，相当于 2pc的测距误差 ，显然它和卫星的测距误差是无法直接比较的！但 如果我们引入 相对误差 的概念，它们的测距误差就 有了可比性。 7实验数据分析方法实验数据分析方法 _Chap.5 n 被测量的绝对误差 与其真值 a之比定义为 这个量的 相对误差 ，并用下式表示： n 当误差较小时，相对误差式中真值 a可用给 定值代替。对于上面的例子，它们测距的 相对误差分别为 1 10和 1 10-1。 即三角视差测量的相对误差反而要比 卫星激光测距的相对误差小！ 8实验数据分析方法实验数据分析方法 _Chap.5 n 由观测的环境因素差异、仪器性能、不同的观 测者等因素造成的按某一确定的规律变化的误差 称为 系统误差 。系统误差的大小和符号在多次重 复观测中几乎相同，通常使观测值往一个方向偏 离。另外，这种误差可以归结为某一因素或某几 个因素的函数，而这种函数通常可以用解析公式 表达出来。人们总是设法找出代表系统误差的解 析表达式，然后在观测结果中 扣除 。 n 由某些难以控制的随机因素造成的，绝对值和 符号的变化时大时小、时正时负，以不可预测的 方式变化的误差称为 随机误差 。虽然就其个体而 言，随机误差没有规律、不可预料，但就其总体 而言，随着观测次数的增加，它又服从某种 统计 规律 。下面我们将从概率论的角度出发讨论随机 误差所满足的统计规律。 5.1.2 系统误差、随机误差和过失误差 9实验数据分析方法实验数据分析方法 _Chap.5 n 古典误差理论认为，随机误差服从 正态分布 ，因 此我们可以用正态分布密度曲线来表征随机误差 ，随机误差的分布密度曲线可表为： 其被称为 高斯误差方程 ，其相应图形也常被称为高斯误差 曲线。式中 称为精密度指数， =x-a, 为随机误差 的均方差。 n 高斯误差方程的一般表达式： 10实验数据分析方法实验数据分析方法 _Chap.5 n 随机误差有下列统计特征，当观测样 本足够大时： (1) 绝对值相等、符号相反的正负误差近于 相等。因此，随机误差的算术平均值随着 观测次数的增加愈来愈小，以零为极限。 (2) 误差的概率与误差的大小有关，绝对值 小的误差出现的概率比绝对值大的误差出 现的概率大，绝对值很大的误差出现的概 率很小。 根据随机误差的这些特征，当不存在系统误差的影响时 ，多次测量结果的平均值将更接近于真值。随机误差产生 的原因很多，观测时环境因素的微小变化，设备中的热噪 声等都是产生随机误差的重要原因。 11实验数据分析方法实验数据分析方法 _Chap.5 n 实际上，系统误差和随机误差之间并没有明显的 界限 有时，我们把一些具有复杂规律但暂末 掌握的系统误差都当作随机误差处理。而随着人 们对误差及其规律的认识的加深，就有可能把这 些以往认识不到因而归之于随机误差的这类误差 确认为系统误差。反之，在一个较短时期内可能 呈现出某种规律，故而归为系统误差，但经过一 段较长时间的观测，发现这种变化规律破坏了， 并呈现出随机性，这就是说，随着时间的推移， 两种不同性质的误差有可能互相转化。 n 过失（人为）误差 是指测量结果与事实明显不符 的一种误差。如观测时对错星或观测过程中望远 镜 /记录仪器的小故障等过失原因造成的结果异常 。这种误差一般比较容易发现，而且只要观测人 员认真细致，基本上是可以避免的。 12实验数据分析方法实验数据分析方法 _Chap.5 n 数据处理中一个很重要的方面是评定一列观测值的 可靠程度。它是指观测结果与真值的一致程度，是观 测结果中系统误差和随机误差大小的综合度量，常用 准确度 这个词来表征。在消除了系统误差之后，观测 的可靠程度由随机误差的大小来衡量。一列观测值精 度高低必须从全列观测值的误差来衡量，而不能只根 据个别值的误差来判断。 n 另外，观测的目的是要从一列观测值中确定 (直接 地或间接地 )被测量的 真值 ，但由于观测手段和观测 次数的限制，真值实际上是测不到的，只能得到它的 一个近似值或估计值。在天文学中通常把最接近于被 测量的真值的一个 近似值 称为它们的 最或然值 ，因此 ，数据处理的又一个重要的问题是给出被测量的最或 然值及其精度。最或然值的精度是衡量观测结果的精 度和处理方法有效性的综合指标。 5.2 观测精度 13实验数据分析方法实验数据分析方法 _Chap.5 n 标准偏差 (又称 均方误差 )是用来衡量一列观测值精 度高低的一个较好指标。 n 设 为被测量的一组观测值， a为被测 量的真值，且 xi 中只包含随机误差，则 称为 xi 的真误差，我们定义 真误差的平方的算 术平均值的平方根 为这列观测值的 标准偏差 或 标准 误差 ，天文上又常称之为中误差，并用 表示，即： 这里定义的标准误差和统计学中从方差的正平方根定义的 标准差是一致的，因为从概率论的角度来说， xi的真值可用 其 数学期望 表示。 5.2.1 精度标准 14实验数据分析方法实验数据分析方法 _Chap.5 n 下面我们来说明标准偏差的大小为什么可以用来 衡量一列观测值的精度高低： 由正态分布的性质可知，观测值 xi在（ a， a+）区间上的概率，或说 i出现在（ ， +） 范围内的概率为 68.3%, 已知 1 2. 则区间（ a1， a+1）小于（ a2 , a+2），也就是说 = 1的观测数据在 a 周围的分布较 密集 ，而 = 的观测值在 a 周围的分布较 分散 ，即标准偏差 的大小可以衡量 一列观测值在真值周围分布的 密度程度 ，而这种密集程度是具有概率含义的， 即误差在（ ， +）内的 置信水平 是 68.3%。 15实验数据分析方法实验数据分析方法 _Chap.5 n 下表列出了一些常用的置信水平误差限： 置信水平 误差限 置信水平 误差限 50.0 % 68.3 % 95.0 % 0.674 1.0 1.96 95.5 % 99.0 % 99.7 % 2 2.58 3 可见，误差落在 3中的概率为 99.7 ,亦即绝对值大于 3的误差仅有 0.3%,这显然是一个小概率事件。所以在有限 次观测中，误差值大于 3的观测值可能含有过失误差，应 考虑舍去该观测值 ; 当然 , 也有可能这个值并不含有过失误 差 , 如舍去它会犯 “ 弃真 ” 错误，但这种误差的最大概率也 只有 0.3%。这种取舍观测值的原则称为 拉依达准则拉依达准则 或简称 为 3 准则准则 。 16实验数据分析方法实验数据分析方法 _Chap.5 高斯函数的性质 17实验数据分析方法实验数据分析方法 _Chap.5 n 在比较两个观测结果时，应在相同的置信水平上比较它们的 误差限，误差限较小的观测较精确，为了说明观测的精度， 通常把观测结果报导为 (置信水平 )。 凡是没有注明置信水平的，一般均指 68.3， 相应的 误差限即为标准误差。 n 在上述各式中，真值 a(或 x)通常是未知的，因此真误差 也是未知的，通常用被测量的最或然值或真值的估计值代替 真值，观测值与其最或然值之差称为观测值的 残差 或 离差 。 标准误差不取决于观测中个别误差的符号，对观测值中较大 误差和较小误差比较灵敏，是表示精度的较好方法。实际应 用中，有时也常用平均误差 离差绝对值的算术平均值来表 示精度；也有时采用概率误差：即绝对值比它大的误差和绝 对值比它小的误差出现的可能性一样大，将误差绝对值按大 小顺序排列，序列的中位数即为概率误差。平均误差和概率 误差只有当 N较大时才较可靠。 n 天体物理中还经常采用 半峰宽度 来表示观测的精度，所谓 半峰宽度 ，即观测值分布曲线在极大值半高度处的全宽（ Full Width at Half Maximum)。 18实验数据分析方法实验数据分析方法 _Chap.5 n 在很多实际问题中 , 待求量往往不能直接观测得到 ，但它们可通过对其它量的观测，再利用它们之间 的函数关系换算求得：这种情况就称为 间接观测 。 间接观测在天文观测中是普遍存在的，例如：在人 造卫星的定轨预报中要测的是卫星在某一历元的轨 道根数，但它们不能直接测得而只能通过测定卫星 的赤经、赤纬换算而得到。 n 对于间接观测的情况，应首先由直接观测量求出间 接观测量的最或然值，然后由直接观测量的精度估 计出间接观测量的精度。 n 通常用下面的式子表示间接观测量 y与 m个直接观 测量 xk (k = 1 m)的关系： 5.2.2 误差传递公式 19实验数据分析方法实验数据分析方法 _Chap.5 n 为了求得间接观测时误差传递的关系，需要对上式 进行线性化处理 如果直接观测量的误差相对于 它们的观测值来说是较小的量，则非线性函数可以 在各个观测值的邻近点上展开成泰勒级数，然后取 误差的一阶项而略去一切高阶误差项： 式中 为观测量 xk的离差，我们把它记为 k。若对 xk(k=1m)各进行了 N次观测，设间接观测量任一次观测的离 差为 y yy0， y0 f(x10， x20， ， xm0), 将 y y+y0, k=xkxk0 代入上式，可得： n 直接观测量 xk的误差以 的形式出现在间接观测量 y的 误差中，或说间接观测量 y的误 差是 m个直接观测量的误差加权和，权重因子 称为 y的 误差传递系数 。 20实验数据分析方法实验数据分析方法 _Chap.5 n 设 m个直接观测量的标准偏差为 ， 根据标准偏差的定义及随机变量方差的运算法则 ，可得间接观测量 y的标准偏差为： 式中 kj为第 k个观测量与第 j个观测量的相关系数。当各 个直接观测量相互独立时，有 kj=0, 则有： 上式通常称为 独立观测量的误差合成定理 。 n 若间接观测量与直接观测量的关系为线性关系时 ，即： 。 则有： 此式即为 线性 情况下的 标准偏差传递公式 。 21实验数据分析方法实验数据分析方法 _Chap.5 例：例： 利用 IRAF进行 测光时，其会根据误差传递以如下 的公式给出测光误差： 根据信噪比的定义： S/N = Flux/Err，故 1/MerrS/N，即 IRAF里给出的测光误差的倒数即为信噪比。除了信噪比会引 起测光误差外，还有很多其他的因素也会带来误差，如减本底 、除平场、减暗流等过程都会带来附加的误差：一般平场的精 度可以达到千分之五左右。 目标源的 测光误差 可以按如下形式给出： Eref为多颗比较星测光误差的平均值， Eobj为目标源测光 误差， Eothers为其他误差，根据不同的情况确定，比如误差 小于千分之五的时候 “ 其他误差 ” 就可能需要包括平场误差 ，再比如比较星的定标误差等。 22实验数据分析方法实验数据分析方法 _Chap.5 n 观测精度的高低是由观测条件决定的，它包括观测的手段 、仪器的精度、观测的次数、观测者技术熟练的程度等 , 因 此我们按观测时的条件把观测分成两大类 : 如果某一列观测 是在完全相同的条件下进行的，则为等精度观测，所得到的 序列称为等精度观测列；如果某一列观测是在不同的条件下 进行的，称为非等精度观测，相应的观测序列为非等精度观 测列。 等精度观测列的标准偏差 对于等精度观测列，可以用全列观测值的标准偏差来衡 量这列观测值的精度。但是，由于观测值的真误差一般是未 知的，为此通常用观测值的残差代替真误差。而对于一列等 精度观测值来说，被测量的最或然值就是这列观测值的算术 平均值 ，则有残差 ，而真误差为 : 5.2.3 等精度观测和非等精度观测 23实验数据分析方法实验数据分析方法 _Chap.5 n 为算术平均值的真误差，对上式两边求平 方和，得： 并有： n 由 线性 情况下的标准偏差传递公式 ，并将 算术平均值的标准偏差 代入上 式则得： n 整理后得到一等精度观测列用残差表示的标准偏 差公式 (这里用高斯符号 表示求和 )： 24实验数据分析方法实验数据分析方法 _Chap.5 权与非等精度观测列 n 处理非等精度观测序列的情况在天文学中是很普遍 的： 例如 利用观测星表编制基本星表就是一个典型的 例子。各种星表中的星位都具有误差；即使是在同一 星表中，它所包含的星位也不都具有相同的标准偏差 。 它们大多数和观测次数的多少有关，故而大多 数星表中有一栏同时列出了各恒星观测的次数，相应 的精度随所用的观测数目的增加而增加。因此，在编 制基本星表时，需根据它们精度的高低区别对待。在 数据处理中，通常用数值 pi表示对某一观测结果 xi的 重视程度，并称之为 权 。观测值精度的高低是和其误 差大小密切相关的：误差越大 , 观测值精度就越低， 对它的重视程度也应相应减小。在观测值只包含随机 误差的情况下，通常定义权与标准偏差的平方成反比 。 25实验数据分析方法实验数据分析方法 _Chap.5 n 设非等精度观测列的标准偏差分别为 1,2, N ,通 常把和最大的标准偏差对应的观测值的权定为 1，设 1 maxi (i = 1N )，则标准偏差为 i的观测值 xi的 权为： 不难看出 p1=1，故 x1被称为单位权观测值。对非等精 度观测序列被测量的最或然值需要加权平均 , 即 : 标准偏差公式为 : n 权只是从 相对意义 上表示一个量的精确程度：我们 同样可以取和最小的 i对应的观测值为单位权观测值 ；这时虽然各个观测值权的数值和原来不同了，但这 些观测值权的比值并未改变。有时为了使所有观测值 的权均为整数，可以根据要求选取单位权观测值。 26实验数据分析方法实验数据分析方法 _Chap.5 n 由于被测量的真值在有限次观测中是无法得到的，数 据处理的任务是通过对被测量的有限次观测求出被测 量的最接近于真值的量，即被测量的最或然值。 5.3 直接观测量的最或然值及其精度 5.3.1 最小二乘准则 n 最小二乘法是求解被测量最或然值的基本方法。按照 最或然值的定义，它是最接近于真值的值。设一组观 测值为 x1, x2, , x N ， 待求的最或然值为 x*，则它们 的残差为 i = xi - x* (i =1 N), 最小二乘准则 就是选 择 x*，使得残差平方和为最小。即 x*必须满足： 27实验数据分析方法实验数据分析方法 _Chap.5 n 对于一列等精度观测列，设由最小二乘准则求出的 最或然值为 x*，由 N个观测值可得 N个残差方程： i = xi - x* (i = 1 N) 根据最小二乘准则，最或然值 x*应满足： 由极值原理，有： 于是得： n 设观测值的标准偏差为 ，则由上式并利用标准偏差的传递 公式得： 5.3.2 等精度观测列的最或然值及精度 多次观测取平多次观测取平 均可以减小观均可以减小观 测结果的随机测结果的随机 误差！误差！ 28实验数据分析方法实验数据分析方法 _Chap.5 n 设 x1, x2, , x N为一非等精度观测列， x*为被测量 的最或然值，由于各个 xi的精度不同，不能像处理 等精度观测列那样直接应用 来求解 x*，而必 须先将它转化为等精度观测列，再利用等精度观测 列的最小二乘准则来求最或然值及其精度。 n 设观测值 xi的权为 pi，可以证明，只要将每个观测值 乘以相应的权的平方根，就可以把原来的非等精度 观测列转化为一等精度观测列 ，与 之对应的残差序列为 。由最小二乘 准则有： 5.3.3 非等精度观测列的最或然值及精度 则 非等精度观测列的加权平均值非等精度观测列的加权平均值 为 29实验数据分析方法实验数据分析方法 _Chap.5 n 非等精度观测列的最或然值的标准偏差为： 由于非等精度观测列中每个观测值的标准偏差可表 示为 ，则上式又可写为： 其中 为单位权标准偏差，它可按等精度观测列的标 准偏差公式计算，但它对应的残差是 ， 最后得： 实例 30实验数据分析方法实验数据分析方法 _Chap.5 n 间接观测中一种较普遍的情况是观测量为待求量的线性函 数。设对直接观测量进行了 N次观测，待求的未知量为 xk (k = 1 m)，则可得 N个观测方程： 如果 li没有误差且各方程是独立的，则由其中 m(m N)个方程可以解 出 m个未知量的 真值 。 n 但实际上观测值总会有误差。如果我们用未知量的最或然 值代入上式 , 则观测量 li与待求量的最或然值的关系可表 示成如下的方程组： 5.4 间接观测量的最或然值及其精度 5.4.1 误差方程 式中 1, 2, N 分别 为 l1,l2, lN 的残差。 31实验数据分析方法实验数据分析方法 _Chap.5 n 通常称以上方程组为 误差方程 或 条件方程 ，在这个 方程组中有 N个方程， m+N个未知量 , 即使不考虑 vi 的影响，也不能找出严格满足所有方程的解，更何 况残差 i必须要考虑，但它又是未知的。因此，要 求出未知量必须要有 附加条件 ，而使用 最小二乘准 则 能得到这个方程圆满的解。 n 根据最小二乘准则，在 等精度观测列 的情况下，未 知量的最或然值是使残差平方和最小的那些值，即 n 由极值原理， xk (k 1 m) 应满足： 5.4.2 正态方程 32实验数据分析方法实验数据分析方法 _Chap.5 即： 经过简单整理并引用高斯符号， 则由此可得到线性 方程组 常称以上方程组为 正态方程 或 法方程 。 33实验数据分析方法实验数据分析方法 _Chap.5 n 间接观测另一种常见的情况是观测值是待求量的 非线性函数 。 n 例如，人造卫星的轨道改正中，观测量是某一历 元卫星的球面坐标，待求量是相应历元的六个轨道 根数，它们之间的关系是很复杂的 非线性关系 ；利 用甚长基线 (VLBI)观测测定地球自转参数，观测量 是来自射电源同一波前到达 VLBI两个测站的钟面时 之差即几何延迟，待求量是地球自转参数，它们之 间的关系也是很复杂的非线性关系；又如，利用食 双星的光变曲线确定其轨道要素是目前测定 食双星 轨道要素的惟一方法，而食双星的光变曲线不仅和 轨道根数有关，还依赖于其它一些因素：包括两颗 子星的大小、光度、形状等 因此，利用光变曲线 得到食双星的轨道要素 (称为食双星的测光轨道解 ) 是一个典型的复杂非线性间接观测问题。 34实验数据分析方法实验数据分析方法 _Chap.5 n 观测量 yi与待求量 xk (k = 1 m)之间的非线性关系 可写为 n 设 x0k为 xk的近似值 (或初值 )，并用 xk表示 xk与其 近似值之差。则由上式可以算出已知待求量近似 值的函数 y0k，并记 yi yi y0i，对上式在 x0k (k 1 m)上进行泰勒展开，并略去 xk的二次及二 次以上的项，这样可得： 其中 (k 1m)当 x0k给定时为己知系数，下面我们 用 bik(k 1m)表示。因为观测值 yi有误差，因此必须考虑 yi中的误差，故而得到误差方程： 35实验数据分析方法实验数据分析方法 _Chap.5 n 利用最小二乘准则，可得到法方程： 解此方程得到 xk(k 1m)，分别加上近似值 x0k(k 1m)，就可得 待求量的最或然值。 n 当 |xk|较大时，可将得到的 xk代替原来的近似值 x0k重新算 出系数 bik 和 yi并解法方程得到新的 xk。这种过程可以反复 迭代，直到最后的 |xk|值小于给定的误差限为止，这时最后 得到的 xk即为所求。这种算法常被称为 高斯 牛顿法 或 泰勒 展开法 ，此法在求解过程中需反复迭代和修正，逐次迭代的 结果将使最后的 xk更接近真解。当 初值 选得较好时，随着迭 代次数的增加，修正值 |k|将越来越小，即为迭代 “ 收敛 ” ；否则称迭代 “ 发散 ” ：迭代得到的新值可能比原来的值更 远离真解，而这种情况在实际应用中时有发生，所以 初值的 选取 是至关重要的。 36实验数据分析方法实验数据分析方法 _Chap.5 n 为了求最或然值的标准偏差，必须要知道它们与观 测值 li的标准偏差之间的关系以及 li的标准偏差； 要求 li的标准偏差，首先要求出 li的残差，而这只 要将从法方程解得的未知量的最或然值代入误差方 程便可得到。 (由残差求标准偏差的公式推导请详 见书中叙述 ) 观测值的标准偏差为： 其中， N m称为 自由度 ，意思是指求解 m个未知量只需 在 m个不同条件下测得 m个观测值；但现有 N m个测得值 ，故而多测了 N m个值。 n 从上面的推导可知，用最小二乘法求解未知量时 ，为了得到较小的标准偏差；通常要求 N - m越大 越好 。 5.4.3 最或然值的标准偏差 = 37实验数据分析方法实验数据分析方法 _Chap.5 n 有了观测值的标准偏差后 , 就可以求 m个最或然值的 标准偏差。设 m个最或然值的标准偏差为 对应的权分别为 ，则由非等精度观测列 的标准偏差公式 可以得到： n 式中 按观测值标准偏差公式计算。 pxk的计算可借 助法方程求得，即只要将法方程右端项 b1l, b2l, bml 改为 1,0,0,0 ，解此法方 程得到的 x1即为 ；若把法方程右端项分别改为 0,1,0 则由 可解得 px2。依次类推 38实验数据分析方法实验数据分析方法 _Chap.5 5.5 最小二乘曲线拟合 n 天文工作中常遇到达样两个问题，其一是： y和 x 是可被观测的天文量 , 且 y是 x的函数，它们的函数 关系由公式 (曲线 )： y f (x， ck) (k 1 m) 给 出， 但式中含有 m个未知参数 ck (k 1 m )。我 们的任务是根据 y和 x的 N组观测值寻求参数 ck的最 佳估计 k，进而得到以上公式 (曲线 )具体形式的最 佳估计； n 另一问题是： y和 x之间的函数形式 未知 ，而需要 利用对 y和 x的观测求出 y和 x之间关系的一个经验公 式 (或经验曲线 )。 39实验数据分析方法实验数据分析方法 _Chap.5 n 由于观测值总含有误差，通常只能用曲线拟合的 方法由 y和 x的 观测值 (yi, xi)i=1N，求得理论曲线 或经验曲线中参数的 估计值 。曲线拟合的特点在于 ，被确定的曲线原则上并不特别要求真正通过给定 的所有观测点，而只要尽可能在绝大多数观测点附 近通过。这对于含有误差的观测来说较之过所有点 的曲线拟合更合理，并有利于减小对未知数据进行 预测时的偏差 *。 n 确定表达式中的参数是曲线拟合中的基本问题。 另外，经验公式的确定又是参数估计的基础，但它 与客观实际联系紧密，必须结合专业知识并依据经 验才能得到较好的解决。 40实验数据分析方法实验数据分析方法 _Chap.5 41实验数据分析方法实验数据分析方法 _Chap.5 n 最小二乘法 （又称 最小平方法 ）是一种 数学 优化 技 术，它通过 最小化 误差 的平方和找到一组数据的最 佳 函数 匹配；其是用最简的方法求得一些绝对不可 知的真值，而令误差平方之和为最小； 最小二乘法 通常用于 曲线拟合 。很多其他的优化问题也可通过 最小化 能量 或最大化 熵 用最小二乘形式表达。 n 1801意大利天文学家 朱赛普 皮亚齐 发现了第一颗小行星 谷 神星 ，在 40天的跟踪观测后，谷神星运行至太阳背后。皮亚 齐失去了谷神星的位置。随后全世界的科学家通过皮亚齐的 观测数据开始了寻找谷神星的行动。但是大多数的计算都没 有结果，只有当时年仅 24岁的高斯成功计算出了谷神星的轨 道，奥地利天文学家 海因里希 奥尔伯斯 在高斯计算出的轨 道上重新发现了谷神星，从此高斯闻名世界。他的这个最小 二乘的方法发表在 1809年的著作 天体运动论 中。法国科 学家 勒让德 也于 1806年独立发明最小二乘法。 1829年，高斯 提供了这个方法较其它方法为优的证明：最小二乘法在很大 方面上优化效果强于其它方法，被称为 高斯 -莫卡夫 定理。 42实验数据分析方法实验数据分析方法 _Chap.5 n 理论曲线 (或经验公式 )中参数的估计问题可用如 下的数学语言描述：若 y是关于自变量 x和待定参 数 ck(k=1m)的形式已知的函数： y f(x,c)。 n 今给出 (x, y)的 N对观测值 (xi, yi) (i=1N)，要确定 参数 ck(k=1m)，使某个目标函数 取极值 (极大值或极小值 )。 n 因此曲线拟合就是对目标函数进行 最优化计算 ，寻求使目标函数 d取极值的一组参数值。目标函 数的具体形式可根据具体问题的要求来选取 , 可 以在非最小二乘意义下确定 c使得 : 5.5.1 目标函数和最优化 43实验数据分析方法实验数据分析方法 _Chap.5 达到极小。也可以在最小二乘意义 下求解 c，即使目标函数： 达到极 小。 n 我们称这种 选取各观测点的残差平方和 作为目标函数的拟合 为 最小二乘曲线拟合 n 最小二乘曲线拟合用拟合的 2量： 作为目标函数。 寻求使 2最小 的参数 c作为参数的 估计值。其中 pi为观测值 yi的权重因子： 5.5.2 最小二乘曲线拟合 44实验数据分析方法实验数据分析方法 _Chap.5 n 满足最小二乘准则的参数值 可由下列方程组解出 ，即由： 解此参数的最小二乘估计 k (k=1m)。 线性情况： 理论曲线是未知参数的线性情况时，它的一般形 式可表示为 对于 N组观测值 (xi, yi)，把线性函数代入上述 方程 组 ，则可得到未知参数 c的线性方程组： 45实验数据分析方法实验数据分析方法 _Chap.5 例如在 m=2且为等精度的情况下（ 2 个未知参数），方程组化为： 已知： 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 的最佳估计值 46实验数据分析方法实验数据分析方法 _Chap.5 n 把参数估计值代入理论关系式，可以得到对应各个 自变量 xi的 y的估计值： n 线性情况最典型的例子是： 这是标准的线性模型，形式简单。但是有些看来 较复杂的模型，常常可以通过变量代换的方法简化 成这样的形式。下面我们给出几个例子： 例 1： 是一个多项式模型，尽管观测值 y对自变量而言是非线性的 ，但它对参数是线性的，因此仍属线性问题。只要作变量 代换： 则多项式即可化为标准的线性形式 47实验数据分析方法实验数据分析方法 _Chap.5 例 2： 观测量 y对自变量 x及参数均为非线性，但通过变量代换仍可 化为线性问题来处理。即对两边取导数，得 令 ，得 例 3： 这是标准的直线模型，解出这是标准的直线模型，解出 C0， Cl后，用逆变换求后，用逆变换求 c0， c1： 式中 Aj， j (j 1， 2)分别为周期函数的振幅和初相位， 它们都是拟合过程中待估计的参数 pj为已知的周期。 48实验数据分析方法实验数据分析方法 _Chap.5 n 这个函数形式是非线性的，但我们亦可以通过变量变换将 其转化为 线性 的： 这是以 c1, c2, c3, c4参数的标准化模型。由线性情况的 最小二乘拟合的参数估计公式解得参数 c1,c2,c3 ,c4 后 可得周期函数的拟合参数 将它们代入周期函数公式中即得周期函数拟合曲线 变量变换的方法可以把看来较复杂的模型化简，且变换既适用于待定参 数也适用于观测量和自变量。这种能通过变量代换的方法化为线性模型的 理论或经验公式称为 广义线性模型广义线性模型 。 49实验数据分析方法实验数据分析方法 _Chap.5 n 思考 ： 若把观测量 y和 x进行调换，最终由上 式得到的最小二乘拟合结果是否不变？ 对 dy 对 dx 对 (dx2+dy2)1/2 50实验数据分析方法实验数据分析方法 _Chap.5 实例： 测定星系中心大质量黑洞的质量 Using RBLR the central mass is: V is the BLR clouds velocity (either from FWHM or LINE) 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.) 51实验数据分析方法实验数据分析方法 _Chap.5 Continuum luminosity vary. BLR respond to the variations (via photoionization). 测定 RBLR： Reverberation Mapping The entire BLR does not respond at the same time. q from the central source and angle RA cloud at a distance to the line of sight will appear to respond after a time: q Line Continuum For a thick shell BLR the response to a continuum flash will be: Time Line flux 52实验数据分析方法实验数据分析方法 _Chap.5 Time Light curves Line Flux Continuum Flux bH Kaspi et al. 2000 53实验数据分析方法实验数据分析方法 _Chap.5 BLR size (RBLR) vs. Luminosity Both are fundamental measured quantities. Peterson et al. (2004) compiled all studies to date. 35 objec