地质统计学(北京科技大学 张树泉)
PPT
汇报提纲
研究现状 目标与内容 过程论述 结论
研究现状
多点地质统计学是目前储层地质随机建模的研 究热点方向。相对于传统的两点地质统计学方 法，多点地质统计学在进行储层地质建模时能 够利用更多空间点数据（n>2）之间的相关性 进行统计计算，因此建立的模型更加满足实际 的情况。 多点地质统计学算法snesim是其中一种基于 概率统计的方法，由于snesim算法的参数设 置较为复杂，因此该算法的各种参数对模型的 控制需要进行研究和分析，从而指导如何使用 该算法建立高质量的储层模型。
硬数据条件化
向网格中加载硬数据：
硬数据条件化模拟结果：
非平稳模拟
旋转角控制河道局部的方向 仿射性控制河道宽度
手绘旋转与仿射性图像：
1.3
0.8 0.5 15°
75°
45°
格式转换：
参数文件复制：
参数设置：
运行得到实现：
涠洲11-2油田4井区岩相模拟
研究区属于辫状河三角洲前缘亚相,微相种类较多,面积较大，但 是钻井较少且分布集中,为了减少没有井控制的地方模拟的不确 定性,本次采用录井解释的砂泥岩相进行岩相模拟，相的种类减 少，有助于减少模拟的不确定性.
利用不同油组的平面沉积微相图作为训练图像：
序贯指示法模拟结果：
snesim算法模拟结果：
用多点法模拟的岩相展布图较为连续且具有一定的真实性
结论
设置不同参数对程序的影响： 目标比例能有效控制模拟结果相比例 目标相比例的设置应与训练图像边际 概率相近，否则要考虑更换训练图像 旋转角控制河道局部的方向，仿射性 控制河道宽度 用多点法模拟的岩相展布图较为连续 且具有一定的真实性
地统计学知识点
地统计学知识点地统计学是一门结合了地质学、统计学和数学等多学科知识的交叉学科，主要用于研究自然现象的空间变异性和相关性。
它在地质学、土壤学、生态学、气象学等众多领域都有着广泛的应用。
接下来，让我们逐步深入了解地统计学的一些重要知识点。
一、地统计学的基本概念首先，我们要明白什么是区域化变量。
区域化变量是指在空间上具有随机性和结构性的变量。
比如，某一地区的土壤湿度、温度、海拔高度等都可以看作是区域化变量。
半方差函数是地统计学中的一个核心概念。
它用于描述区域化变量的空间自相关性。
简单来说，就是衡量两个位置的变量值之间的相似程度随着距离的增加是如何变化的。
克里金估计则是基于区域化变量的空间结构和已知数据，对未观测点进行无偏最优估计的一种方法。
二、数据的采集与预处理在进行地统计分析之前，数据的采集至关重要。
为了获得可靠的数据，需要精心设计采样方案，确保样本具有代表性。
采样的方法有很多种，比如简单随机采样、系统采样、分层采样等。
采集到的数据往往需要进行预处理，包括数据的筛选、清理和转换。
例如，去除异常值、填补缺失值、对数据进行标准化或归一化处理，以满足后续分析的要求。
三、空间自相关性分析空间自相关性反映了一个变量在不同位置上的相似程度。
正的空间自相关性表示相似的值倾向于聚集在一起，负的空间自相关性则表示不同的值倾向于聚集。
常用的分析空间自相关性的方法除了前面提到的半方差函数外，还有 Moran's I 指数等。
通过这些方法，可以了解变量的空间分布模式，为进一步的分析和建模提供基础。
四、克里金方法克里金方法有多种类型，如普通克里金、简单克里金、泛克里金等。
每种方法都有其适用的条件和场景。
在使用克里金方法进行估计时，需要确定合适的半方差函数模型。
常见的模型有球状模型、指数模型、高斯模型等。
通过拟合实际数据，选择最优的模型参数，以提高估计的准确性。
五、地统计模拟地统计模拟是通过随机过程生成多个可能的实现，以反映区域化变量的不确定性。
319数理学院 - 中国地质大学(北京)
中国地质大学（北京）硕士研究生《统计学》考试大纲科目名称：统计学科目代码：432一、考试性质《统计学》是为我校招收应用统计硕士专业学位研究生设置的入学资格考试科目。
其目的是科学、公平、有效地测试考生是否具备攻读应用统计硕士专业学位所必须具备的基本素质、应用能力和培养潜能，以利选拔优秀人才入学, 为国家的经济建设培养具有优良的职业道德、法制观念、国际视野、及较强分析与解决实际问题能力的高层次、应用型、复合型统计专业人才。
考试的目标是测试考生是否掌握数据收集、整理、分析和描述等一些基本统计方法。
具体考试要求是：1. 掌握基本的概率论知识、原理和方法；2. 掌握数据收集、统计分析、统计处理的基本原理和方法；3. 具有统计建模及用统计方法解释数据的基本能力；4. 具有使用常用统计软件分析统计数据、做出统计推断的能力。
二、考试形式与试卷结构1. 试卷总分：150分。
2. 考试时间：180分钟（3小时）。
3. 试卷形式与结构：概率论60分，由以下三种题型构成：单项选择题：10题，每小题2分，共计20分；填空题：10题，每小题2分，共计20分；计算与证明题：2题，每小题10分，共计20分。
统计学90分，由以下三种题型构成：单项选择题：15题，每小题2分，共计30分；填空题：15题，每小题2分，共计30分；计算与分析题：3题，每小题10分，共计30分。
4. 答题方式：闭卷、笔试。
5. 注意事项：允许使用仅具备四则运算和开方运算功能的简单计算器，但不得使用带有公式和文本存储等复杂功能的高级计算器。
三、考试内容1. 概率论⑴. 随机试验、样本空间、随机事件；⑵. 事件的关系、运算及运算性质；⑶. 概率的定义与性质、加法公式；⑷. 条件概率、乘法公式、全概率公式、贝叶斯公式；⑸. 随机变量、分布函数的概念；⑹. 离散型随机变量及其分布、常用离散分布：两点分布，二项分布、泊松分布、负二项分布、几何分布、超几何分布；⑺. 连续型随机变量及其分布、常用连续分布：均匀分布、指数分布、正态分布；⑻. 多维向量的分布、边缘分布、条件分布；常用多维离散分布：多项分布；常用多维连续分布：均匀分布、正态分布；⑼. 数字特征：期望、方差、协方差、相关系数、矩；⑽. 大数定律，中心极限定理。
CAP方法反演震源机制的误差分析_以胶东半岛两次显著中等地震为例
第５８卷第２期２０１５年２月地 球 物 理 学 报ＣＨＩＮＥＳＥ ＪＯＵＲＮＡＬ ＯＦ ＧＥＯＰＨＹＳＩＣＳＶｏｌ．５８，Ｎｏ．２Ｆｅｂ．，２０１５郑建常，林眉，王鹏等．２０１５．ＣＡＰ方法反演震源机制的误差分析：以胶东半岛两次显著中等地震为例．地球物理学报，５８（２）：４５３－４６２，ｄｏｉ：１０．６０３８／ｃｊｇ２０１５０２０９．Ｚｈｅｎｇ Ｊ Ｃ，Ｌｉｎ Ｍ，Ｗａｎｇ Ｐ，ｅｔ ａｌ．２０１５．Ｅｒｒｏｒ ａｎａｌｙｓｉｓ ｆｏｒ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ｆｒｏｍ ＣＡＰ ｍｅｔｈｏｄ ｉｎｖｅｒｓｉｏｎ：Ａｎ ｅｘａｍｐｌｅ ｏｆ ２ｍｏｄｅｒａｔｅ ｅａｒｔｈｑｕａｋｅｓ ｉｎ Ｊｉａｏｄｏｎｇ Ｐｅｎｉｎｓｕｌａ．Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．（ｉｎ Ｃｈｉｎｅｓｅ），５８（２）：４５３－４６２，ｄｏｉ：１０．６０３８／ｃｊｇ２０１５０２０９．ＣＡＰ方法反演震源机制的误差分析：以胶东半岛两次显著中等地震为例郑建常，林眉，王鹏，徐长朋山东省地震局，济南 ２５００１４摘要 利用区域波形数据使用ＣＡＰ方法反演中强地震的震源机制正逐渐得到广泛应用．本文以胶东半岛近期发生的两次显著中等地震为例，讨论了使用ＣＡＰ方法反演震源机制时的误差估计，展示了反演结果的不确定性分析过程．２０１３年１１月２３日和２０１４年１月７日在山东莱州和乳山分别发生了Ｍ４．６和Ｍ４．３级中等地震，两次事件均造成了较大影响．我们基于ＣＡＰ方法，使用自助抽样（ｂｏｏｔｓｔｒａｐ）技术多次重复反演过程，得到大样本量的震源机制解数据；基于这些数据，使用粒子群算法和聚类分析技术给出了优化解，估计了震源机制解的误差范围，并利用震源机制解的Ｐ、Ｔ轴给出了震源球上的概率密度分布．关键词 莱州地震；乳山地震；波形反演；聚类分析；不确定性ｄｏｉ：１０．６０３８／ｃｊｇ２０１５０２０９中图分类号 Ｐ３１５收稿日期２０１４－０２－２７，２０１４－０７－２４收修定稿基金项目 国家科技支撑计划项目（２０１２ＢＡＫ１９Ｂ０４－０１－０５）、山东省自然科学基金（ＺＲ２０１２ＤＱ００６）及中国地震局监测预报司震情跟踪工作专项（２０１４０２０１０３）资助．作者简介 郑建常，１９７８年生，山东临清人，副研，２０１１年于中国地震局地球物理研究所获博士学位，主要从事地震活动性及数字地震学研究．Ｅ－ｍａｉｌ：ｚｊｃｍａｉｌ＠ｙｅａｈ．ｎｅｔＥｒｒｏｒ ａｎａｌｙｓｉｓ ｆｏｒ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ｆｒｏｍ ＣＡＰ ｍｅｔｈｏｄ ｉｎｖｅｒｓｉｏｎ：Ａｎ ｅｘａｍｐｌｅ ｏｆ ２ｍｏｄｅｒａｔｅ ｅａｒｔｈｑｕａｋｅｓ ｉｎ Ｊｉａｏｄｏｎｇ ＰｅｎｉｎｓｕｌａＺＨＥＮＧ Ｊｉａｎ－Ｃｈａｎｇ，ＬＩＮ Ｍｅｉ，ＷＡＮＧ Ｐｅｎｇ，ＸＵ Ｃｈａｎｇ－ＰｅｎｇＥａｒｔｈｑｕａｋｅ Ａｄｍｉｎｉｓｔｒａｔｉｏｎ ｏｆ Ｓｈａｎｄｏｎｇ Ｐｒｏｖｉｎｃｅ，Ｊｉ′ｎａｎ ２５００１４，ＣｈｉｎａＡｂｓｔｒａｃｔ Ａｓ ａｎ ｅｆｆｅｃｔｉｖｅ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍ ｉｎｖｅｒｓｉｏｎ ｍｅｔｈｏｄ ｆｏｒ ｒｅｇｉｏｎａｌ ｅａｒｔｈｑｕａｋｅｓ，ＣＡＰ（Ｃｕｔａｎｄ Ｐａｓｔｅ）ｉｓ ｗｉｄｅｌｙ ｕｓｅｄ ｉｎ Ｃｈｉｎａ ｉｎ ｒｅｃｅｎｔ ｙｅａｒｓ．Ｉｔｓ ｑｕａｌｉｔｙ ａｎｄ ｅｒｒｏｒ ｌｅｖｅｌ ｎｅｅｄ ｔｏ ｂｅ ｅｖａｌｕａｔｅｄｆｏｒ ｓｕｃｈ ｓｏｌｕｔｉｏｎｓ ａｒｅ ｉｎｃｒｅａｓｉｎｇｌｙ ｒｅｔｒｉｅｖｅｄ．Ｏｎ ｔｈｅ ｏｔｈｅｒ ｈａｎｄ，ｗｈｅｎ ｕｓｉｎｇ ｔｈｅ ＣＡＰ ｍｅｔｈｏｄ，ｔｈｏｓｅ ｐｈａｓｅｓ ｗｈｉｃｈ ｆｉｔｔｅｄ ｗｅｌｌ ａｒｅ ｕｓｕａｌｌｙ ｃｈｏｓｅｎ ｆｏｒ ｉｎｖｅｒｓｉｏｎ，ａｎｄ ｔｈｅ ｏｔｈｅｒ ｐｈａｓｅｓ ｗｈｉｃｈ ａｒｅｔｈｏｕｇｈｔ‘ｂａｄ’ｏｒ ｆｉｔｔｅｄ ｎｏｔ ｓｏ ｇｏｏｄ ａｒｅ ｉｇｎｏｒｅｄ．Ｉｔ ｈａｓ ｂｅｅｎ ｆｏｕｎｄ ｔｈａｔ ｄｉｆｆｅｒｅｎｔ ｓｔａｔｉｏｎｃｏｍｂｉｎａｔｉｏｎｓ ｗｉｌｌ ｌｅａｄ ｔｏ ｖａｒｉｅｄ ｒｅｓｕｌｔｓ ｗｉｔｈ ｕｎｎｅｇｌｅｃｔａｂｌｅ ｄｉｓｃｒｅｐａｎｃｉｅｓ．Ｏｂｊｅｃｔｉｖｅｌｙ ｓｐｅａｋｉｎｇ，ｉｎａ ｓｃｉｅｎｔｉｆｉｃ ｐｅｒｓｐｅｃｔｉｖｅ，ｔｈｉｓ ａｒｔｉｆｉｃｉａｌｌｙ ｓｅｌｅｃｔｅｄ ｐｒｏｃｅｓｓ ｗｉｌｌ ｉｎｃｒｅａｓｅ ｕｎｃｅｒｔａｉｎｔｉｅｓ ｉｎ ｆｉｎａｌｉｎｖｅｒｓｉｏｎ ｒｅｓｕｌｔｓ，ｅｓｐｅｃｉａｌｌｙ ｕｎｄｅｒ ｔｈｅ ｐｒｅｓｅｎｔ ｉｎｓｔｒｕｍｅｎｔ ｓｔａｔｕｓ ｉｎ Ｃｈｉｎａ．Ｆｕｒｔｈｅｒｍｏｒｅ，ｂｅｃａｕｓｅｔｈｅ ｇｒｉｄ ｓｅａｒｃｈ ｓｃｈｅｍｅ ｕｓｅｄ ｉｎ ｔｈｅ ＣＡＰ ｍｅｔｈｏｄ ｗｈｉｃｈ ｉｓ ｎｏｔ ｅｖｅｎｌｙ ｄｉｓｔｒｉｂｕｔｅｄ ｏｎ ｆｏｃａｌ ｓｐｈｅｒｅ，ｗｅ ｃａｎ ｎｏｔ ｇｉｖｅ ａ ｃｏｎｖｉｎｃｉｎｇ ｐｒｏｏｆ ｔｏ ｉｌｌｕｓｔｒａｔｅ ｔｈａｔ ｗｈｅｔｈｅｒ ｔｈｅ ｃａｕｓｅ ｏｆ ｂａｄｌｙ－ｆｉｔｔｅｄ ｗａｖｅｆｏｒｍｓｃｏｍｅｓ ｆｒｏｍ ｄａｔａ ｅｒｒｏｒ ｏｒ ｆｒｏｍ ｕｎ－ｓｕｆｆｉｃｉｅｎｔ ｓｅａｒｃｈｅｄ ｓｏｌｕｔｉｏｎｓ．Ｔｗｏ ｅａｒｔｈｑｕａｋｅｓ ｏｆ Ｍ４．６ａｎｄ Ｍ４．３ｏｃｃｕｒｒｅｄ ｉｎ Ｓｈａｎｄｏｎｇ Ｐｅｎｉｎｓｕｌａ ｏｎ ２３Ｎｏｖｅｍｂｅｒ ２０１３地球物理学报（Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．）５８卷 ａｎｄ ７Ｊａｎｕａｒｙ ２０１４，ｒｅｓｐｅｃｔｉｖｅｌｙ．Ｔｈｅ ｆｏｒｍｅｒ ｉｓ ｔｈｅ ｌａｒｇｅｓｔ ｅｖｅｎｔ ｉｎ ｔｈｅ Ｓｈａｎｄｏｎｇ ａｒｅａ ｓｉｎｃｅ１９９５．Ｆｏｒ ｃｏｎｖｅｎｉｅｎｃｅ，ｗｅ ｌａｂｅｌ ｔｈｅ ｆｏｒｍｅｒ ｅｖｅｎｔ ａｓ ｔｈｅ Ｌａｉｚｈｏｕ ｅａｒｔｈｑｕａｋｅ，ａｎｄ ｔｈｅ ｌａｔｔｅｒ ａｓ ｔｈｅＲｕｓｈａｎ ｅａｒｔｈｑｕａｋｅ，ａｃｃｏｒｄｉｎｇ ｔｏ ｔｈｅｉｒ ｅｐｉｃｅｎｔｅｒｓ．Ｔａｋｉｎｇ ｔｈｅｓｅ ｔｗｏ ｅｖｅｎｔｓ ｆｏｒ ｅｘａｍｐｌｅ，ｔｈｉｓｐａｐｅｒ ｄｉｓｃｕｓｓｅｓ ｅｒｒｏｒ ｅｓｔｉｍａｔｉｏｎ ｆｏｒ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍ ｉｎｖｅｒｓｉｏｎ ｕｓｉｎｇ ｔｈｅ ＣＡＰ ｍｅｔｈｏｄ．Ｔｈｅ ｐａｐｅｒａｌｓｏ ｐｒｅｓｅｎｔｓ ａｎ ｕｎｃｅｒｔａｉｎｔｙ ａｎａｌｙｓｉｓ ｐｒｏｃｅｓｓ ｆｏｒ ｔｈｅ ｉｎｖｅｒｓｉｏｎ ｒｅｓｕｌｔｓ．Ｂｒｉｅｆｌｙ，ａ ｂｏｏｔｓｔｒａｐ ｔｅｃｈｎｉｑｕｅ ｉｓ ａｄｏｐｔｅｄ，ｗａｖｅｆｏｒｍｓ ａｒｅ ｒａｎｄｏｍｌｙ ｓａｍｐｌｅｄ ｗｉｔｈ ｅｑｕａｌｐｒｏｂａｂｉｌｉｔｙ ｆｒｏｍ ｏｒｉｇｉｎ ｄａｔａｓｅｔ，ａｎｄ ｔｈｅｎ ｕｓｅｄ ａｓ ｄａｔａ ｆｏｒ ａ ｒｅｐｅａｔｉｎｇ ｉｎｖｅｒｓｉｏｎ ｐｒｏｃｅｄｕｒｅ．Ａｆｔｅｒ ａｌａｒｇｅ ｎｕｍｂｅｒ ｏｆ ｉｎｖｅｒｓｉｏｎｓ，ｅ．ｇ．，１０００ｔｉｍｅｓ，ｗｅ ｆｉｎａｌｌｙ ｇｅｔ ｂｏｏｔｓｔｒａｐ ｒｅｓｕｌｔｓ ｃｏｎｓｉｓｔｉｎｇ ｏｆ １０００ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ．Ｂａｓｅｄ ｏｎ ｔｈｅｓｅ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ，ｗｅ ｃｏｎｄｕｃｔ ｔｈｅ ｆｏｌｌｏｗｉｎｇ ｗｏｒｋ．（１）Ｗｅ ｅｍｐｌｏｙ ａ ＰＳＯ（ｐａｒｔｉｃｌｅ ｓｗａｒｍ ｏｐｔｉｍｉｚａｔｉｏｎ）ａｌｇｏｒｉｔｈｍ ｔｏ ｓｅａｒｃｈ ａ ｓｏｌｕｔｉｏｎ，ｏｆ ｗｈｉｃｈｔｈｅ Ｋａｇａｎ ａｎｇｌｅ ｉｓ ｍｉｎｉｍａｌ ｔｏ ａｌｌ ｔｈｅ ｄｏｕｂｌｅ ｃｏｕｐｌｅ ｍｏｄｅｌｓ，ａｎｄ ｕｓｅ ｉｔｓ ｓｔａｎｄａｒｄ ｄｅｖｉａｔｉｏｎ ａｓ ｔｈｅｉｎｖｅｒｓｉｏｎ ｒｅｓｕｌｔｓ ｕｎｃｅｒｔａｉｎｔｙ ｒａｎｇｅ．Ｔｈｅ ｕｎｃｅｒｔａｉｎｔｙ ｏｆ ｓｏｌｕｔｉｏｎｓ ｆｏｒ ｔｈｅ Ｌａｉｚｈｏｕ ｅｖｅｎｔ ｉｓ±２３．７°；ｆｏｒＲｕｓｈａｎ ｅｖｅｎｔ，ｉｓ±６．４°．Ｔｈｅｓｅ ｔｗｏ ｓｏｌｕｔｉｏｎｓ ｃａｎ ｂｅ ｅｖａｌｕａｔｅｄ ａｓ‘Ａ’ｌｅｖｅｌ ａｃｃｏｒｄｉｎｇ ｔｏ ｔｈｅＨａｒｄｅｂｅｃｋ′ｓ ｉｎｄｉｃａｔｏｒ ｏｆ ｍｅｃｈａｎｉｓｍｓ ｑｕａｌｉｔｙ．（２）Ａ ｃｌｕｓｔｅｒｉｎｇ ａｎａｌｙｓｉｓ ｉｓ ｕｓｅｄ ｆｏｒ ｂｏｏｔｓｔｒａｐ ｒｅｓｕｌｔｓ．Ｆｏｒ ｔｈｅ Ｒｕｓｈａｎ ｅｖｅｎｔ，ｔｈｅ ｃｌｕｓｔｅｒｉｎｇｃｅｎｔｅｒ ｉｓ ｃｏｉｎｃｉｄｅｎｔ ｔｏ ｔｈｅ ＰＳＯ ｏｐｔｉｍｉｚｅｄ ｓｏｌｕｔｉｏｎ；ｆｏｒ ｔｈｅ Ｌａｉｚｈｏｕ ｅｖｅｎｔ，ｓｅｖｅｒａｌ ｃｅｎｔｅｒｓ ａｒｅｆｏｕｎｄ ｉｎ ｃｌｕｓｔｅｒｉｎｇ，ｄｅｓｐｉｔｅ ｉｓｏｌａｔｅ ｓｏｌｕｔｉｏｎｓ．Ｔｈｅｒｅ ａｒｅ ｔｗｏ ｃｌｕｓｔｅｒｉｎｇ ｃｅｎｔｅｒｓ，ｏｆ ｗｈｉｃｈ ｔｈｅｃｏｒｒｅｓｐｏｎｄｉｎｇ ｄａｔａ ｐｒｏｐｏｒｔｉｏｎ ｉｓ ａｂｏｕｔ ９８．７ｐｅｒｃｅｎｔ．（３）Ｐｒｏｊｅｃｔｉｎｇ Ｐａｎｄ Ｔａｘｅｓ ｏｆ ｂｏｏｔｓｔｒａｐ ｒｅｓｕｌｔｓ ｏｎｔｏ ａ ｆｏｃａｌ ｓｐｈｅｒｅ，ｃａｌｃｕｌａｔｉｎｇ ｉｔｓｐｒｏｂａｂｉｌｉｔｙ ｄｅｎｓｉｔｙ，ｗｅ ｇｅｔ ｔｈｅ ｐｒｏｂａｂｉｌｉｔｙ ｄｅｎｓｉｔｙ ｄｉｓｔｒｉｂｕｔｉｏｎ ｏｆ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ｏｎ ｔｈｅ ｆｏｃａｌｓｐｈｅｒｅ．Ｔｈｅｎ ｗｅ ｃａｎ ｇｉｖｅ ｔｈｅ ｃｏｎｆｉｄｅｎｃｅ ｉｎｔｅｒｖａｌ ｏｎ ｄｉｆｆｅｒｅｎｔ ｌｅｖｅｌｓ ｆｏｒ ｍｅｃｈａｎｉｓｍ ｓｏｌｕｔｉｏｎｓ．Ｔｈｅ ｍｅｔｈｏｄ ｄｅｍｏｎｓｔｒａｔｅｄ ｉｎ ｔｈｉｓ ｐａｐｅｒ ｉｓ ｎｏｔ ｃｏｎｆｉｎｅｄ ｔｏ ａｃｈｉｅｖｉｎｇ ｍｏｒｅ ａｃｃｕｒａｔｅ ｆｏｃａｌｍｅｃｈａｎｉｓｍ ａｎｄ ｏｂｔａｉｎｉｎｇ ｒａｔｉｏｎａｌ ｉｎｖｅｒｓｉｏｎ ｅｒｒｏｒ，ｗｈｉｌｅ ｉｔ ｃａｎ ａｌｓｏ ｂｅ ｕｓｅｄ ｔｏ ｅｘｃｌｕｄｅ ｉｓｏｌａｔｅ ａｎｄｉｎｃｏｒｒｅｃｔ ｓｏｌｕｔｉｏｎｓ ｅｆｆｅｃｔｉｖｅｌｙ，ａｎｄ ａｖｏｉｄ ｔｈｅ ｅｆｆｅｃｔ ｏｆ ｄａｔａ ｆｒｏｍ ｓｔａｔｉｏｎｓ ｗｉｔｈ ｌａｒｇｅｒ ｄｉｓｔｕｒｂａｎｃｅｓ．Ｔｈｅｒｅｆｏｒｅ，ｔｈｉｓ ｍｅｔｈｏｄ ｃａｎ ｂｅ ｕｓｅｄ ｔｏ ｉｎｖｅｒｔ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ａｕｔｏｍａｔｉｃａｌｌｙ ｉｍｍｅｄｉａｔｅｌｙ ａｆｔｅｒｍｏｄｅｒａｔｅ ｅａｒｔｈｑｕａｋｅｓ ｏｃｃｕｒ．Ｓｕｓｔａｉｎｉｎｇ ｂｙ ｐｏｗｅｒｆｕｌ ｃｏｍｐｕｔａｔｉｏｎａｌ ｃａｐａｂｉｌｉｔｉｅｓ，ｗｅ ｃａｎ ｇｅｔ ｍｏｒｅａｃｃｕｒａｔｅ ａｎｄ ｒｅｌｉａｂｌｅ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍ ｒｅｓｕｌｔｓ ｗｉｔｈｏｕｔ ｍａｎｕａｌ ｗｏｒｋ．Ｋｅｙｗｏｒｄｓ Ｌａｉｚｈｏｕ ｅａｒｔｈｑｕａｋｅ；Ｒｕｓｈａｎ ｅａｒｔｈｑｕａｋｅ；Ｗａｖｅｆｏｒｍ ｉｎｖｅｒｓｉｏｎ；Ｃｌｕｓｔｅｒｉｎｇ ａｎａｌｙｓｉｓ；Ｕｎｃｅｒｔａｉｎｔｙ１ 引言据山东台网测定，２０１３年１１月２３日１３时４４分在山东省莱州市（３７．１０°Ｎ，１２０．０２°Ｅ）发生Ｍ４．６级地震，这次地震是山东陆地地区自１９９５年苍山５．２级地震后发生的最大地震，影响范围广，山东东部市地普遍有感；２０１４年１月７日２２时２４分在山东乳山（３６．８０°Ｎ，１２１．７０°Ｅ）发生Ｍ４．３级地震，这次地震也造成胶东地区大面积有感．这两次事件是１９７０年以来胶东半岛陆地及近海地区发生的最强烈的地震活动，其中莱州地震震中区在１９７０年以来的小震目录上属于典型的少震、弱震区，活动水平不高，很少有ＭＬ≥３．０地震发生，仅在１９９１年２月以及２０１２年７月分别发生最大ＭＬ３．８级小震序列各一次；乳山地震震中区历史上曾发生公元１０４６年岠嵎山５１／２级和１９３９年乳山下初５１／２级地震．虽然这两次地震的震中区历史上没有强烈地震活动，但胶东半岛北部近海曾发生多次６、７级强震，如１５４８年渤海海峡７．０、１９４８年威海近海６．０以及１９６９年渤海７．４级等．因此确定这两次显著中等地震的震源机制对于研究区域地质构造的活动特征，以及研判该地区的地震危险性等具有重要的科学价值．我们使用近年来在国内得到广泛使用的ＣＡＰ（Ｃｕｔ ａｎｄ Ｐａｓｔｅ）方法反演这两次地震的震源机制．为了得到更准确的解，并且合理地估计反演结果的４５４ ２期郑建常等：ＣＡＰ方法反演震源机制的误差分析：以胶东半岛两次显著中等地震为例不确定性，我们使用自助抽样法（ｂｏｏｔｓｔｒａｐ）对反演过程随机重复，在大样本量的反演结果基础上，使用粒子群算法搜索优化解，利用动态聚类技术对结果进行聚类分析，从而得到了更加稳定可靠的断层面解，给出了可能的误差范围，并进一步给出了震源球上Ｐ、Ｔ轴的概率密度分布．２ 理论与方法２．１ ＣＡＰ方法反演震源机制震源机制和传播效应决定了观测波形的变化．如果地壳模型已知，可以准确地计算波形传播过程中的效应，因此我们可以通过理论波形ｓ（ｔ）和观测波形ｕ（ｔ）的拟合来估计震源的断层面参数．双力偶震源产生的理论位移ｓ（ｔ）可以表示为（Ｚｈｕ ａｎｄ Ｈｅｌｍｂｅｒｇｅｒ，１９９６）：ｓ（ｔ）＝Ｍ０∑３ｉ＝１Ａｉ（φ－θ，δ，λ）Ｇｉ（ｔ），（１）其中，ｉ＝１，２，３对应三种基本断层响应，即：垂直走滑、垂直倾滑以及倾角为４５°的倾滑；Ｇｉ为格林函数，Ａｉ是辐射系数，φ是台站方位角，Ｍ０为标量地震矩．θ，δ，λ分别为断层的走向、倾角、滑动角．系数Ａｉ由６个矩张量分量和台站方位角表示如下（Ｊｏｓｔａｎｄ Ｈｅｒｒｍａｎｎ，１９８９）：Ａ１＝１２（Ｍｙｙ－Ｍｘｘ）ｃｏｓ２φ＋Ｍ１２ｓｉｎ２φ，Ａ２＝Ｍｘｚｃｏｓφ＋Ｍｙｚｓｉｎφ，Ａ３＝－１２（Ｍｘｘ＋Ｍｙｙ），Ａ４＝１２（Ｍｘｘ－Ｍｙｙ）ｓｉｎ２φ－Ｍｘｙｃｏｓ２φ，Ａ５＝－Ｍｙｚｃｏｓφ＋Ｍｘｚｓｉｎφ．（２） 走向θ、倾角δ、滑动角λ，以及标量地震矩Ｍ０等可以通过求解以下方程进行估计：ｕ（ｔ）＝ｓ（ｔ）．（３） 波形反演可以使用全波形数据，也可以单独使用体波或面波震相进行拟合．ＣＡＰ方法是一种联合使用体波和面波进行反演的方法，近年来在国内得到了广泛的应用（吕坚等，２００８；黄建平等，２００９；郑勇等，２００９；龙锋等，２０１０；韩立波等，２０１２），由于该方法分别截取波形的Ｐｎｌ部分和面波部分分别拟合（Ｚｈａｏ ａｎｄ Ｈｅｌｍｂｅｒｇｅｒ，１９９４；Ｚｈｕ ａｎｄ Ｈｅｌｍｂｅｒｇｅｒ，１９９６），并在反演的过程中允许它们在适当的时间变化范围内相对移动，在一定程度上避免了因为地壳模型不准确而引起的震相到时的误差因素，对速度模型和地壳横向变化的依赖性较小，因此在实际的区域地震震源机制求解中有明显的优势．ＣＡＰ方法使用频率Ｆ－波数Ｋ法（Ｚｈｕ ａｎｄ Ｒｉｖｅｒａ，２００２）计算格林函数，使用网格搜索方法搜寻最优震源机制参数和震源深度．考虑到波形随震中距的衰减，方法定义误差函数如下：ｅ＝‖ｒｒ（）０ｐ‖·‖ｕ－ｓ‖，（４）式中，ｒ为台站震中距，ｒ０为选定的参考震中距，ｐ为指数因子．参考有关研究，对体波ｐ＝１，面波ｐ＝０．５（韩立波和蒋长胜，２０１２）．２．２ ＣＡＰ方法的优化解及其不确定性估计地球内部的任意震源可以表示为６个独立分量的矩张量，由于ＣＡＰ方法限制震源为双力偶模型，并且无需发震时刻的对齐，因此只需对震源模型的三个角度，即走向θ、倾角δ、滑动角λ，以及标量地震矩Ｍ０进行搜索，理论上而言，仅需要２个台站的三分向波形就可以求解；虽然研究显示，对于大多数３个三分向台的组合，使用波形反演就可以得到相对准确的震源机制，但实际情况也显示，不同的台站组合波形反演得到的解之间仍然存在一定的差异（Ｇｏｄａｎｏ ｅｔ ａｌ．，２００９；郑建常和陈运泰，２０１２）．目前国内台网密度已经达到相当水平，在东部地区，一个中等地震通常有数十甚至上百个宽频带台能够记录到清晰的波形，以此次莱州地震为例，通过对原始波形进行去均值、去趋势、积分等简单变换后，根据直观的观察，震中距３００ｋｍ范围内，采样率１００Ｈｚ的宽频带三分向波形有近４０个台站的资料可用．在使用ＣＡＰ方法求解震源机制时，一般的做法是选择部分波形拟合较好的台进行反演，删去拟合不好的台或者震相；有些情况下，甚至仅使用面波部分而删除体波震相，需知面波尤其是径向和切向分量，很容易受到台站下方浅层地壳结构的影响．由于ＣＡＰ方法是采用网格搜索的方法，因而这种人为的选择，必然会为反演结果增加主观的不确定性因素．我们无法令人信服地说明，拟合不好的波形究竟是数据本身确实存在干扰，还是说搜索到的解无法满足该条数据．另外，ＣＡＰ方法虽然可以在最后的输出结果中给出断层面参数的不确定性，但该估计值只是面向所使用的台站数据的结果，在上述的人为选择下，该不确定性估计能够在多大程度上客观地反映最终解的整体不确定性，是无法说明的．为了求得更加稳定可靠的解并且合理客观地给出解的误差估计，我们在相对丰富的观测数据基础５５４地球物理学报（Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．）５８卷 上，采用自助抽样统计方法进行分析．具体方法是在可用的观测台站中可重复地随机抽取一定数量的台站组成新的台站组合，使用该台站组合的观测数据重复反演过程．在大量的重复计算后（例如，超过１０００次），可以有效地排除观测质量不高或存在较大干扰误差的数据的影响，从而得到更加接近真实解的结果，并且可以有效地给出解的不确定性．另外，由于ＣＡＰ方法在搜索断层面解时采用的是网格搜索的方法，然后通过插值计算误差函数ｅ的最小值，并且由于固定步长的走向、倾角、滑动角的尝试位置在震源球上的分布是不均匀的（许向彤等，１９９５），因此在最终解中可能会有空缺（ｇａｐ）的存在．为了求解优化解，我们进一步使用Ｋａｇａｎ（１９９１）定义的双力偶模型最小空间旋转角，对上面自助抽样得到的大量满足条件的震源机制结果进行分析，定义与所有解的空间偏转角度和为目标函数，使用粒子群非线性优化方法搜索该目标函数最小的结果，视为最优解．２．３ 聚类分析在震源机制求解中，常见的情况是在震源球上存在几簇相对集中分布的解，对这些可能的解直接取数学平均是不甚合理的，并且在数据存在较大误差或干扰的情况下，满足条件的可能解的分布范围也许会相当大．因此针对这一现象，刁桂苓等（１９９２）、俞春泉等（２００９）分别使用系统聚类和动态聚类技术，对所有的可能解进行聚类分析，求取聚类中心作为反演的优化解，数值试验和实际应用都有很好的效果．聚类分析可以很好地排除孤立解和错误解，从而在大量的数据中获取更加接近真实解的结果．本文在使用不同台站组合重复进行波形反演后，同样得到了大量的震源机制解数据，受台站布局和数据误差的影响，这些解或多或少存在差别，因此对这些结果进行聚类分析是很有必要的．３ 数据与资料本文使用了山东台网提供的波形资料，其中还包括了邻省如辽宁、河北、江苏等省交换资料的部分台站．图１给出了本项研究使用的台站分布，其中个别台如ＪＩＭ、ＺＳＬ、ＨＵＤ等为短周期台，在波形反演中没有使用．本文研究中，首先由观测记录直接读取初动符号，用于约束波形反演；然后将观测数据扣除仪器响应，经过去均值、零漂等预处理后积分至位移记录，旋转到Ｚ－Ｒ－Ｔ坐标系，对观测波形和理论波图１ 本文研究的两次地震震中及山东台网台站分布图Ｆｉｇ．１ Ｍａｐ ｏｆ ｓｔａｔｉｏｎｓ ｉｎ Ｓｈａｎｄｏｎｇ Ｎｅｔｗｏｒｋａｎｄ ｔｗｏ ｅａｒｔｈｑｕａｋｅｓ ｓｔｕｄｉｅｄ ｉｎ ｔｈｉｓ ｐａｐｅｒ．Ｒｅｄ ｃｉｒｃｌｅｓ ｄｅｎｏｔｅ ｅｐｉｃｅｎｔｅｒ，ｔｒｉａｎｇｌｅｓ ａｒｅｓｔａｔｉｏｎｓ，ａｎｄ ｓｏｌｉｄ ｂｌａｃｋ ｌｉｎｅｓ ａｒｅ ｆａｕｌｔｓ．形同样进行带通滤波，然后用于反演．使用Ｃｈａｎｇ等（２００６）给出的朝鲜半岛南部至黄海地区的中上地壳速度结构模型用于本文的震源机制反演．相关地质资料显示，胶东半岛、南黄海以及朝鲜半岛南部在大地构造分区上都属于下扬子地块，地质构造属性相对较为一致（Ｒｅｅ ｅｔ ａｌ．，１９９６）．４ 结果与分析４．１ 乳山Ｍ４．３震源机制选择震中距在２５０ｋｍ以内的１５个台站的宽频带波形记录进行反演，Ｐｎｌ和面波的反演波段分别选择０．０５～０．１５Ｈｚ和０．０３３～０．０６７Ｈｚ频段．图２给出了不同深度的最佳双力偶解，及拟合误差随不同深度变化的关系，由图可见，震源深度在４ｋｍ时观测波形和理论波形的错配值最小，说明事件的震源深度较浅．由ＣＡＰ方法反演得到的最佳震源机制：节面Ａ的参数为：走向２０２°、倾角７５°、滑动角１５３°；节面Ｂ的参数为：走向２９９．５°、倾角６４°、滑动角１６．７°；参考乳山序列的双差定位结果（李冬梅和郑建常，２０１４）分析认为，节面Ｂ可能是乳山地震的发震断层；震源机制显示为左旋走滑型，反演得到此次地震的矩震级ＭＷ＝４．２．６５４ ２期郑建常等：ＣＡＰ方法反演震源机制的误差分析：以胶东半岛两次显著中等地震为例图２ ２０１４年１月７日乳山Ｍ４．３级地震不同震源深度的波形拟合误差及最佳震源机制解Ｆｉｇ．２ Ｗａｖｅｆｏｒｍ ｆｉｔ ｅｒｒｏｒｓ ａｎｄ ｂｅｓｔ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ａｓｆｕｎｃｔｉｏｎ ｏｆ ｄｅｐｔｈ ｆｏｒ Ｊａｎ．７，２０１４Ｒｕｓｈａｎ Ｍ４．３ｅｖｅｎｔ 图３给出了对应最优解的理论波形和观测波形的拟合情况．１５个台一共７５个震相，其中理论波形与观测波形相关系数大于０．９的有３９个，超过５０％；相关系数大于０．６（相关性较好）的有６７个，约占８９．３％；最佳解的方差减少（ｖａｒｉａｎｃｅ ｒｅｄｕｃｔｉｏｎ）为７０．３％，说明理论波形很好地拟合了观测波形，反演结果是可靠的．个别台（如ＷＥＨ）平均相关系数较差，可能与台站位于震源机制解的节面线附近，振幅相对较小所致；另外如ＣＨＤ台的拟合程度不好，可能与该台处于海域、噪声干扰较大有关．使用自助抽样的统计方法，对乳山地震震源机制解的不确定性进行估计．选用震中距在３００ｋｍ以内的采样率为１００Ｈｚ的２２个三分向宽频带台的观测波形组成原始数据集，为了保证用于反演的数据的样本量，设用于反演的台站数为２０个，对原始数据集进行每个台站等概率、可重复地随机抽取，抽取出的台站波形组成新的数据集，然后用于ＣＡＰ方法的波形反演．对上述的抽取台站反演过程重复１０００次，将反演得到的震源机制的断层节面解和Ｐ、Ｔ轴绘制在一个震源球上，见图４．可以看出，反演中除去个别反演过程的断层面解出现一定程度的偏离外，其余结果集中分布，均显示为近走滑的机制；图４中的Ｐ、Ｔ轴位置和断层节面线集中成丛，大致显示出断层面解的误差范围．使用粒子群非线性优化方法，以与自助抽样给出的１０００个机制解（图４）的Ｋａｇａｎ角之和为目标函数，搜索最优解．结果显示最优解为，节面Ａ：走向２０８．４°，倾角８９．７°，滑动角１５４．３°；节面Ｂ：走向２９８．５°，倾角６４．３°，滑动角０．３°；最优解与图４所示表１ 波形反演乳山Ｍ４．３地震震源机制解结果Ｔａｂｌｅ １ Ｐａｒａｍｅｔｅｒｓ ｏｆ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍ ｒｅｓｕｌｔｓ ｆｒｏｍ ｗａｖｅｆｏｒｍｉｎｖｅｒｓｉｏｎ ｆｏｒ Ｊａｎ．７，２０１４Ｒｕｓｈａｎ Ｍ４．３ｅｖｅｎｔ方法节面Ａ节面Ｂ走向（°）倾角（°）滑动角（°）走向（°）倾角（°）滑动角（°）ＣＡＰ ２０２ ７５ １５３ ２９９．５ ６４．０ １６．７粒子群优化２０８．４ ８９．７ １５４．３ ２９８．５ ６４．３ ０．３聚类分析２０８．０ ８９．３ １５４．１ ２９８．４ ６４．１ １．８１０００个解的平均夹角４．３７°，以其与所有解Ｋａｇａｎ角的２倍标准差为震源机制解的误差范围，结果显示不确定性为６．４４°（图５）．对自助抽样结果进行动态聚类分析，结果显示最优解为：节面Ａ：走向２０８．０°，倾角８９．３°，滑动角１５４．１°，节面Ｂ：走向２９８．４°，倾角６４．１°，滑动角１．８°，与粒子群优化解非常一致（见表１）．将自助抽样结果中的Ｐ、Ｔ轴投影到震源球上（图４），对其进行概率密度统计分析，结果见图６．４．２ 莱州Ｍ４．６震源机制使用ＣＡＰ方法对２０１３年１１月２３日莱州Ｍ４．６地震进行反演（郑建常等，２０１５），同样进行ＣＡＰ反演情况的自助抽样统计分析．选用震中距在２７０ｋｍ以内的采样率为１００Ｈｚ的２２个三分向宽频带台的观测波形组成原始数据集，采用全样本随机抽取方法，自助抽样反演１０００次，图７给出了反演得到的震源机制的断层节面解和Ｐ、Ｔ轴在震源球上的分布情况．结果显示，莱州地震的自助抽样结果同样很好地显示出了反演得到震源机制解的误差范围，相对于乳山地震，出现了极个别反演过程的结果偏离较大的情况．使用粒子群非线性优化方法，搜索与自助抽样结果的旋转角最小的解．结果显示最优解为，节面Ａ：走向２３６．９°，倾角７６．２°，滑动角－１６９．３°；节面Ｂ：走向１４４．３°，倾角７９．６°，滑动角－１４．０°；与所有自助抽样解的平均偏转角１７．４°（图８），以其与所有解Ｋａｇａｎ角的２倍标准差（图８红色虚线所示）为误差范围，结果显示震源机制解的不确定性为２３．７°．从自助抽样得到的所有机制解在震源球上的分布情况（图７）可以直观地看出，断层节面线尤其是北西向节面呈现出两组集中．由于我们定义的粒子群优化的目标函数是搜索与所有自助抽样解的空间旋转最小，因此从图７可以看出，最优解的节面位置处于其中一组的边缘位置，在此情况下，对自助抽样结果进行聚类分析是有意义的．７５４地球物理学报（Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．）５８卷图３ ２０１４年１月７日乳山Ｍ４．３地震最优解的理论波形（红）与观测波形（黑）波形图下方第一行数字为各段理论地震波形相对实际观测波形的移动时间，正值表示理论波形相对观测波形超前．第二行数字为理论波形与观测波形的相关系数（百分比）．波形图左侧字母为台站，其下数字分别为台站震中距（ｋｍ）和方位角（°）．图左侧的震源球上红色区域代表压缩区，白色代表拉张区，震源球采用下半球投影．震源球上标注的“＋”和“－”表示反演使用台站的Ｐ波初动．Ｆｉｇ．３ Ｃｏｍｐａｒｉｓｏｎ ｂｅｔｗｅｅｎ ｓｙｎｔｈｅｔｉｃｓ（ｒｅｄ）ａｎｄ ｏｂｓｅｒｖｅｄ（ｂｌａｃｋ）ｓｅｉｓｍｏｇｒａｍｓ ｏｆ Ｊａｎ．７，２０１４Ｒｕｓｈａｎ Ｍ４．３ｅｖｅｎｔＴｈｅ ｎｕｍｂｅｒｓ ｏｎ ｔｈｅ ｌｏｗｅｒ ｌｅｆｔ ｓｉｄｅ ｏｆ ｅａｃｈ ｓｅｉｓｍｏｇｒａｍ ａｒｅ ｔｈｅ ｔｉｍｅ ｓｈｉｆｔｓ（ｕｐｐｅｒ）ａｎｄ ｃｒｏｓｓ－ｃｏｒｒｅｌａｔｉｏｎ ｃｏｅｆｆｉｃｉｅｎｔ ｉｎ ｐｅｒｃｅｎｔ（ｌｏｗｅｒ）．Ｐｏｓｉｔｉｖｅ ｔｉｍｅ ｓｈｉｆｔｓ ｍｅａｎ ｔｈａｔ ｔｈｅ ｏｂｓｅｒｖｅｄ ｄａｔａ ｈａｖｅ ｂｅｅｎ ｄｅｌａｙｅｄ．Ｔｈｅ ｌｅｔｔｅｒｓ ｏｎ ｔｈｅ ｌｅｆｔ ｓｉｄｅ ａｒｅ ｓｔａｔｉｏｎｓ，ｔｈｅ ｎｕｍｂｅｒｓ ｂｅｌｏｗ ｉｔ ａｒｅｅｐｉｃｅｎｔｒａｌ ｄｉｓｔａｎｃｅ（ｉｎ ｋｍ）ａｎｄ ａｚｉｍｕｔｈ（ｄｅｇｒｅｅ）．Ｔｈｅ ｒｅｄ ｃｏｌｏｒ ｉｎ ｂｅａｃｈ－ｂａｌｌ ｄｅｎｏｔｅｓ ｃｏｍｐｒｅｓｓｉｏｎ ａｒｅａ，ｗｈｉｌｅ ｗｈｉｔｅ ｉｓ ｅｘｔｅｎｓｉｏｎ．Ｔｈｅ‘＋’ａｎｄ‘－’ｓｉｇｎｓ ｏｎ ｂｅａｃｈ－ｂａｌｌ ｉｎｄｉｃａｔｅ ｐｏｌａｒｉｔｉｅｓ ｏｎ ｉｎｖｅｒｓｉｏｎ ｕｓｅｄ ｓｔａｔｉｏｎｓ．Ｌｏｗｅｒ ｈｅｍｉｓｐｈｅｒｅ ｐｒｏｊｅｃｔｉｏｎ ｉｓ ｕｓｅｄ．８５４ ２期郑建常等：ＣＡＰ方法反演震源机制的误差分析：以胶东半岛两次显著中等地震为例图４ 自助抽样得到的１０００次乳山地震震源机制解及粒子群最优解（下半球极射投影）震源球上黑色细线条表示自助抽样结果的断层节面线，红色线条表示粒子群优化解的节面线．Ｆｉｇ．４ １０００ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ｏｆ Ｒｕｓｈａｎ Ｍ４．３ｅｖｅｎｔｒｅｔｒｉｅｖｅｄ ｂｙ ａ ｂｏｏｔｓｔｒａｐｐｉｎｇ ｐｒｏｃｅｓｓ，ａｌｌ ｎｏｄａｌ ｌｉｎｅｓ（ｂｌａｃｋ）ａｎｄ Ｐ，Ｔａｘｅｓ（ｂｌｕｅ ａｎｄ ｒｅｄ ｐｏｉｎｔｓ，ｒｅｓｐｅｃｔｉｖｅｌｙ）ａｒｅｐｌｏｔｔｅｄ ｏｎ ｏｎｅ ｂｅａｃｈ－ｂａｌｌ．Ｔｈｅ ｒｅｄ ｌｉｎｅｓ ｏｎ ｔｈｅ ｂｅａｃｈ－ｂａｌｌｓｈｏｗ ｔｈｅ ｏｐｔｉｍｉｚｅｄ ｓｏｌｕｔｉｏｎ ｇｉｖｅｎ ｂｙ ａ Ｐａｒｔｉｃｌｅ ＳｗａｒｍＯｐｔｉｍｉｚａｔｉｏｎ ｍｅｔｈｏｄ．Ｉｔｓ ｃｏｒｒｅｓｐｏｎｄｉｎｇ Ｐ，Ｔａｘｅｓ ａｒｅ ａｌｓｏｄｉｓｐｌａｙｅｄ ｏｎ ｔｈｅ ｂｅａｃｈ－ｂａｌｌ（ｙｅｌｌｏｗ ａｎｄ ｇｒｅｅｎ ｐｏｉｎｔ，ｒｅｓｐｅｃｔｉｖｅｌｙ）．Ｌｏｗｅｒ ｈｅｍｉｓｐｈｅｒｅ ｐｒｏｊｅｃｔｉｏｎ ｉｓ ｕｓｅｄ．图５ 乳山地震自助抽样结果与粒子群优化解的Ｋａｇａｎ角分布Ｆｉｇ．５ Ｋａｇａｎ ａｎｇｌｅｓ ｏｆ ｂｏｏｔｓｔｒａｐ ｒｅｓｕｌｔｓ ｔｏｔｈｅ ＰＳＯ ｓｏｌｕｔｉｏｎ ｆｏｒ Ｒｕｓｈａｎ ｅｖｅｎｔ 以两个震源机制解之间的最小空间旋转角（即Ｋａｇａｎ角）为距离的定义，对１０００次自助抽样结果进行聚类分析，图９给出了聚类谱系图（由于完整的聚类树过于密集和庞大，因此我们只显示了Ｋａｇａｎ＞７°的部分），以５０°为阈值，可以将结果分为５类．图１０给出了聚类分析的结果，属于Ⅰ类的数据占３２．１％，Ⅱ类６６．６％，其余三类数据合计仅有１．３％．由图６ 乳山地震自助抽样结果的震源球概率密度分布俯视图（未进行极射投影）色标中正值表示Ｔ轴的概率密度分布，负值表示Ｐ轴的概率密度分布．Ｆｉｇ．６ Ｐｒｏｂａｂｉｌｉｔｙ ｄｅｎｓｉｔｙ ｄｉｓｔｒｉｂｕｔｉｏｎ ｏｆ ｓｏｌｕｔｉｏｎｓ ｏｎｂｅａｃｈ－ｂａｌｌ（ｔｏｐ ｖｉｅｗ ｏｆ ｌｏｗｅｒ ｈｅｍｉｓｐｈｅｒｅ，ｗｉｔｈｏｕｔ ｐｒｏｊｅｃｔｉｏｎ）Ｐｏｓｉｔｉｖｅ ｖａｌｕｅｓ ｏｎ ｔｈｅ ｃｏｌｏｒ ｓｃａｌｅ（ｃｏｒｒｅｓｐｏｎｄｉｎｇ ｔｏ ｔｈｅ ｒｅｄ ａｒｅａｏｎ ｂｅａｃｈ－ｂａｌｌ）ｉｎｄｉｃａｔｅ ｐｒｏｂａｂｉｌｉｔｙ ｏｆ Ｔ ａｘｉｓ，ｗｈｉｌｅ ｎｅｇａｔｉｖｅｖａｌｕｅｓ（ｃｏｒｒｅｓｐｏｎｄｉｎｇ ｔｏ ｂｌｕｅ ａｒｅａ）ｍｅａｎ ｐｒｏｂａｂｉｌｉｔｙ ｏｆ Ｐａｘｉｓ．图７ 自助抽样得到的莱州地震震源机制解及粒子群最优解（下半球极射投影）黑色节面线为自助抽样得到的震源机制解；红色节面线为粒子群最优解．Ｆｉｇ．７ １０００Ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ｏｆ Ｌａｉｚｈｏｕ Ｍ４．６ｅｖｅｎｔｒｅｔｒｉｅｖｅｄ ｂｙ ａ ｂｏｏｔｓｔｒａｐｐｉｎｇ ｐｒｏｃｅｓｓ，ａｌｌ ｎｏｄａｌ ｌｉｎｅｓ（ｂｌａｃｋ）ａｎｄ Ｐ，Ｔａｘｅｓ（ｂｌｕｅ ａｎｄ ｒｅｄ ｐｏｉｎｔｓ，ｒｅｓｐｅｃｔｉｖｅｌｙ）ａｒｅ ｐｌｏｔｔｅｄ ｏｎ ｏｎｅ ｂｅａｃｈ－ｂａｌｌ．Ｔｈｅ ｒｅｄ ｌｉｎｅｓ ｏｎ ｔｈｅ ｂｅａｃｈ－ｂａｌｌ ｓｈｏｗ ｔｈｅ ｏｐｔｉｍｉｚｅｄ ｓｏｌｕｔｉｏｎ ｇｉｖｅｎ ｂｙ ａ ＰＳＯ ｍｅｔｈｏｄ．Ｉｔｓ ｃｏｒｒｅｓｐｏｎｄｉｎｇ Ｐ，Ｔａｘｅｓ ａｒｅ ａｌｓｏ ｄｉｓｐｌａｙｅｄ ｏｎ ｔｈｅｂｅａｃｈ－ｂａｌｌ（ｃｙａｎ ａｎｄ ｙｅｌｌｏｗ ｐｏｉｎｔ，ｒｅｓｐｅｃｔｉｖｅｌｙ）．Ｌｏｗｅｒｈｅｍｉｓｐｈｅｒｅ ｐｒｏｊｅｃｔｉｏｎ ｉｓ ｕｓｅｄ．图１０可以看出，其余三类的断层节面线和Ｐ、Ｔ轴位置明显偏离集中区且机制解类型与绝大部分结果（走滑型）不一致，是典型的孤立解．孤立解（或错误解）的出现，可能说明我们使用的数据中个别台站（或分向）存在较大干扰．使用俞春泉等（２００９）的方法求取了四类解的聚类中心，其中Ｉ类解的聚类中９５４地球物理学报（Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．）５８卷图８ 莱州地震自助抽样结果与粒子群优化解的Ｋａｇａｎ角分布Ｆｉｇ．８ Ｋａｇａｎ ａｎｇｌｅｓ ｏｆ ｂｏｏｔｓｔｒａｐ ｒｅｓｕｌｔｓｔｏ ｔｈｅ ＰＳＯ ｓｏｌｕｔｉｏｎ ｆｏｒ Ｌａｉｚｈｏｕ ｅｖｅｎｔ图９ 莱州地震自助抽样结果的聚类谱系图Ｆｉｇ．９ Ｄｅｎｄｒｏｇｒａｍ ｐｌｏｔ ｏｆ ｔｈｅ ｈｉｅｒａｒｃｈｉｃａｌｂｉｎａｒｙ ｃｌｕｓｔｅｒ ｔｒｅｅ ｆｏｒ Ｌａｉｚｈｏｕ ｅｖｅｎｔ图１０ 莱州地震自助抽样结果的聚类分析（ａ）断层节面线的分类显示；（ｂ）机制解Ｐ（＋）、Ｔ（⊙）轴位置的分类显示．断层节面线和Ｐ、Ｔ轴颜色表示分类，与图９分类颜色一致．Ｆｉｇ．１０ Ｃｌｕｓｔｅｒｉｎｇ ｒｅｓｕｌｔｓ ｏｆ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ｆｒｏｍ ａ ｂｏｏｔｓｔｒａｐ ｐｒｏｃｅｓｓ ｆｏｒ Ｌａｉｚｈｏｕ ｅｖｅｎｔＴｈｅ ｄｉｆｆｅｒｅｎｔ ｃｏｌｏｒｓ ｏｆ ｎｏｄａｌ ｌｉｎｅｓ ａｎｄ Ｐ（“＋”ｓｉｇｎ ｉｎ ｒｉｇｈｔ ｐａｎｅｌ），Ｔ（“⊙”ｉｎ ｒｉｇｈｔ ｐａｎｅｌ）ａｘｅｓｄｅｎｏｔｅ ｄｉｆｆｅｒｅｎｔ ｃｌａｓｓｅｓ，ｗｈｉｃｈ ａｒｅ ｃｏｒｒｅｓｐｏｎｄｉｎｇ ｔｏ ｃｏｌｏｒｓ ｓｈｏｗｎ ｉｎ Ｆｉｇ．９心的断层面参数（设北东向节面为发震断层面）为：走向２３１．６°，倾角８８．７°，滑动角－１６８．２°；Ⅱ类解的聚类中心为：走向２３８．５°，倾角７４．１°，滑动角－１６４．８°．将自助抽样结果中的Ｐ、Ｔ轴投影到震源球上，对其进行概率密度统计分析，结果见图１１．可见，同聚类分析的结果一致，Ｐ轴位置的概率密度在震源球上出现了两个极值区，分别对应Ⅰ类和Ⅱ类两个聚类中心．５ 讨论与结论基于山东省宽频带数字地震波形资料，本文首先使用ＣＡＰ方法反演了近期胶东半岛地区发生的两次显著中等地震活动的震源机制，讨论了如何合理地估计ＣＡＰ方法反演震源机制的误差范围以及如何确定优化解的问题．我们首先使用自助抽样方法，对原始数据进行等概率随机抽样，多次重复波形反演过程，排除了人为选择数据的干扰，得到大样本量的震源机制解数据；在此基础上，我们（１）使用了粒子群优化算法从中搜索与这些机制解空间偏转角最小的解当作优化解，以Ｋａｇａｎ角的二倍标准差作为反演结果的不确定性范围，结果显示：乳山地震的粒子群优化解为：走向２９８．５°，倾角６４．３°，滑动角０．３°，不确定性为±６．４°；莱州地震的优化解为：走向２３６．９°、倾角７６．２°、滑动角－１６９．３°，不确定性为±２３．７°．（２）对自助抽样结果进行聚类分析，其中：乳山地震结果的聚类中心与粒子群优化解基本一致；莱州地震结果存在多个聚类，排除孤立解后，有两个聚类中心，其对应两类数据合计占结果的９８．７％，说０６４ ２期郑建常等：ＣＡＰ方法反演震源机制的误差分析：以胶东半岛两次显著中等地震为例图１１ 莱州地震自助抽样结果的震源球概率密度分布（ａ）下半球俯视图（未进行极射投影）；（ｂ）东南方向，４５°三维侧视图．色标说明同图６．Ｆｉｇ．１１ Ｐｒｏｂａｂｉｌｉｔｙ ｄｅｎｓｉｔｙ ｄｉｓｔｒｉｂｕｔｉｏｎ ｏｆ ｓｏｌｕｔｉｏｎｓ ｏｎ ｂｅａｃｈ－ｂａｌｌ（ａ）Ｔｏｐ ｖｉｅｗ ｏｆ ｌｏｗｅｒ ｈｅｍｉｓｐｈｅｒｅ，ｗｉｔｈｏｕｔ ａｎｙ ｐｒｏｊｅｃｔｉｏｎ；（ｂ）Ｓｉｄｅ ｖｉｅｗ ｏｆ ｗｈｏｌｅ ｂｅａｃｈ－ｂａｌｌ，ｆｒｏｍ ｓｏｕｔｈ－ｅａｓｔ ４５°ｐｏｓｉｔｉｏｎ．Ｐｏｓｉｔｉｖｅ ｖａｌｕｅｓｏｎ ｔｈｅ ｃｏｌｏｒ ｓｃａｌｅ（ｃｏｒｒｅｓｐｏｎｄｉｎｇ ｔｏ ｔｈｅ ｒｅｄ ａｒｅａ ｏｎ ｂｅａｃｈ－ｂａｌｌ）ｉｎｄｉｃａｔｅ ｐｒｏｂａｂｉｌｉｔｙ ｏｆ Ｔａｘｉｓ，ｗｈｉｌｅ ｎｅｇａｔｉｖｅ ｖａｌｕｅｓ（ｃｏｒｒｅｓｐｏｎｄｉｎｇｔｏ ｂｌｕｅ ａｒｅａ）ｍｅａｎ ｐｒｏｂａｂｉｌｉｔｙ ｏｆ Ｐａｘｉｓ．明此次地震的真实解在这两类数据范围内．（３）将自助抽样结果中的Ｐ、Ｔ轴投影到震源球，对其进行概率密度统计，给出了机制解在震源球上的概率密度分布图．本文方法不单可以得到更准确的震源机制优化解、给出科学合理的误差估计，而且可以有效地排除孤立解和错误解，克服存在较大干扰台站的数据的影响，因此在震后应急的震源机制自动化求解中可以发挥作用．在强大计算能力的支持下，无须人工干预即可得到准确可靠的震源机制结果，从而为震害评估、趋势分析等提供重要的科学依据．自助抽样结果显示乳山地震震源机制解的不确定性要小于莱州地震，笔者推测可能有莱州地震使用的台站中个别台的干扰较大的原因，另外也无法排除莱州地震的震源破裂过程可能更加复杂的可能．ＣＡＰ方法中用于计算格林函数的Ｆ－Ｋ方法使用狄拉克－Ｄｅｌｔａ函数作为震源时间函数（Ｚｈｕ ａｎｄＲｉｖｅｒａ，２００２），对于小震级的事件该简化方案更为适用，莱州地震（Ｍ４．６）与乳山地震（Ｍ４．３）震级相差不大，但莱州地震的自助抽样结果出现了两个概率较高的聚类中心，在使用大部分相同台站的情况下，这可能意味着莱州地震的震源破裂随时间的变化可能与狄拉克－Ｄｅｌｔａ函数存在一定的偏离．Ｒｏｄｒíｇｕｅｚ－Ｌｏｚｏｙａ等（２００８）的研究显示，区域中等地震也可能有复杂的震源破裂过程，在该问题上的深入研究需要更进一步的工作．致谢 感谢两位匿名审稿专家提出的宝贵意见．聚类分析中使用了俞春泉等（２００９）提供的部分开放代码，粒子群搜索使用了Ｓ Ｃｈｅｎ给出的Ｍａｔｌａｂ软件包，在此一并表示感谢！ＲｅｆｅｒｅｎｃｅｓＣｈａｎｇ Ｓ Ｊ，Ｂａａｇ Ｃ Ｅ．２００６．Ｃｒｕｓｔａｌ ｓｔｒｕｃｔｕｒｅ ｉｎ Ｓｏｕｔｈｅｒｎ Ｋｏｒｅａｆｒｏｍ ｊｏｉｎｔ ａｎａｌｙｓｉｓ ｏｆ ｒｅｇｉｏｎａｌ ｂｒｏａｄｂａｎｄ ｗａｖｅｆｏｒｍｓ ａｎｄ ｔｒａｖｅｌｔｉｍｅｓ．Ｂｕｌｌ．Ｓｅｉｓｍ．Ｓｏｃ．Ａｍｅｒ．，９６（３）：８５６－８７０．Ｄｉａｏ Ｇ Ｌ，Ｙｕ Ｌ Ｍ，Ｌｉ Ｑ Ｚ．１９９２．Ｈｉｅｒａｒｃｈｉｃａｌ ｃｌｕｓｔｅｒｉｎｇ ａｎａｌｙｓｉｓｏｆ ｔｈｅ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍ ｓｏｌｕｔｉｏｎ－Ｔａｋｉｎｇ ｔｈｅ Ｈａｉｃｈｅｎｇ ＥａｒｔｈｑｕａｋｅＳｅｑｕｅｎｃｅｓ ｆｏｒ ｅｘａｍｐｌｅ．Ｅａｒｔｈｑｕａｋｅ Ｒｅｓｅａｒｃｈ ｉｎ Ｃｈｉｎａ（ｉｎＣｈｉｎｅｓｅ），８（３）：８６－９２．Ｇｏｄａｎｏ Ｍ，Ｒｅｇｎｉｅｒ Ｍ，Ｄｅｓｃｈａｍｐｓ Ａ，ｅｔ ａｌ．２００９．Ｆｏｃａｌｍｅｃｈａｎｉｓｍｓ ｆｒｏｍ ｓｐａｒｓｅ ｏｂｓｅｒｖａｔｉｏｎｓ ｂｙ ｎｏｎｌｉｎｅａｒ ｉｎｖｅｒｓｉｏｎ ｏｆａｍｐｌｉｔｕｄｅｓ：ｍｅｔｈｏｄ ａｎｄ ｔｅｓｔｓ ｏｎ ｓｙｎｔｈｅｔｉｃ ａｎｄ ｒｅａｌ ｄａｔａ．Ｂｕｌｌ．Ｓｅｉｓｍ．Ｓｏｃ．Ａｍｅｒ．，９９（４）：２２４３－２２６４．Ｈａｎ Ｌ Ｂ，Ｊｉａｎｇ Ｃ Ｓ，Ｂａｏ Ｆ．２０１０．Ｓｏｕｒｃｅ ｐａｒａｍｅｔｅｒ ｄｅｔｅｒｍｉｎａｔｉｏｎｏｆ ２０１０Ｔａｉｋａｎｇ ＭＳ４．６ｅａｒｔｈｑｕａｋｅ ｓｅｑｕｅｎｃｅｓ．Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．（ｉｎ Ｃｈｉｎｅｓｅ），５５（９）：２９７３－２９８１，ｄｏｉ：１０．６０３８／ｊ．ｉｓｓｎ．０００１－５７３３．２０１２．０９．０１６．Ｈａｎ Ｌ Ｂ，Ｊｉａｎｇ Ｃ Ｓ．２０１２．Ｆｏｃａｌ ｍｅｃｈａｎｉｓｍ ｉｎｖｅｒｓｉｏｎ ｏｆ ８Ｊｕｎ ２０１１Ｔｏｋｓｕｎ ＭＳ５．３ｅａｒｔｈｑｕａｋｅ．Ａｃｔａ Ｓｅｉｓｍｏｌｏｇｉｃａ Ｓｉｎｉｃａ（ｉｎ Ｃｈｉｎｅｓｅ），３４（３）：４１５－４２２．Ｈｕａｎｇ Ｊ Ｐ，Ｎｉ Ｓ Ｄ，Ｆｕ Ｒ Ｓ，ｅｔ ａｌ．２００９．Ｓｏｕｒｃｅ ｍｅｃｈａｎｉｓｍ ｏｆ ｔｈｅ２００６ ＭＷ５．１Ｗｅｎ′ａｎ ｅａｒｔｈｑｕａｋｅ ｄｅｔｅｒｍｉｎｅｄ ｆｒｏｍ ａ ｊｏｉｎｔ ｉｎｖｅｒｓｉｏｎ ｏｆ１６４地球物理学报（Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．）５８卷 ｌｏｃａｌ ａｎｄ ｔｅｌｅｓｅｉｓｍｉｃ ｂｒｏａｄｂａｎｄ ｗａｖｅｆｏｒｍ ｄａｔａ．Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．（ｉｎ Ｃｈｉｎｅｓｅ），５２（１）：１２０－１３０．Ｊｏｓｔ Ｍ Ｌ，Ｈｅｒｍａｎｎ Ｒ Ｂ．１９８９．Ａ ｓｔｕｄｅｎｔ′ｓ ｇｕｉｄｅ ｔｏ ａｎｄ ｒｅｖｉｅｗ ｏｆｍｏｍｅｎｔ ｔｅｎｓｏｒｓ．Ｓｅｉｓｍ．Ｒｅｓ．Ｌｅｔｔ．，６０（２）：３７－５７．Ｋａｇａｎ Ｙ Ｙ．１９９１．３－Ｄ ｒｏｔａｔｉｏｎ ｏｆ ｄｏｕｂｌｅ－ｃｏｕｐｌｅ ｅａｒｔｈｑｕａｋｅｓｏｕｒｃｅｓ．Ｇｅｏｐｈｙｓ．Ｊ．Ｉｎｔ．，１０６（３）：７０９－７１６．Ｌｉ Ｄ Ｍ，Ｚｈｅｎｇ Ｊ Ｃ．２０１４．Ｒｅｌｏｃａｔｉｏｎ ａｎａｌｙｓｉｓ ｏｆ Ｏｃｔ．１，２０１３Ｒｕｓｈａｎ ｓｗａｒｍｓ．Ｑｉｌｕ Ｅａｒｔｈｑｕａｋｅ Ｓｃｉｅｎｃｅｓ（ｉｎ Ｃｈｉｎｅｓｅ），１０（１）：２６－２８．Ｌｏｎｇ Ｆ，Ｚｈａｎｇ Ｙ Ｊ，Ｗｅｎ Ｘ Ｚ，ｅｔ ａｌ．２０１０．Ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓｏｌｕｔｉｏｎｓ ｏｆ ＭＬ≥４．０ｅｖｅｎｔｓ ｉｎ ｔｈｅ ＭＳ６．１Ｐａｎｚｈｉｈｕａ－Ｈｕｉｌｉｅａｒｔｈｑｕａｋｅ ｓｅｑｕｅｎｃｅ ｏｆ Ａｕｇ ３０，２００８．Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．（ｉｎ Ｃｈｉｎｅｓｅ），５３（１２）：２８５２－２８６０，ｄｏｉ：１０．３９６９／ｊ．ｉｓｓｎ．０００１－５７３３．２０１０．１２．００８．ＬüＪ，Ｚｈｅｎｇ Ｙ，Ｎｉ Ｓ Ｄ，ｅｔ ａｌ．２００８．Ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ａｎｄｓｅｉｓｍｏｇｅｎｉｃ ｓｔｒｕｃｔｕｒｅｓ ｏｆ ｔｈｅ ＭＳ５．７ａｎｄ ＭＳ４．８Ｊｉｕｊｉａｎｇ－Ｒｕｉｃｈａｎｇ ｅａｒｔｈｑｕａｋｅｓ ｏｆ Ｎｏｖ．２６，２００５．Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．（ｉｎ Ｃｈｉｎｅｓｅ），５１（１）：１５８－１６４，ｄｏｉ：１０．３３２１／ｊ．ｉｓｓｎ：０００１－５７３３．２００８．０１．０２０．Ｒｅｅ Ｊ－Ｈ，Ｃｈｏ Ｍ，Ｋｗｏｎ Ｓ－Ｔ，ｅｔ ａｌ．１９９６．Ｐｏｓｓｉｂｌｅ ｅａｓｔｗａｒｄｅｘｔｅｎｓｉｏｎ ｏｆ Ｃｈｉｎｅｓｅ ｃｏｌｌｉｓｉｏｎ ｂｅｌｔ ｉｎ Ｓｏｕｔｈ Ｋｏｒｅａ：ＴｈｅＩｍｊｉｎｇａｎｇ ｂｅｌｔ．Ｇｅｏｌｏｇｙ，２４（１２）：１０７１－１０７４．Ｒｏｄｒíｇｕｅｚ－Ｌｏｚｏｙａ Ｈ Ｅ，Ｑｕｉｎｔａｎａｒ Ｌ，Ｏｒｔｅｇａ Ｒ，ｅｔ ａｌ．２００８．Ｒｕｐｔｕｒｅ ｐｒｏｃｅｓｓ ｏｆ ｆｏｕｒ ｍｅｄｉｕｍ－ｓｉｚｅｄ ｅａｒｔｈｑｕａｋｅｓ ｔｈａｔｏｃｃｕｒｒｅｄ ｉｎ ｔｈｅ Ｇｕｌｆ ｏｆ Ｃａｌｉｆｏｒｎｉａ．Ｊ．Ｇｅｏｐｈｙｓ．Ｒｅｓ．，１１３（Ｂ１０３０１），ｄｏｉ：１０．１０２９／２００７ＪＢ００５３２３．Ｘｕ Ｘ Ｔ，Ｘｕ Ｚ Ｈ，Ｚｈａｎｇ Ｄ Ｎ．１９９５．Ａ ｐｒｏｂａｂｉｌｉｓｔｉｃ ｇｒｉｄ ｔｅｓｔｍｅｔｈｏｄ ｏｆ ｄｅｔｅｒｍｉｎｉｎｇ ｅａｒｔｈｑｕａｋｅ ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ｕｓｉｎｇ Ｐ－ｗａｖｅ ｏｎｓｅｔ ｐｏｌａｒｉｔｙ ｄａｔａ．Ｓｅｉｓｍｏｌｏｇｉｃａｌ ａｎｄ ＧｅｏｍａｇｎｅｔｉｃＯｂｓｅｒｖａｔｉｏｎ ａｎｄ Ｒｅｓｅａｒｃｈ（ｉｎ Ｃｈｉｎｅｓｅ），１６（４）：３４－４２．Ｙｕ Ｃ Ｑ，Ｔａｏ Ｋ，Ｃｕｉ Ｘ Ｆ，ｅｔ ａｌ．２００９．Ｐ－ｗａｖｅ ｆｉｒｓｔ－ｍｏｔｉｏｎ ｆｏｃａｌｍｅｃｈａｎｉｓｍ ｓｏｌｕｔｉｏｎｓ ａｎｄ ｔｈｅｉｒ ｑｕａｌｉｔｙ ｅｖａｌｕａｔｉｏｎ．Ｃｈｉｎｅｓｅ Ｊ．Ｇｅｏｐｈｙｓ．（ｉｎ Ｃｈｉｎｅｓｅ），５２（５）：１４０２－１４１１，ｄｏｉ：１０．３９６９／ｊ．ｉｓｓｎ．０００１－５７３３．２００９．０５．０３０．Ｚｈａｏ Ｌ Ｓ，Ｈｅｌｍｂｅｒｇｅｒ Ｄ Ｖ．１９９４．Ｓｏｕｒｃｅ ｅｓｔｉｍａｔｉｏｎ ｆｒｏｍ ｂｒｏａｄｂａｎｄｒｅｇｉｏｎａｌ ｓｅｉｓｍｏｇｒａｍｓ．Ｂｕｌｌ．Ｓｅｉｓｍ．Ｓｏｃ．Ａｍｅｒ．，８４（１）：９１－１０４．Ｚｈｅｎｇ Ｊ Ｃ，Ｌｉ Ｄ Ｍ，Ｗａｎｇ Ｐ，ｅｔ ａｌ．２０１５．Ｆｏｃａｌ ｍｅｃｈａｎｉｓｍｓ ａｎｄｓｅｉｓｍｉｃ ｔｅｃｔｏｎｉｃ ｆｅａｔｕｒｅｓ ｏｆ ｔｈｅ ２０１３Ｌａｉｚｈｏｕ Ｍ４．６ｅａｒｔｈｑｕａｋｅｓｅｑｕｅｎｃｅ．Ｓｅｉｓｍｏｌｏｇｙ ａｎｄ Ｇｅｏｌｏｇｙ（ｉｎ Ｃｈｉｎｅｓｅ），ｉｎ ｐｒｅｓｓ．Ｚｈｅｎｇ Ｊ Ｃ，Ｃｈｅｎ Ｙ Ｔ．２０１２．Ｓｔａｂｉｌｉｔｙ ｏｆ ｓｐａｒｓｅ ｓｔａｔｉｏｎ ｄａｔａｉｎｖｅｒｓｉｏｎ ｆｏｒ ｄｅｖｉａｔｏｒｉｃ ｍｏｍｅｎｔ ｔｅｎｓｏｒ ｓｏｌｕｔｉｏｎ ｏｆ ｒｅｇｉｏｎａｌｅａｒｔｈｑｕａｋｅｓ．Ａｃｔａ Ｓｅｉｓｍｏｌｏｇｉｃａ Ｓｉｎｉｃａ（ｉｎ Ｃｈｉｎｅｓｅ），３４（１）：３１－４３．Ｚｈｅｎｇ Ｙ，Ｍａ Ｈ Ｓ，ＬüＪ，ｅｔ ａｌ．２００９．Ｓｏｕｒｃｅ ｍｅｃｈａｎｉｓｍ ｏｆ ｓｔｒｏｎｇａｆｔｅｒｓｈｏｃｋｓ（Ｍｓ≥５．６）ｏｆ ｔｈｅ ２００８／０５／１２Ｗｅｎｃｈｕａｎ ｅａｒｔｈｑｕａｋｅ ａｎｄｔｈｅ ｉｍｐｌｉｃａｔｉｏｎ ｆｏｒ ｓｅｉｓｍｏｔｅｃｔｏｎｉｃｓ．ＳｃｉｅｎｃｅｉｎＣｈｉｎａＳｅｒｉｅｓＤ：Ｅａｒｔｈ Ｓｃｉｅｎｃｅｓ，５２（６）：７３９－７５３．Ｚｈｕ Ｌ Ｐ，Ｈｅｌｍｂｅｒｇｅｒ Ｄ Ｖ．１９９６．Ａｄｖａｎｃｅｍｅｎｔ ｉｎ ｓｏｕｒｃｅ ｅｓｔｉｍａｔｉｏｎｔｅｃｈｎｉｑｕｅｓ ｕｓｉｎｇ ｂｒｏａｄｂａｎｄ ｒｅｇｉｏｎａｌ ｓｅｉｓｍｏｇｒａｍｓ．Ｂｕｌｌ．Ｓｅｉｓｍ．Ｓｏｃ．Ａｍｅｒ．，８６（５）：１６３４－１６４１．Ｚｈｕ Ｌ Ｐ，Ｒｉｖｅｒａ Ｌ Ａ．２００２．Ａ ｎｏｔｅ ｏｎ ｔｈｅ ｄｙｎａｍｉｃ ａｎｄ ｓｔａｔｉｃｄｉｓｐｌａｃｅｍｅｎｔｓ ｆｒｏｍ ａ ｐｏｉｎｔ ｓｏｕｒｃｅ ｉｎ ｍｕｌｔｉｌａｙｅｒｅｄ ｍｅｄｉａ．Ｇｅｏｐｈｙｓ．Ｊ．Ｉｎｔ．，１４８（３）：６１９－６２７．附中文参考文献刁桂苓，于利民，李钦祖．１９９２．震源机制解的系统聚类分析———以海城地震序列为例．中国地震，８（３）：８６－９２．韩立波，蒋长胜，包丰．２０１２．２０１０年河南太康ＭＳ４．６地震序列震源参数的精确确定．地球物理学报，５５（９）：２９７３－２９８１，ｄｏｉ：１０．６０３８／ｊ．ｉｓｓｎ．０００１－５７３３．２０１２．０９．０１６．韩立波，蒋长胜．２０１２．２０１１年６月８日新疆托克逊ＭＳ５．３地震震源机制解反演．地震学报，３４（３）：４１５－４２２．黄建平，倪四道，傅容珊等．２００９．综合近震及远震波形反演２００６文安地震（ＭＷ５．１）的震源机制解．地球物理学报，５２（１）：１２０－１３０．李冬梅，郑建常．２０１４．２０１３年１０月１日乳山震群重新定位结果分析．齐鲁地震科学，１０（１）：２６－２８．龙锋，张永久，闻学泽等．２０１０．２００８年８月３０日攀枝花—会理６．１级地震序列ＭＬ≥４．０事件的震源机制解．地球物理学报，５３（１２）：２８５２－２８６０，ｄｏｉ：１０．３９６９／ｊ．ｉｓｓｎ．０００１－５７３３．２０１０．１２．００８．吕坚，郑勇，倪四道等．２００８．２００５年１１月２６日九江－瑞昌ＭＳ５．７，ＭＳ４．８地震的震源机制解与发震构造研究．地球物理学报，５１（１）：１５８－１６４，ｄｏｉ：１０．３３２１／ｊ．ｉｓｓｎ：０００１－５７３３．２００８．０１．０２０．许向彤，许忠淮，张东宁．１９９５．求震源机制Ｐ波初动解的格点尝试概率法．地震地磁观测与研究，１６（４）：３４－４２．俞春泉，陶开，崔效锋等．２００９．用格点尝试法求解Ｐ波初动震源机制解及解的质量评价．地球物理学报，５２（５）：１４０２－１４１１，ｄｏｉ：１０．３９６９／ｊ．ｉｓｓｎ．０００１－５７３３．２００９．０５．０３０．郑建常，王鹏，李冬梅等．２０１５．２０１３年莱州Ｍ４．６地震序列震源机制与发震构造初探．地震地质，待刊．郑建常，陈运泰．２０１２．稀疏台网反演区域地震偏量矩张量解的稳定性．地震学报，３４（１）：３１－４３．郑勇，马宏生，吕坚等．２００９．汶川地震强余震（Ｍｓ≥５．６）的震源机制解及其与发震构造的关系．中国科学Ｄ辑：地球科学，３９（４）：４１３－４２６．（本文编辑 何燕）２６４。
地质统计学反演技术专题之三：基于MCMC的StatMod及RockMod
地质统计学反演技术专题之三：基于MCMC的StatMod及RockMod概述在前二期的文章中，我们介绍了HampsonRussell工具包中的地质统计学地震储层表征模块GeoSI（地质统计学随机反演）。
在这篇文章中，我们将介绍CGG GeoSoftware的另外两个地质统计学地震储层表征工具，即Jason地学软件工具包中的StatMod和RockMod。
如同GeoSI，这两个工具也是在地层地质模型框架内精细地整合地质信息和地震数据，生成多个高精度的实现。
这些实现都是对储层的预测，也可用于不确定性分析和风险评估。
StatMod和RockMod紧密融合地球科学领域内所有相关的不同种类数据，从而产生多个储层模型。
这种高度跨学科一致性的特点，确保了油藏模型是符合实际生产开发情况的，并最大限度地挖掘了测量数据和推断信息的价值。
其成果是一系列精确的深度域的储层模型，可用于预测油田储量、流体流动样式和产量估算。
同时，这些模型也为定量估计不确定性提供了可靠的依据；不确定性评估与先验信息、专家知识、井资料和地震具有直接关系。
这些储层模型在远离井点处具有较高预测性，而这恰恰是传统地质统计学建模的痛点。
StatMod使用一个（全）叠加地震数据，所以只在纵波阻抗足以区分岩相时使用。
RockMod同时使用多个（部分）叠加（AVO/AVA）地震数据，因此在需要多个弹性参数组合（如纵波阻抗、横波阻抗和密度）才能区分岩相时使用。
利用地震数据，以及不同领域来源的数据，RockMod可同时得到岩相体、弹性参数体和油藏工程属性体。
图1. RockMod同时反演出弹性属性、岩相和油藏工程属性简要研发历史最早的Jason地质统计学储层表征方法可追溯至1996年。
它始于基于褶积模型的叠后地震反演。
模拟方法包含简单的序贯高斯模拟（SGS）、序贯指示模拟（SIS）和SGS同时模拟、带趋势的序贯指示模拟（SISTR）和序贯高斯协模拟（SGCS）。
光谱地壳计划探索新一代矿产勘查技术
"修回日期 收稿日期 $ % ! $ % # $ = $ % ! $ % > ! " ? $ $ 基金项目国家自然科学基金项目# 和地质调查局项目# 资助% " ! % @ $ $ " = ! $ ! $ % ! ! ! $ % $ $ $
作者简介 王晋年# 男! 研究员! 长期致力于航空遥感& 高光谱遥感! 以及空间信息在地质找矿与生态环境中的应用 ! ! A & & $
3, 2 5 航空成像光谱数据分析和提取了铀矿信息" 国土资源航空物探遥感中心和澳大利亚合作! 利用
射计 , 欧空局自制小卫星 K 5 4 0 *! * ( R , 搭载的 紧密型高分辨率成像仪 6 和 N* 2 5 0 + 7 2 5 , 4 卫星 搭载的中分辨率成像光谱仪 30 都得到了广 ! * 2 5 泛的应用%另外! 澳大利亚& 加拿大& 德国& 法国& 芬 兰& 日本等国也已相继发射或即将发射搭载高光谱 成像仪的卫星%故此! 许多国家在矿产勘查和资源 利用等方面取得了显著的效果% 我国的一些相关部门和研究单位也利用国内 外的高光谱卫星数据开展了许多岩矿地质方面的 应用研究%国土资源航空物探遥感中心应用 0 ( 1 ! 的高光谱传感器 N 对西藏驱龙地区 Q H D 8 < 数据! O P
6 8 B I L 8 J S金矿及其相邻地区进行了机载成像光谱 技术的矿物填图% 美国 ) $ % % @ 年! 5 / 5用 N 3 C O P
%
高光谱遥感技术在油气探测& 固体矿产资源普 查与勘查等方面发挥了重要作用% 它为成矿系统 和成矿过程提供完备& 规模化的可靠稳定的岩矿光 谱信息源! 满足地质勘查业务化的要求% 美国& 澳 大利亚等国家的矿业部门已大量采用高光谱遥感 技术& 高光谱岩芯探测编录技术用于精细采矿& 选 矿等工作! 大大提高了工作效率% 诸如! 利用地面 成像光谱仪获取地表岩石矿物的精确高光谱信息" 岩芯光谱扫描系统获取地下岩石矿物的高光谱数 据! 以开展矿物的 地 表 和 地 下 的 全 方 位 立 体 填 图 成为了目前世界各国研究的热点% 等! 在航天高光谱遥感方面! 美国国家航空航天局
地质统计学在固体矿山中的应用
ｓ ｅ ｒ ｖ ｅ ｓ ｏ ｆ ｔ ｈ ｅ ｍｉ ｎ ｅ ｒ ａ ｌ ｒ ｅ ｓ ｏ ｕ ｒ ｃ ｅ ｓ ，ｔ ｈ ｅ ｄ ｉ ｓ ｔ ｒ ｉ ｂ ｕ ｔ ｉ ｏ ｎ ａ ｌ ｃ ｈ ａ ｒ ａ ｃ ｔ ｅ ｒ ｉ ｓ ｔ ｉ ｃ ｓ ｏ ｆ ｔ ｈ ｅ ｍｉ ｎ ｅ ｒ ａ ｌ ｓ ，ｔ ｈ ｅ ｃ ｌ ａｓ ｓ ｉ ｆ ｉ ｃ ａ ｔ ｉ ｏ ｎ ｏ ｆ ｔ ｈ ｅ ｒ ｅ ｓ ｅ ｒ ｖ ｅ ｓ ， ｔ ｈ ｅ ｏ ｐ ｔ ｉ — ｍｉ ｚ ａ ｔ ｉ ｏ ｎ ｏ ｆ ｔ ｈ ｅ ｄ ｅ ｔ ｅ ｒ ｍｉ ｎ ｅ ｅ ｘ ｐ ｌ ｏ ｒ ａ ｔ ｏ ｒ ｙ ｇ ｒ ｉ ｄ ａ ｎ ｄ ｔ ｈ ｅ ｉ ｎ ｖ ｅ ｓ ｔ ｉ ｇ ａ ｔ ｉ ｏ ｎ ｏ ｆ ｔ ｈ ｅ ｍｉ ｎ ｅ ｒ ａ ｌ ｓ ． Ｋｅ ｙ ｗｏ ｒ ｄ ｓ ：ｇ ｅ ｏ ｓ ｔ ａ ｔ ｉ ｓ ｔ ｉ ｃ ｓ ， ｋ ｒ ｉ ｇ ｉ ｎ ｇ， ｍｉ ｎ ｅ ｒ ａ ｌ ｒ ｅ ｓ ｏ ｕ ｒ ｃ ｅ ｓ ， ｍｉ ｎ ｅ ｒ ａ ｌ ｅ ｘ ｐ ｌ ｏ ｒ ａ ｔ ｉ ｏ ｎ
摘 要 ： 地质 统计 学作 为一 门新 兴 的 交叉 学科 ， 在近 ５ ０年 的研 究 和 实践 中得 到 了很 大的 发展 ， 近 年
来 又被 称 为 空 间信 息统计 学 。 国 内外 的生产 实践 表 明 ， 地 质统 计 学除 了在 地 学科研 方 面具有 明 显 的优越
性， 在 固体矿 山中的应 用也 越 来越 广 泛 。本 文主 要 介 绍 了地 质 统 计 学在 矿 产 资源储 量计 算 、 矿 产 分 布特
２ ．Ｆａ ｃ ｕ ｌ ｔ ｙ ｏ ｆ Ｍｅ ｔ ａ ｌ ｌ ｕ ｒ ｇ ｉ ｃ ａ ｌ ａ ｎ ｄ Ｅｎ ｅ ｒ ｇ ｙ Ｅｎ ｇ ｉ ｎ ｅ ｅ ｒ ｉ ｎ ｇ，Ｋｕ ｎ ｍｉ ｎ ｇ Ｕｎ ｉ ｖ ｅ ｒ ｓ ｉ ｔ ｙ ｏ ｆ Ｓ ｃ ｉ ｅ ｎ ｃ ｅ ａ ｎ ｄ Ｔｅ ｃ ｈ ｎ ｏ ｌ ｏ ｇ ｙ ， Ｋｕ ｎ ｍｉ ｎ ｇ ６ ５ ０ ０ ９ ３ ，Ｃｈ ｉ ｎ ａ ）
2010各专业最受欢迎十大教授榜
所在学校
所在院系
丁为祥
陕西师范大学
政治经济学院
金林南
河海大学
公共管理学院
林美茂
中国人民大学
哲学院
刘宏九
辽宁大学
哲学与公共管理学院
聂敏里
中国人民大学
哲学院
徐英瑾
复旦大学
哲学学院
臧峰宇
中国人民大学
哲学院
赵林
武汉大学
哲学学院
周濂
中国人民大学
哲学院
朱良志
北京大学
哲学系
2010政治学与公共管理专业最受欢迎十大教授榜(211院校类)
姚淦铭
江南大学
文学院
易平
南昌大学
人文学院
声明：搜狐教育、评师网联合发布，转载请注明出处！
（注：可编辑下载，若有不当之处，请指正，谢谢!）
2010计算机专业最受欢迎十大教授榜(211院校类)
老师姓名
所在学校
所在院系
董洁
北京科技大学
信息工程学院
韩智
南开大学
软件学院
李胜睿
厦门大学
信息科学与技术学院
李雄飞
吉林大学
计算机科学与技术学院
卢红星
郑州大学
信息工程学院
卢燕宁
西北大学
信息科学与技术学院
孙效里
大连理工大学
电子与信息工程学院
王宇颖
哈尔滨工业大学
计算机科学与技术学院
张树泉
北京科技大学
土木与环境工程学院
朱仲涛
清华大学
智能技术与系统国家重点实验室
2010经济及管理专业最受欢迎十大教授榜(211院校类)
老师姓名
所在学校
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
-1.0
-0.5
0
+0.5
正相关程度增加
+1.0
r
负相关程度增加
• 变异函数(Variogram)
• 协方差函数(Covariance)
• 平稳假设(Stationary assumption)
• 内蕴假设(intrinsic assumption)
• 估计方差(Estimation variance) 1
• 变异函数及变异曲线
• 变异函数及变异曲线
• Z(x)和Z(x+h)的相关与变异函数 的关系
•变异函数与协方差之间的关系
• 存在趋势的变异函数
•具有空穴效应的变异函数
•几何异向性和带状异向性比较
• 变异函数在原点处的性状
• 变异函数的理论模型
C0-块金常数 a—变程 C0+C—基台
p=(i-0.5)/n 0.025 0.125 0.175
标准正态分位数 -1.9600 -1.1503 -0.9346
5
7 10 12 14 15 17 18 19 20
39
43 49 51 56 59 61 74 75 99
0.225
0.325 0.475 0.575 0.675 0.725 0.825 0.875 0.925 0.975
1、区域化变量理论
• 区域化变量 G.马特隆定义区域化变量是：一种在空间上具有数值 的实函数，它在空间的每一个点取一个确定的数值，即当 由一个点移到下一个点时，函数值是变化的。 从地质及矿业角度来看，区域化变量具有如下性质： （1）空间局限性：即它被限制在一个特定的空间（如一 个矿体内）；该空间称为区域化的几何域；区域化变量是 按几何支撑定义的。 （2）连续性：不同的区域化变量具有不同的连续性，这 种连续性是通过相邻样品之间的变异函数来描述的。 （3）异向性：当区域化变量在各个方向上具有相同的性 质时称各向同性，否则称各向异性。 （4）相关性：一定范围内、一定程度上的空间相关性， 当超出这一范围后相关性减弱以至消失。 （5）对于任一区域化变量而言，特殊的变异性可以叠加 在一般规律之上。
泛克立格法方程组举例
c( x1 , x1 ) c( x1 , x2 ) c( x1 , xn ) c( x2 , x1 ) c( x2 , x2 ) c( x2 , xn ) c( xn , x1 ) c( xn , x2 ) c( xn , xn ) 1 1 1 x1 x2 xn y1 y2 yn 2 2 2 x x x 1 2 n 2 2 y12 y2 yn x1 y1 x2 y2 xn yn 1 x1 y1 1 x2 1 xn 0 0 0 0 0 0 0 0 0 0 0 0 x12 y12 x1 y1 1 c( x1 , x0 ) 2 2 c ( x , x ) y2 x2 y2 x2 y2 2 1 0 2 2 c ( x , x ) y n xn y n x n y n n n 0 1 0 0 0 0 0 0 0 0 0 1 x0 0 0 0 0 2 y0 2 x 0 0 0 0 3 0 2 y0 0 0 0 0 4 0 0 0 0 5 x0 y0
• 地质统计学诞生过程
上世纪40年代后期，当南非统计学家H.S西奇 尔（Sichel）判明南非各金矿的样品品位呈对数 正态分布以后，才真正确立了地质统计学的开端。 1951年，南非的矿山工程D.G.克立格Daniel Krige）在H.S西奇尔研究的基础上提出一个论点： “可以预计，一个矿山总体中的金品位的相对变 化要大于该矿山某一部分中的金品位的相对变 化”。换句话说，以较近距离采集的样品很可能 比以较远距离采集的样品具有更近似的品位。这 一论点是描述在多维空间内定义的数值特征的空 间统计学据以建立的基础。
对特高品位进行处理。 计算各个方向的变异函数并确定采用球状或其 他模型，同时确定结构参数。
确定椭球体的长轴、次轴和最小轴的方向及比 例，通过交叉验证。 确定块模型的尺寸及相关的属性。 应用普通克立格等方法对矿块进行估值。 对储量进行分类统计，做品位-吨位曲线。
5. 建立反应空间变异性的 椭球体 6. 创建一个空的块模型 7. 块模型估值 8. 矿块储量的基本统计
负线性相关
不相关
正线性相关
相关系数
(计算公式)
• 样本相关系数的计算公式
XY
r
Cov ( X,Y) D( X ) D( Y )
( x x )( y y ) (x x) ( y y)
2
2
相关系数
(取值及其意义)
条件模拟：计算正态分位数计算
条件模拟：产生随机数
• 产生一个介于0和1之间的一个随机数从CDF （累计分布曲线）中读取一个模拟值
•协同克里格法（Co-Kriging）
矿业软件Surpac的应用
•基本统计分析 •空间变异性分析 •选用不同的方法估值
应用Surpac估算储量的步骤
1. 样品组合 2. 基本统计分析 地质统计学要求数据支撑大小相同 研究数据的分布特征。如数据不符合正态分布 则进行变换。
3. 特异值的处理 4. 建立变异函数模型
• 变异函数的理论模型
•三种有基台值标准模型比较
•具有空穴效应的变异函数
• 无基台值标准模型
• 变异函数结构分析
• 几何异向性结构的套合
•带状异向性结构的套合
•比例效应
• 相对变异函数改正
• 变异函数的套合
• 普通克立格法(Ordinary Kriging)
• 普通克立格方程组
地质统计学是在1977年由美国福禄尔采 矿金属有限公司（Flour Mining & Meta Incorporation）H.M.Parker博士随美中贸易 全国委员会矿业代表团来华访问，传入我 国，继而得到进一步的发展。1989年11月 召开的全国第一届地质统计学学术讨论会， 地质统计学在我国的发展进入了一个新的 阶段，理论研究更加深入，涉及的方法原 理更加广泛。地质统计学已经被广泛地承 认是矿床评价的必要部分，在我国已经认 可用地质统计学对矿床进行评价的地质报 告。
地质统计学
张树泉
北京科技大学
一、地质统计学的发展历史和现状
• 什么是地质统计学？ 地质统计学(Geostatistics)包含经典统计学 与空间统计学，按其基本原理可定义为： 地质统计学是以区域化变量理论为基础， 以变异函数为主要工具，研究那些在空间 分布上既有随机性，又有结构性的自然现 象的科学 。
应用Surpac估算储量的步骤
1. 样品组合 2. 基本统计分析 地质统计学要求数据支撑大小相同 研究数据的分布特征。如数据不符合正态分布 则进行变换。建立组合样文件。
3. 特异值的处理 4. 建立变异函数模型
对特高品位进行处理。 计算各个方向的变异函数并确定采用球状或其 他模型，同时确定结构参数。
确定椭球体的长轴、次轴和最小轴的方向及比 例，通过交叉验证。 确定块模型的尺寸及相关的属性。 应用普通克立格等方法对矿块进行估值。 对储量进行分类统计，做品位-吨位曲线。
5. 建立反应空间变异性的 椭球体 6. 创建一个空的块模型 7. 块模型估值 8. 矿块储量的基本统计
i 1 3 4
数据 20 23 37
到上世纪60年代，才认识到需要把样品值之间的相似 性作为样品间距离的函数来加以模拟，并且得出了半变异 函数。法国概率统计学家马特隆（Matheron）创立了一个 理论框架，为克立格作出的经验论点提供了精确而简明的 数学阐释。马特隆创造了一个新名词“克立格法” （Kriging），藉以表彰克立格在矿床的地质统计学评价工 作中所起到的先驱作用。即1962年，马特隆在克立格和西 奇尔研究的基础上，将他们的成果理论化、系统化，并首 先提出了区域化变量（Regionalized variable）的概念， 为了更好地研究具有随机性及结构性的自然现象，提出了 地质统计学（Geostatistics）一词，发表了《应用地质统 计学》，该著作的出版标志着地质统计学作为一门新兴边 缘学科而诞生。地质统计学开始进入了学术界。在法国枫 丹白露成立了地质统计学中心（Centre de Geostatistiques），培养了一大批学员，不仅为地质统计 学的研究而且为它的传播起到了巨大的作用。
• 估计方差(Estimation variance) 2
• 估计方差(Estimation variance) 3
• 离差方差 (Dispersion variance)
2、变异函数及结构分析
为表征一个矿床金属品位等特征量的变 化，经典统计学通常采用均值、方差等一 类参数，这些统计量只能概括该矿床中金 属品位等特征量的全貌，却无法反映局部 范围和特定方向上地质特征的变化。地质 统计学引入变异函数这一工具，它能够反 映区域化变量的空间变化特征——相关性 和随机性，特别是透过随机性反映区域化 变量的结构性，故变异函数又称结构函数。
f ( x, y) a0 a1x a2 y a3 x2 a4 y 2 a5 xy
• 指示克立格法(Indicator Kriging)
• 指示克立格法(Indicator Kriging)
•协同克里格法（Co-Kriging）
• 协同克里格法（Co-Kriging）
•协同克里格法（Co-Kriging）
正态分布的误差图示
x z 2 x
- 2.58x -1.65 x
x
+1.65x + 2.58x
x
-1.96 x
+1.96x