弓长岭露天矿采空区顶板位移超前预测算法

李相熙，朱万成，任 敏

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

露天矿边坡下层采空区不可避免的会因为爆破、机械凿岩、地表行车、地下水侵蚀等外部扰动发生顶板崩落、整体失稳，最终影响露天矿山的生产作业安全[1-2]。为了保证露天矿山人员和设备安全，提高生产效率和降低开采难度，有必要对露天矿山边坡下层采空区状态进行超前预测。

随着采空区稳定性监测技术的不断提高，已基本实现对采空区位移场、应力场等的自动化实时监测[3-4]。由于受静电、气温变化和机械振动等因素的干扰，导致采集的数据携带大量的误差数据，给后期数据的分析和应用带来巨大困扰。常见的时间序列数据超前预测方法对数据质量均有严格要求，直接使用时间序列原始数据使预测结果失真，进而无法提供长期、稳定、可信的预测结果。这就需要对所采集数据中的误差数据进行科学剔除，以筛选出有效的数据集。目前剔除时间序列数据中误差数据的方法主要基于数理统计[5]、距离计算、密度分析、聚类分析或小波分析等[6-8]中的一种算法来实现。其中，数理统计方法适用于对随机误差数据的剔除；基于距离计算、密度分析和聚类分析的方法是对“离群”较大的杂乱数据进行剔除。但是，露天矿山长期监测到的采空区顶板位移数据中随机分布着大量杂乱信号和随机误差，单一的误差数据识别剔除方法不足以对其充分剔除。

在数据超前预测算法方面，目前对时间序列数据长期超前预测技术的研究还处于初步探索阶段。徐峰[9]等提出将三峡库区八字门滑坡体的位移分解为趋势项位移和周期项位移，分别采用灰色模型及AR自回归模型进行位移预测；LIU Zhongqiang[10]等通过建立长短期记忆法、随机森林法等多种神经网络模型，实现对三峡库区“阶跃”型边坡位移的预测，预测准确度较高；赵静波[11]等利用灰色理论分阶段预测了露天矿边坡变形指标，结果表明预测值与实测值发展趋势一致；陶慧[12]等基于遗传算法对多个地压监测变量进行重构，通过广义回归神经网络( GRNN )实现对冲击地压的预测；蔡舒凌[13]等通过萤火虫模型( FA模型 )与非线性自回归模型( NAR模型 )的结合，提高单一自回归模型预测的准确度，实现了对隧道围岩变形的预测；曹文贵[14]等提出动态改变时间序列步长，再进行灰色预测的方法，提高了基桩峰值载荷预测的准确度；李恒凯[15]等建立由多种修改的灰色模型组成的预测模型库，应用灰色评价模型以挑选最佳预测模型，并在盘古山钨矿变形数据预测中取得了理想的效果。

针对露天矿边坡下层采空区顶板位移超前预测的问题，本文以鞍钢集团弓长岭铁矿独木采区采空区顶板位移-时间序列为基础，对K-means聚类算法进行改进，采用聚类-数理统计的组合模型实现对位移-时间序列数据的诸多误差数据进行剔除，以保证监测结果的准确性，降低误差数据对采空区顶板位移超前预测算法的干扰。在此基础上，基于灰色系统理论实现对位移-时间序列数据进行超前预测，以提前预测采空区顶板位移量的变化，为矿山的安全生产和高效开采提供保障。

1 采空区顶板位移超前预测算法

首先应用层次迭代聚类算法实现杂乱数据的粗略检测及剔除，使位移-时间序列数据能够满足正态分布；然后基于拉依达准则精细判断随机误差数据，采用插值替补法替换随机误差数据；最后使用单一变量的一阶灰色预测模型实现采空区顶板位移-时间序列数据的超前预测，并采用残余误差检验法、级比偏差检验法对预测模型进行检验，检验结果作为预测值的可信程度参考指标。技术路线如图1所示。

图1 技术路线 Fig. 1 Technology roadmap

1.1 误差数据剔除组合模型

首先将误差数据分为杂乱数据和随机误差数据两大类，在此基础上提出层次迭代聚类算法，实现误差数据粗略检测及剔除，然后基于拉依达准则实现随机误差数据的精细判别，并使用线性内插法替换判别出的随机误差数据。

1.1.1 层次迭代聚类算法

在对采空区顶板位移的实时监测中，不可避免地采集到诸多杂乱数据，杂乱数据一般具备以下2点较为明显的特征：① 杂乱数据一般偏差较大；

通过对弓长岭露天边坡下层采空区顶板多点位移计NA2测点的已有历史数据进行统计分析，将杂乱数据按照值大小粗略分为3级，见表1。

表1 杂乱数据分级标准 Table 1 Clutter data grading criteria

依据表1所示分级标准，对2019年12月31日采集的位移-时间序列进行分级，结果如图2所示。

图2 原始数据分级现象示意 Fig. 2 Schematic diagram of original data classification phenomenon

从图2可以看出位移原始数据大部分集中于( 0，5 )的正常值范围内，少量杂乱数据在表1所划分的3级区间内均有分布，故可以采用聚类的方法分步剔除这些杂乱数据。由于杂乱数据中少量第3级数据“离群”过远，传统基于距离的聚类方法，如Kmeans聚类算法，会使部分“离群”不太远的杂乱数据被计入正常值，导致杂乱数据剔除效果有限。为改善位移-时间序列数据中对杂乱数据的聚类效果，本文对K-means聚类算法进行改进，用于检测并剔除杂乱数据，改进后的算法称为层次迭代聚类算法。该算法的基本流程如下：

( 1 ) 对位移原始数据序列D0取每条数据的绝对值，得到绝对值序列D(0)。部分原始数据中杂乱数据既“离群”很远，又有正有负，直接进行聚类会出现3个聚类中心，为了简化聚类算法，不直接基于原始时间序列数据聚类，而是基于原始时间序列的绝对值序列进行聚类，这样聚类算法始终都只有2个聚类中心。

( 2 ) 获取绝对值序列D(0)的最大和最小值dmax，dmin，分别作为绝对值序列D(0)的2个簇中心，2个簇分别记为簇D1，簇D2。

( 3 ) 计算绝对值序列D(0)中每条数据d(i0)到2个簇的欧氏距离分别为Δ1，Δ2，并比较Δ1，Δ2的大小。Δ1计算公式(Δ2计算方法相同 )为

式中，Δ1为数据到簇D1的欧氏距离，mm；di(0)为第i个数据点的值，mm；dmin为簇D1的簇中心，mm。

( 4 ) 遍历完绝对值序列D(0)并将所有数据进行聚类后，将簇D1中的成员记为数据序列D(1)。

( 5 ) 计算数据序列D(1)所有成员的均值，以及数据序列D(0)的均值，如果有，则结束聚类迭代，数据序列D(1)为剔除杂乱数据后的初步清洗序列，否则令D(0)＝D(1)，D(1)清空后循环执行步骤( 1 )～( 4 )。

1.1.2 拉依达准则

采空区顶板位移监测过程中，因电压不稳等偶然因素产生的随机误差数据波动较小，层次分级不明显，导致聚类时随机误差数据和正常值始终在同一个簇。经层次迭代聚类得到的数据序列D(1)中不含杂乱数据，满足正态分布，故可以使用数理统计方法进行随机误差数据的检测。检测出的随机误差数据采用线性内插法进行插值替换。

应用拉依达准则( 也称做3σ准则 )对数据序列D(1)中可能含有的随机误差数据进行精细判定。拉依达准则认为：在一个数据序列中，如果有数据点的残余误差大于数据序列标准差的3倍，则认为该数据点含有偶然误差[16]。

拉依达准则精细判定随机误差数据的过程如下：

( 1 ) 计算数据序列D(1)的均值dabs：

式中，dabs为数据序列D(1)的均值，mm；di为数据序列D(1)中第i个数据点。

( 2 ) 计算数据序列D(1)的标准差σ：

( 3 ) 计算数据序列D(1)中每个数据点的残差，并做3σ判定：

根据标准正态分布函数，满足式( 4 )的数据点占数据序列D(1)的0.973；不满足式( 4 )的数据点占数据序列D(1)的0.027[17]。对于不满足拉依达准则的数据点，采用线性内插法进行替补，线性内插法公式为

1.1.3 误差数据剔除效果检验

通过计算时间序列光滑度来检验误差数据剔除效果。较为理想的误差数据剔除效果是：上述2种随机分布的误差数据被剔除或替补后得到的时间序列中，任意相邻2个数据点的值不会出现较大波动，时间序列曲线是平滑的。本文采用级比检验法检验位移-时间序列数据的光滑度。假设剔除误差数据后的位移-时间序列为D0＝{d1，d2，…，dn}，则时间序列D0的级比检验公式为

式中，λi(0)为时间序列第i个数据的级比；di为第i个数据点，mm。

如果时间序列D0中任意相邻2个数据点的级比满足式( 7 )，那么可以判定位移-时间序列D0误差数据剔除效果较好。

式中，n为时间序列长度。

1.2 位移-时间序列灰色预测模型

开采扰动、地质条件、地下水等诸多因素都会对采空区顶板位移产生作用，但是在对顶板位移量的预测研究中，不可能实现所有因素的量化，有些未知因素甚至无法进行量化描述，在这种信息量十分有限的情况下，可以将采空区顶板位移发展看作一个灰色系统。灰色系统是由我国学者邓聚龙教授在20世纪80年代提出并不断进行完善的一种惯性抽象系统[18]。灰色系统认为自然界有些事件在发展过程中，如果没有受到外界环境突变影响，那么这一事件将会沿着一条不易改变的迹线产生惯性演变。灰色模型就是要基于有限的已知数据找出这条迹线。

本文对位移-时间序列使用单一变量的一阶灰色模型进行超前预测，按照灰色理论，这种灰色模型被称为GM( 1，1 )模型[19-20]。对位移-时间序列数据进行预测的基本原理：通过对剔除误差数据后的位移-时间序列数据进行生成处理( 一般使用累加生成、指数生成等方法 )，使其形成呈指数发展趋势的递增序列。这种累加序列可实现原始离散数据的规律化，弱化灰色系统中的不确定性因素，在后续的递减还原中又可还原原始数据中的特征数据。然后针对递增序列建立一阶差分方程，由一阶差分方程变换得到参数未知的灰色预测模型。基于最小二乘原理近似计算出灰色模型计算参数，完成灰色模型建模。最后灰色预测模型通过递减还原处理，得到模型值的重预测值及未来一段时间的预测值。

1.2.1 建立迭代灰色预测模型

位移-时间序列数据以( 时间，位移值 )键值对的形式存在，假设通过光滑度检验的位移绝对值时间序列为G(0)，做一次累加后得到的累加序列为G(1)，累加公式为

式中，gi(0)为第i个数据点，mm；gi(1)为第i个累加值，mm。

针对累加序列G(1)建立一阶差分方程：

式中，a为灰色模型发展系数；b为灰色模型灰作用量。

基于最小二乘法基本原理，求解未知参数a，b。构造一个将a，b作为自变量的函数：

上式分别对a，b求偏导，求出使ω最小、拟合效果最优的a，b的值，代入式( 10 )整理得：

设通过最小二乘法得到的累加序列G(1)的拟合曲线为G(2)，则有灰色模型：

式( 12 )是累加序列G(1)的拟合曲线( 散点 )，对此拟合曲线做递减还原，得到的时间序列即为原始时间序列G(0)的重预测值及超前预测值组成的预测时间序列G(3)：

式中，gi(3)为第i个预测值，mm。

通过以上过程可实现灰色模型建模，但灰色模型参数a，b是固定值，这就使得灰色模型不可应用于任意长度的时间序列。因为随着时间序列步长的增加，参数a，b所确定的模型预测精度将会逐渐降低，直至模型失真。本文选择时间序列的步长为10，即n=10，确保现在值与过去值的比例适中，再基于“先进先出”的原则，用现在值替补过去值，计算a，b的动态值，在迭代中及时修正灰色模型。迭代的灰色预测可以提高预测准确度，提供长期可信的稳定预测结果。

1.2.2 计算预测准确度指标

预测结果是否可信需要有直观的评价指标作为参考。本文使用残余误差( 相对误差 )检验法和级比偏差检验法对灰色预测模型的可信度进行评价。

残余误差检验法使用重预测值的平均相对误差做为预测结果可信度评价指标。级比偏差可以判断一组重预测值的光滑度。灰色模型能够成立的假设是：在稳定环境下原始时间序列能够按照一条平滑迹线发展。那么预测时间序列无论预测值是否准确，都应该首先满足：预测时间序列也是一条光滑曲线。级比偏差检验法计算公式为

式中，λi(0)为原始时间序列的级比；λi(3)为重预测时间序列的级比；为平均级比偏差。

2 位移-时间序列超前预测微服务

基于Java语言完成算法程序的设计与开发，并创建Maven工程以方便项目管理。项目工程中使用Mybatis插件实现数据库连接与数据读取，使用Quartz插件实现算法程序的定时执行。最后使用Maven工具将项目打包成war包，直接运行于Tomcat服务器。

弓长岭矿独木采区因早期开采技术落后、盗采等原因形成的采空区一直未能得到有效治理。为确保采区生产安全，在已探明的采空区顶板布置多点位移计监测系统，对采空区顶板位移进行监测。每日采集的位移-时间序列数据最多达3 000条，其中随机分布着诸多杂乱数据、随机误差。该位移监测系统自2018年6月开始运行，在运行初期阶段，多点位移计运行状况较好，监测值从0开始逐渐平稳增大，后在1～1.4 mm范围内波动，监测数据质量较好。但由于长期受地表行车、降雨等影响，监测数据中逐渐出现“离群”较远的误差数据，位移监测值出现较大波动，导致监测数据序列产生整体偏移。多点位移计监测系统中NA2测点距离地表约13.61 m，距离采空区顶板约7.2 m，监测状态稳定，能够反映采空区顶板位移变化趋势。故选取NA2测点的位移-时间序列数据，自2019年9月开始进行长期预测。

将存储NA2测点位移原始数据的数据表接入到数据预测微服务程序中，借助Java Quartz定时任务插件，于每日0时读取前10 d的位移监测值，通过误差数据剔除程序及灰色预测程序得到新一天的采空区顶板位移预测值。

2.1 NA2测点剔除误差数据

NA2测点平均每日可采集230条数据，若不进行误差数据的识别与剔除，直接对原始数据求平均值，得到的结果如图3所示。

图3 NA2测点原始数据 Fig. 3 Original data of NA2 measuring point

由图3可以看出，2019-09—2020-02期间，NA2测点处顶板位移在1.5～1.7 mm范围内波动，由于受到误差数据的干扰，监测数据中夹杂有较多“离群”较远的散点，导致监测数据无法反映采空区顶板真实的位移变化。直接使用未经处理的原始数据进行预测，将会使预测结果失真。本文使用迭代聚类+数理统计的组合模型对误差数据进行处理，其程序设计如图4所示。

图4 误差剔除程序设计 Fig. 4 Eliminate error programming

在针对NA2测点位移-时间序列的微服务程序中，程序自动剔除NA2测点的误差数据，得到NA2测点剔除误差数据后的位移-时间序列，如图5所示。

图5 NA2测点位移-时间序列误差剔除效果 Fig. 5 Error elimination effect of NA2 measuring point displacement-time series

在剔除误差数据后，NA2测点位移监测值在1.4～1.6 mm范围内波动，整体呈现增大趋势。对该曲线进行光滑度检验，根据式( 6 )计算得到剔除误差数据之后位移-时间序列的级比如图6所示，符合由式( 7 )计算得到的级比区间( 0.833 8，1.199 4 )，说明误差数据剔除效果较好。

图6 光滑度检验结果 Fig. 6 Inspection result of smoothness

同时该曲线符合监测点处于采空区顶板中较密实位置，受外界环境影响相对较小，采空区顶板应力处于相对平衡状态，故在很长一段监测时间内，测点处累积位移量仅有轻微增加。说明本文通过迭代聚类+数理统计的组合模型，能有效剔除位移-时间序列中的杂乱数据、随机误差数据，剔除误差数据后的位移-时间序列可以较准确地反映采空区顶板实际监测位移量。综合考量以上定量计算及定性分析结果，可以认为误差数据剔除算法合理有效。

2.2 NA2测点灰色预测

使用剔除误差数据后的位移-时间序列执行灰色预测程序进行超前预测。预测程序首先得到位移-时间序列的累加序列。引入累加序列可以很好地将对波动较大的原始位移-时间序列曲线的拟合转化为对累加序列曲线的拟合，即将弱规律或无规律的复杂曲线拟合问题转化为单调递增的简单曲线拟合，拟合难度大大降低。然后基于最小二乘原理得到累加序列的有参拟合函数，节选部分位移-累加序列如图7所示。

由图7可以看出，原本复杂、波动的位移-时间序列得到了一个可以使用简单函数进行精准拟合的累加序列。基于此累加序列的拟合结果进行递减还原，得到采空区顶板历史位移数据的重预测值及未来预测值。预测程序在2019-09-01—2020-02-10期间运行得到的预测结果如图8所示。

图7 节选部分位移累加序列 Fig. 7 Extract part of the displacement accumulation sequence

从图8可以看出剔除误差数据后，灰色预测模型预测效果较好，对采空区顶板NA2测点位移自2019年9月—2020年2月的预测中，预测值与实测值基本吻合。得到的预测结果能够反映采空区顶板位移量的变化。分别采用残余误差检验法和级比偏差检验法对预测结果的可信度进行检验，NA2测点位移-时间序列超前预测微服务自运行以来，预测值的残余误差率均值为0.22%，最大不超过0.96%，预测值的级比偏差均值为0.21%，最大不超过0.89%。

图8 位移-时间序列预测结果 Fig. 8 Predictive results of displacement-time series

3 结 论

( 1 ) 结合位移-时间序列数据特征，在时间序列数据中含有误差数据的情况下，提出层次迭代聚类+拉依达准则组合算法，实现误差数据的自动检测及剔除。通过改进K-means聚类算法得到的层次迭代聚类算法可以有效剔除位移-时间序列数据中“离群”较远的杂乱数据，为确保超前预测程序的长期稳定运行提供了保障。

( 2 ) 基于灰色系统理论，对采空区顶板位移-时间序列数据进行长期超前预测分析。采空区顶板位移的历史预测结果与对应的实测值基本吻合，预测值的残余误差和级比偏差均在1%以内，能够正确预测弓长岭露天边坡下层采空区顶板位移量的变化趋势，预测结果与实际采集值具有对比意义。

( 3 ) 从近6个月的位移预测结果及实测值可以看出，弓长岭露天边坡下层采空区顶板位移不明显，位移量呈现出极缓慢的增加趋势( 6个月内采空区顶板仅产生0.2 mm左右的位移 )，短期内没有发生失稳破坏的风险。