龙门山断裂带发震时刻间隔的SARIMA模型分析预测

叶秋吟　薛源　曾文骏　杨静　饶敏

摘要： 基于龙门山断裂带2012年1月—2021年9月MS2.5及以上地震目录数据，按震级分组建立发震时刻间隔序列，然后对各序列进行平稳白噪声检验，自相关、偏相关性分析，使用SARIMA模型对其进行短、中、长周期拟合及预测。通过分析模型拟合效果，得到不同震级序列的最优模型参数及相应周期数值。其中，序列MS≥2.5及序列MS≥3.0各模型调整R2均在0.86以上，最高达0.911;两序列对应模型的短时预测表现良好，预测RMSE分别为10.686及8.800。模型预测结果表明，龙门山断裂带后续发震时刻间隔总体趋势平稳，序列MS≥3.0预测结果趋势有微弱增长，一段时间内龙门山断裂带MS≥3.0地震发震次数将略微下降，地震活动性降低。该分析结果可为地震活动研究工作提供科学依据，其分析方法及过程为地震发震时间的分析预测提供了有效可行的途径。

关键词： 龙门山断裂带; 发震时刻间隔; SARIMA模型; 时间序列分析

中图分类号： P315;O21文献标志码：A 文章编号： 1000－0844（2023）04-0991-10

DOI：10.20000/j.1000－0844.20220421002

Analysis and prediction of the occurrence time interval of earthquakes along the Longmenshan fault zone using the SARIMA model

YE QiuyinXUE YuanZENG WenjunYANG JingRAO Min1，2

Abstract：  Different sequences of earthquake occurrence time intervals were first established in accordance with the magnitude based on the catalog data of MS≥2.5 earthquakes along the Longmenshan fault zone from January 2012 to September 2021. Testing of stationary white noise and the analysis of autocorrelation and partial correlation were then conducted， and short－， medium－， and long－period fitting and prediction were performed for each sequence by using the SARIMA model. Finally， the optimal model parameters and corresponding period values of different magnitude sequences were obtained by analyzing the model fitting effect. Among them， the model－adjusted R2of MS≥2.5 and MS≥3.0 sequences are more than 0.86， with a maximum value of 0.911; the short－term prediction effects of corresponding models are good and the values of predicted RMSE are 10.686 and 8.800， respectively. The prediction results of the models show that the overall trend of the occurrence time interval of subsequent earthquakes along the Longmenshan fault zone is stable. The trend of the MS≥3.0 sequence shows a slight increase; that is， the number of MS≥3.0 earthquakes along the Longmenshan fault zone will decrease slightly in a given time period， and the seismicity will decrease. Analysis results can provide a scientific basis for the study of seismic activity， and the analysis method and process can introduce an effective and feasible way for the analysis and prediction of earthquake occurrence times.

Keywords： Longmenshan fault zone; occurrence time interval; SARIMA model; time series analysis

0 引言

龍门山断裂带位于青藏高原东缘、四川盆地西北边缘，以北东－南西向展布，北端抵达秦岭断裂带，南端止于鲜水河—小江断裂带，长约500 km，宽约50 km，主要由后山断裂、中央断裂、前山断裂和山前隐伏断裂4条主干断裂及其控制的逆冲推覆体组成［1－2］。该断裂带范围内分布有都江堰等数十个重要城镇，主山前边界大断裂东南方向靠近人口密集的成都平原地区及双城经济圈，地理位置特殊。探究龙门山断裂带的发震时间规律，有助于为该断裂带地震活动波及地区的地震预测工作提供科学的地震管理决策依据。

如今地震研究是跨学科的综合性领域，对于地震时间序列数据，很多研究者利用统计、概率的方法和模型对地震发生时间、位置及危险性等重要信息进行预测。部分学者着眼于对震级时间序列进行探究［3－4］，如和宏伟［3］利用门限自回归等统计方法处理震级序列。另外也不乏对发震时刻间隔序列的研究，其中多为分析地震时间间隔的概率分布，如王璇［5］，Sadeghian［6］，Liu［7］，Khan等［8］及Sil等［9］将其与统计学中的分布相比较。也有对时间间隔开展预测研究工作，如郭安宁等［10］利用倍周期叠加黄金分割法对西北地区的MS8.0地震进行的验证性预测，并发现有很好的吻合性。Mariani等［11］提出将独立O－Uh过程叠加的随机微分方程应用于地震研究中，并用该模型拟合了南美MS8.0以下的地震序列。此外，还有学者拓展预测途径，如陈翰林等［12］，游佳等［13］，朱海宁［14］及Asencio－Cortes等［15］运用SVM、LSTM及回归算法对地震时间进行预测。季节性自回归综合移动平均（Seasonal Autoregressive Moving Average，SARMA）模型既能用于预测自然现象如降雨量［16－17］、河流流量［18］等，也可以对社会经济现象如交通运输客流量、贸易量等方面进行预测，如学者于琛等［19］，易靖韬等［20］，孙越等［21］及Milenkovic等［22］所做研究。在地震预测领域，有学者如张继艳等［23］，李磊等［24］及Saqib等［25］运用求和自回归移动平均（Autoregressive Integrated Moving Average，ARIMA）模型对地震前兆异常进行识别，通过探测震前电离层异常预报地震。但鲜有将ARIMA或季节求和自回归移动平均（Seasonal Autoregressive Integrated Moving Average，SARIMA）模型用以拟合发震时刻的数据。

本文选取中国四川地区龙门山断裂带2012年至2021年9月的地震事件作为研究对象，从发震时刻间隔的角度，利用SARIMA模型分别按震级对序列进行拟合，以此预测发震时刻间隔的发展趋势及下一次发震的时间。

1 SARIMA模型

1.1 模型理论

ARIMA模型于20世纪70年代被Box及Jenkins首次提出，全称为差分自回归移动平均模型，记作ARIMA（p，d，q），其中p和q分别表示自回归和移动平均阶数;d为趋势差分阶数，基本模型结构为：

ARIMA模型可对无季节效应的非平稳时间序列建模，而现实生活中很多时间序列存在着一定周期性，对于既含季节效应又含长期趋势效应且两者存在复杂交互影响关系的序列，可采用SARIMA模型。SARIMA模型又称乘积季节模型，最初由Wang等［27］提出通用表达式，记为ARIMA（p，d，q）×（P，D，Q）s，其中P和Q分别为季节性自回归和季节性移动平均阶数;D为季节差分的阶数;s为季节周期数。作为ARIMA模型的扩展，该模型在原模型基础上对序列进行了季节差分提取季节效应信息，以达到对季节性时间序列建模的目的。

1.2 预测模型构建

针对小样本数据的预测，传统时间序列预测模型可以较充分地提取数据所含信息，同时也不易出现过拟合致使预测结果不理想的情况。本文在分析数据结构和特点后，选择利用SARIMA模型对地震时间序列进行建模。作为传统统计预测模型，SARIMA模型对小样本数据的拟合程度较好、预测精度高、训练速度快且能够能较准确反映数据动态变化，在捕捉序列长期趋势变化的同时提取其周期性波动的信息［16］。

对按震级分组后的序列进行差分，得到发震时刻间隔序列，得到不同参数的SARIMA模型：模型拟合值序列 Y′t及预测值序列Yt，并以此估计下一次发震时刻。预测时选取了建模效果最好的长中短周期，表示为sa、sb、sc，以保证模型对短期与长期时序信息的有效提取，提高各期预测的精度。模型分析流程图如下：

2 龙门山断裂带发震时间间隔研究

2.1 数据与样本

本文采用2012年1月9日—2021年9月24日中国地震台网（http：//www.ceic.ac.cn）地震目录数据。地震事件序列有发震时间、经纬度、震级等因素，具有波动性、非线性的特点，因此通过数据预处理将其确定为特定空间和一定震级下的数据，即龙门山断裂带经纬度范围内（29.5°～33.5°N，102°～107°E） MS2.4及以上的时间序列。对数据进一步筛选，考虑其发震断裂位置得实验样本数据437条，发震时刻2012 年1月9日—2021年8月，MS∈（2.5，7）。

对选出样本的发震时刻差分处理得发震时刻间隔序列，该样本具有如下特征：（1）样本量相对较小。由于设备和技术等原因，早期监测数据不够完备，由此计算时间间隔数值准确性将显著降低，因此选取近10年地震序列，有效避免该问题，但同时也减少了数据量;（2）时间跨度大。数据时间从2012年1月—2021年8月，包含龙门山断裂带近10年间地震事件;（3）非平稳、强波动性。数据的影响因素较多且可能有复杂交互影响，且地震序列本身就有复杂的周期性特点。

另外，序列中记录的第58次发震为2013年雅安芦山MS7.0大地震，随后伴有非常紧密的、大小不一的数次余震，由此造成监测数据中第58次至203次发震前后时间间隔短，一段时间内序列趋近零值。以Xt表示原序列，序列如图2所示。

2.2 模型确定

分析数据所含特征，對序列分解，提取其趋势性、季节性及残差序列如图3所示，表明该序列有一定趋势性和季节性，初步推测最小周期为5天。

SARIMA模型需要序列满足平稳或差分后平稳的前提条件。对Xt进行平稳性检验，确定差分次数d=2，差分后序列记为Xd。观察Xd序列图（图4），该序列似乎不具有明显趋势性。进一步采用ADF检验进行判断，检验原理即通过考察序列中是否含有单位根来判断该序列是否平稳：若序列平稳，则不存在单位根。结果得Xd的ADF检验统计量为-15.463 4，p值小于0.05，表明该序列已平稳。

按9∶1划分Xd数据为模型训练集、测试集，依据数据的ACF和PACF为模型定阶。ACF用于衡量该序列的当前观测值与其滞后项间相关程度;PACF用于衡量除去之前滞后已解释的影响后，当前观测值与其滞后项间相关程度。考察两者的截尾、拖尾性质，可对 SARIMA模型季节性及非季节性的自相关、偏自相关阶数给出参考。由图5中ACF及PACF图初步确定模型非季节性阶数。PACF图可看做拖尾或9、10阶截尾。

结合模型调整的R2选择模型最优阶数，不同阶ARIMA模型拟合效果如表1所列，确定模型p=9，q=1。

再考察对Xd进行1阶s步差分后，滞后对应周期整数倍阶的ACF及PACF。以s=22为例，如图6，滞后22阶的ACF及PACF系数均在2倍标准差范围外，且随后逐渐收敛。初步假设P=1，Q=1，其他周期情况大致相同。进一步对该假设模型进行参数检验和LB检验。参数检验用于检验模型是否成立，即模型所有参数是否均显著不为0;LB检验则用于检验模型拟合后的残差序列是否为白噪声序列，若为白噪声（p≥0.05），则表示残差序列中无更多可提取信息，即模型已充分提取原数据蕴含的信息。

检验结果如表2所列，模型中含不显著参数（1阶季节自回归系数）。经比较，确定模型最佳季节性阶数P=0，Q=1。

2.3 模型评价

研究所用地震序列有大小周期复杂交叠的情况，在分析数据周期性时，针对短、中、长期的预测各选取不同季节期数。对序列Xd做不同跨度的中心移动平均，所得序列记为Xn，n∈（4，45）为跨度。比较各Xn与Xd的差异，观察各周期对季节性消除程度，拟合s∈（4，45）的模型。为了更好地比较不同模型间预测效果，本文采用常见评价指标中的均方根误差（RMSE）来代表模型的拟合误差：

（3）

式中：xt为真实值;yt为拟合值。对于通过参数检验和LB检验的模型，可分析其调整R2及拟合RMSE来判定其拟合效果。AIC准则（最小信息化准则）也常用于衡量统计模型复杂度及拟合优良性，在鼓励数据拟合的优良性的同时尽量避免过度拟合，其AIC值越小表明模型在该准则下性能越优。本文模型选取的原则为：优先采用调整R2数值高的模型，在该值差距不显著的情况下，考察各模型拟合RMSE，同时以模型AIC值作为辅助参考，确定最优模型。表3提供了部分已通过参数检验及LB检验的不同周期模型拟合情况。

由此可以确定地震数据短、中、长期预测模型分别为ARIMA（9，2，1）×（0，1，1）6，ARIMA（9，2，1） ×（0，1，1）22，ARIMA（9，2，1）×（0，1，1）42，模型拟合效果如图7所示。

3 震级分组

对原序列按震级分组，得序列MS≥3.0共390条数据，序列MS≥4.0共58条数据。对其差分后得时刻间隔序列，两序列分别以Xd3、Xd4表示。按上述步骤对各序列进行分析建模，考虑对Xd3拟合模型ARIMA（7，2，1）×（0，1，1）s，对Xd4拟合模型ARIMA（3，2，1）×（0，1，0）s。各组部分不同周期的模型对比如表4及表5所列。

按上述原则，综合各模型指标，确定以s=7，s=16，s=32对Xd3进行短、中及长期预测，预测模型为ARIMA（7，2，1）×（0，1，1）7，ARIMA（7，2，1）×（0，1，1）16，ARIMA（7，2，1）×（0，1，1）32。对于Xd4，由于MS4.0及以上地震序列数据量较少，难以捕捉长周期变化趋势，相应的长周期模型均参数不显著或拟合效果不佳;同时，短周期模型拟合效果亦较差，序列Xd4仅表现出较明显的中周期性质。因此以s=11对Xd4进行中期预测，预测模型为ARIMA（3，2，0）×（0，1，0）11。各序列不同周期拟合效果如图8、图9所示。

利用最优模型预测各序列测试集数据，将预测结果与真实值对比，如图10、11所示。对于Xd及Xd3，其测试集中数据较多，属于长时预测，对应模型能捕捉其周期性变化，但预测值整体较真实值偏高;对Xd4进行短时预测，ARIMA（3，2，0）×（0，1，0）11能很好地预测近几期该序列发展趋势，且未表现出预测结果偏高的现象。各模型预测RMSE如表6所列。

进一步计算序列Xd及Xd3真实值及预测值的平均水平，两序列对应各周期模型的平均误差最低

为17.578天，最高为32.967天，参考目前发震时刻预测的普遍准确度，认为预测值平均水平的偏高程度在能接受的范围内。修正序列的预测值，对序列Xd及Xd3模型预测进行修正，以序列Xd周期为22及Xd3周期为32对应模型为例，预测值减去对应平均误差后效果如图12，表明模型趋势及周期性预测性能较好。

由于序列Xd和Xd3的预测步长大于40，可能会降低预测精度，因此我们考虑对序列Xd进行15步预测，对数据量相对较小的序列Xd3进行5步预测。选取上述两模型为例，预测结果如图13所示，预测RMSE分别为10.686 0、8.800 9。图中序列Xd所预测发震间隔无明显趋势，而Xd3的预测序列则稍有增长趋势。

补充截止2022年3月6日龙门山断裂带MS≥4.0的地震序列新增数据共2条于数据集中，利用ARIMA（3，2，0）×（0，1，0）11模型对最近发生MS≥4.0的地震时刻间隔进行外推预测。對比预测值与真实值，如图14所示。该模型对最新数据的预测RSME为55.711 2，发现其仍很好地捕捉到其趋势，可在此基础上对模型精度继续改进。

4 结论

本文提出以分析龙门山断裂带发震时刻间隔的角度，通过时间序列分析方法，利用SARIMA模型对数据进行建模，分析地震时间序列隐含的信息，并由此预测下一次的发震时间。通过对模型拟合度及拟合效果的分析，为该时刻间隔序列选取适合的模型，并针对不同周期预测的需求，选择适合的模型。根据结果，有以下几点结论：

（1）SARIMA模型适用于地震时间间隔序列的分析预测，各序列最优模型调整R2均在0.86以上，最高达0.911。

（2）长时预测时，序列Xd及Xd3对应各周期模型均有预测值较真实值偏高的情况，对偏低数值（时间间隔趋近0）的预测性能较为欠缺。

（3）短时预测时，序列Xd及Xd3最优模型预测效果均良好，能较好预测出序列周期性。序列MS≥3.0预测结果趋势有微弱增长，一段时间内龙门山断裂带MS≥3.0地震发震次数将略微下降，地震活动性降低;序列Xd4中周期性质较明显，对应模型在短时预测中能很好地捕捉其发展变化，预测趋势与真实情况贴合，模型预测精度可进一步提升。

参考文献（References）

［1］许康生，王维欢，李英，等.2016—2017年龙门山地区微震活动的演化：基于对微震目录的模糊聚类分析［J］.地球物理学进展，2022，37（1）：69－77.

XU Kangsheng，WANG Weihuan，LI Ying，et al.Evolution of microseismic activity in Longmenshan area from 2016 to 2017：based on fuzzy clustering analysis of microseismic catalogue［J］.Progress in Geophysics，2022，37（1）：69－77.

［2］杨宜海，张雪梅，花茜，等.龙门山断裂带的分段性特征：来自密集震源机制解的约束［J］.地球物理学报，2021，64（4）：1181－1205.

YANG Yihai，ZHANG Xuemei，HUA Qian，et al.Segmentation characteristics of the Longmenshan fault：constrained from dense focal mechanism data［J］.Chinese Journal of Geophysics，2021，64（4）：1181－1205.

［3］和宏伟.门限自回归模型在地震预报中的应用初探［J］.地震研究，1992，15（2）：154－161.

HE Hongwei.Study of the threshold autoregressive model and its application in earthquake prediction［J］.Journal of Seismological Research，1992，15（2）：154－161.

［4］KUNDU S，OPRIS A，YUKUTAKE Y，et al.Extracting correlations in earthquake time series using visibility graph analysis［J］.Frontiers in Physics，2021，9：656310.

［5］王璇.四川地区M5.0及以上地震时间间隔统计分析与概率预测［J］.四川地震，2021（4）：40－45.

WANG Xuan.Statistical analysis and probability prediction of earthquake time interval for earthquakes greater than M5.0 in Sichuan［J］.Earthquake Research in Sichuan，2021（4）：40－45.

［6］SADEGHIAN R.Finding a probability density function for time intervals between earthquake occurrences［J］.Arabian Journal of Geosciences，2018，11（17）：1－10.

［7］LIU G L.A new physical model for earthquake time interval distribution［J］.Physica A：Statistical Mechanics and Its Applications，2017，465：62－65.

［8］KHAN M Y，IQBAL T，IQBAL T，et al.Probabilistic modeling of earthquake interevent times in different regions of Pakistan［J］.Pure and Applied Geophysics，2020，177（12）：5673－5694.

［9］SIL A，SITHARAM T G，HAIDER S T.Probabilistic models for forecasting earthquakes in the northeast region of India［J］.Bulletin of the Seismological Society of America，2015，105（6）：2910－2927.

［10］郭安寧，徐爱信.用倍周期叠加黄金分割法对西北地区八级地震的时间间隔的分析及预测［J］.西北地震学报，2005，27（增刊1）：138－139.

GUO Anning，XU Aixin.Analyzing the time interval of great earthquakes （M≥8） in northwest China and making prediction by the method of multiplied period with golden section［J］.Northwestern Seismological Journal，2005，27（Suppl01）：138－139.

［11］MARIANI M C，TWENEBOAH O K，GONZALEZ－HUIZAR H，et al.Stochastic differential equation of earthquakes series［J］.Pure and Applied Geophysics，2016，173（7）：2357－2364.

［12］陈翰林，蔡晋安，黄伟.基于LSTM方法的新疆地震时间预测试验［C］//第三届地球物理信息前沿技术研讨会论文摘要集. 贵阳，2021：81.

CHEN Hanlin，CAI Jin'an，HUANG Wei.Earthquake time prediction experiment in Xinjiang based on LSTM method［C］//Abstracts of the Third Symposium on Geophysical Information Frontier Technology.Guiyang，2021：81.

［13］游佳，楊鑫，刘议丹，等.基于灰色预测和SVR算法的地震预测及可视化［J］.闽江学院学报，2021，42（5）：1－7.

YOU Jia，YANG Xin，LIU Yidan，et al.Earthquake prediction and visualization based on grey prediction and SVR algorithm［J］.Journal of Minjiang University，2021，42（5）：1－7.

［14］朱海宁.基于改进支持向量机回归的地震预测方法研究［D］.合肥：安徽大学，2016.

ZHU Haining.Research on earthquake prediction method based on improved support vector machine regression［D］.Hefei：Anhui University，2016.

［15］ASENCIO－CORTS G，MORALES－ESTEBAN A，SHANG X，et al.Earthquake prediction in California using regression algorithms and cloud－based big data infrastructure［J］.Computers & Geosciences，2018，115：198－210.

［16］王坤，阮金梅，邓妮.基于SARIMA模型的曲靖市空气质量指数预测［J］.曲靖师范学院学报，2018，37（3）：25－29.

WANG Kun，RUAN Jinmei，DENG Ni.Prediction of air quality index in Qujing based on SARIMA model［J］.Journal of Qujing Normal University，2018，37（3）：25－29.

［17］RAY S，DAS S S，MISHRA P，et al.Time series SARIMA modelling and forecasting of monthly rainfall and temperature in the South Asian countries［J］.Earth Systems and Environment，2021，5（3）：531－546.

［18］MOEENI H，BONAKDARI H，EBTEHAJ I.Monthly Reservoir inflow forecasting using a new hybrid SARIMA genetic programming approach［J］.Journal of Earth System Science，2017，126（2）：18.

［19］于琛，付玉慧，张逸飞，等.基于ARIMA－BIGRU的船舶航迹预测［J］.船海工程，2021，50（6）：147－152.

YU Chen，FU Yuhui，ZHANG Yifei，et al.The track prediction method based on ARIMA－BIGRU neural network［J］.Ship & Ocean Engineering，2021，50（6）：147－152.

［20］易靖韬，严欢.基于小波分解和ARIMA－GRU混合模型的外贸风险预测预警研究［J］.中国管理科学，2023，31（6）：100－110.

YI Jingtao，YAN Huan.Early prediction and warning of international trade risks based on wavelet decomposition and ARIMA－GRU hybrid model［J］.Chinese Journal of Management Science，2023，31（6）：100－110.

［21］孙越，宋晓宇，金莉婷，等.基于ARMA－LSTM组合模型的铁路客流量预测［J］.计算机应用与软件，2021，38（12）：262－267，273.

SUN Yue，SONG Xiaoyu，JIN Liting，et al.Railway passenger flow forecast based on ARMA－LSTM combined model［J］.Computer Applications and Software，2021，38（12）：262－267，273.

［22］MILENKOVI C′ M，VADLENKA L，MELICHAR V，et al.SARIMA modelling approach for railway passenger flow forecasting［J］.Transport，2016：1－8.

［23］张繼艳，王新安，雍珊珊，等.基于ARIMA模型的九寨沟7.0级地震前兆异常检测［J］.华北地震科学，2019，37（1）：28－33.

ZHANG Jiyan，WANG Xin'an，YONG Shanshan，et al.Precursory anomaly detection of Jiuzhaigou M7.0 earthquake based on ARIMA model［J］.North China Earthquake Sciences，2019，37（1）：28－33.

［24］李磊，张宁，尹淑慧，等.基于残差修正的ARMA模型探测门源MS6.4地震前电离层TEC异常［J］.大地测量与地球动力学，2021，41（4）：382－386.

LI Lei，ZHANG Ning，YIN Shuhui，et al.ARMA residual correction model for detecting ionospheric TEC anomalies before Menyuan MS6.4 earthquake［J］.Journal of Geodesy and Geodynamics，2021，41（4）：382－386.

［25］SAQIB M，RK E，SAHU S A，et al.Comparisons of autoregressive integrated moving average （ARIMA） and long short term memory （LSTM） network models for ionospheric anomalies detection：a study on Haiti （MW=7.0） earthquake［J］.Acta Geodaetica et Geophysica，2022，57（1）：195－213.

［26］易丹辉，王燕.应用时间序列分析［M］.5版.北京：中国人民大学出版社，2019.

YI Danhui，WANG Yan.Applied time series analysis［M］.5th ed.Beijing：China Renmin University Press，2019.

［27］WANG W Q，GUO Y.Air pollution PM2.5 data analysis in los angeles long beach with seasonal ARIMA model［C］//International Conference on Energy and Environment Technology.Guilin，China.IEEE，2009：7－10.

（本文编辑：贾源源）

收稿日期：2022-04-21

基金项目：四川省科技计划项目《四川地区主要断裂系统地震活动性统计预测与技术应用研究》

第一作者简介：叶秋吟（1999-），女，硕士研究生，主要从事应用统计、数学方法应用、地学数据统计方面研究。E－mail：352299099@qq.com。