你还记得那年的东日本大地震吗第一波新闻里写的是8.8级第二天多家机构修定为9.0级。这不是简单的四舍五入而是随着更多台站数据汇入震级被重新测定了。这个我们从小就听熟的“几级”其实是地球物理学里最核心的概念之一magnitude震级。我一直觉得震级是地震科学里最容易被误解、也最值得亲手算一遍的东西。这篇博文想把震级从一个新闻名词拆成一套你可以自己复现的数据分析流程——从公开地震目录下载、Python数据处理到计算b值、画震级-频度曲线。想真正把“震级”弄明白的人不管是做数据处理的学生还是对地震科普有硬核兴趣的爱好者这篇内容应该能帮你在拿到目录后不再一头雾水。1. 震级到底是什么——先从那个“几级”说起1.1 为什么必须用对数标尺震级的定义绕不开对数。这不是数学家为了显得高深才这么干的而是地震释放能量的跨度实在太大。人感觉不到的微震能量大概只有几万焦耳1960年智利Mw 9.5大地震释放的能量大约在10的18次方焦耳量级两者相差十几个数量级。如果按“差多少焦耳”来画图小震全都挤在零点附近什么也看不出来。这就像你要比较一个人的存款是10万还是1000亿你会说“差5个数量级”而不是关心具体差多少块。里克特在1935年定义里氏震级时采用的办法就是取对数把距离震中100公里处、标准伍德-安德森地震仪记录到的最大振幅单位微米取以10为底的对数再减去一个校准常数。换句话说振幅每放大10倍震级就加1。这么做的直接好处是地震目录里的震级数字变得极其直观。6级、7级、8级数字相差不大但它们对应的振幅差了100倍能量差了大约31600倍。正是因为取了对数人类才能在一张纸上同时写下微震和大震的观测记录。这个数学选择是后面所有统计分析和历史目录比对的前提。1.2 ML、Ms、Mb、Mw同一场地震为什么有好几个数字实际操作中你会发现一场地震往往有多个震级数字这不是新闻乱报而是不同测量方式算出来的结果。里氏震级ML是最早的定义但它依赖的是南加州浅源地方震的观测所以后来人们又根据不同的地震波震相发展出了几种新标度。体波震级mb使用1秒左右短周期P波深源地震也能记到早期全球速报常用但它和地震真正释放能量之间的关系不够直接。面波震级Ms使用20秒左右周期的面波对浅源大震更稳定但震级超过一定范围后会“饱和”。矩震级Mw由地震矩M0推导而来M0等于剪切模量乘以断层破裂面积再乘以平均滑动量。公式是Mw (2/3)log10(M0) - 6.07M0单位为N·m。为什么现在新闻里的大震级别都以Mw为准因为矩震级的物理意义明确它直接对应断层破裂的规模和位错量而且不会因为地震波高频部分饱和而封顶。1960年智利地震早期用面波震级测定只有8.3级后来用矩震级测定为9.5级差别就是这么来的。震级类型基于的观测优点局限ML短周期地方震振幅定义早、资料多仅适合浅源地方震易饱和mb1秒周期体波P波可测深源地震与释放能量相关性较弱Ms20秒周期面波适合远震浅源超过8.0到8.5会饱和Mw地震矩M0物理意义明确、不饱和需要反演震源机制计算慢1.3 一份目录里到底存的是哪种震级我见过太多人从USGS下载一份CSV就开始统计画图算b值结果画出来的曲线在某处莫名其妙拐弯。问了一句“你用的哪个震级字段”对方才意识到问题CSV里有一个mag字段还有一个magType字段后者标明了这一条记录的震级标度。同一份目录里可能是mb、ML、Ms、Mw混着来的。所以拿到任意一份地震目录后的第一件事不是急着画图而是先看元数据。具体到这个csv里magType都有哪些取值各占多少比例。如果这个步骤不做后面所有分析都建立在沙子上。2. 数据获取拉一份真实的地震目录2.1 公开数据源怎么选全球公开的地震目录有好几套最常用的是USGS的ANSS综合目录覆盖全球API开放字段完整适合大多数分析场景。如果你做更底层的科研ISC国际地震中心目录数据质量和完整性更好但下载需要注册和邮件确认门槛略高。GCMT项目专门提供矩心矩张量解只有Mw但事件数量有限不适合做小震级统计分析。对于初学者来说USGS是最合适的起点。我平时做快速验证也直接用USGS它的FDSN事件查询接口可以直接用URL拼参数返回CSV用pandas读进来就行不需要任何特殊权限。如果想了解中国及周边区域的地震活动中国地震台网中心也提供公开目录但字段规范和数据格式跟USGS不太一样后期清洗会稍微多花点时间。2.2 用Python拉取USGS目录USGS的接口格式很友好核心路径是earthquake.usgs.gov/fdsnws/event/1/query配合formatcsv和起止时间、震级范围就能拼出下载URL。我用pandas的read_csv直接读远程URL几秒钟就能把多年目录拉下来。import pandas as pd start 2015-01-01 end 2023-12-31 url ( https://earthquake.usgs.gov/fdsnws/event/1/query f?formatcsvstarttime{start}endtime{end} minmagnitude5.0 ) df pd.read_csv(url) print(df.shape) print(df[magType].value_counts())初次运行可能会因为网络或限流失败建议加一层重试或者按年份循环下载再用pd.concat合并。另一个容易忽略的点是时间USGS接口里的时间是UTC如果你按当地时间划分研究窗口务必先加上时区处理否则对比历史报告时会错位。2.3 字段解读与数据清洗拿到表格后先花点时间理解字段。time、latitude、longitude、depth这些都是基本定位信息。depth单位是千米偶尔会出现负值表示震源深度在海平面以上多见于海洋区域的特殊算法结果。mag是震级数值magType紧跟其后magSource表示测定机构。nst是参与定位的台站数gap是台站方位角空隙这两个字段直接反映定位质量。nst小于3的记录定位可靠性很低建议删除。gap大于180说明台站几乎只分布在震中一侧震中位置可能偏差很大。mag为空或者为NaN的记录直接删掉别含糊。df df.dropna(subset[mag, latitude, longitude]) df df[df[nst] 3] if nst in df.columns else df df df[df[gap] 180] if gap in df.columns else df另外USGS目录偶尔会有重复事件我遇到过同一条地震同时被几个机构上报的情况。简单去重可以按time、lat、lon、mag四列做drop_duplicates更稳妥的办法是根据事件ID字段去重。清洗完成后最好把结果存成parquet或者csv后面反复读取会快很多。3. Python实操震级分布与b值计算3.1 Gutenberg-Richter关系震级-频度的核心规律地震不是随机事件它有一个著名的统计规律叫古登堡-里克特关系log10(N) a - b·M。这里的N是震级大于等于M的地震数。意思是说震级每降低1级地震数量大约变成原来的10倍。b值通常接近1但不同区域的b值有差异反映了地壳应力状态和断层类型的区别。我第一次亲手拟合这条关系时才真正理解为什么有人说地震具有“自相似性”。从M2的微震到M8的大震只要数据足够完整它们之间的数量比例关系基本稳定。这也是为什么在地震危险性评估中人们常用b值反推未来大震的发生概率。3.2 计算b值Aki最大似然估计与Bootstrap误差计算b值最常用的方法是Aki在1965年提出的最大似然估计公式非常简单b log10(e) / (M_mean - M_c ΔM/2)。其中M_mean是样本平均震级M_c是目录完整性震级ΔM是震级分档宽度。麻烦的是M_c到底取多少这决定了计算结果的可靠性。实操中我会先画累积震级-频度图找到偏离直线的那一段把M_c取在曲线开始变直线的位置。然后对M_c以上的数据进行Aki估计再用Bootstrap重采样1000次估计误差。import numpy as np mag df[mag].values mag mag[mag M_c] b np.log10(np.e) / (mag.mean() - M_c 0.05) print(fb value: {b:.3f}) rng np.random.default_rng(42) boot_b [] for _ in range(1000): sample rng.choice(mag, sizelen(mag), replaceTrue) boot_b.append(np.log10(np.e) / (sample.mean() - M_c 0.05)) ci np.percentile(boot_b, [2.5, 97.5]) print(f95% CI: {ci[0]:.3f} - {ci[1]:.3f})注意M_c取高一点虽然会损失样本量但能保证进入统计的都是被台网完整记录的事件。如果把不完整的低震级数据也算进来b值会被明显拉低。3.3 能量对比1级之差差了31倍还是32倍这是科普里最常见的误区。震级每差1级振幅差10倍能量差大约31.6倍而不是人们常说的“10倍”。能量和震级的经验关系是log10(E) 1.5·M 4.8E单位是焦耳。这意味着M从6.0涨到7.0能量约增大10^1.5倍也就是31.62倍。把这个关系落到具体数字上会更直观一次6级地震释放的能量约等于6.3×10^13焦耳7级则是2.0×10^15焦耳左右。8级地震释放的能量相当于大约1000次6级地震或者31.6次7级地震。公众听到“8级是6级的100倍”大概率是指这个量级关系。我在目录里随便抽了几个大型地震事件算了算映象最深的是大震与小震之间的能量差远超直觉。这也是为什么一次7级地震造成的破坏可能比几十次6级地震加在一起还要大因为能量是按指数增长的。3.4 可视化画一张震级-频度图数据有了b值算了最后用图把所有内容串起来。画直方图和累积频度曲线时竖轴用对数刻度这样G-R关系的直线形态才看得清楚。import matplotlib.pyplot as plt bins np.arange(mag.min(), mag.max() 0.1, 0.1) counts, edges np.histogram(mag, binsbins) centers (edges[:-1] edges[1:]) / 2 plt.figure(figsize(8, 5)) plt.bar(centers, counts, width0.1, alpha0.6, labelhistogram) plt.yscale(log) plt.xlabel(Magnitude) plt.ylabel(Number of events) plt.title(Magnitude-Frequency Distribution) plt.legend() plt.grid(alpha0.3) plt.show()画完之后理想情况下你会看到尾部拖得很长的高震级区以及直方图的整体包络基本是一条直线。如果这条直线在中段有转折先别急着解释地质原因回头检查是不是M_c取错或者震级标度混用。4. 踩过的坑与排查技巧4.1 震级标度混用造成的假象最典型的坑就是不分青红皂白地把mb、ML、Ms、Mw当成同一个量纲去统计。我试过把某一年全球目录里的所有震级直接拿来画直方图结果在M5附近出现了一个明显的“鼓包”一开始还以为是区域地震活动异常后来查了magType才发现那堆数据大多是mb而其他标度在这个区间的分布并不一样。建议拿到数据后先跑一句value_counts看分布。如果混用情况不严重可以只保留一种标度分析如果必须合并要给出明确的转换关系并且说明这种转换本身有误差。比如用经验公式把mb转换成Mw这只是一个近似不适合用于精度要求极高的研究。4.2 震级饱和为什么大震被低估了震级饱和是一个很隐蔽的问题。里氏震级ML依赖短周期波振幅当断层破裂尺度非常大时短周期波的振幅增长会放缓导致ML在大约7级以后“封顶”。面波震级Ms在8.0到8.5之间也会出现饱和。这就是为什么一些历史大地震在旧目录里震级偏低。处理办法是涉及大震研究时只用Mw历史目录中如果只有Ms要查阅当时的测定方法必要时查阅已发表的震级转换文献。我踩过这个坑之后现在只要看到目录里有一堆8.0以上的Ms记录第一反应就是这批数据需要谨慎处理。4.3 震级完整性与小样本问题b值估计的另一个致命问题是样本量和完整性。Aki公式依赖平均值而平均值对样本量有最低要求。超过一定震级的事件本来就少如果M_c取得太低混入大量检测不完备的小震b值会偏低如果M_c取得太高样本量不足误差又太大。实操经验是样本少于50个时基本别谈b值。我个人的习惯是至少取100个事件再算这样Bootstrap置信区间才勉强能看。另外很多目录里的震级只保留了一位小数这相当于人为引入了分档宽度ΔMAki公式里补上ΔM/2这一项就是这个原因。4.4 数据清洗的其他细节深度为负值的记录虽然在USGS里不多但会影响定位筛选。主震后的余震会让短时窗内的震级-频度关系严重偏离背景值研究背景地震活动时需要做余震去丛。同一区域、不同时间段的数据完整性可能不同比如早期全球台网稀疏M5以下事件漏记严重。这些都是我会在正式分析前跑一遍检查的步骤。虽然看起来琐碎但少了任何一步都可能让最后的统计结果出现偏差。最后再分享一个我后来做数据检查时很依赖的小技巧。拿到一份新地震目录后先用震级对nst台站数画一张散点图通常震级低到一定程度时台站数会明显坠落。那个坠落位置就是该台网在这片区域的完整检测边界。用这个办法你甚至不需要专门计算M_c就能对目录完整性有一个直观把握。我每次拿到新数据都会先画这张图再决定从哪里开始算b值。这个习惯帮我少走了很多弯路。
震级、b值与Python：从地震目录到Gutenberg-Richter分析的完整指南
你还记得那年的东日本大地震吗第一波新闻里写的是8.8级第二天多家机构修定为9.0级。这不是简单的四舍五入而是随着更多台站数据汇入震级被重新测定了。这个我们从小就听熟的“几级”其实是地球物理学里最核心的概念之一magnitude震级。我一直觉得震级是地震科学里最容易被误解、也最值得亲手算一遍的东西。这篇博文想把震级从一个新闻名词拆成一套你可以自己复现的数据分析流程——从公开地震目录下载、Python数据处理到计算b值、画震级-频度曲线。想真正把“震级”弄明白的人不管是做数据处理的学生还是对地震科普有硬核兴趣的爱好者这篇内容应该能帮你在拿到目录后不再一头雾水。1. 震级到底是什么——先从那个“几级”说起1.1 为什么必须用对数标尺震级的定义绕不开对数。这不是数学家为了显得高深才这么干的而是地震释放能量的跨度实在太大。人感觉不到的微震能量大概只有几万焦耳1960年智利Mw 9.5大地震释放的能量大约在10的18次方焦耳量级两者相差十几个数量级。如果按“差多少焦耳”来画图小震全都挤在零点附近什么也看不出来。这就像你要比较一个人的存款是10万还是1000亿你会说“差5个数量级”而不是关心具体差多少块。里克特在1935年定义里氏震级时采用的办法就是取对数把距离震中100公里处、标准伍德-安德森地震仪记录到的最大振幅单位微米取以10为底的对数再减去一个校准常数。换句话说振幅每放大10倍震级就加1。这么做的直接好处是地震目录里的震级数字变得极其直观。6级、7级、8级数字相差不大但它们对应的振幅差了100倍能量差了大约31600倍。正是因为取了对数人类才能在一张纸上同时写下微震和大震的观测记录。这个数学选择是后面所有统计分析和历史目录比对的前提。1.2 ML、Ms、Mb、Mw同一场地震为什么有好几个数字实际操作中你会发现一场地震往往有多个震级数字这不是新闻乱报而是不同测量方式算出来的结果。里氏震级ML是最早的定义但它依赖的是南加州浅源地方震的观测所以后来人们又根据不同的地震波震相发展出了几种新标度。体波震级mb使用1秒左右短周期P波深源地震也能记到早期全球速报常用但它和地震真正释放能量之间的关系不够直接。面波震级Ms使用20秒左右周期的面波对浅源大震更稳定但震级超过一定范围后会“饱和”。矩震级Mw由地震矩M0推导而来M0等于剪切模量乘以断层破裂面积再乘以平均滑动量。公式是Mw (2/3)log10(M0) - 6.07M0单位为N·m。为什么现在新闻里的大震级别都以Mw为准因为矩震级的物理意义明确它直接对应断层破裂的规模和位错量而且不会因为地震波高频部分饱和而封顶。1960年智利地震早期用面波震级测定只有8.3级后来用矩震级测定为9.5级差别就是这么来的。震级类型基于的观测优点局限ML短周期地方震振幅定义早、资料多仅适合浅源地方震易饱和mb1秒周期体波P波可测深源地震与释放能量相关性较弱Ms20秒周期面波适合远震浅源超过8.0到8.5会饱和Mw地震矩M0物理意义明确、不饱和需要反演震源机制计算慢1.3 一份目录里到底存的是哪种震级我见过太多人从USGS下载一份CSV就开始统计画图算b值结果画出来的曲线在某处莫名其妙拐弯。问了一句“你用的哪个震级字段”对方才意识到问题CSV里有一个mag字段还有一个magType字段后者标明了这一条记录的震级标度。同一份目录里可能是mb、ML、Ms、Mw混着来的。所以拿到任意一份地震目录后的第一件事不是急着画图而是先看元数据。具体到这个csv里magType都有哪些取值各占多少比例。如果这个步骤不做后面所有分析都建立在沙子上。2. 数据获取拉一份真实的地震目录2.1 公开数据源怎么选全球公开的地震目录有好几套最常用的是USGS的ANSS综合目录覆盖全球API开放字段完整适合大多数分析场景。如果你做更底层的科研ISC国际地震中心目录数据质量和完整性更好但下载需要注册和邮件确认门槛略高。GCMT项目专门提供矩心矩张量解只有Mw但事件数量有限不适合做小震级统计分析。对于初学者来说USGS是最合适的起点。我平时做快速验证也直接用USGS它的FDSN事件查询接口可以直接用URL拼参数返回CSV用pandas读进来就行不需要任何特殊权限。如果想了解中国及周边区域的地震活动中国地震台网中心也提供公开目录但字段规范和数据格式跟USGS不太一样后期清洗会稍微多花点时间。2.2 用Python拉取USGS目录USGS的接口格式很友好核心路径是earthquake.usgs.gov/fdsnws/event/1/query配合formatcsv和起止时间、震级范围就能拼出下载URL。我用pandas的read_csv直接读远程URL几秒钟就能把多年目录拉下来。import pandas as pd start 2015-01-01 end 2023-12-31 url ( https://earthquake.usgs.gov/fdsnws/event/1/query f?formatcsvstarttime{start}endtime{end} minmagnitude5.0 ) df pd.read_csv(url) print(df.shape) print(df[magType].value_counts())初次运行可能会因为网络或限流失败建议加一层重试或者按年份循环下载再用pd.concat合并。另一个容易忽略的点是时间USGS接口里的时间是UTC如果你按当地时间划分研究窗口务必先加上时区处理否则对比历史报告时会错位。2.3 字段解读与数据清洗拿到表格后先花点时间理解字段。time、latitude、longitude、depth这些都是基本定位信息。depth单位是千米偶尔会出现负值表示震源深度在海平面以上多见于海洋区域的特殊算法结果。mag是震级数值magType紧跟其后magSource表示测定机构。nst是参与定位的台站数gap是台站方位角空隙这两个字段直接反映定位质量。nst小于3的记录定位可靠性很低建议删除。gap大于180说明台站几乎只分布在震中一侧震中位置可能偏差很大。mag为空或者为NaN的记录直接删掉别含糊。df df.dropna(subset[mag, latitude, longitude]) df df[df[nst] 3] if nst in df.columns else df df df[df[gap] 180] if gap in df.columns else df另外USGS目录偶尔会有重复事件我遇到过同一条地震同时被几个机构上报的情况。简单去重可以按time、lat、lon、mag四列做drop_duplicates更稳妥的办法是根据事件ID字段去重。清洗完成后最好把结果存成parquet或者csv后面反复读取会快很多。3. Python实操震级分布与b值计算3.1 Gutenberg-Richter关系震级-频度的核心规律地震不是随机事件它有一个著名的统计规律叫古登堡-里克特关系log10(N) a - b·M。这里的N是震级大于等于M的地震数。意思是说震级每降低1级地震数量大约变成原来的10倍。b值通常接近1但不同区域的b值有差异反映了地壳应力状态和断层类型的区别。我第一次亲手拟合这条关系时才真正理解为什么有人说地震具有“自相似性”。从M2的微震到M8的大震只要数据足够完整它们之间的数量比例关系基本稳定。这也是为什么在地震危险性评估中人们常用b值反推未来大震的发生概率。3.2 计算b值Aki最大似然估计与Bootstrap误差计算b值最常用的方法是Aki在1965年提出的最大似然估计公式非常简单b log10(e) / (M_mean - M_c ΔM/2)。其中M_mean是样本平均震级M_c是目录完整性震级ΔM是震级分档宽度。麻烦的是M_c到底取多少这决定了计算结果的可靠性。实操中我会先画累积震级-频度图找到偏离直线的那一段把M_c取在曲线开始变直线的位置。然后对M_c以上的数据进行Aki估计再用Bootstrap重采样1000次估计误差。import numpy as np mag df[mag].values mag mag[mag M_c] b np.log10(np.e) / (mag.mean() - M_c 0.05) print(fb value: {b:.3f}) rng np.random.default_rng(42) boot_b [] for _ in range(1000): sample rng.choice(mag, sizelen(mag), replaceTrue) boot_b.append(np.log10(np.e) / (sample.mean() - M_c 0.05)) ci np.percentile(boot_b, [2.5, 97.5]) print(f95% CI: {ci[0]:.3f} - {ci[1]:.3f})注意M_c取高一点虽然会损失样本量但能保证进入统计的都是被台网完整记录的事件。如果把不完整的低震级数据也算进来b值会被明显拉低。3.3 能量对比1级之差差了31倍还是32倍这是科普里最常见的误区。震级每差1级振幅差10倍能量差大约31.6倍而不是人们常说的“10倍”。能量和震级的经验关系是log10(E) 1.5·M 4.8E单位是焦耳。这意味着M从6.0涨到7.0能量约增大10^1.5倍也就是31.62倍。把这个关系落到具体数字上会更直观一次6级地震释放的能量约等于6.3×10^13焦耳7级则是2.0×10^15焦耳左右。8级地震释放的能量相当于大约1000次6级地震或者31.6次7级地震。公众听到“8级是6级的100倍”大概率是指这个量级关系。我在目录里随便抽了几个大型地震事件算了算映象最深的是大震与小震之间的能量差远超直觉。这也是为什么一次7级地震造成的破坏可能比几十次6级地震加在一起还要大因为能量是按指数增长的。3.4 可视化画一张震级-频度图数据有了b值算了最后用图把所有内容串起来。画直方图和累积频度曲线时竖轴用对数刻度这样G-R关系的直线形态才看得清楚。import matplotlib.pyplot as plt bins np.arange(mag.min(), mag.max() 0.1, 0.1) counts, edges np.histogram(mag, binsbins) centers (edges[:-1] edges[1:]) / 2 plt.figure(figsize(8, 5)) plt.bar(centers, counts, width0.1, alpha0.6, labelhistogram) plt.yscale(log) plt.xlabel(Magnitude) plt.ylabel(Number of events) plt.title(Magnitude-Frequency Distribution) plt.legend() plt.grid(alpha0.3) plt.show()画完之后理想情况下你会看到尾部拖得很长的高震级区以及直方图的整体包络基本是一条直线。如果这条直线在中段有转折先别急着解释地质原因回头检查是不是M_c取错或者震级标度混用。4. 踩过的坑与排查技巧4.1 震级标度混用造成的假象最典型的坑就是不分青红皂白地把mb、ML、Ms、Mw当成同一个量纲去统计。我试过把某一年全球目录里的所有震级直接拿来画直方图结果在M5附近出现了一个明显的“鼓包”一开始还以为是区域地震活动异常后来查了magType才发现那堆数据大多是mb而其他标度在这个区间的分布并不一样。建议拿到数据后先跑一句value_counts看分布。如果混用情况不严重可以只保留一种标度分析如果必须合并要给出明确的转换关系并且说明这种转换本身有误差。比如用经验公式把mb转换成Mw这只是一个近似不适合用于精度要求极高的研究。4.2 震级饱和为什么大震被低估了震级饱和是一个很隐蔽的问题。里氏震级ML依赖短周期波振幅当断层破裂尺度非常大时短周期波的振幅增长会放缓导致ML在大约7级以后“封顶”。面波震级Ms在8.0到8.5之间也会出现饱和。这就是为什么一些历史大地震在旧目录里震级偏低。处理办法是涉及大震研究时只用Mw历史目录中如果只有Ms要查阅当时的测定方法必要时查阅已发表的震级转换文献。我踩过这个坑之后现在只要看到目录里有一堆8.0以上的Ms记录第一反应就是这批数据需要谨慎处理。4.3 震级完整性与小样本问题b值估计的另一个致命问题是样本量和完整性。Aki公式依赖平均值而平均值对样本量有最低要求。超过一定震级的事件本来就少如果M_c取得太低混入大量检测不完备的小震b值会偏低如果M_c取得太高样本量不足误差又太大。实操经验是样本少于50个时基本别谈b值。我个人的习惯是至少取100个事件再算这样Bootstrap置信区间才勉强能看。另外很多目录里的震级只保留了一位小数这相当于人为引入了分档宽度ΔMAki公式里补上ΔM/2这一项就是这个原因。4.4 数据清洗的其他细节深度为负值的记录虽然在USGS里不多但会影响定位筛选。主震后的余震会让短时窗内的震级-频度关系严重偏离背景值研究背景地震活动时需要做余震去丛。同一区域、不同时间段的数据完整性可能不同比如早期全球台网稀疏M5以下事件漏记严重。这些都是我会在正式分析前跑一遍检查的步骤。虽然看起来琐碎但少了任何一步都可能让最后的统计结果出现偏差。最后再分享一个我后来做数据检查时很依赖的小技巧。拿到一份新地震目录后先用震级对nst台站数画一张散点图通常震级低到一定程度时台站数会明显坠落。那个坠落位置就是该台网在这片区域的完整检测边界。用这个办法你甚至不需要专门计算M_c就能对目录完整性有一个直观把握。我每次拿到新数据都会先画这张图再决定从哪里开始算b值。这个习惯帮我少走了很多弯路。