土壤水分是连接水循环、能量循环和碳循环的关键水文状态变量，在水文、气象和农业等学科领域中扮演着重要的角色^[1,2]。获取大范围土壤水分的时空动态变化对于作物干旱监测^[3,4]、水资源管理^[5] 以及区域和全球气候变化研究^[6,7]具有十分重要的意义。目前土壤水分监测方法主要有地面站点观测和遥感观测。地面站点观测法只能实现对点尺度的土壤水分观测，并且需要耗费大量的人力和物力，然而土壤水分在不同时空尺度上具有高度的异质性^[8,9]，其大范围的实时监测更多需要依靠遥感技术手段。

微波遥感具有全天时、全天候的观测能力，并且凭借其对土壤介电常数的敏感性，在土壤水分监测中发挥着重要的作用^[10]。微波遥感观测土壤水分主要是基于微波信号与土壤介电性质之间的紧密联系，随着土壤水分的增加，从干土到纯水，介电常数发生3 ~80 F·m^-1的变化^[11]。微波遥感按照观测方式可以分为主动微波和被动微波。被动微波基于微波辐射计探测地物辐射亮温来实现土壤水分的反演，常用的被动微波产品有AMSR-E、AMSE2、SMOS、SMAP(主被动)和SSM/I等。主动微波传感器则通过传感器发射电磁波并接收目标地物回波信号实现。常用的搭载主动微波传感器的卫星有ENVISAT/ASAR、ALOS/PALSAR、RadarSat、ASCAT、Sentinel-1及我国的高分三号(GF-3)等。当前，已发布的土壤水分产品空间分辨率较低，如ASCAT (产品分辨率约为25 km)，AMSR-E (分辨率约为25 km)，这限制了土壤水分微波产品在田间尺度上的应用。Sentinel-1卫星提供高分辨率的合成孔径雷达(SAR) 观测数据(空间重采样可达10 m)。自2016年起，其A/B双星的配置使得部分区域重访时间缩短至6 d。Sentinel-1卫星的这些特性为需要连续时相观测数据的变化检测算法反演高分辨率的土壤水分提供有力的保障。

SAR数据在反演土壤水分过程中需要考虑植被覆盖及地表粗糙度等地表情况对反演精度的影响。通常根据实际情况选用不同的散射模型。对于裸露地表，现有的散射模型主要有经验/半经验模型，如Oh模型^[12]和Dubois模型^[13]；理论模型如积分方程模型 (IEM)和高级积分方程模型 (AIEM)^{[14,15,16,17,18]}。植被覆盖区常用的散射模型主要有水云模型^[17]和MIMICS模型^[20]。当前，大部分方法通过耦合裸土散射模型和植被散射模型来消除地表粗糙度和植被覆盖的影响。这一过程中，通过一系列条件假设对模型进行简化是土壤水分反演的不确定性主要来源。

变化检测算法可以避免微波散射模型中对实际的复杂散射过程进行简化所带来的不确定性。其核心思想是充分利用高时间分辨率的微波数据，这样在较短的重返周期内，连续时段的影像中地表粗糙度和植被覆盖度不会发生明显的变化，可以视作常量。由此判定后向散射系数的变化主要是由土壤水分波动带来的介电常数变化所引起，从而避免直接求解复杂的散射方程。Wagner等^[21]在1999年针对ERS散射计发展了基于变化检测的土壤水分反演方法。该方法更简洁有效地避免了复杂参数和辅助数据的输入，充分利用了卫星提供长时间序列数据的优势。多项研究对此算法进行了验证和改进。万幼川等^[22]使用SCAT和ASCAT卫星数据，采用移动时间窗口函数替代原有的经验函数描述σ'(θ)(后向散射关于入射角θ的一阶导数)的季节性变化，摆脱了经验函数的局限性，使算法具备良好的精确度和可移植性。Zribi 等^[23]使用高分辨率的ENVISAT ASAR数据将变化检测方法应用于半干旱区域，发现土壤相对湿度指数的均方根误差（RMSE）为0.13(约为0.035 m³/m³)。Pathe等^[24]以美国俄克拉荷马州为研究区，证实了ENVISAT ASAR数据在变化检测方法中的适用性。Gao等^[25]使用近期发射的Sentinel-1卫星数据反演土壤水分，并对此变化检测方法进行简化，使用两个连续时相的土壤水分变化值表征当前与最干燥情况下的后向散射系数差，最后得到较好的反演效果(RMSE为0.059 m³/m³)。

已有的研究清晰地展示了将变化检测算法运用于高分辨率的SAR数据进行田间尺度土壤水分反演的良好潜力。但在实际应用中，该方法仍有几个关键问题需要解答。①后向散射时序差值-植被指数 (Δσ-VI)空间中所选取的采样点比例是否会影响反演的精度？②选择不同植被指数对变化检测方法构建的Δσ-VI经验关系式是否有影响？③变化检测方法作为一种递推算法，具有误差累积的特性，其估算精度对于关键参数，如土壤水分初始值和土壤水分尺度缩放因子 (即两个连续时相土壤水分变化的最大值ΔMsmax) 的敏感性如何？这些问题的解答有助于研究者选择合适的植被指数产品，并选择合适的参数 (单个站点或区域所有站点的土壤水分观测值，或低分辨率的微波土壤水分反演值)，以构建变化检测模型从而高效反演土壤水分。本文围绕上述问题，采用Sentinel-1卫星SAR数据，结合Landsat 8和中分辨率成像光谱仪（MODIS）两种传感器的归一化植被指数（NDVI）产品，将变化检测算法应用于黑河中游地区，详细探讨变化检测算法中具体参数设置对估算精度的影响。

黑河流域生态—水文过程综合遥感观测联合试验(简称“黑河生态水文遥感试验”)，是围绕黑河流域生态—水文集成研究中的核心科学目标，于2012~2015年在黑河流域开展的卫星和航空遥感及地面观测互相配合的多尺度综合观测试验^[26,27]。

研究区位于黑河生态水文遥感试验区的张掖南部中游人工绿洲试验区(图1)，该试验区为地形开阔的盆地，海拔在1 450~1 680 m之间，地势自南向北逐渐降低。研究区气候类型为温带干旱荒漠气候，据中国气象数据网的统计结果显示，该区域多年年平均降水量为132.6 mm，年降水量变化范围71.6 ~ 216.3 mm；多年年平均气温为7.8 ℃，年均温变化范围0.9 ~ 16.2 ℃。试验区北部为湿地，南部为荒漠，西部为戈壁，东部为沙漠，中部为人工绿洲，绿洲中有农田、防护林、村庄、道路及渠系等多种地类，主要作物为制种玉米和大田玉米等^[28]。本文选取试验区中部45 km × 44 km的人工绿洲作为研究区，介于38.7°~39.2° N，100.2°~100.7° E之间。研究区站点信息如表1所示。

Sentinel-1卫星是欧洲空间局 (ESA) 在“哥白尼”计划中发射的卫星，包括Sentinel-1A和Sentinel-1B，分别于2014年和2016年发射。Sentinel-1卫星搭载C波段 (频率为5.405 GHz，波长为5.6 cm)合成孔径雷达，干涉宽幅模式(IW)空间分辨率为5 m × 20 m，单颗卫星重访时间为12 d，A/B双星重访时间为6 d，采用VV和VH双极化方式。诸多研究表明，VH交叉极化对植被变化更为敏感，而VV极化对于植被和土壤水分均具有敏感性，后者在土壤水分反演中具有更大的潜力^[29,30,31] 。因此本研究只选用VV极化后向散射信号。

哨兵数据产品可以从欧空局Sentinel-1数据网站上下载 (https://scihub.copernicus.eu/)。本文使用了覆盖研究区的GRD产品影像，成像模式选取IW模式 ，共选取2017年1月至12月间的46景影像。

数据具体处理步骤包括：①应用精密轨道文件，提高影像配准及地理编码的精度；②热噪声去除，消除雷达系统内部的热噪声；③辐射定标，减少辐射偏差；④使用Sigma Lee滤波方法进行斑点滤波，消除相干斑噪声对SAR图像的干扰；⑤地形校正，采用30 m 的SRTM DEM改正地物在雷达影像中发生的畸变，如透视收缩、叠掩等。

数据预处理全部完成后, 从影像中提取出线性化至dB单位的后向散射信号。

Landsat 8卫星于2013年2月发射，携带热红外传感器 (TIRS) 和陆地成像光谱仪 (OLI)。OLI有11个波段， 其中全色波段分辨率为15 m，多光谱波段分辨率为30 m，重访时间为16天。

本文选取了2017年1月至12月期间包含研究区的20景影像 (https://earthexplorer.usgs.gov/)，使用Savitzky–Golay (SG) 滤波平滑信号，消除噪声^[32]。 Chen等^[33]通过大量试验表明，对NDVI而言，采用窗口大小为9，滤波阶数为6的SG滤波效果最佳，故在此研究中，选用该SG滤波参数对原始数据进行滤波处理。由于Landsat 8卫星与Sentinel-1卫星过境时间不相同，本研究利用双线性插值法将Landsat 8数据插值为日数据，并根据Sentinel-1卫星过境时间选取对应的数据。本研究使用Landsat 8 红光波段和近红外波段计算NDVI。具体计算公式如下：

其中:Red和NIR分别为红光和近红外波段的反射率。

搭载于Terra卫星上的MODIS 传感器获取的NDVI数据由美国国家航空航天局NASA提供 (https://search.earthdata.nasa.gov/search)。该数据集来源于MOD13Q1产品，空间分辨率为250 m，时间分辨率为16 d。数据时间覆盖度为2017年1月至12月共23景影像。同样使用SG滤波进行异常点滤除及插值，滤波参数及插值方法与2.2.2节中Landsat 8 NDVI 相同。

本研究地面实测数据集来源于黑河生态水文遥感试验中的流域水文气象观测^[26,34]，研究站点为大满超级站及黑河遥感自动气象站，站点信息见表1^[35,36,37]。

站点土壤水分探头埋设在2、4、10、20、40、80、120和160 cm处，观测频率为每10 min一次，本研究中使用4 cm处的土壤水分观测值作为验证数据。

本文所用的变化检测方法的核心思想在于将两个连续时相的雷达后向散射系数差归结为土壤水分的变化，然后通过递推的方式求得土壤水分。研究时段内Sentinel-1卫星数据重访时间为6 d，假设地表粗糙度在重访周期内保持不变，植被变化在两个连续时相内变化也不明显，则两者在短时间内均可视为常量。因此获取两个连续时相的雷达后向散射系数之差可以消除地表粗糙度的影响，即雷达信号变化表示为土壤水分变化的函数。两个连续时相t和t+1的雷达后向散射系数差异可以表示为：

其中，在特定的植被指数VI下，σVI0t+1是在时间t+1的后向散射系数，σVI0t是在时间t的后向散射系数。

两个连续时相的雷达后向散射系数的差值（变化）分为正向变化和负向变化，两种变化的绝对值都随着VI增大而减小，所以雷达后向散射系数时间变化的敏感度随VI增大而下降。理论上，正向和负向变化的趋势是一种对称、线性模式^[25]。

为了获取VI和连续时相后向散射系数差的线性关系，取每个VI对应的一定比例的后向散射系数观测作为样本点，用于模拟VI与对应的后向散射系数变化最大值之间的关系。雷达后向散射系数变化的最大值，可以用一元线性函数来描述：

即：

其中：Δσmaxbare是两个连续时相裸土的后向散射系数差，对应于最大的土壤水分变化。bVI是指植被对后向散射系数的影响，当VI为0时，则ΔσmaxVI等于Δσmaxbare，b是经验函数f的斜率。

如前所述，变化检测算法发展的重要假设之一是两个连续时相的雷达后向散射系数的差可以表示为土壤水分差异的函数，因此土壤水分变化可以表示为：

其中：ΔMsmax是观测期间内两个连续时相之间的土壤水分的最大变化，即土壤水分缩放因子，该值可以从站点的实测数据中得到，也可以通过卫星土壤水分产品获取其近似估计值。函数S(t,t+1)是指时间t和t+1之间的土壤水分变化。最后通过递推的方法，从连续时序的第一天递推得到其后的土壤水分，公式如下：

……

其中：Ms (t) 是时间t的土壤水分反演值，Ms (1) 为土壤水分的起始值，理想情况下可以通过地面实际观测得到，也可以通过微波产品获取其近似值。

反演土壤水分的过程中，在后向散射系数差-植被指数 (Δσ-VI)空间中选取不同的采样点比例进行线性建模对反演精度具有一定的影响。为了量化这一影响，研究中设置20个不同等级的采样点比例，从1%的起始比例按1%的步长递增至20%。将Landsat 8 NDVI和MODIS NDVI这两种植被产品，分别按照这20个等级的采样点比例进行建模，并进行结果的验证。反演误差与采样点比例的关系如图2所示。

从图中可知，对于MODIS NDVI和Landsat 8 NDVI两种植被指数，随着采样点比例的增高，土壤水分反演结果RMSE都是先减小后增大，分别在2%和4% 时RMSE达到最小值，其后RMSE随着Δσ-VI采样点比例的增大而增大。值得注意的是，RMSE达到最小值时对应的最优采样点比例并非固定的，而是随ΔMsmax的变化而略有变化(但都低于5%)。

本文选取了Landsat 8 NDVI和MODIS NDVI两种植被产品，根据公式(4) 构建变化检测的经验模型，Δσ-VI采样点比例由4.1节的最优结果确定。由于雷达后向散射系数变化有正向变化和负向变化，所以拟合的一元线性函数分为上函数（负向变化值的拟合结果）和下函数（正向变化值拟合结果）。图3展示的是两种指数的VI和雷达后向散射系数差值Δσ的经验函数。

从图3可知，雷达后向散射系数变化随着植被覆盖的增加而不断衰减，拟合的上函数和下函数都呈近似对称的线性关系。此外，两种植被指数构建的Δσ-VI拟合关系式都通过p = 0.01的显著性检验，但两套拟合关系之间也有所差异。首先，从拟合模型来看， MODIS NDVI 构建的上函数和下函数具有明显的对称性，拟合的直线斜率绝对值近似相等(图3（b）),而Landsat 8 NDVI的拟合直线对称效果稍差(图3（a）)。其次，由于MODIS植被指数变化范围较Landsat 8更大(研究区MODIS NDVI值范围约为0~0.9，而Landsat 8 NDVI变化范围为0~0.6)，故MODIS NDVI拟合的线性函数斜率较小。

基于4.2节采用的两种植被指数 (Landsat 8 NDVI、MODIS NDVI) 所构建的后向散射系数随植被指数变化的线性方程，分别进行2017年研究区土壤水分反演，其估算值和实测值的散点图如图4所示。总体上，Landsat 8 NDVI 反演结果略微优于MODIS NDVI，两者均方根误差RMSE分别为0.040 m³/m³和0.044 m³/m³,相关系数R分别为0.86和0.83。这是由于Landsat 8 NDVI的空间分辨率较MODIS NDVI更接近Sentinel-1，故其构建的变化检测模型更能反映真实的植被-后向散射系数的变化关系。但由于Landsat 8数据集在研究时段内的完整性较MODIS更差，故在本方法进行VI插值过程中将引入更大的误差。因此，根据研究区的特点选取合适的植被指数数据集，有助于得到更为准确的土壤水分反演结果^[38]。值得注意的是，根据两种植被指数构建的变化检测模型，土壤水分估算偏差 (Bias)均为0.009 m³/m³，在水分的高值区和低值区都没有明显的系统偏差。相比Gao等^[25]在西班牙的研究，本研究选取了更大的研究范围(45 km×44 km)，并且得到了更优的结果，证实了本方法利用Sentinel-1数据反演高精度土壤水分信息的巨大潜力。

变化检测方法中，土壤水分反演是一个递推的过程，其中有两个关键的辅助参数，分别是土壤水分反演的起始值Ms1，以及两个连续时相土壤水分变化的最大值，或称为土壤水分缩放因子ΔMsmax。在最理想的情况下，Ms1和ΔMsmax都来源于站点观测的“真实值”，将这种理想情况称为方案一。而实际的土壤水分反演流程中，两个参数通常无法同时获取，考虑到这一实际情况，研究中设置另外两个对比方案。方案二中，假设ΔMsmax无法由点尺度的观测信息获取，只能通过低分辨率的微波土壤水分产品(3 km Level-2 SMAP/Sentinel-1L2_SM_SP数据) 进行估算。方案三，假设初始时刻的土壤水分真实值Ms1无法准确获取，同样以初始时刻的SMAP/Sentinel-1 土壤水分作为变化检测算法的起始值。以上这3种方案都是基于4.1~4.3节中选出的Landsat 8 NDVI在4%的最优(Δσ-VI)空间采样点比例下构建的经验模型。

图 5 是3种方案中土壤水分实测值与反演值的散点图，与预期一致，方案一的精度最高 (RMSE为0.04 m³/m³)，方案三均方根误差次之，方案二误差最大。由此可见，土壤水分缩放因子ΔMsmax对反演精度的影响大于土壤水分初始值Ms1的影响。方案二是在方案一的基础上，将土壤水分缩放因子的理想“真实”值替换为3公里Level-2 SMAP/Sentinel-1土壤水分产品，除了尺度差异带来的不确定性，还引入了SMAP产品约为0.04 m³/m³的误差^[39]，最终方案二RMSE约为 (0.04 + 0.04) m³/m³= 0.08 m³/m³，相对最优的方案一反演效果降低最为明显。

图6展示了大满站和黑河遥感站采用3种方案反演得到的土壤水分时间序列。首先通过对比方案一和方案三来研究引入初始值误差的影响，由于大满站实测土壤水分的初始值与SMAP/Sentinel-1初始值大致相等，故大满站的方案一与方案三反演得到的时间序列基本重叠。而黑河遥感站中SMAP/Sentinel-1获取的土壤水分初始值比实测值高0.021 m³/m³，故黑河遥感站的方案三相对于方案一有稳定的系统偏差，该偏差的大小明显取决于土壤水分初始值的误差，且通过起始值这一参数引入的误差不会造成土壤水分反演误差的累积。对于缩放因子引入的误差 (对比方案一和方案二)，两个站点的土壤水分反演误差时间序列都具有较大的波动，且在水分值较高时具有较大的误差。此外，反演误差还随着缩放因子误差的增大而明显增大，大满和黑河遥感两个站点的缩放因子误差分别为0.11 和0.06 m³/m³，由此造成的两个站点的土壤水分反演误差则为0.058和0.025 m³/m³。此外，反演的土壤水分时间序列也进一步说明，土壤水分缩放因子对反演结果的影响明显大于土壤水分初始值对反演结果的影响。

采用Landsat 8 NDVI植被指数构建线性方程，并从2017年3月20日3公里Level-2 SMAP/Sentinel-1数据获取土壤水分起始值Ms1和缩放因子ΔMsmax，应用于公式(5)~(7)进行2017年6月29日黑河中游土壤水分反演(共进行17次递推)。最后通过双线性法将反演的5 m×20 m土壤水分重采样得到研究区范围内100 m×100 m土壤水分空间分布(图7(a))。作为对比，图7(b)展示同一天的3公里Level-2 SMAP/Sentinel-1土壤水分空间分布。该日期内研究区土壤水含量介于0.02 ~0.35 m³/m³之间，符合半干旱区的地表水分状态特征。

总体而言，黑河中游人工绿洲试验区土壤水分空间异质性明显，土壤水分含量呈中 部高、四周低的趋势，高值区主要分布在中部和东部的农业地带，土壤含水量介于0.25~0.35 m³/m³之间。这是由于中部农田西侧有黑河的支流经过，给中部农田提供充足的水分，且研究区中部植被覆盖率较高，使土壤汇流和截流能力高于四周地区。低值区则主要位于四周的荒漠戈壁，东北部与西南部的荒漠区土壤含水量总体上低于0.15 m³/m³。相较于3 km的SMAP/Sentinel-1产品，研究区高分辨率的土壤水分反演结果能更加精细地刻画土地利用类型及人类活动对该区域土壤水分空间异质性的影响。

本文基于欧空局Sentinel-1数据，以及相同时段的MODIS、Landsat 8 NDVI植被产品构建变化检测模型，估算黑河中游土壤水分，探讨变化检测算法中具体参数设置对土壤水分估算精度的影响，得出以下主要结论：

（1)在构建后向散射系数差-植被指数的线性模型过程中，MODIS NDVI和Landsat 8 NDVI两种植被指数所建模型的反演精度都会受到所选取的采样点比例的影响，但都在采样点比例低于5%时取得最高精度。

（2)总体上，Landsat 8 NDVI构建的变化检测模型略优于MODIS NDVI构建的模型，这主要是由于Landsat 8 NDVI的空间分辨率更接近Sentinel-1，故其构建的后向散射系数差-植被指数空间更准确。实际应用中应根据研究区的特点和研究目标选择合适的植被指数。

（3)对于变化检测算法的两个关键参数——土壤水分初始值和土壤水分缩放因子，缩放因子这一参数的误差将带来土壤水分估算精度更为显著的降低。