物理指标与数据特征融合驱动的冲击地压时序预测方法

曹安业 ，刘耀琪，杨 旭，李 森，王常彬，白贤栖，刘亚鹏

(1.中国矿业大学 矿业工程学院,江苏 徐州 221116；2.中国矿业大学 江苏省矿山地震监测工程实验室,江苏 徐州 221116；3.徐州物硕信息技术有限公司,江苏 徐州 221116；4.中国矿业大学 计算机科学与技术学院,江苏 徐州 221116)

冲击地压灾害由采掘空间煤岩体集聚弹性能的瞬间释放所引起，具有突发、急剧以及猛烈的特点，易造成巷道严重损毁和人员重大伤亡[1-2]。究其根源，冲击地压的诱发机理不清、冲击灾变信息判识不明以及监测技术精度不足是主要原因[3]。

相较于冲击地压诱发机理研究的突破可有效指导针对性制定防控措施，冲击地压的监测预警更注重冲击灾变信息的全方位捕捉与灾变前兆信息的有效辨识。冲击灾变信息的捕捉方法主要可分为：①岩石力学方法，主要包括钻屑法[4-5]、煤岩体变形观测法、煤岩体应力测量法[6-7]、煤炮及矿震统计法、顶板来压预测法、工作面“见方”预测法等；②地球物理方法，主要包括微震法[8-10]、电磁辐射法[11]、声发射(地音)法[12-14]、电荷感应法[15]等，通过连续记录煤岩体内出现的动力现象，实时监测冲击危险状态。其中，微震监测方法是目前煤矿开采过程中监测范围最大、信息量最多的方法，我国冲击地压矿井几乎都安装了该类监测系统，已成为煤矿冲击地压区域性监测预警的最主要手段。齐庆新等[16-17]较早开展了冲击地压的微震监测理论及应用研究，姜福兴等[9]研发并在国内多个矿井建立了矿山微震监测系统，潘一山等[18]研制了矿震监测定位系统并进行了现场应用，窦林名等[19]在冲击地压微震监测和防治方面进行了大量研究，取得了显著成果。此外，冯夏庭、唐春安等[20-22]采用微震监测技术对金属矿山的岩爆进行监测，取得了较好的效果。

冲击地压危险的预测具体可分为冲击危险性的时序预测与冲击危险区域的空间预测，两者各有侧重，均具有广阔的研究前景。就冲击危险时序预测而言，诸多学者基于微震监测数据提出基于物理指标的时序预测方法，该方法的特点是致力于寻找冲击地压发生的前兆信息。如田向辉等[23]提出了一种基于微震能量、频次的定量-趋势冲击危险预警方法，李宏艳等[24]初步建立了以响应能量异常系数和无响应时间异常系数为指标的冲击危险性微震评价方法。同时，已有学者发现单一物理指标仅代表冲击孕育过程中某些物理特征的变化趋势，导致使用单一物理指标预警准确率难以进一步提升，据此提出采用多参量综合预警的思路，在现场应用获得较好效果[25-26]。然而，面对日益复杂的采掘工程环境，如“三高一扰动”等复杂环境下的冲击地压机理模型有待进一步完善，而基于物理指标的预测方法往往基于一定的物理假设与力学模型，导致其预测模型的适用性可能大打折扣。此外，在煤矿智能化建设的大背景下，如何高效利用海量的监测数据，充分挖掘其内在的灾变信息特征，以及预警模型的快速更新迭代是制约煤矿冲击地压智能化预警的新的瓶颈。

面对此困境，已有学者展开积极探索。国家重点研发计划“煤矿典型动力灾害风险判识及监控预警技术研究”已布局开展基于人工智能的冲击地压多源多场智能监测预警的阶段性研究[27]，魏立科等[28]提出了冲击地压灾害风险监察智能分析系统关键技术架构，陈结等[29]提出了冲击地压“双驱动”智能预警架构与方法。此外，由于岩爆与冲击地压均为矿山典型动力灾害，具有一定的相似性，在监测预警理念与方法上可相互借鉴。就机器学习算法在岩爆预测的研究方面，已有学者结合神经网络、支持向量机、随机森林、深度学习等人工智能方法对岩爆灾害预测开展初步研究。如孙臣生等[30]以非线性科学理论为指导，结合BP 神经网络对隧道岩爆预测模型进行改进，葛启发、冯夏庭[31]结合集成学习算法和BP 人工神经网络算法进行了岩爆等级预测模型研究，张航[32]开展了基于深度学习的隧道微震信号处理及岩爆智能预警研究，李明亮等[33]结合相关系数与支持向量机算法开展了岩爆烈度等级预测研究，田睿等[34]开展了基于深度残差神经网络的岩爆预测方法研究，DONG 等[35]结合随机森林算法进行了岩爆等级预测方法研究。岩爆预测领域的相关研究进展可为煤矿冲击地压智能预警研究提供一定借鉴。然而，基于机器学习所构建的模型具有数据集构建成本高、可解释性差、稳定性低等弊端。物理指标驱动和基于数据驱动的预警方法各具优势，如何将两者优势互补，做到物理指标与数据驱动的有效融合，是未来冲击地压预测预警的主要发展方向。

笔者以陕西彬长矿区某强冲击危险工作面为背景，统计分析了多次大能量事件发生前物理指标的变化特征并剖析了仅使用物理指标驱动的冲击地压危险预测方法的短板与不足，在此基础上提出采用物理指标与数据特征融合驱动的冲击地压时序预测方法。该预测模型包括数据预处理、特征提取以及预警模型构建3 个模块，数据预处理将原有微震监测数据处理为具有特定时间窗的前兆模式序列；特征提取主要包括基于物理指标的显式特征以及基于卷积神经网络的数据隐式特征提取；因此，提出基于注意力机制的显式特征和隐式特征的深度融合方法，并通过全连接网络实现预测模型分类，实现对不同冲击危险等级的大能量事件进行预测，并进行了方法应用及实践。

1 工程背景

彬长矿区某矿主采4 煤层，平均厚度9.43 m，埋深介于800～1 000 m，直接顶为厚度1～4 m 的泥岩，基本顶为厚度约10 m 的细粒砂岩，上覆多层厚度不等的粗粒砂岩。4 煤层经鉴定具有强冲击倾向性，其顶、底板分别具有弱冲击倾向性。

该矿生产期间大能量事件频发，并诱发过多起冲击地压显现。如2021-09-29 监测到一起1.6×105J 的微震事件，事件发生时地面有明显震感，井下冒顶范围可达55 m，冒顶高度约2 m，顶板锚杆索脱落、崩断，顶板钢带破断；底板与巷道两帮无明显变形。现场典型显现场景如图1 所示。

图1 彬长某矿典型冲击地压显现场景Fig.1 Typical coal burst scene of a coal mine in Binchang mine area

本文研究对象302 工作面为该矿三盘区首采工作面，平均埋深约980 m，顶板上方存在多层坚硬砂岩顶板。工作面区域煤层起伏较大，其中部与终采线附近分布有较大幅度的向斜，构造应力水平较高。仅2020 年11 月—2021 年2 月回采期间就发生大能量事件(能量大于105J)19 起，大能量震源分布及工作面布置如图2 所示，综合柱状如图3 所示。

图2 302 工作面布置Fig.2 Layout of longwall 302

图3 302 工作面综合柱状Fig.3 Synthesis column around lonwall 302

302 工作面自2020 年6 月—2021 年5 月回采期间监测到多次微震事件，各能级微震事件频次与能量占比如图4 所示。工作面回采过程中能量大于5×104J的事件频次占比仅为3.81%，而能量占比却可达32.67%，可见工作面开采应力调整剧烈，大能量事件对工作面的冲击危险性起到重要影响。同时，工作面监测到的大量微震信息可为本文后续研究提供数据基础。

图4 302 工作面回采过程中微震事件频次及能量占比Fig.4 Microseismic event frequency and energy ratio during the longwall 302 mining

302 工作面回采前采用综合指数法及多因素耦合法对冲击危险进行了预评估。冲击危险综合指数为0.86，具有强冲击危险性；受大采深(采深接近千米)、断层以及褶曲等地质构造影响，工作面以中等危险区域为主，并多处具有强冲击危险性。为实时掌握工作面冲击危险程度，指导现场卸压解危工作，有必要进一步开展工作面危险预警方法研究，尤其应提高大能量事件的预警准确率。

2 物理指标驱动的冲击危险预警实例分析

为实现对大能量事件的精准预警，本节首先分析工作面大能量事件发生前常见物理预警指标，如矿震时间信息熵Qt、矿震活动度S、a、b、A(b)、Zmap、矿震活动标度 ΔF、缺震以及断层总面积A(t)等[25]指标的变化规律，进而对物理指标驱动的冲击危险预警的有效性进行分析(震活动度S结合了矿震频次、矿震最大能量、矿震平均能量和震源分布集中度等常用的指标，并与震源分布集中度的记忆效应相结合，其值的大小描述了矿震活动的强度。A(b)综合考虑某一时间段某一区域内矿震事件的能量、频次和活动性等各种因素，直观定量地对某一区域矿震活动程度为“增强”或“平静”做出评价。a的大小可以代表矿震活动水平的高低，a越大说明矿震发生的频次越高，矿震活动越频繁。b的大小代表了矿震活动强度的高低，b越小，发生强矿震的可能性就越大，发生冲击地压的可能性也越大)。

如图5 所示，选取2020 年11 月—2021 年2 月期间的微震数据进行分析，分别计算了该时间段19 次大能量事件发生前上述9 个指标的变化情况(计算时间窗口与滑移步长为1 d，大能量震源分布如图2 所示)。

图5 大能量事件发生前后各物理指标时序变化特征Fig.5 Temporal variation characteristics of physical indexes before and after the occurrence of high-energy events

利用物理指标变化特征通常从2 方面进行预警(此处分析的时间窗为3 d，滑移步长为1 d)：①趋势预警，如图5 中各指标存在急剧上升(3 d 内指标连续上升)、急剧下降(3 d 内指标连续下降)、先升后降(3 d 内指标先升后降)、先降后升(3 d 内指标先降后升)等情况；实际应用过程中根据大能量事件发生前各指标“异常变化趋势”开展分析预警。②阈值预警，通过大能量事件发生前各指标所达到的阈值进行预警，如图5中各指标存在局部极值的情况，实际应用过程中可通过各指标值的“局部极值”是否超过阈值进行预警。

通过及时总结各物理指标的变化规律，可有效改善冲击地压预警的准确率。但是仅凭人工进行庞大监测数据统计，难免工程量巨大。此外，基于人工经验对变化趋势、阈值进行标定时，难免由于个人经验等因素产生非必要的误差，例如对指标变化趋势中“上升”还是“下降”的判断，若选择不同的统计窗口可能会得出完全相左的结果；此外，冲击危险程度不同，大能量事件发生前各指标的变化速率也显著不同，凭人工进行判断势必会产生误差，而通过统计方法计算时不同的时间窗口也会产生差异。对采用“阈值”进行预警而言，同样存在采用统计方法与统计窗口的不同，可产生误差甚至导致预警失败。

为进一步分析采用“趋势”法进行大能量事件预警的有效性。对图5 中19 次大能量事件发生前各物理指标的变化率、当日预警值、前3 d 最大值以及变化趋势进行了统计(缺震只包含1 和0 两种结果，此处不做统计)，各指标变化趋势的统计分布如图6 右轴与上轴所示。

图6 大能量事件发生前后各物理指标变化统计特征Fig.6 Statistical characteristics of physical index changes before and after high-energy events

如图6 所示，各物理指标的4 类变化特征在大能量事件发生之前存在明显的规律性。以图6(a)为例，大能量事件发生前Zmap的4 类指标变化特征的概率分布存在较为明显的集中分布特征，如当日预警值与变化率的概率峰值集中在0.3 与0.4 左右，前3 d 最大值与变化趋势的概率峰值均集中在0.4 附近；通过统计其他7 类指标(图6(b)～(h))也均可获得类似的统计特征。

理论上，通过统计各类指标的变化特征可对大能量事件进行预警。但其实不然，主要原因有：①各指标变化特征的概率分布虽具有明显的集中分布特征，但是集中度并不十分突出，如图6(c)中S的当日预警值与变化率概率分布曲线较为平缓，在实际判断过程中可能难以得出准确的预警值；②由于各物理指标往往基于一定的物理含义，不同地质与采掘环境下指标的敏感性可能存在差异，而井下采掘过程中条件复杂，孕灾环境多变，这会导致物理指标预警的敏感性显著不同，以本文302 工作面各物理指标统计特征为例(图6(a)～(h))，在本工作面开采背景下各指标的指标统计特征存在明显的概率分布差异可证明这一点；③井下采掘过程中微震数据时刻产生，积累了大量的数据，如何自动、高效的获取物理指标的前兆信息也面临挑战。

因此，笔者后续尝试结合大数据与机器学习方法对预警指标的前兆信息特征进行提取，并提出一种可大规模处理与实时分析的冲击地压时序预测方法。

3 物理指标与数据特征融合驱动的预测方法

3.1 预测模型架构

如图7 所示，本文提出的预测方法主要包括3 个部分：数据预处理、特征提取以及预测模型。

图7 物理指标与数据融合驱动的冲击地压预测模型架构Fig.7 Coal burst prediction model driven by combined physical index and data characteristics

数据预处理将微震原始数据包括发生时间、能量以及震源坐标等转化为用于模型输入的前兆模式序列数据。具体过程为首先将原始数据进行统计分析，计算得到每日的最大能量和平均能量，生成以日为最小单位的时间序列数据。指定前兆模式序列长度，生成多前兆模式序列及其标签。

特征提取包括知识驱动的显式特征提取和数据驱动的隐式特征提取。显式特征是指利用物理指标所计算的相关指标，并基于子集搜索方法确定最优的显式特征指标组合。隐式特征是指利用深度学习方法从海量数据中所挖掘的隐藏规律信息，本文利用卷积神经网络实现隐式特征提取。

预测模型包括特征融合和分类网络。在特征融合部分，提出了基于注意力机制的显式特征和隐式特征的深度融合方法。在分类网络部分，通过全连接网络拟合实现分类，进而构建冲击地压大能量事件预测模型。

3.2 数据预处理过程

如图8 所示，为使数据适用于模型的训练和预测，对微震原始数据包括时间、能量以及震源坐标等进行预处理。

图8 数据预处理过程Fig.8 Process of data preprocessing

首先，以固定时间窗统计数据，构建时间序列数据集。假设第i个时间窗计算得到的数据记录为mi，其可以表示为

其中，id为时间窗编号；emax为时间窗内微震最大能量；emean为时间窗内的微震平均能量；f为时间窗内微震频次。当数据被划分为n个时间窗并对数据进行遍历可得时间数据序列数据集M=(m0,m1,m2,···,mn-1)，进一步根据时间序列数据集M构建前兆模式序列，假设其第i个前兆模式序列为si，可表示为

其中，p为前兆模式序列长度；j为采样步长。因此，在n个时间窗口且n≫p的前提下，生成的前兆模式序列集S可表示为

式中，D为前兆模式序列样本在预测时间范围为N小时情况下的前兆模式序列数量；N为统计时长。

本文的目的是预测大能量事件发生的概率，在利用深度学习模型训练之前，需对不同冲击危险等级的微震能量进行定义及标记，具体对应关系见表1。此外，便于后续模型训练，需建立前兆模式序列集合对应的标签集合T表示为

表1 标记的微震能量及其与冲击危险的对应关系Table 1 Labeled microseismic energy and its relation to coal burst risk

其中，ti为前兆模式序列si的标签。将前兆模式序列si未来N小时内最大能量表示为ei，取值范围见表1。

3.3 特征提取方法

3.3.1 物理指标的显式特征提取

显式特征是基于物理指标提取的特征指标，本文将第2 节所述9 个指标作为候选集，采用反馈式选择方法进行最优指标筛选，其主要思想是动态添加指标，具体利用已有微震数据计算不同指标添加后的模型精度，根据精度判断该指标集合是否最优，步骤如下：给定微震指标集合 {a1,a2,···,am}，将每个指标看作一个候选子集，对m个候选指标子集分别进行评价，假定第1 轮最优候选集是{a2}，则将其作为第1 轮选定集；后续在该选定集中加入一个指标，构成了包含2个指标的候选子集，假定这m-1 个候选子集中{a2,a5}最优，则将其作为第2 轮选定集；直到r+1 轮，最优候选子集精度低于上一轮，则将第r轮指标集合做为最终的特征选择结果。在本研究中，一个特征子集优劣的评价标准是特征子集在不同数据集和不同采样方法上的模型平均准确率。反馈式微震指标选择机制解决了微震指标选择问题，一定程度上提高了模型预测准确率。

依据“反馈式选择”的思路，笔者选择预测时长为3 d，选用样本训练集和测试集样本数量比为1 919∶822，选用模型参数见表5，依次对不同物理指标组合预测的F1(计算公式见式(9)，试验结果见表2)进行了比较，笔者发现使用b、a、A(b)、缺震组成的指标组合预警准确率较高(F1=0.854)，若在此基础上再增添指标预警准确率上升有限，却会大幅增加运算时间，对实现该方法的实时预警不利，因此本文将上述4 个指标作为优选的显式特征指标，其物理含义及计算方法见表3。

表2 不同物理指标组合预测的预测F1Table 2 F1 values predicted by different physical indicator couples

表3 优选的物理指标及其含义Table 3 Physical indicators and their meanings

3.3.2 数据驱动的隐式特征提取

本文提出基于深度学习的隐式特征提取方法，该方法利用深度学习中的卷积神经网络实现隐式特征提取。卷积神经网络其隐含层内的卷积核参数共享和层间连接的稀疏性使得卷积神经网络能够以较小的计算量实现隐式特征提取。本文采用的卷积神经网络为3 层卷积，总共包含112 个卷积核，卷积核大小均为3×3，步长为1，填充为1。经过3 层卷积后进行扁平化操作，再通过一个线性层，即可输出1 000 维度的隐式特征向量，该向量即为本文所要提取的隐式特征。

3.4 预测模型构建

基于3.2 节创建的训练数据集，可建立目标函数，利用反向传播算法进行网络权重的更新，使得目标函数损失最小化，从而生成预测模型。训练完的模型用于大能量事件预测，即输入前兆模式序列，用于预测未来一段时间内发生大能量事件的等级。本文任务属于时间序列分类任务，其输入是前Nd 的矿震数据，输出为(N+1)d 的冲击地压等级。模型整体网络架构如图9 所示，包括特征提取、特征融合、分类网络3 部分。

图9 预测模型网络架构Fig.9 Network architecture for prediction model

3.4.1 特征融合

在预测模型中，需将显式特征和隐式特征进行深度融合。由于显示特征和隐式特征的复杂性和异构性，导致简单的加权特征融合方法并不适用，因此本文利用深度学习中的注意力机制法实现特征融合，并实现显式特征和隐式特征内每个维度的加权。

如图10 所示，首先，将显式特征Fe和隐式特征Fi进行合并，从而得到初始特征向量Fs，有

图10 预测模型中显式与隐式特征融合示意Fig.10 Schematic of explicit and implicit feature fusion in prediction model

式中，⊙为哈达玛积，融合后的特征用于后续分类网络输入。

3.4.2 分类网络

预测模型中的分类网络包括一个全连接层和一个激活函数，在深度学习中，全连接层起到分类器的作用，其可以将学到的分布式特征表示映射到样本标记空间中，本模型中全连接层由2 000 个神经元组成。然后，利用激活函数进行归一化处理，得到每个类别的概率，概率最大的类别即为输出类别，即发生该能量区间(表1)的大能量事件概率最大。

3.4.3 模型训练

考虑到实际数据中大能量事件发生的次数较少，导致标签为小能量事件的数据远多于大能量事件，进一步导致训练数据集样本不平衡。若使用传统的深度学习模型训练方法，则会导致模型在分类时偏向学习样本较多的一类。针对该问题，本文采用“再放缩”的思想，在训练过程中，模型根据批样本和总体样本的分布情况，动态调整各类别在计算损失的权重。在网络训练中，使用加权交叉熵损失函数计算模型的损失，进而不断更新神经网络模型中的参数，使得模型在训练数据集上的损失最小化。此外，在冲击地压预测过程中，本文在损失函数中添加了各类别的学习权重z0和z1，通过调整大能量事件的学习权重，使得该模型可以更偏向大能量样本的预测，以降低大能量事件的漏报率。本文的加权交叉熵损失L可表示为

其中，Li为第i个前兆模式序列的损失；N为前兆模式序列数量；M为能量等级类别数量；zc为类别c的学习权重；wc为类别c的样本分布权重，若第i个前兆模式序列的标签为c，则yic=1，否则yic=0；pic为观测样本i为类别c的预测概率；W为所有类别样本的权重之和；yic为真实样本类别概率分布，即若第i个前兆模式序列的标签为c，则yic=1，否则yic=0。上述方法有效改善了微震数据不平衡的问题，加快了模型收敛速度，并提高了模型预测准确率。此外，通过深度学习中常用的“dropout”方法来避免模型过拟合问题。

4 模型训练及预测实验

4.1 评价指标

本实验利用混淆矩阵记录模型预测结果(表4)。对于每个类别，若实际情况为真，且预测为真，则记作TP(True Positive)；若实际情况为真，且预测为假，则记为FN(False Negative)。若实际情况为假，且预测为假，则记TN(True Negative)；若实际情况为假，且预测为真，则记为FP(False Positive)。

表4 记录预测结果的混淆矩阵Table 4 Confusion matrix of the predicted results

在模型评估中，本文使用F1作为评估指标。F1能综合考虑查准率(Precision)和召回率(Recall),计算公式为

式中，Pr和Re分别为查准率和召回率，计算方法为

4.2 模型训练及预测

图11 为302 工作面2020-11-01—2021-03-01回采期间微震日最大能量分布情况，受复杂采掘环境等因素影响，工作面多数情况下处于中等冲击危险状态，甚至局部时间段频繁发生105J 以上的大能量事件。一方面说明在302 工作面实施大能量事件预测的必要性，另一方面丰富的数据可满足预警模型训练的需要。以该时间段共131 d 的微震数据为样本进行训练和测试，具体将数据以7∶3 的比例进行随机划分，生成训练集和测试集。在预测未来1、2、3 d 情况下，生成的训练集和测试集样本数量比分别为1 952∶836、1 936∶829、1 919∶822。后续在预警效能分析时采用训练集数据完成模型训练后，使用预警模型对测试集中的大能量事件进行预测以检验预测模型的有效性。数据预处理参数设置见表5。设置模型训练时的批大小为16，最大训练轮次为130，学习率为0.000 1。不同预测实验中训练集与测试集的目标函数损失变化如图12 所示，实验过程中当目标函数损失连续20 轮不降低或者F1连续20 轮不上升即停止训练。

表5 数据预处理参数取值Table 5 Parameter value of warning model

图11 302 工作面2020-11-01—2021-03-01 回采期间日最大能量分布Fig.11 Maximum daily energy distribution during longwall 302 retreat from 2020-11-01 to 2021-03-01

冲击地压危险实时预测有助于保障人员及时撤出，灾情及时处理及指导卸压防控措施等。因此，预测冲击地压或大能量事件何时、何地以及能级大小便是关键所在。针对预测冲击地压危险“何时”发生，为验证本文提出的融合驱动的冲击地压预测方法的预测效能，采用训练完成的预测模型对未来一定时间内(预测1、2 及3 d)的大能量事件进行预测，并采用混淆矩阵分析其F1，预测结果的混淆矩阵以及F1见表6～8。

表6 预测未来1 d 的预测结果混淆矩阵Table 6 Predicted results in the future one day

表8 预测未来3 d 的预测结果混淆矩阵Table 8 Predicted results in the future three day

由表6～8 可知，预测实验F1最高可达0.956，最低为0.854，可见本模型对冲击危险的预测效果较好。因此，现场选用本模型进行预测时，可根据现场危险情况(采用综合指数法等可提前预估采掘工作面整体冲击危险等级)确定预测时长：若提前预测采掘工作面冲击危险等级较高，例如本文实验的302 工作面，大能量事件频发(图11)，冲击危险性较高，需实施大量的卸压工程、耗费大量时间，因此需尽量提前预测危险，此时选用预测时长3 d 为宜；若预测现场预估的整体冲击危险等级较低或者预估的冲击危险等级较高同时需采取的卸压工程量较小时，选用预测时长1 d 为宜。

通过本实验不难发现，在工作面采掘过程中可保证微震数据正常采集情况下，使用本模型可有效预测冲击地压危险；未来在此基础上，笔者将进一步开展实验，以探求本模型在预测时长更短时能否精准预测危险的可能性。在此过程中，需保证模型在微震数据不足的情况下获得较好的训练结果，拟解决的途径主要有：①试验本模型在小样本情况下的适应性，适时调整参数，以期获得较好的训练结果；②通过云平台等建立冲击地压事件数据库，采用大批量的模型调整模型参数，使预测模型具有更强的泛化能力；③研究基于迁移学习的模型泛化方法，使模型具有较强的自适应及调整能力，实现模型参数的及时更新迭代与调整。

4.3 工程应用性能测试

为进一步评估模型总体性能，本节进行工程应用测试。实验设计思路如下：①利用302 工作面前110 d的数据训练初始模型，然后预测第111～113 天内发生各能量等级事件的概率；②模拟时间推移1 d，对模型进行训练并更新，具体采用前111 d 的数据重新训练，预测第112～114 天内发生各能量等级事件的概率；③再模拟时间推移1 d，对模型再进行更新，利用前112 d 的数据重新训练，预测第113～115 天内的冲击危险性；以此类推，共进行了50 次模拟预测。由表9 可知，在2021-02-27—2021-04-17 这50 d 内，302 共发生6 起105J 以上的大能量事件。对比实际与预测危险等级不难发现，本模型通过每日训练更新迭代后准确预测了6 起大能量事件中的5 起，将另外1 起具有强冲击危险性的大能量事件预测为中等冲击危险。为进一步对模型的工程应用效果进行分析，统计了50 d 模拟预测结果中误差分布(图13)，不难发现50 次预测中26 次预测结果准确，20 次预测结果为误差1 级(例如实际危险等级为弱而预测等级为中等，两者相差1 级)，仅有4 次预测结果误差等级可达2 级。工程应用结果表明本融合驱动预测模型可对强冲击危险尤其是大能量事件提前预测，进而为针对性防控措施的制定提供重要参考。

表9 工程应用中实际与预测危险等级对比Table 9 Comparison of actual and predicted risk levels in engineering applications

图13 工程应用中预测结果误差分布Fig.13 Error distribution of prediction results in engineering applications

5 结论及展望

(1)通过实例分析了物理指标驱动用于冲击地压危险预测的不足与短板。结果表明物理指标主要利用冲击危险来临前各物理指标的异常变化进行预警，如典型的“趋势法”与“阈值法”，而2 者的概率分布特征分析表明难以得出明显的统计特征；各物理指标往往仅适用于特定采掘条件，模型泛化与跨场景预测效率较差；物理指标往往通过特定统计数据得出，难以结合大数据分析技术进行模型更新及实时分析。

(2)提出物理指标与数据特征融合驱动的预测模型。预测模型包括数据预处理、特征提取以及预测模型融合。数据预处理将原有微震监测数据处理为具有特定时间窗的前兆模式序列；特征提取主要包括基于物理指标的显式特征以及基于卷积神经网络的数据隐式特征提取；提出基于注意力机制的显式特征和隐式特征的深度融合方法，并通过全连接网络实现预测模型分类，实现对不同冲击危险等级的大能量事件进行预测。

(3)采用混淆矩阵对预测结果进行分析，预测实验结果表明预测时长为未来1 d、未来2 d 以及未来3 d时预测F1分别可达0.956、0.950、0.854，现场可根据需求选用预测时长；工程应用结果表明融合驱动模型可对大能量事件准确预测，误差分析结果表明模型预测准确率较高，可满足现场需求。

煤岩动力灾害的智能化预测是煤矿智能化建设的重要组成部分，而机器学习与大数据技术是实现智能化预测与提高预测准确率的关键所在。本文提出的基于物理指标与数据特征融合驱动的预测方法仅是将机器学习技术用于冲击地压危险预测的初步尝试，经笔者粗略总结，未来仍有面向不同煤矿预测模型迁移、小样本数据情况下的机器学习模型构建以及多源异构数据融合等关键技术问题需进行攻关。