结合SBAS-InSAR与支持向量回归的开采沉陷监测与预测

师芸，李杰，吕杰，马东晖

(1.西安科技大学 测绘科学技术学院，西安 710054；2.自然资源部煤炭资源勘查与综合利用重点实验室，西安 710021)

0 引言

煤炭是我国的主要消耗能源，开采后采空区周围岩体失去原来的平衡状态而发生移动，诱发滑坡、塌陷、地裂缝等地质灾害[1]。煤炭开采引起地面沉降、塌陷等矿区地质灾害问题亟待当地政府和企业解决。分析预测开采沉陷的可能性对矿区安全开采具有积极的意义[2]。

合成孔径雷达差分干涉测量技术(differential interferometric synthetic aperture radar，D-InSAR)对矿区进行沉降监测已被众多国内外学者证实具有可行性，其监测精度可达到cm甚至mm级[3-4]。然而，传统D-InSAR技术观测单一，易受时空基线的限制导致失相干现象[5]，从而无法获得矿区内时序沉降量。为了克服上述问题，时序InSAR技术[6-7]逐渐发展了起来，以小基线集InSAR(small baseline subset InSAR，SBAS-InSAR)为代表的多主影像时序InSAR技术能够进行大尺度长时间序列的地表形变监测，广泛应用于城市地表形变[8-9]、地质灾害监测[10]、矿区地表形变等领域[11]。

目前，常用的矿区沉降预测方法主要有概率积分法[12]、神经网络(neural networks，NN)[13]、灰色模型(grey models，GM)[14]和支持向量机(support vector machine，SVM)[15]等。支持向量机是由Vapnik提出的一种机器学习算法，具有估算精度高、可解决小样本和非线性问题、泛化能力强等特点。然而，SVM的性能受其参数影响较大，模型参数选取存在一定的主观性与随意性使其效率低下。遗传算法(genetic algorithm，GA)[16]是一种通过模拟自然界遗传机制和生物进化论来进行随机优化的算法，它具有较强的全局搜索能力，不依赖特定的求解模型，但基本遗传算法(sample genetic algorithm，SGA)存在欺骗问题，即进化过程缓慢，易早熟收敛。针对这一问题，Srinivas等[17]提出了一种根据种群适应度自动调整的自适应遗传算法(adaptive genetic algorithm，AGA)，对SVM模型中的参数进行寻优，可以避免选参的盲目性从而得到全局最优解，同时提高预测效率和准确度。

本文提出将AGA-SVR模型应用于开采沉陷预测，利用17景Sentinel-1A影像，采用SBAS-InSAR技术实现对矿区的长时间序列的地表沉降监测，将获得的沉降结果作为支持向量回归(support vector regression，SVR)的学习训练样本，采用AGA对SVR模型中的参数进行寻优并得到预测结果。

1 矿区沉降监测与预测原理

1.1 SBAS-InSAR基本原理

SBAS-InSAR技术是在传统InSAR技术的基础上发展起来的，其主要思想是采用多主影像的方式，通过设置合理的时间与空间基线阈值将SAR影像组合成若干个小基线集合从而获取干涉对，其次在每个子集内利用最小二乘法(least squares，LS)求解对应的地面沉降值，对不同子集间采用奇异值分解法(singular value decomposition，SVD)进行联合求解，可以获取观测时间段内矿区的沉降时间序列。

假设在时间t0，t1，…，tn内获取了N+1幅同一研究区域的SAR影像，选择一幅作为主影像，将其余影像统一配准至该影像上，设置合适的时空基线后生成M幅干涉图，M满足式(1)。

(N+1)/2≤M≤N(N+1)/2

(1)

对于由ta和tb(ta

(2)

式中：λ代表雷达中心波长；dtb、dta分别代表以t0为初始时刻，对应于tb与ta时刻在雷达视线方向(line of sight，LOS)的累积形变量；Δφtop，j、Δφatm，j、Δφnoise，j分别代表残余的地形相位差、大气相位差、噪声相位差。去除形变量以外的相位后，干涉相位简化由式(3)至式(4)表示。

(3)

dtb-dta=νi(tb-ta)

(4)

式中：νi代表ta至tb时间段雷达视线方向平均形变速率。由此解缠后的差分干涉图相位可由矩阵表示为式(5)。

Aν=δφ

(5)

式中：A是一个M×N的秩亏矩阵，需利用SVD法进行求解，从而对得到的各时间段内的平均速率进行积分最终获取时间序列形变量。考虑到矿区地表以下沉为主，且水平移动对雷达视线方向的贡献远小于下沉值，所以可将LOS向形变dLOS直接转换为垂直沉降值W(式(6))。

W=dLOS/cosθ

(6)

式中：θ为雷达入射角[18]。

1.2 SVR算法原理

SVR是建立在SVM基础上的回归算法，在预测回归方面具有良好的泛化能力。其基本思想是通过核函数定义的非线性特征，将低维样本数据映射至高维空间中，并在高维空间中对样本进行线性回归，获取最优决策函数f(x)(式(7))。其中，输入样本为(x1，y1)，(x2，y2)，…，(xn，yn)。

f(x)=ωφ(x)+b

(7)

式中：ω为权值矢量；x为输入变量；b为阈值。其优化问题用式(8)表示。

(8)

满足：

(9)

引入拉格朗日函数，并将式(9)转换为对偶形式，同时选取径向基函数作为SVR的核函数，满足KKT(Karush-Kuhn-Tucker，KKT)约束条件后，可求解出SVR回归函数模型(式(10))。

(10)

(11)

1.3 自适应遗传算法参数寻优

在SVR模型中，C、ε、g3个参数的选取对预测精度有很大的影响，对SVR进行参数寻优能有效避免算法陷入局部最优解。首先，建立SVR目标函数与AGA适应度函数的关系，在全局范围内求解目标参数的最优解；然后，随机创建初始种群并开始最优解迭代搜索，计算种群个体的适应度；最后，以最大适应度为标准自动调整交叉和变异的概率并生成下一代种群。循环迭代这一过程直到满足优化终止条件，得到最优个体。算法流程如图1所示。

图1 自适应遗传算法优化流程图

适应度是描述种群个体的主要指标，根据适应度的大小对个体进行优胜劣汰。适应度函数是根据目标函数确定的判别种群中个体优劣性的评价标准。选择SVR目标函数的倒数作为适应度函数，适应度函数表示为式(12)。

(12)

式中：e为避免分母为零所加的一个常数。

遗传算法的主要操作包括选择、交叉、变异等。选择操作对种群中的个体按适应度大小进行排序，根据个体的适应度从上一代种群中选择优良个体复制到下一代。交叉是遗传算法的核心操作，即模拟遗传学的杂交原理，随机结合种群中的2个个体，交换部分染色体从而形成新的个体。变异操作以变异概率改变个体的基因组，从而维持种群的多样性，避免算法过早收敛。交叉概率Pc和变异概率Pm分别如式(13)、式(14)所示。

(13)

(14)

式中：fmax为每代种群的最大适应度；favg为每代种群的平均适应度；f′为要交叉的2个个体中较大的适应度；f为变异个体的适应度。

自适应遗传算法的核心思想是在算法前期保持较大的交叉和变异概率，以丰富种群的多样性；在算法后期降低交叉和变异的概率，保证优良的个体不受破坏。同时，对于适应度高于种群平均适应度的个体，减小Pc和Pm，使得优良个体能够进入下一代；对于低于平均适应度的个体，增大Pc和Pm，淘汰劣质个体。

2 实验结果与分析

2.1 研究区概况和实验数据

李家壕煤矿位于内蒙古自治区鄂尔多斯中南部，地处鄂尔多斯高原腹地。井田东西长约8.0 km，南北约9.0 km，于2011年5月开始试生产。井田构造为向西南倾斜的单斜构造，倾向220°～260°，煤层倾角小于5°，采用综合机械化采煤方法进行综采，其中12113工作面宽度300 m，采高2.3 m，埋深170 m，2017年可推进长度630 m，采煤量可达52万吨。本文数据来源为哨兵1号卫星。哨兵卫星是欧洲航天局哥白尼计划中的地球观测卫星，能够获取C波段的雷达影像，其中Sentinel-1A重访周期为12 d。选用2017年10月2日至2018年4月24日共17景Sentinel-1A影像，影像模式为干涉宽幅模式IW(interferometric wide)，分辨率为5 m×20 m，条带宽度240 km，极化方式为HH。同时，使用精密轨道文件去除平地效应和大气相位的影响。

2.2 时序InSAR监测及结果分析

时序InSAR通过对多时相的SAR影像进行参数估计，可以有效减少平地效应和大气延迟的影响。SBAS技术是一种时序技术，按照空间基线的大小将所有数据组成若干个集合，每个集合的时间序列形变量可以通过最小二乘法求得，使用奇异值分解将多个小基线集进行联合求解就可以获得整个观测周期形变信息，最终利用入射角将视线向形变量转换为垂直沉降信息。

图2表示2017年10月至2018年4月时间段该矿区年平均沉降速率图，该工作面是自东向西开采。由图2可以看出，整个工作面自东向西呈现逐渐下沉的趋势，这是因为东侧要比西侧先开采，同时比西侧更早进行回填。最大沉降速率为360 mm/a。

图3表示以2017年10月2日为时间参考基准的累计沉降量时间序列图。由图3可知，随着时间变化，该工作面的沉降量不断增大，同时沉降范围沿着工作面开采方向由东向西不断扩大，形成明显的沉降漏斗，截至2018年4月24日，最大的累计沉降量为163 mm。

图2 年平均沉降速率图

图3 研究区时序累计沉降量

为了量化分析该工作面的时序累计沉降值，在该工作面中心做剖面提取若干像点进行分析，图4为走向时序累计沉降值折线图。由图4可知，该工作面自东向西进行开采并随着时间的累计沉降值不断变大，出现一个快速下沉盆地，开采中心沉降值从开始的20 mm以360 mm/a平均速率增加到163 mm，说明该工作面还未稳定，随着时间推移沉降值还会不断增大且沉降中心持续向西移动。结合矿区资料可知，自2018年1月起对该工作面开采结束的部分进行回填，导致沉降的区域有所抬升。故图4中 35至55像元处出现明显的抬升区域。

图4 走向时序累计沉降

为了验证SBAS-InSAR技术在矿区形变监测的可靠性，将GPS获得的三维形变投影至雷达视线方向并于与SBAS-InSAR监测结果进行对照。由于GPS的观测数据表示一个点的形变，而InSAR的监测结果表示一个像元分辨率大小的形变，当地表变形不剧烈时，可认为GPS和InSAR获得的沉降值近似相等[19]。图5、图6为SBAS-InSAR技术和GPS技术在观测点的形变序列结果，G1点GPS与InSAR获得的沉降量最大绝对误差、均方误差分别为4.7 mm和3.4 mm；G2点最大绝对误差、均方误差分别为4.8 mm和3.6 mm，符合矿区沉降监测的要求[20]。

图5 G1点雷达视线向SBAS-InSAR与GPS形变结果对比

图6 G2点雷达视线向SBAS-InSAR与GPS形变结果对比

2.3 矿区沉降预测及结果分析

本文取17景哨兵1号影像，共获取16期沉降量，将前12期数据作为训练样本，后4期数据作为测试样本进行预测。采用AGA对SVR模型中C、g、ε3个参数进行优化，初始种群个数为30，最大进化次数为100，C的寻优范围为[0，100]，g的寻优范围为[0，0.01]，ε的寻优范围为[0，0.1]。迭代100次后，AGA优化后的参数C=80，g=0.001，ε=0.05。为了验证本文预测模型的可靠性，与GM、SVR、GA-SVR 3种算法进行对比，预测结果如图7所示。

图7 各模型预测结果

由图7可知，在前6期的观测中观测点位沉降较为缓慢，此时该工作面刚开始开采，从第7期开始，该点位随着工作面的掘进沉降速率逐渐增大。4种预测模型对于前12期观测数据拟合较为理想，但从13期开始，GM和SVR模型在预测时出现了较大的偏差，沉降速率变化趋势与原始值不相符，而AGA-SVR模型与原始值反映的沉降速率基本一致，能够更好地预测矿区地表沉降。各模型预测精度如表1所示。由表1可知，GM预测结果较差，13期和15期2期预测结果的绝对误差均大于10 mm，而SVR模型的15期预测结果绝对误差超过15 mm，与SBAS-InSAR结果偏差较大。与未经参数优化的SVR算法相比，GA-SVR模型预测效果相对较好，但仍存在过早收敛等问题，采用自适应遗传算法进行参数寻优的SVR算法取得了良好的预测结果，4期数据的绝对误差均在5 mm以内，与GA-SVR算法相比，绝对误差分别减小了37.56%、26.31%、24.73%、31.40%。

表1 各预测模型精度对比 mm

分别计算各预测模型的平均绝对误差(mean absolute error，MAE)、均方根误差(root mean square error，RMSE)和运行时间，如表2所示。GM和SVR模型虽然运算时间短但精度较差，MAE和RMSE均大于11 mm，GA-SVR模型的精度有所提升，但运算时间也大幅增加，可能与陷入局部最优解有关。相比其他3种模型，本文提出的AGA-SVR模型表现出更优的精度，MAE和RMSE均小于5 mm，同时运算时间相比 GA-SVR 算法有所下降，避免了遗传算法中存在的欺骗问题，提高了参数寻优效率。

表2 各个预测模型精度对比

为了进一步验证AGA-SVR算法的有效性，在图4剖面自沉降边缘至沉降中心等间距选取3个点A1、A2、A3进行预测。剖面上3点预测结果如表3所示。由表3可知，不同沉降区域的3个点位均取得了较好的预测结果，绝对误差均在5 mm以内。

表3 剖面上所选点位预测结果 mm

从时间距离和预测精度的关系考虑，3个点的预测精度均与时间距离成反比，即随着观测期数的增加，各点预测精度逐渐降低。本文对第13～16期数据采用了滚动预测的方法，即先用第1～12期的观测数据预测第13期的数据，然后将第13期的预测结果编入训练样本，用第2～13期的数据对下一期数据进行预测，以此类推完成所有的预测。由于第1～12期的数据均为独立观测获得，而预测第14～16期数据所用的训练样本中均含有预测值，且越往后各个训练样本间的相关性越强，导致预测模型训练结果较差，因此时间距离越短，预测的精度相对越高。

3个点的预测精度如表4所示。由于A3点的累计沉降量最大，因此A3点的相对误差最小。各个点位的MAE和RMSE均在5 mm以内，取得了较好的预测结果，进一步验证了AGA-SVR算法是一种较为可靠的预测方法。

表4 各个点位预测精度

3 结束语

针对传统D-InSAR技术存在的失相干等问题，利用SBAS-InSAR技术对内蒙古李家壕煤矿某工作面进行开采沉降监测，获得了矿区连续地表形变特征，在此基础上使用AGA-SVR模型对矿区沉降值进行预测。监测结果表明，自2017年10月该工作面开始开采，该地区存在连续下沉现象，截至2018年4月，该工作面最大累计垂直沉降量为163 mm，最大沉降速率为360 mm/a。SBAS-InSAR技术监测矿区沉降的最大绝对误差为4.8 mm。与其他预测模型相比，AGA-SVR模型能够避免算法陷入局部最优解，显著提高了预测精度，最大绝对误差为5.0 mm，最大相对误差为4.08%，MAE和RMSE均优于4.5 mm，监测精度和预测精度均满足工程实践需要。

在矿区地表形变监测中，利用InSAR技术获取开采沉陷的影响范围与发展趋势，对于形变剧烈的区域，应考虑在实地布设GPS监测点并定期进行水准观测。同时，建立预测模型，将InSAR技术提取的形变量作为训练样本进行开采沉陷动态预测，当预测结果的形变速率显著增大时，应结合其他监测手段进行分析并对矿区开采活动发出预警，以保障矿区生产活动的平稳进行。下一步的工作将考虑增加数据量并采用不同极化方式的数据以提取矿区的三维形变信息。