综述：结合InSAR技术模拟矿区地表形变速率衰减时间

《Canadian Journal of Remote Sensing》：Combination of InSAR Techniques for Modeling The Decay Time of Surface Deformation Rates in Mining-Impacted Areas

本文综述了利用合成孔径雷达干涉测量（InSAR）技术，结合Knothe时间函数，预测矿区在停止开采后地表形变（主要为沉降）衰减时间的新方法。研究在波兰的Legnica-Glogow铜矿带开展，通过小基线子集（SBAS）方法处理Sentinel-1 SAR数据，获得了高分辨率的形变时间序列。结果表明，地表变形在采矿停止后最迟可持续约8年，且形变衰减时间与开采期间测得的地表沉降值成反比。该研究为后矿区土地利用规划提供了更精确的空间信息，相较于传统大地测量方法，能更精细地分析地表变化的动力学过程。

　　

地下矿产开采活动总会对地表地形产生影响，导致岩石层迁移、地下形成空洞，并在地表形成沉降槽，对建筑基础设施造成破坏性影响。这种地表变形通常用沉降（w）、倾斜（T）、曲率（K）和水平应变（ε）等指标来描述。研究表明，约80%的直接影响变化在开采过程中已在地表显现，其余20%的变化则是后续过程的结果，传统大地测量估计这些过程发生在停止开采后的二到五年内，之后是长期的残余变形过程。然而，传统的精密大地测量方法（如水准测量和GNSS测量）在矿山关闭后往往停止，限制了我们对后矿区地表变化动力学的理解。

近年来，干涉合成孔径雷达（InSAR）方法的应用日益广泛。这项技术能够大范围、高精度地监测地表位移，即使是在矿山关闭之后。已有研究利用InSAR方法确认了采矿变形在矿场关闭后仍会持续发生，例如在比利时、德国、荷兰、法国、英国和波兰等多个矿区都观测到了显著的地表位移。

地表形变预测与Knothe函数

然而，采矿区和后矿区的安全在很大程度上依赖于掌握变形过程随时间推移的动力学信息，即对这些变形的预测。这类预测基于描述采矿变形动力学的理论模型，例如Knothe模型。在该模型中，采矿活动停止后形变发生的时间是使用Knothe时间函数计算的，其关键参数是时间系数c。这个系数描述了采矿影响在岩体中的迁移速率，因而反映了岩体的特性。传统上，系数c是基于开采期间进行的大地测量，通过最小二乘法确定的，但这种方法信息有限，且通常在空间上分辨率不足。

近期，Apanowicz等人提出了一种创新方法，将永久散射体干涉测量（PSI）时间序列方法与Knothe时间函数相结合来确定时间系数c。他们的研究表明，基于SAR数据可以比经典大地测量更精确地确定系数c，从而能更准确地估算矿场关闭后的形变衰减时间。本文的研究正是基于这一方法论的扩展和改进。

研究区域与方法

本研究在波兰的Legnica-Glogow铜矿带（LGCB）进行。该地区铜矿资源丰富，目前由三个矿山在深度超过1200米处开采。采矿作业采用房柱式后退开采法，并在需要保护地面基础设施的区域进行水力充填。该地区同时存在采矿相关变形、矿震和岩体排水引起的变形，使其成为研究采矿活动对地表变形影响的宝贵区域。

研究采用了改良的方法来确定Knothe函数的时间系数c。首先，利用小基线子集（SBAS）方法处理了2020年1月至2022年10月期间获取的Sentinel-1 SAR数据，计算了视线向（LOS）位移。数据处理使用了GMTSAR软件包，并进行了二维分解，以获取垂直方向和东西方向的形变分量。由于Sentinel-1卫星近极地轨道几何的限制，雷达测量对南北方向的位移敏感性有限，因此主要关注垂直位移。SBAS处理直接产生了整个研究期间的累积地表位移时间序列，随后根据改良的方法，使用最小二乘法从SAR数据结果中确定了Knothe函数的时间系数c。

地表形变衰减时间的空间分析

基于从SBAS处理的SAR图像中得出的地表沉降值，并利用Knothe时间函数，对研究区域的理论形变衰减时间进行了分析。时间系数c及其标准差是从沉降值的幂律趋势模型中确定的。分析聚焦于垂直沉降超过0.2米的区域，这对应于该研究地点采矿区域分类中的1类区（理论上的最低类别，认为地表变形不可忽略）。研究区域面积约为1.05平方公里，共包含853个像元。

估算出的时间系数c值范围在0.55到2.04之间，绝大多数情况下的标准差在系数值的10%到30%之间。假设99%的最大地表沉降已经发生（这是沉降建模中常用的理论渐近线），计算出采矿停止后的地表形变衰减时间从2年3个月到8年4个月不等。研究发现，在沉降最大（约40厘米）的位置，衰减时间最短（约两到三年）。随着测得的沉降值减小，沉降过程持续增加的时间则变长。超过5年的最长沉降持续时间，出现在沉降槽东部边缘、北部东西边缘、西部边缘以及盆地周边个别位置的小区域。这表明，形变衰减时间与开采期间测得的地表沉降值成反比。

时间系数c的分布分析：χ²检验与Kolmogorov-Smirnov检验

Knothe理论基于对采矿引起的垂直地表沉降采用高斯分布的假设。因此，有理由推测Knothe时间函数的参数c也近似遵循高斯分布。对确定的853个系数进行的统计分析显示，其分布与高斯分布一致，大多数值聚集在均值（和中位数）1.29附近。为验证这一假设，在显著性水平α=0.05下进行了两项统计检验。

首先进行了卡方（χ²）独立性检验，以验证时间系数c的经验分布与理论高斯分布的独立性。检验结果显示，χ²统计量之和为33，大于自由度为5、α=0.05时的临界值11.07。因此，拒绝原假设（认为分布独立），接受备择假设，即时间系数c的分布是高斯分布。

随后进行了Kolmogorov-Smirnov检验，以验证经验分布是否就是高斯分布。计算出的检验参数D为0.033，经验参数D√n小于临界值λ_k。因此，不能拒绝原假设。基于Kolmogorov-Smirnov检验，时间系数c的分布是高斯分布。

两项统计检验的结果表明，所分析的时间系数c统计样本遵循高斯分布。此外，还计算了Cramér’s V系数以评估经验分布与理论正态分布的拟合强度，其值为0.20，表明拟合强度较弱，这在分布图中也可观察到。

结果验证

研究结果与位于研究区域内的一条由13个测量点组成的水准线进行了比较。首先，根据每个点实测的沉降值，确定了时间系数c及其相应的标准差。接着，在同一位置，通过SBAS方法处理SAR图像计算了沉降值。从SBAS结果中，分别基于SBAS推导的沉降值、对数趋势模型和幂律趋势模型，三次计算了系数c。计算再次使用了最小二乘法和Knothe函数。

基于水准测量确定的系数c值在0.79到1.97之间，但部分点的标准差过大（超过系数值的50%），这些点被排除在比较分析之外。从原始SBAS时间序列数据确定的系数c，其标准差值比系数值大数倍，因此结果不可靠。从对数趋势模型确定的系数c在所有测量点均为0.31，标准差为0，这表明变异性无法量化，结果同样不可靠。从幂律趋势模型确定的系数c范围在0.77到2.31之间，标准差在0.24到0.30之间，相对偏差（占系数值的百分比）大多数情况下不超过28%，因此被认为是三者中最可靠的，并用于比较分析。

从大地测量结果和SBAS时间序列结果确定的Knothe函数时间系数c值之间的差异，从点6127的0.02到点6128的0.19不等。沿线测量点的均方根误差（RMSE）为0.23。所有点的时间系数c平均值，大地测量数据为1.29，SBAS数据为1.24。基于这些平均值，并假设99%的采矿沉降已经发生，估算出的研究区域形变衰减时间分别为约3年7个月（基于大地测量）和约3年8个月（基于SBAS），两者仅相差约一个月。这种微小差异表明了结果的高度一致性，并证实了所应用模型的适用性，在形变时间预测的背景下是一个非常理想的结果。

结论与讨论

本文展示了SAR数据和Knothe时间函数在采矿作业停止后理论形变衰减时间空间分析中的应用。研究采用了基于PSI方法处理SAR图像的既定方法，但为了适应SBAS方法进行了修改。这些修改包括处理两个相对轨道（一个升轨，一个降轨）的SAR图像，并分析地面变形的垂直分量以便与大地测量结果进行更有潜力的比较。同时，由于SBAS方法直接产生位移时间序列，因此省略了原方法中与PSI方法特点相关的步骤。此外，SBAS方法能够分析更大区域，这优于主要局限于城市区域的PSI方法。

对垂直地表沉降超过0.2米区域的形变衰减时间分析表明，估算的衰减时间范围从2年3个月到8年4个月，这从理论上证明了变形在采矿活动结束最多约8年后才可能显现。以往基于大地测量的研究通常提及的是5年周期，差异可能源于大地测量在空间上的局限性，以及其通常在地表变形已经显著显现后才开始监测。本研究结果表明，结合InSAR方法与Knothe时间函数，可以比使用传统大地测量方法更精确地分析地表变化的动力学。这对于后矿区土地规划、风险评估和基础设施建设具有重要的现实意义。