引用本文: 陈蒙, 王华. 2022. 地震动强度参数估计的可解释性与不确定度机器学习模型. 地球物理学报, 65(9): 3386-3404, doi: 10.6038/cjg2022P0428 通讯作者: 王华, 男, 1982年生, 教授, 国家特聘青年专家.2005年本科毕业于中国石油大学(华东)勘查技术与工程专业, 2012年博士毕业于中国石油大学(北京), 为中国石油大学(北京)与美国麻省理工学院联合培养博士, 美国麻省理工学院地球大气行星科学系博士后.主要从事诱发地震、测井、人工智能等研究.E-mail: huawang@uestc.edu.cn

中图分类号:

准确预测地震动强度参数(峰值加速度PGA、峰值速度PGV等)对于震后应急和地震危险性概率分析至关重要.作为地震动强度参数预测的新手段, 机器学习算法具有优势, 但也存在可解释性差和难给出预测结果不确定度的问题.本文提出采用自然梯度提升(NGBoost)算法在预测结果的同时提供其不确定度, 并结合SHAP值解释机器学习模型.基于NGA-WEST2强震动数据库, 本文训练出了适合预测活跃构造区地壳地震的PGA和PGV概率密度分布的机器学习模型.测试集数据PGA和PGV的预测值与真实值的相关系数可达0.972和0.984, 并可给出预测结果的合理概率密度分布.通过SHAP值, 我们从数据角度弄清了各输入特征(矩震级 M _W 、Joyner-Boore断层距 R _jb 、地下30 m平均S波速度 V _S30 、滑动角Rake、断层倾角Dip、断层顶部深度 Z _TOR 和 V _S 达到2.5 km·s ^-1 时的深度 Z _2.5 )对机器学习模型预测结果的影响机理.SHAP值显示, 基于NGBoost算法的机器学习模型的预测方式基本与物理原理相符, 说明了机器学习模型的合理性.SHAP值还揭示出一些以往研究忽视的现象: (1)对于活跃构造区地壳地震, 破裂深度较浅( Z _TOR ＜~5 km)时, Z _TOR 的SHAP值低于破裂深度较深( Z _TOR ＞~5 km)时的值, 表明浅部破裂可能主要受速度强化控制, 地震动强度较弱.并且 Z _TOR 的SHAP值随 Z _TOR 值增大而减小, 表明地震动强度可能还受破裂深度变化引起的几何衰减变化影响; (2)破裂深度较深时, Z _TOR 的SHAP值随 Z _TOR 值增大而增大, 表明深部破裂的地震动强度可能受和破裂深度变化相关的应力降或品质因子 Q 的变化影响; (3) Z _2.5 较小( Z _2.5 ＜~1 km)时, Z _2.5 的SHAP值的变化规律对于PGA和PGV预测是相反的, 表明加速度和速度频率不同, 受浅层沉积物厚度变化引起的共振频率变化影响不同.

强地面运动 可解释机器学习 自然梯度提升 Abstract:

Ground motion parameters prediction (peak ground acceleration, PGA and peak ground velocity, PGV) is of the essence in rescue efforts aftermath of earthquakes and seismic hazard analysis. As new developed approaches for predicting ground motion parameters, machine learning algorithms do have some advantages, but also have difficulties in estimating predictive uncertainties and interpreting machine learning models. In this study, we use the natural gradient boosting (NGBoost) algorithm to evaluate predictive uncertainties, and use the SHAP values to interpret trained machine learning models. Based on NGA-WEST2 database, we trained machine learning models which are suitable for predicting PGA and PGV in active tectonic regions. The correlation coefficients between the predicted PGA and PGV and observations in testing dataset reach up to 0.972 and 0.984, respectively. The trained machine learning models also provide reasonable probability distributions of predicted values. With the computed SHAP values, we figured out the influence of the input features (moment magnitude, M _W ; Joyner-Boore distance, R _jb ; V _S over top 30 m, V _S30 ; rake angle, Rake; dip angle, Dip; depth to the top of fault, Z _TOR ; and depth to V _S =2.5 km·s ^-1 , Z _2.5 ) on the outputs of machine learning models. According to the SHAP values of input parameters, we find that the predicting mechanisms of trained machine learning models make sense in physics which illustrates the machine learning models are reasonable. In addition, SHAP values also revealed some facts which are ignored in previous studies: (1) The SHAP values of Z _TOR in general are low when the depths of rupture planes are shallow ( Z _TOR < ~5 km), indicating that the ground motions from ruptures in the shallow part of crust may be controlled by velocity strengthening and are systematically weaker.The SHAP values of Z _TOR decrease with Z _TOR , which indicate ground motions from ruptures in the shallow part of crust may also be affected by depth-varying geometrical attenuation; (2) When depths of ruptures are large ( Z _TOR > ~5 km), the SHAP values of Z _TOR increase with Z _TOR , which indicate ground motions from ruptures in the deep part of crust may highly be impacted by depth-varying stress drops or quality factors ( Q ); (3) The variations of SHAP values of Z _2.5 are different for predictions of PGA and PGV when Z _2.5 are low ( Z _2.5 < ~1 km), which indicate impacts of differences in resonance frequencies of sediments caused by variations of Z _2.5 on PGA and PGV are different, since the frequencies of velocity and acceleration are different.

Key words: Strong ground motion Explainable machine learning Natural gradient boosting Geological hazard Figure 3.

For prediction of PGA and PGV, the variations of average scores of validation dataset with the number of estimators for different learning rate (0.1 and 0.01) and max depth of estimator (3~10). All the input features ( M _W + R _jb + V _S30 + Z _TOR +Rake+Dip+ Z _2.5 ) are used when training these models for prediction of PGA and PGV

Campbell K W, Bozorgnia Y. 2013. NGA-West2 Campbell-Bozorgnia ground motion model for the horizontal components of PGA, PGV, and 5%-damped elastic pseudo-acceleration response spectra for periods ranging from 0.01 to 10 sec. Berkeley, CA: Pacific Earthquake Engineering Research Center, University of California. Figure 1.

Data distribution of strong ground motion used for machine learning study. The data is filtered from NGA-WEST2 database

Figure 2.

Block diagram of NGBoost algorithm for ground motion parameters prediction

Figure 3.

For prediction of PGA and PGV, the variations of average scores of validation dataset with the number of estimators for different learning rate (0.1 and 0.01) and max depth of estimator (3~10). All the input features ( M _W + R _jb + V _S30 + Z _TOR +Rake+Dip+ Z _2.5 ) are used when training these models for prediction of PGA and PGV

Figure 4.

Comparisons between the machine learning model for predicted and observed PGA and PGV. The training (80%) and testing (20%) dataset are randomly selected from NGA-WEST2 strong motion dataset

Figure 5.

Histograms of ratios of machine learning model for predicted and observed PGA and PGV. The training (80%) and testing (20%) dataset are randomly selected from NGA-WEST2 strong motion dataset

Figure 6.

The predicted probability density distributions by the proposed NGBoost model. Dataset are the 2004 Parkfield M _W 6.0 earthquake, 2009 San Bernardino M _W 4.45 earthquake and 2007 Chuetsu-oki M _W 6.8 earthquake

Figure 7.

The computed SHAP values of different input features for three strong motion records (the 2007 Chuetsu-oki M _W 6.8 earthquake with distances at 20.4, 107.1 and 279.3 km)

Figure 8.

For prediction of PGA and PGV, the SHAP summary plot of all input features

Figure 9.

For PGA prediction, the SHAP dependence plots of different input features which show how the outputs of model depend on input features

Figure 10.

For PGV prediction, the SHAP dependence plots of different input features, which show how the outputs of model depend on input features

Figure A1.

For different times (1933—2000, 2001—2005 and 2006—2011) and moment magnitudes, the maximum Joyner-Boore epicentral distances used for selecting seismic data

Figure A2.

For the 2004 Parkfield M _W 6.0 earthquake, 2009 SanBernardino M _W 4.45 earthquake and 2007 Chuetsu-oki M _W 6.8 earthquake, the comparisons of prediction results using M 5′ model trees and classification and regression trees as base learners

Figure A3.

For nearfield (< 2 km) and farfield (> 10 km) PGA and PGV prediction, the SHAP dependence plots of M _W , which show how predicted PGA and PGV depend on M _W