广西盘龙铅锌矿涌水量时间序列变化特征分析及ARIMA预测

吴卫忠，邓 忠，陈余道，卢丹美，陈 盟，夏 源，邹志坤，陆仁骞

(1.桂林理工大学 岩溶地区水污染控制与用水安全保障协同创新中心，广西 桂林 541006；2.广西壮族自治区水文地质工程地质队，广西 柳州 545006；3.自然资源部南方石山地区矿山地质环境修复工程技术创新中心，南宁 530031)

0 引 言

矿山涌水的发生和演化受多重因素影响, 具有典型的非线性动力特征, 能够导致突发性涌水事件[1-2]。矿山突水是矿井事故的重要诱因[3], 开展矿山涌水量时间序列分析及预测研究, 对防治矿山突水、 保障安全生产具有重要意义[4-6]。

近年来, 时间序列的确定性分析法和非线性分析法在多领域的推广应用, 为发展矿山涌水量预测方法提供了新的视角[4, 7]。除了统计学方法外, 重标极差法(R/S分析)和自回归移动平均模型(ARIMA)较为典型。前者能够用来分析时间序列的长程依赖性(LRD), 刻画出序列变化的主要控制因素, 能较好地解释非线性参数特征[8-9]、 降水量序列特征[10-11]和径流特征[12]等; 后者作为时间序列预测的一种模型, 能较准确地开展短期趋势预测[13-15]。这两种方法在矿山涌水量序列分析方面均有应用的案例, 比如山东郓城煤矿山工作面涌水量预测等[14]、 青龙煤矿月最大涌水量预测[16]、 东欢坨煤矿矿井涌水量预测[17]、 赵家寨煤矿矿井涌水量预测[18]等。这些应用多集中在我国北方地区且以煤矿为主, 取得了很好的预测效果。对于南方降水充沛、 地表水丰富且与地下水联系密切的岩溶区有色金属矿山, 吴松明[19]采用水均衡法、 大井法和水文地质比拟法等3种方法, 预测出铝土矿区的矿坑涌水量, 但应用R/S分析和ARIMA模型分析涌水量时间序列的案例未见报道。为此, 对广西大水量岩溶矿山——盘龙铅锌矿开展了涌水量时间序列变化特征分析和ARIMA预测研究。

盘龙铅锌矿位于广西武宣县黔江河畔。自2008年开采以来, 矿业安全一直面临地表水和地下水的侵害影响。根据矿山规划, 未来将向深部开采, 开采标高-440～-1 100 m; 加上2019年大藤峡水利枢纽黔江河段蓄水, 矿坑涌水也将成为影响安全生产的重要因素之一。

本文根据2010—2021年盘龙铅锌矿矿坑涌水量监测数据, 结合矿区水文地质条件, 在分析矿坑涌水量时间序列结构特征基础上, 采用R/S方法分析涌水量序列LRD特征, 并应用ARIMA模型进行涌水量分段预测, 以期为矿山水安全监测和水灾害防御提供参考。

1 矿区水文地质概况

盘龙铅锌矿所在区域属亚热带温湿气候区, 年平均气温21.1 ℃, 多年平均降水量1 370 mm, 降水多集中在5～9月份, 占全年降水量的69%。黔江是矿区东侧毗邻的河流, 多年平均流量4 240 m3/s, 常年水位38～45 m。地貌类型主要为构造侵蚀、 溶蚀峰丛谷地及残丘洼地(平原)。

盘龙铅锌矿是广西大瑶山西侧的大型沉积岩型铅锌矿床, 大地构造上位于桂中凹陷带与大瑶山隆起的结合部位, 经历了多期构造作用[20]。矿区出露的地层有第四系、 石炭系、 泥盆系、寒武系(详见图1),主要由白云岩、 石灰岩、 泥灰岩、硅质岩等组成。其中, 上伦组上段(D1sl2)白云岩为主要赋矿层位,在矿区出露面积最大(厚度约970 m), 是构成大岭矿段矿坑顶底板直接充水的含水层; 上伦组下段(D1sl1)为泥质灰岩与灰岩互层,富水性较弱。矿区发育多条断层, 其中北东向逆断层F1构成了矿区西侧边界; 近南北向平移断层F2则将矿区分为崩山矿段和大岭矿段(图1), 前者已在2009年停止开采, 后者2008年以来处于开采状态。

矿区水文地质条件复杂, 以溶蚀裂隙发育为主, 富水性由强到弱不等, 可分为4个地下水子系统: 大岭矿地下水子系统(Ⅱ-1)、 大坪岭地下水子系统(Ⅱ-2)、 崩山矿地下水子系统(Ⅱ-3)与六沙地下水子系统(Ⅱ-4)。各子系统相互之间有一定的水力联系。大岭矿地下水子系统(Ⅱ-1), 南侧主要由下泥盆统上伦组下段、 郁江组(D1y)、 那高岭组(D1n)、 莲花山组(D1l)泥灰岩、浅变质砂岩和泥页岩等构成,地层倾角约50°～70°;北侧由下泥盆统二塘组(D1e)泥灰岩构成;东侧毗邻黔江河流;西侧与F2断层相接(图1)。该子系统边界条件可概括为:南侧和北侧为隔水边界,东部为水头边界,西部为弱透水边界[21]。

2 矿坑涌水来源分析

根据矿区地下水系统划分(图1)及开拓中段涌水量监测, 大岭矿段坑道涌水主要有垂向入渗和东、 西侧补给3个方向的水源:

(1)大气降水、 低洼地表水入渗补给: 矿区近5年(2017—2021年)年均降水量1 403.24 mm, 积水区主要位于重晶石采坑、 水塘、 积水坑和沟渠等, 多为常年积水, 加大了地表水渗入补给地下水的水量。

(2)东部黔江河水侧向反补给: 大岭矿段距离黔江河岸约500 m, 近河岸的地下水水位常年低于黔江(平水期水位约38 m), 矿坑疏干水位为-115.00 m, 地下水水力坡度达到34%。矿坑地下水水位随着黔江水位升降而升降, 水力联系较密切。黔江河床底部为裸露的碳酸盐岩, 河水可以通过溶洞、 溶蚀裂隙带补给矿坑疏干区。2019年大藤峡水利枢纽蓄水后黔江水位提高20～30 m, 矿山在东部建成了挡水帷幕(图1), 帷幕底部高程-125 m, 削弱了黔江对矿区的侧向补给。

(3)西部崩山矿地下水子系统的侧向补给: 西部崩山矿地下水子系统与大岭矿地下水子系统被弱透水F2断层分隔。大岭矿段疏干排水时, 崩山矿地下水子系统的地下水能穿越F2断层渗流补给大岭矿坑, 补给强度受裂隙空间大小、 连通性、 水头高度等因素影响, 补给量有限[22]。

盘龙铅锌矿坑有-70、 -120、 -170、 -220、 -270、 -320、 -380和-440 m共计8个标高的开采中段, 各中段出水量成为矿坑涌水量主要来源, 2010年1月—2021年1月期间的监测结果如图2所示。据矿区勘察, -120 m标高以上中段岩层平均渗透系数K=4.85 m/d。溶蚀裂隙、 溶洞发育, 钻孔遇溶洞率48.24%, 巷道遇溶洞率74.07%, 岩溶中等发育。-120 m标高以下中段岩层的平均渗透系数为0.40 m/d, 钻孔遇溶洞率9.52%, 巷道遇溶洞率18.52%, 岩溶弱发育。随着开采中段的纵深开拓, 涌水点数量增加, 但涌水量主要来源于-120 m以上中段。开采初期(2014年之前),-70、-120和-170 m中段平均涌水量分别为550.5、93.1和126.2 m3/h,以-70 m中段涌水为主; 2014年之后, -70 m中段涌水量平均为433.3 m3/h(在2021年9月异常增大, 达775.1 m3/h), -120～-270 m中段涌水量较小, -320 m中段涌水量平均为317.4 m3/h(最大值为2020年7月的597 m3/h), 2021年增加了-380和-440 m涌水点, 平均涌水量分别为104.1和79.5 m3/h。因此, 总体上历年涌水量以-70 m和-320 m中段为主。

图2 大岭矿段各中段2010—2021年月均涌水量监测结果Fig.2 Monitoring results of water inflow in each middle section of Daling ore block from 2010 to 2021

3 涌水量时间序列结构特征分析

根据盘龙铅锌矿2010—2021年的月均涌水量时间序列(图3), 年内涌水量随时间上下波动, 与降水特征相似, 受黔江水位影响, 具有季节性, 其中6—8月份涌水量较高, 具有中间高两侧低的分布特征。

图3 盘龙铅锌矿2010—2021年月均矿坑涌水量时间序列Fig.3 Time series of mine water inflow of Panlong lead-zinc deposit from 2010 to 2021

3.1 涌水量时间序列结构参数分析

对涌水量时间序列进行统计, 得到平均值、 标准差、 偏度、 峰度和增强迪基-福勒(ADF)检验值(表1)。2010—2021年涌水量多年平均值为792.7 m3/h, 最小年均值为707.3 m3/h, 最大年均值为982.7 m3/h。2019、 2020和2021年年均涌水量有了明显提高, 分别达865.8、 942.7和982.7 m3/h, 这不仅与矿坑纵深掘进涌水量增加有关, 而且可能与黔江蓄水增加侧向补给有关[21]。

开采初期, 2010和2012年标准差较大, 与浅部含水层富水性差异大有关, 导致涌水量波动幅度大; 随着纵深开拓, 深部富水性差异逐渐变小, 地下涌水量序列趋于稳定; 2019年1月大藤峡水利枢纽蓄水后又出现了涌水量大幅度增大。总体上, 2010—2021年涌水量序列的偏度为0.842, 具有右偏分布趋势及其长尾特征, 说明后续随着矿区的开采出现涌水量超常值可能性相对比较大。由于2021年矿区涌水量异常增长, 也使得2010—2021年涌水量时间序列呈现尖峰态特征; ADF检验统计量t值和概率P值分别为0.133和0.723, 涌水量序列为非平稳时间序列[23]。

3.2 涌水量时间序列的长程依赖性分析

3.2.1 R/S分析方法 设不同时间t1，t2，…，tn所对应的矿区涌水量分别为x1，x2，…，xn，记为

时间序列。在τ=tn-t1时间段内, 其矿区涌水量平均值为

(1)

在tj(tj-t1)这一时期的矿区涌水量对于平均值的积偏差为

(2)

不同的τ对应不同的x(t,τ)序列, 将每一个τ值对应的x(t,τ)序列中的极差定义为区间, 记为R(τ)

R(τ)=maxx(t,τ)-minx(t,τ)。

(3)

同一个τ值下, 与R(τ)区间相对应的矿区涌水量标准偏差为

(4)

英国科学家赫斯特(Hurst)等认为, 时间序列时段τ内的极差R(τ)和标准偏差S(τ)的比值与τ/2之间属于一种幂律关系[24]。引入无量纲比值R(τ)/S(τ), 对R(τ)重新进行标度, 有

(5)

时间序列时段τ内, 极差和标准偏差的比值与τ/2的关系

(6)

式中:τ/2的指数H被称为Hurst指数。由不同时间序列得出的R/S值建立双对数坐标系ln(τ/2)-ln(R/S), 并应用最小二乘法进行线性拟合, 得出的斜率即为Hurst指数(H)。Hurst指数是反映动力系统标度不变性定律的重要指标, 可用来反映时间序列变化趋势的持久性和非持久性, 即长程依赖性, 通常在0～1, 无量纲。若0

3.2.2 涌水量序列的长程依赖性 根据2010—2021年盘龙铅锌矿矿坑涌水量时间序列所作的R/S分析结果如表2所示。各年涌水量序列的Hurst指数范围为0.780 1～0.994 9(平均0.918 6),表明矿区涌水量在年内的变化趋势十分相似, 持续性强,其原因与年内降水量和黔江水位呈季节性变化的规律有关。然而,2020和2021年Hurst指数分别低至0.839 6和0.780 1,相对于之前的年份序列持续性下降明显。结合表1中的序列峰度, 可以认为: 2020和2021年这两年涌水量异常升高、 呈现尖峰态特征是Hurst指数下降的主要原因。

表2 涌水量时间序列R/S分析

对2010—2021年际间涌水量序列开展R/S分析, 可得到Hurst指数为0.700 5(范围0.5～1), 说明其具有较弱的长程依赖性, 未来序列变化趋势可能会相对过去趋势出现偏移。推测其原因, 可能与水利工程枢纽蓄水(2019年1月开始蓄水)导致地下水系统水动力场异常变化有密切关系。

4 涌水量序列预测

4.1 时间序列ARIMA模型建立

对于非平稳的涌水量时间序列, 经差分转换为平稳时间序列, 可利用ARIMA模型对该序列进行预测, 其模型结构为[13, 25]

yt=c+φ1Xt-1+φ2Xt-2+…+φpXt-p+et-

θ1et-1-θ2et-2-…-θqet-q,

(7)

式中:yt为t时序数据;c为常数项;et为白噪声序列;p、q为模型的阶数;φp、θq为自回归和移动平均系数。

根据矿区坑道涌水点分布及2010—2021年涌水量组成(图2), 按-120 m标高以上段(-70和-120 m中段)和以下段(-170～-400 m中段), 利用2010—2019年月涌水量分别开展短期2020年和长期2021年各月涌水量预测, 通过对比预测值和实测值反映ARIMA模型的适用性。

首先对盘龙铅锌矿2010—2019年-120 m标高上段涌水量序列进行ADF检验, 统计量t值为-0.951(大于1%水平临界值),概率P为0.303(大于10%概率),可判别为非平稳序列;对数据进行一阶差分,再进行ADF检验,t值为-7.209(小于10%水平临界值),概率P为0.000(小于1%概率),表明涌水量序列在一阶差分后变成了平稳序列;进一步对平稳序列进行自相关性(AC)和偏自相关性(PAC)检验, 表明可应用ARIMA模型[26]进行分析预测。

通过赤池信息准则(AIC)与施瓦茨准则(SC)确定ARIMA模型的阶数p和q, 当p和q均等于2时, ARIMA模型的AIC与SC数值均最小, 说明该模型与原始矿坑涌水量序列拟合效果最佳。由此得到ARIMA模型表达式

yt=-2.648+1.521yt-1-0.722yt-2+et+

1.371et-1-0.371et-2。

(8)

对ARIMA模型进行误差检验, 表明该ARIMA模型符合要求可进行涌水量预测。同理, 对盘龙铅锌矿2011—2019年-120 m标高下段涌水量序列建立ARIMA模型, 可得到模型表达式

yt=4.419+1.070yt-1-0.649yt-2+et+

0.867et-1-0.161et-2。

(9)

4.2 预测结果与讨论

利用上述ARIMA模型对2020和2021年盘龙铅锌矿-120 m标高上段、 下段涌水量序列进行拟合预测, 结果见表3, 模型相关的平均绝对误差(MAE)、 均方根误差(RMSE)、 平均绝对百分比误差(MAPE)与全年矿坑总涌水量误差见表4。

表3 -120 m标高上、下段涌水量分段预测结果

表4 -120 m标高上、下段涌水量预测评价结果

2020年-120 m标高上段和下段矿坑涌水量预测值分别为4 323.7和6 274.0 m3/h, 与实测值误差分别为2.49%和3.14%; 2021年预测值分别为 000.4和6 837.0 m3/h, 与实测值误差分别为26.38%和9.64%。

2020年分段预测的MAPE均小于10%, 说明ARIMA模型预测准确, 并且-120 m标高上段预测值的MAE、RMSE与全年矿坑涌水量总误差表现优于-120 m标高下段, 均能够满足实际工程短期预报的需求。模型具有可靠性, 可用于盘龙铅锌矿矿坑突水预防工作。对于2021年涌水量预测结果, 预测时间较远, 且矿区开拓层不断加深, 其精度比2020年结果差, 主要原因是模型模拟未能考虑2020和2021年的年度异常涌水。因此, 用ARIMA模型作短期预报更为合适。

5 结 论

(1)2010—2021年盘龙铅锌矿矿坑涌水量年内变化受当地季节性降雨量变化的影响, -70和-320 m标高中段是主要的涌水段; 涌水量时间序列总体上呈现右偏长尾、 尖峰-平峰态和非平稳特征。

(2)年内涌水量序列具有明显的长程依赖性, 未来发展趋势与历年趋势之间呈现出较强的持续性, 与气象水文要素的季节性变化有关。年际涌水量序列(2010—2021年)具有较弱的长程依赖性, 过去的趋势存在偏移, 可能与矿山开采深度加大、 水文地质条件变化和黔江蓄水活动有关。

(3)ARIMA模型可用于盘龙铅锌矿矿坑涌水量的短期预测, 按-120 m标高上下分段预测较为合理, 能满足矿坑突水防御的需要。该模型预测需及时更新监测数据, 校正模型, 以提高模型精度。