浅埋采空区顶板应力变形与地表沉降的动态关系分析

任 敏 朱万成 贾瀚文 程关文 李斓堃 代 风

（东北大学资源与土木工程学院，辽宁沈阳110819）

目前，地下采空区引起的地表沉降逐渐成为了一种常见地质灾害。以鞍钢弓长岭露天铁矿为例，早在日伪时期就进行过地下富矿开采，造成地下多处采空区。随着露天采场开采范围拓展与延伸及深部铁矿资源的开采和利用，地表沉降严重威胁着矿山大型设备的安全，给正常采矿生产带来了极大困难。为了能够实时掌握弓长岭露天铁矿浅埋采空区顶板岩体的稳定性和地表沉降的发展演化，研究人员在采空区顶板和露天坑底的地表建立了集锚杆应力计、多点位移计、液压水准仪为一体的联合监测系统。因此，基于该联合监测系统的监测信息，分析浅埋采空区顶板岩体内应力、位移的变化与地表沉降之间的影响关系与响应机制，对于多传感器监测信息有效融合和精准的预测预警模型构建意义重大。

为避免地表沉降等灾害的发生，越来越多的监测技术与监测设备应运而生，为地表沉降监测提供了多种多样的监测手段和监测方法。目前普遍常用的地表沉降监测手段主要有全站仪［1］、水准仪［2］、GPS站点［3］等定点式监测方法，以及三维激光扫描［4］、雷达扫描［5］等区域式监测方法。在获取现场实时监测数据的前提下，为避免地表沉降变形引起灾害及人员、财产损失，开展地表沉降预测研究十分必要。近年来，众多研究者致力于利用不同的方法手段实现这一目标，主要包括人工智能法、灰色预测法、数学函数法等。如李勇发等［6］基于PS-InSAR 雷达干涉测量技术对地表沉降进行了监测，对获取的实时沉降数据利用遗传算法优化的反馈神经网络算法进行了地表沉降预测，获得了良好的预测效果；张涛［7］将DInSAR 技术与小基线集技术进行有机结合对地表沉降进行了监测，并使用粒子群算法优化支持向量机模型对沉降数据进行了拟合与预测；党星海等［8］采用卡尔曼滤波方法对三维激光扫描沉降监测获取的初始监测数据进行了去噪处理，而后利用GM（1，1）灰色预测模型进行了沉降变形预测，得到了优于传统灰色预测模型的预测结果；李篷［9］采用GPS方法进行了煤矿地表非采动沉降监测，并利用GM（1，1）灰色预测模型对煤矿主井、副井、煤仓附近区域开展了地表沉降预测；Nie 等［10-11］基于采空区上覆岩体不同破裂带的“S”型变形特征，以及塌陷坑GPS 监测获取的“S”型地表沉降曲线，提出了采空区的反正切函数预测模型；袁兴明等［12］利用数字水准仪进行了沉降观测，并分别利用BP 神经网络算法、改进BP 神经网络算法对沉降数据进行了预测。就目前研究来看，基于单一监测信息的地表沉降预测存在一定的缺陷，是因为数据源本身既无法充分反映地表沉降的多种诱导因素，也无法全面体现沉降过程中的多种响应信息，从而导致预测效果不稳定且与实际存在一定的偏差。

为了获取更全面、精准的监测数据，从多种信息融合角度判断监测对象的发展状态，采用多种监测手段针对地表沉降变形建立联合监测系统，开展多手段联合监测已经得到广泛研究。如何倩等［13］采用三维激光扫描与D-InSAR 联合监测矿区地表动态沉降；臧巍［14］将测量机器人与电子水准仪相结合应用于隧道监测工程中；郭利民等［15］利用GPS 和PS-In-SAR 联合监测西秦岭断裂北缘地表运动；李杨杨等［16］将GPS-RTK联合水准仪应用于充填工作面地表高精度动态监测中。目前，学术界对于多手段联合监测数据的处理及应用方面存在两种研究方向：一是将多种监测技术进行有机结合，得到精度较高的地表变形监测结果，如余礼仁［17］在对淮北袁店二矿煤矿区地表沉降进行监测时，结合了D-InSAR 和GIS技术得到监测数据，而后采用三次样条插值与灰色预测组合模型GM（1，1）对该矿7255工作面中心区域的下沉趋势进行了拟合并预测；邢学敏［18］将CR-In-SAR 与PS-InSAR 相结合，对矿区时序地表变形进行了监测，提高了矿区地表变形监测的精度和可靠性；二是根据模糊系统理论，将多个传感器获取的信息进行融合，进而对地表沉降状态进行综合估计与评价，如Wu 等［19］采用野外调查、地球物理勘探、钻井勘探和D-InSAR 技术等联合监测手段，应用模糊综合评价法对地表沉降进行了风险分区。由此可见，对于监测数据而言，虽然多传感器联合监测总体上提高了数据样本的多样性和精度，但是在基于多传感器联合监测结果的预测方法研究与应用中，对多传感器监测获得的地表变形时序数据未能进行有效融合，对于提升预测预警效果作用不显著。

程关文等［20-21］认为深部岩体变形是形成地表变形区域的原因，并采用数值模拟方法分析了随着不同开挖步的增加，深部岩体变形产生地表沉降变形的整个过程。宋长根等［22］对深部岩体变形对地表沉降变形的影响进行了研究，认为当深部岩体发生破坏后，地表会进入加速变形阶段，地表与构筑物不久便出现裂缝，严重影响了矿区安全生产与构筑物正常使用。本研究以弓长岭露天铁矿某空区顶板为例，根据锚杆应力计、多点位移计、液压水准仪对空区顶板岩体及地表沉降的联合监测结果，采用自回归分布滞后模型（ARDL）分析浅埋采空区顶板岩体内的应力、位移变化对地表沉降的影响关系和响应机制，确定该空区顶板岩体内发生应力、位移变化影响到地表沉降的滞后时间，这种发展关系的探究不仅可以为多传感器的信息融合提供有效途径，还可以为地表沉降灾害的预测预警提供支撑。

1 弓长岭铁矿某采空区概况

近年来，弓长岭露采分公司采用地球物理探测方法研究采空区的空间分布和大小情况，取得了显著成效，基本查清了一定深度空区的分布、大小等信息。但是，仅仅研究采空区的分布、大小、深度远远无法满足矿山安全生产的需要。因此，结合弓长岭铁矿的工程实际情况，建立了一套先进的监测系统，包括锚杆测力计、多点位移计、液体静力水准仪等，实现对地表沉降变形和采空区顶板岩体的实时在线监测，掌握空区顶板岩体内应力、位移变化及地表沉降的发展，以确定机械设备、车辆能否安全通过，对于露天矿安全开采至关重要。系统测点布置如图1所示。

本研究选择2018-06-11—2018-07-01获取的时间序列监测数据进行分析，测点为研究区内（图1）地表沉降液压水准仪测点S4、A 孔多点位移计位于135.2 m 水平的测点uA4和位于129.2 m 水平的测点uA1，A 孔锚杆测力计位于138.7 m 水平的测点σA3和位于131.7 m 水平的测点σA1。选择该类测点进行分析，原因在于：①3 种传感器测点基本处于同一地质环境，上覆岩体中积聚产生的应力和变形会以较单一的形式影响到地表；②各个测点空间距离较近，且在垂直方向上多点位移计测点与锚杆测力计测点在采空区顶板岩体中交错布置，基本能够监测到岩体内发生的应力和变形进而影响地表产生沉降的全过程。因此，可根据该类测点的数据，研究上覆岩体内的应力、位移变化对地表沉降的影响机制。

2 模型构建及实证分析

2.1 ARDL模型设定

本研究采用自回归分布滞后模型（ARDL）分析由锚杆应力计、多点位移计监测的浅埋采空区顶板岩体内的应力、位移对液压水准仪监测获得的地表沉降之间的影响关系、响应机制和滞后时间。Pesaran等［23］认为ADRL 模型的优势在于不需要变量同阶单整，也可以检验变量之间的长期关系，当解释变量为内生变量时，ARDL 模型的估计也不会受到影响，并且在小样本情况下ARDL的估计结果也更加稳定。基于此，本研究构建的地表沉降与顶板岩体应力、位移之间的ARDL模型为

式中，S 为地表沉降，为被解释变量；u 和σ 分别为采空区顶板岩体的位移和应力，为解释变量；q1、q2、q3分别为各变量对地表沉降长期影响的滞后期数；α1、α2、α3为各变量及其滞后期数对被解释变量的影响程度，若解释变量的滞后期t 值显著，则该变量的滞后期对地表沉降有显著影响；β0为常数项；ut为残差项。

2.2 基于ARDL模型实证分析

2.2.1 单位根检验

为了避免在实证检验过程中发生伪回归，需要首先对数据进行单位根检验。尽管ARDL 模型不要求变量具有同阶单整，但应用该模型的前提条件需确保变量单整阶数不能超过1。本研究采用最常用的ADF 检验统计量对模型变量进行单位根检验，结果如表1所示。

注：单位根检验形式（C，T，L）中，C为截距项，T为趋势项，L为滞后阶数。

由表1可知，在5%的显著性水平下，包括采空区顶板岩体的位移、应力以及地表沉降在内的所有变量均为非平稳序列。对这3 个变量进行一阶差分后的ADF 检验，在5%的显著性水平上都是平稳序列，说明这3 个变量均是一阶单整（即I(1)）序列，符合ARDL模型的应用要求。

2.2.2 协整检验

ARDL 模型的第一步是检验变量的协整关系。由于协整关系对滞后阶数的选择非常敏感，因此，需要估计向量自回归模型（Vector Autoregression，VAR），确定滞后阶数。通过信息准则和似然比（LR）检验法得到VAR 模型滞后阶数的确定结果如表2 所示。在滞后阶数为1 时，最终预测误差FPE 和3 个信息准则AIC、SC、HQ 均最小，由4 个参考值确定1 阶滞后为最优，因此可得滞后阶数为1。

注：LogL 表示对数似然函数；表中用“*”表示从每一列标准中选的滞后阶数；

滞后阶数确定后，为方便起见，本研究采用Johansen 协整检验方法进行协整检验，结果如表3 所示。

?

由表3可知：地表沉降S4与采空区顶板岩体内的位移uA1、uA4和应力σA1、σA3之间迹统计量和最大特征值统计量的检验结果均存在协整关系，因而可得出地表沉降与浅埋采空区顶板岩层中的位移、应力之间存在长期关系。

2.2.3 长期均衡系数估计

确定各变量之间存在协整关系后，可以用ARDL模型进行长期关系和短期关系的估计。在估计长期关系的变量系数之前需要确定模型中的滞后阶数，本研究使用SC准则确定模型中各变量的最优滞后阶数。估计结果显示：地表沉降S4与位移uA1、uA4，应力σA1、σA3以ARDL( 1,1 ,1,1,1 )模型最为合适，长期系数如表4所示。

注：D 表示差分，是因为单位根检验已得到的各变量均为一阶差分平稳，故所有变量均采用一阶差分形式；表示一阶差分后的各变量由ARDL 模型选择的滞后一阶形式表示一阶差分的原变量；C 为截距项；@TREND 为时间趋势项；AIC 为赤池AIC信息准则（Akaike Information Criterion），SC为施瓦茨SIC信息准侧（Schwarz Criterion），两者均是衡量回归模型优良性的标准，数值较小意味着滞后阶数较为合适；R2 统计量用于衡量模型的估计值对实际值拟合的好坏，数值越大，表明模型对因变量的拟合效果越好。

由表4 可知：浅埋采空区顶板岩体内的应力、位移滞后量对地表沉降影响显著，说明上覆岩体内的应力、变形影响传递到地表需要时间，以使得滞后量与地表沉降最相关。变量系数有正有负，变量系数为正，说明采空区顶板岩体内产生的应力、位移使得地表沉降加剧；变量系数为负则较为复杂，存在水平地应力增加的情形，使得上覆岩体内的应力、变形对地表呈现出与塌陷相反方向变形的影响。估计得到的地表沉降S4与位移uA1、uA4，应力σA1、σA3的ARDL( 1,1 ,1,1,1 )模型中，滞后一期的位移uA4、应力σA1和σA3对地表沉降的影响均显著，即前1 d 采空区顶板内产生的应力、位移与当前发生的地表沉降最相关，其中应力σA1对地表沉降S4的影响最大，其系数为0.529 29。

此外，还进行了数值模拟分析，即保证模拟的弓长岭露天铁矿浅埋采空区顶板岩体所受的竖向应力基本控制在监测结果范围内，测点所受应力为0.38～0.74 MPa 时，提取地表测点处的沉降值为0.005 97～0.009 39 m（图2）。由数值模拟结果分析可知：σA1的改变对地表沉降的影响较大，系数达0.572 81 左右，与ARDL 模型的估计结果较吻合。故可将浅埋采空区顶板岩体内的应力、位移作为地表沉降的长期影响指标，度量地表沉降。

2.2.4 短期动态调整关系分析

在估计得到长期关系系数后，可进一步求得基于ARDL模型的误差修正项（ECM）的系数，分析各个解释变量的短期变动对被解释变量的动态调整关系。ARDL的误差修正模型具体估计结果见表5。

注：EC M( -1 )为滞后一期的误差修正项。

由表5 可知：在10%的显著性水平下，滞后一期的地表沉降与滞后一期的位移uA1以及当期的应力σA1对地表沉降均存在显著影响，其中应力的影响最为显著，该变量的正向变动会导致地表沉降加剧。由误差修正项ECM 估计结果表明存在反向调整机制，意味着当受到短期冲击后，系统会较快恢复至长期均衡状态。

因此，基于上述协整检验、ARDL 模型的长期均衡系数估计和短期动态调整关系的分析，证实了地下采空区顶板岩体内产生的应力、位移与地表沉降之间存在长期和短期影响。

3 结论与建议

在地下开采等多因素的影响下，浅埋采空区顶板岩体会产生应力、变形，作为直接影响因素导致地表沉降，体现在传感器上能够捕捉到的信息即为顶板岩体内的应力、位移以及地表沉降。为了给地表沉降多传感器监测信息的有效融合及预测预警模型构建提供支撑，本研究利用ARDL模型系统分析了浅埋采空区顶板岩体内应力、位移的变化和地表沉降之间的影响关系。研究表明：地表沉降与浅埋采空区顶板岩体内的应力、位移之间存在长期和短期的影响；利用SC 准则确定的最优模型ARDL( )1,1,1,1,1分析表明，作为地表沉降解释变量的顶板岩体内的应力、位移滞后期均为1，即当前产生的地表沉降与前1 d 顶板岩体内产生的应力、位移有关，且滞后的解释变量的影响均很显著。

针对地表沉降灾害锚杆应力计、多点位移计、液压水准仪的联合监测，本研究提出如下建议：已证实浅埋采空区顶板岩体内的应力、位移对地表沉降具有长期影响，且滞后的应力、位移影响显著，故而可以采用各种预测模型（如人工智能算法等）将地表沉降时间序列数据作为输出值，采空区顶板岩体内滞后的应力、位移时间序列数据作为输入值，根据监测数据建立变量之间的线性、非线性关系进行预测。基于多传感器监测数据融合构建的预测模型，可以有效体现出地表沉降过程中浅埋采空区顶板岩体的响应机理，优于具有固定理论模型的数值模拟方法，因而具有进一步深入研究的价值和应用前景。