第三章广义反演法精选PPT.ppt

第三章广义反演法第1页，本讲稿共58页内容内容1、广义逆矩阵的概念、广义逆矩阵的概念;2、奇异值分解（、奇异值分解（SVD）和自然逆）和自然逆;3、广义反演法；、广义反演法；4、数据分辨矩阵；、数据分辨矩阵；5、参数分辨矩阵；、参数分辨矩阵；6、特征值的应用；、特征值的应用；7、分辨力高低和方差大小的测度；、分辨力高低和方差大小的测度；8、最佳折衷解；、最佳折衷解；第2页，本讲稿共58页1、广义逆矩阵的概念、广义逆矩阵的概念前面我们讨论了解线性反演问题的长度法，无疑还可以定前面我们讨论了解线性反演问题的长度法，无疑还可以定义其他各式各样的长度，比如义其他各式各样的长度，比如 范数等。但是，由于其范数等。但是，由于其他范数解的应用并非如此广泛，因而，没有必要在这里进他范数解的应用并非如此广泛，因而，没有必要在这里进一步加以论述了。一步加以论述了。这里我们将从另一个角度，即广义逆矩阵的角度讨论线性这里我们将从另一个角度，即广义逆矩阵的角度讨论线性反演问题，并称基于广义逆矩阵建立起来的线性反演法叫反演问题，并称基于广义逆矩阵建立起来的线性反演法叫广义反演法（广义反演法（Gener-alized Inversion），或广义线性反），或广义线性反演法（演法（Generalized Linear Inversion，缩写为，缩写为GLI）。）。第3页，本讲稿共58页设线性反演问题：如果把G看成一个映射算子，那么正演问题就是将模型空间 中的m模型通过算子G映射到数据空间 中的观测数据d通过映射 到模型空间中的模型m的一种运算。第4页，本讲稿共58页第5页，本讲稿共58页 由矩阵理论可知，若G是非奇异矩阵，那么 。这里 是G的逆矩阵，且有：第6页，本讲稿共58页在G是奇异矩阵的情况下，G的逆 并不存在，故我们称 为矩阵G的广义逆。所谓广义逆是矩阵G在常规意义下的逆之推广。普通逆矩阵只是广义逆矩阵的一种特殊形式。显然，在奇异矩阵情况下，：第7页，本讲稿共58页2、奇异值分解（、奇异值分解（SVD）和自然逆）和自然逆为了更好地了解在线性反演中应用相当普为了更好地了解在线性反演中应用相当普遍的奇异矩阵的奇异值分解（遍的奇异矩阵的奇异值分解（Singular Value Decomposition,缩写为缩写为SVD），），我们先从矩阵分解讲起。我们先从矩阵分解讲起。第8页，本讲稿共58页-实对称矩阵的正交分解;任何一个实对称矩阵G均可分解为三个矩阵之连乘积，第一和第三个矩阵分别为G的特征向量矩阵U和它的转置，而第二个矩阵则是G的特征值构成的对角线矩阵 。第9页，本讲稿共58页-非奇异且非对称矩阵的分解;第10页，本讲稿共58页-Lanczos的奇异值分解;（3.25）任何一个MxN阶的矩阵G，均可分解为(3.25）式，即可分解为三个矩阵之乘积。第11页，本讲稿共58页取：（3.29）为矩阵G的逆算子，它被Lanczos称为“自然逆”（natural inverse）。Jackson又称它为Lanczos逆。尔后，大多数学者（如Aki），包括Penros在内都把它称为广义逆。而把基于（3.29）式建立起来的解线性反演问题的方法统称为广义反演法.第12页，本讲稿共58页因而Gm=d的解为：可以证明（3.29）式定义的自然逆满足Penros给出的四个条件。第13页，本讲稿共58页3、广义反演法；、广义反演法；在这节里，我们只涉及基于在这节里，我们只涉及基于Lanczos自然逆而建立起来的自然逆而建立起来的广义反演法，而不讨论基于一般广义逆（即不全部满足广义反演法，而不讨论基于一般广义逆（即不全部满足Penros定义的四个条件的逆）的所谓广义反演法。定义的四个条件的逆）的所谓广义反演法。第14页，本讲稿共58页3、广义反演法；、广义反演法；在这节里，我们只涉及基于在这节里，我们只涉及基于Lanczos自然逆而建立起来的自然逆而建立起来的广义反演法，而不讨论基于一般广义逆（即不全部满足广义反演法，而不讨论基于一般广义逆（即不全部满足Penros定义的四个条件的逆）的所谓广义反演法。定义的四个条件的逆）的所谓广义反演法。第15页，本讲稿共58页设线性反演问题为：Gm=d根据自然逆的定义，有：下面，我们分如下四种情况分别讨论。第16页，本讲稿共58页（1）当M=N=r时，和 均不存在，即 和 均不存在，即 和 都是标准的正交矩阵且：因此，第17页，本讲稿共58页（2）当 时，Gm=d是超定方程。不复存在，但 存在，此时 是正交矩阵，即：而U，是半正交矩阵，即：第18页，本讲稿共58页因此，在这种情况下，广义反演法的解为：第19页，本讲稿共58页（3）当 时，Gm=d是欠定方程。此时，不复存在，而 存在。是正交矩阵，且：而 是半正交矩阵，即：第20页，本讲稿共58页因此，广义反演法的解为：这就是欠定问题的最小长度解，而且解是惟一的。第21页，本讲稿共58页（4）当 时，和 都存在。因此，可以把广义反演解看成是同时在U空间极小 和在V空间极小 的结果。为了帮助大家理解奇异值分解和广义逆的意义，现在分析两个简单的例子。例1:例2:第22页，本讲稿共58页4、数据分辨矩阵；、数据分辨矩阵；用广义反演法解线性反演问题，不但可以求得一个拟合观用广义反演法解线性反演问题，不但可以求得一个拟合观测数据的模型测数据的模型m，而且可以获得一些与观测数据，而且可以获得一些与观测数据d和模型和模型参数参数m有关的辅助信息，例如，数据分辨矩阵（有关的辅助信息，例如，数据分辨矩阵（data resolution matrix）等。）等。第23页，本讲稿共58页假定已经求得模型，即：这里，用 表示用广义反演法构制的模型，以示和真实模型m之区别。试问，能拟合观测数据吗？也就是说，把代入线性方程D=Gm能获得与d相同的重建数据吗？第24页，本讲稿共58页若用 表示重建数据，则有：式中：是 阶方阵，叫数据分辨矩阵（data resolution matrix）或信息密度矩阵（information density matrix）。它是拟合观测数据好坏程度的标志，第25页，本讲稿共58页如图所示，矩阵F的第i行中诸要素 越接近于1，则 越接近于 ，即分辨力越高，因为：第26页，本讲稿共58页由于数据分辨矩阵F主对角线要素 表明 接近 的程度，因此又定义F的对角线矩阵 ，即：为重要性（importance）矩阵。第27页，本讲稿共58页5、参数分辨矩阵；、参数分辨矩阵；由广义反演法构制出来的模型由广义反演法构制出来的模型是真正的模型是真正的模型m吗？为回答这一问题，可先将吗？为回答这一问题，可先将 代入上式，则得：代入上式，则得：式中：式中：R为阶方阵，为阶方阵，R称之为参数分辨矩阵（称之为参数分辨矩阵（parameter resolution matrix）或模型分辨矩阵（）或模型分辨矩阵（model resolution matrix）。它是用广义反演法构制的模型和真正地球物理模型）。它是用广义反演法构制的模型和真正地球物理模型m接近程度的一种重要标志。接近程度的一种重要标志。第28页，本讲稿共58页当 时，。当 时，即在纯超定情况下，才有 。这时R的分辨力最高。当 存在时，。所以 。的每一个要素 ,均可视为m各要素加权的结果，这是因为：第29页，本讲稿共58页如果 ，虽有峰值，但变化比较缓慢，或者其峰值不在R的主对角线上，则R的分辨力不高。分辨力越低，说明模型参数之间越存在相关。和数据分辨矩阵相似，参数分辨矩阵也只是数据核G和反演时所加先验信息的函数，而与观测数据d无关，因此，R矩阵也是实验设计的重要依据。同样，可以定义：为分辨核。越小，R矩阵的分辨能力越高。一般取其倒数 作为分辨能力的（欲称分辨力）的定量度量。第30页，本讲稿共58页6、特征值的应用；、特征值的应用；在这一节中，将要论述利用广义反演提供了一些有用的、在这一节中，将要论述利用广义反演提供了一些有用的、重要的辅助信息重要的辅助信息从特征值所获取的辅助信息。从特征值所获取的辅助信息。第31页，本讲稿共58页1.特征值对观测数据和模型参数的影响将数据核矩阵G作奇异值分解，并代入数据方程得：如用求和形式书写，则有：特征值越大，其对重建观测数据的贡献越大；相反越小，则对的贡献也越小。当反演中大小特征值相差非常悬殊时，小特征值对重建观测数据几乎毫无作用，甚至将它们去掉也不会影响观测数据的重建精度。第32页，本讲稿共58页另一方面，有：其求和形式为：其结论和完全相反，即特征值越小，它对构制的模型参数 影响越大。第33页，本讲稿共58页2.解的方差如果，观测数据具有误差 ，当然用广义反演法所得的结果也有误差 ，且满足：因此，解的协方差矩阵：第34页，本讲稿共58页如果观测数据是统计且独立的，并有相同的方差 ，则：故单位协方差矩阵为：第35页，本讲稿共58页例 解如下联立方程：第36页，本讲稿共58页显然，这是三个未知数三个方程式联立方程。其中：只与第一式有关，各有无限多解，而 与第二、第三式有关，它们却是矛盾方程，无一般意义下的解，现在用广义反演法解之，并分析由此而获得的一些辅助信息，如数据分辨矩阵、参数分辩矩阵和解的方差等。若将上式写成矩阵。即：第37页，本讲稿共58页因此有：第38页，本讲稿共58页由于 和 的特征值完全相同，不难求得它们分别为2，2。显然，此时 是一个混定问题。因此，矩阵G之奇异值分别是：第39页，本讲稿共58页与对称矩阵 和 相对应的特征向量U和V分别是：第40页，本讲稿共58页显然，这是一个混定问题，即 。此时 ，由于矩阵 和 的秩都是2。据此：第41页，本讲稿共58页根据奇异值分解，可得：第42页，本讲稿共58页因而有：能拟合观测数据吗？将 代入数据方程，可得：第43页，本讲稿共58页显然，与d并不完全相同，这是因为数据分辨矩阵：即F不是单位矩阵，所以由：第44页，本讲稿共58页表明，与真值相差甚远。与真实的模型相差多大？可将d=Gm代入：中进行分析，式中：第45页，本讲稿共58页可得：这就说明，而 。由此看来，根据广义反演法，可以惟一地确定 ，而不能惟一地确定 和 的数值大小，只能求得它们的平均值。第46页，本讲稿共58页至于m的协方差：第47页，本讲稿共58页7、分辨力高低和方差大小的测度；、分辨力高低和方差大小的测度；前面讨论了利用观测数据的方差和模型的长度为最小这一前面讨论了利用观测数据的方差和模型的长度为最小这一原则求取线性反演问题的长度解，下面将定义一种利用数原则求取线性反演问题的长度解，下面将定义一种利用数据分辨矩阵据分辨矩阵F，参数分辨矩阵，参数分辨矩阵R和协方差矩阵计算模型参和协方差矩阵计算模型参数的办法。数的办法。第48页，本讲稿共58页由于分辩矩阵（F，R）接近单位矩阵时，说明其分辨力最高，因此一种最好办法是利用非对角线元素之大小（或其展伸情况）来描述分辩力之高低。现以英文Spread表示展伸系数，则有：第49页，本讲稿共58页这种展伸准则有时也叫狄里西莱准则。第50页，本讲稿共58页把目标函数写为：并极小之，可以得到模型的最佳解，式中，为相应项的加权系数。第51页，本讲稿共58页8、最佳折衷解、最佳折衷解在大多数地球物理反演问题中，矩阵在大多数地球物理反演问题中，矩阵G的条件数都很差，的条件数都很差，最大与最小奇异值有时相差几十个级次。我们知道，小的最大与最小奇异值有时相差几十个级次。我们知道，小的奇异值会引起模型参数的最大误差，却能保证模型参数的奇异值会引起模型参数的最大误差，却能保证模型参数的高分辨能力。分辨率和方差是一对矛盾，分辨率高必然方高分辨能力。分辨率和方差是一对矛盾，分辨率高必然方差大；反之，分辨率低、方差也小，二者不可兼得，只能差大；反之，分辨率低、方差也小，二者不可兼得，只能取其折衷。或者以牺牲一些分辨率为代价换取较低的方取其折衷。或者以牺牲一些分辨率为代价换取较低的方差；或者以较大的方差为代价，获得较高的分辨率。差；或者以较大的方差为代价，获得较高的分辨率。第52页，本讲稿共58页Wiggins和Jockson（1972）建议，用广义反演法求解时，设一个最大允许方差t，使即可截断或摒弃小于的特征值。这里t为“方差门槛”值。若特征值按大小顺序排列，即其中仅保留k个大特征值，而截断 个小特征值。显然，应按以下方法计算观测数据的有效自由度q，即有：第53页，本讲稿共58页式中：为矩阵V之要素：去掉 个奇异值相当于把矩阵U和V中最后 个向量用零向量代替。因而，相对应的数据分辨矩阵F和参数分辨矩阵R都发生了变化，设为和 ，则有：第54页，本讲稿共58页此时，相应的广义逆变为：由于将小特征值截断的结果，使定义的分辨力降低了，而方差却大大降低。左图是第k个参数分辨力和方差的示意图。分辨力随k的增加而提高（即 降低），方差则随k值的增大而增大；反演时，令k值从小到大变化，计算 和 ，以 和 为纵、横坐标，作出右图所示的不同k值的折衷曲线，并由此选择最佳折衷k值。第55页，本讲稿共58页第56页，本讲稿共58页在实际地球物理资料反演中，还可以经常见到下图所示的折衷曲线的其他形态。曲线a说明，在最佳折衷点，最佳k值很好确定；当奇异值变化缓慢时，会出现b型折衷曲线或其他类型的折衷曲线，无疑，这都不是我们所期望的。第57页，本讲稿共58页总结：1、广义逆矩阵的概念、广义逆矩阵的概念2、奇异值分解（、奇异值分解（SVD）和自然逆）和自然逆3、广义反演法；、广义反演法；4、数据分辨矩阵；、数据分辨矩阵；5、参数分辨矩阵；、参数分辨矩阵；6、特征值的应用；、特征值的应用；7、分辨力高低和方差大小的测度；、分辨力高低和方差大小的测度；8、最佳折衷解；、最佳折衷解；第58页，本讲稿共58页