基于陆基L波段辐射计的青藏高原地表冻融反演研究

Based on a Relative Frost Factor Algorithm Retrieving Landscape Freeze/Thaw State by Ground-Based L-Band Radiometer

摘要: 青藏高原作为中低纬度地区面积最大的季节性冻土分布区，其地表冻融在气候变化及反馈中扮演了关键性角色。L波段被动微波辐射计观测作为最具潜力的全球尺度冻融反演技术还未被评估其在青藏高原的适用性，因此，本文基于高频时次的陆基辐射计观测配合详尽的站点土壤监测，尝试以相对冻结因子阈值判断算法为基础算法，并详细探讨了应用多种微波极化指数在青藏高原典型季节性冻土区的适用性。结果表明：针对青藏高原季节性冻土，最适采用垂直极化指数作为相对冻结因子阈值判断算法的输入指数，其总精度能够在60˚入射观测角度时达到0.92，而极化差异指数NPR作为最广泛应用的判别指数，并不适合青藏高原地区。本研究为青藏高原地区开展L波段星载微波辐射计反演土壤冻融提供了客观的极化指数选取依据。

Abstract: Soil freezing dramatically alters the soil hydraulic and thermal properties that in turn affect surface energy and water budgets, net ecosystem exchange of carbon and the surface runoff over grassland. Passive microwave observations at L-band (1 - 2 GHz) benefits from a relatively high contrast in the dielectric constant between free water and ice with lower scattering by surface vegetation. Thus, this study uses a relative frost fac-tor algorithm with four polarization ratios to detect temporal changes of soil freeze/thaw (F/T) state based on a soil temperature and freeze/thaw synthesis experiments conducted over the north-eastern part of the Tibetan Plateau. The L-band brightness temperature and the in-situ soil temperature are measured simultaneously for the purpose of evaluating performances of the relative frost factor algorithm on different incident angles. Results of estimated soil F/T state indicate that the algorithm corresponding to Vertical Polarization index (Vpol) achieves the highest preci-sion (0.92) with incident angle of 60˚. NPR index, as wide and successful applications for frost factor algorithm, is still needed to be optimized to decrease its degree of confusion in discriminating F/T state. The evaluation of the performances in this study is useful for selecting optimum polarization index for spaceborne L-band radiometers retrieving soil freeze/thaw state in the seasonally frozen Tibetan grassland.

1. 引言

土壤冻融过程涉及土壤中液态水和固态冰之间的相态交替作用，在陆面过程中显得极为独特 [1]。青藏高原作为地球第三极，由于其独特的高海拔气候条件，成为中低纬度地带面积最广、厚度最大的冻土分布区，其中季节性冻土面积约为1.2 × 106平方千米 [2]。与北美、俄罗斯及极地地区的冻土相比，青藏高原的冻土温度较高且厚度较薄，在气候变化中更为敏感 [3]。土壤的冻融过程通过影响陆面能量与水分收支和碳氮循环等过程，在气候变化反馈中扮演了关键性角色 [4]。更重要的是，青藏高原冻融过程的深度和广度，通过陆面能量再分配和陆面碳循环过程，对地气系统各圈层产生了极其深远的反馈作用 [5]。土壤冻融过程的动态变化与气候变化并非孤立存在，其间存在复杂的相互反馈作用 [6]，土壤冻融的动态变化特征不仅是气候变化的指示器，更为关键的是：土壤冻融过程将通过影响湿地、积雪和植被等地表覆被特征直接作用于流域的产汇流过程。因此，理解流域土壤冻融状态随气候变化的敏感度至关重要 [7]。然而，现今可信度较高且时空连续的土壤冻融监测数据十分缺乏。因此，如何确定表征高原土壤冻融过程的关键物理参数，有效地监测青藏高原土壤冻融过程的时空特征就显得尤为迫切 [8]。

尽管国内外的土壤分析监测网络提供了非常有价值的土壤冻融监测数据 [9] 可以用于模式及卫星遥感反演产品的校准评估 [10]，但下垫面的空间异质性使得观测网络仍难以构建准确的区域尺度的土壤冻融过程的时空分布 [11]，现有的观测台站和网络无法满足对流域及区域尺度的土壤冻融过程监测。而随着卫星遥感对地观测系统的建立和应用，大大提高了地球表层系统科学研究的效率，特别是由于星载微波辐射计的发射和长期在轨运行，在全球和区域尺度上，星基微波遥感技术已经具有精确监测和分析陆面水分循环关键参量的潜力 [12]。而土壤中的水分变化和冻融过程对微波辐射低频波段(L波段)的响应更为敏感 [13]，特别是当土壤冻结时，发射率的显著跃升和下降。而随着土壤湿度和海洋盐度卫星(Soil Moisture and Ocean Salinity, SMOS)和土壤湿度主动/被动卫星(Soil Moisture Active Passive, SMAP)的发射，被动微波遥感观测已成为最有潜力的全球尺度土壤冻融观测手段 [14]。

然而，实际应用被动微波辐射遥感反演土壤冻融循环关键变量将受到地表客观条件的制约，因此，基于陆基微波遥感(L波段)冻融的野外观测试验已陆续开展起来，最具代表性的是在芬兰北方森林索丹屈莱地区 [15] 和加拿大副极地草原萨斯喀彻温地区的土壤冻融观测试验 [16]。陆基遥感观测的优势在于：通过详细的辅助野外观测试验获取典型下垫面卫星遥感印痕的关键参数的优化值，区分下垫面的微波信号贡献和误差来源分析 [15] [16]，进行比星载微波辐射计观测时间分辨率更高频率的长期观测(时间分辨率为15 min)。这对深入理解青藏高原季节性冻土和多年冻土活动层冻融过程的微波辐射信号的响应状态至关重要。

综上所述，本研究将依托中国科学院寒旱区陆面过程与气候变化重点实验室隶属的黄河源陆面过程与遥感观测试验平台，基于L波段陆基被动微波辐计亮温观测数据，着重探讨现有的L波段辐射亮温的土壤冻融指数判别方法在青藏高原地区的适用性。通过对多种冻融指数的分析，阐述不同冻融指数对地表冻融的响应特征，并获取最优参数适应L波段星载微波辐射计的地表冻融判别方法。

2. 研究区域及观测数据源

2.1. 研究区域概况

本研究开展观测试验的区域位于青藏高原东北部的玛曲地区(33˚54′N；102˚09′E)，该地区作为黄河上游最重要的产流区补给了黄河上游水量的45.0%，其水资源和生态环境问题，将威胁到下游数亿人口的供水安全，而研究该区域的土壤冻融过程对认识黄河源区乃至青藏高原的水文和生态过程具有重要的现实意义。该区域植被以高寒草甸和高寒湿地草甸为主，观测站点处海拔3450 m，多年平均降水量约为604.0 mm，年平均气温约为1.8℃，属典型高原亚寒带半湿润大陆性季风气候。该区域典型表层土壤类型为粉砂壤土和砂质壤土，2.5与5.0 cm深度处的土壤电导率实地测量表明，全年最大电导率仅为0.095和0.107 S/m，表明观测站处表层土壤为典型的低盐含量土壤，因此，0.0℃可以作为土壤冻结的合理阈值 [17]。

2.2. 辐射亮温

本研究使用的L波段被动微波辐射计ELBARA-III是瑞士联邦森林-积雪研究所为欧空局研发的一款狄克型L波段微波辐射计(f = 1.41 GHz, λ = 21 cm)。ELBARA-III采用了最优的双极化(水平极化/垂直极化)锥形天线设计(直径1.4 m，天线长度2.7 m，−3dB波束宽度为12˚)，包含两个同步信道(1.40 GHz~1.418 GHz, 1.409 GHz~1.427 GHz)增强频带内的噪声信号识别，自带冷源和热源标定，能将信号干扰和观测误差降至最低 [19]。该型辐射计地表亮温观测的绝对误差小于1.0 K，可响应约0.1 K的地表亮温变幅，辐射计观测以30 min的时间间隔和5˚角度步长开展40˚至65˚观测角的扫描，其中50˚观测角对应的印痕大致范围如图1(c)所示。在本研究中，选取2017年10月1日~2018年6月1日为研究时段，并提取该时段内的双极化亮温观测。

2.3. 土壤温度

为配合微波辐射计的观测时段，获取同期的地表冻融循环期间的土壤温度，本研究同期开展了土壤温度剖面观测，用以获取不同季节的土壤温度和湿度廓线。土壤温度观测廓线是通过美国decagon公司生产的5TM传感器获取，传感器精度为0.1℃，数据输出时间间隔与ELBARA-III相同为30 min，观测温度廓线垂直分布深度为2.5 cm、5.0 cm、10.0 cm、15.0 cm、20.0 cm、25.0 cm、35.0 cm、45.0 cm，如图1(b)所示。需要指出的是，已有的研究表明，无论是高土壤水分含量或低含量条件下 [20]，甚至是冻结期间 [14]，L波段的有效反演深度不不超过2.5 cm，因此本研究选取2.5 cm深度处的土壤温度观测作为研究期的表层土壤温度的有效值。

Figure 1. Overview of field data and L-band passive microwave remote sensing instrument used in this study [18]

图1. ELBARA-III L波段被动微波遥感观测实验概况 [18]

2.4. 星载微波辐射计亮温观测和土壤温度观测网

为了评估冻融反演算法的适用性，我们将利用SMAP (NASA Soil Moisture Active Passive mission)任务卫星，依据陆基微波辐射计观测标定的冻融算法的参数和阈值，开展星载微波辐射计空间尺度的冻融反演。SMAP任务卫星是太阳同步轨道卫星，轨道高度为685 km，升交点和降交点时间分别为6 PM和6 AM。SMAP卫星搭载了L波段(1.41 GHz)被动微波辐射计，它能够以40˚的入射观测角度，对~40 km像元内的双极化亮温(垂直/水平极化)进行观测。本研究直接应用Level-3亮温产品(由国家冰雪数据中心National Snow and Ice Data Center，NSIDC发布，发布网址为https://nsidc.org/data/smap/smap-data.html)。该产品基于EASE-2格点进行投影，空间分辨率为36 × 36 km，研究选取的SMAP的观测时段为2016年9月1日至2017年6月31日。

为客观的评估SMAP亮温观测反演土壤冻融的效果，星载微波辐射计同期观测的土壤温度也将被提取作为反演土壤冻融的基准值。对应星载微波辐射计观测空间尺度的土壤温度数据来源于玛曲区域尺度土壤温湿度监测网(SMST monitoring networks)。该观测网包含20个土壤温湿度廓线观测，观测深度为5，10，20，40，以及80 cm，其中有11个站点位于选定的SMAP像元内。本研究仅选取5.0 cm深度土壤温度数据为反演土壤冻融的基准数据。

3. 研究方法

3.1. 冻融反演方法

本研究所采用的冻融反演方法是基于相对冻结因子算法(relative Frost Factor, FF_rel)，该方法属经验性极化指数阈值判别方法，该方法基于土壤液态水分冻结造成的土壤介电常数的明确响应为前提，利用双极化亮温的变化特征来区分土壤冻融状态 [21] (Rautiainen et al. 2012)。该方法以最优化阈值统计分析得出的阈值∆作为判据，判断 大于∆时，地表为冻结状态，而当FF_rel小于∆时，则为地表冻结，

$F{F}_{rel}\langle \Delta \to freezeandF{F}_{rel}\rangle \Delta \to thaw$ (1)

其中，FF_rel则通过2式计算获得，

$F{F}_{rel}=\frac{NPR\left(t\right)-NP{R}_{fr}}{NP{R}_{th}-NP{R}_{fr}}$ (2)

其中，NPR(t)代表t时刻的冻融指数计算值，NPR_th和NPR_fr代表该指数条件下的冻结和融化基准值，需要指出的是，在本研究中采用的NPR_th和NPR_fr依据五月和一月升轨和降轨的NPR均值获取。

3.2. NPR指数的选取

本研究选取了四种代表性的冻融指数，选取的原则是尽可能包含双极化亮温观测中不相关的信息用以提取冻融特征。选取的指数包含双极化归一化指数NPR_NPR，垂直极化指数NPR_Vpol，标准差异冻融指数NPR_STI和垂直极化差异权重指数NPR_CombV。需要指出的是，总计有9种指数源于先前基于陆基L波段被动微波辐射计的土壤冻融观测研究试验 [14] [15] [21]，但是其它形式的指数和本研究选取的某一种指数间必定存在高的线性相关关系(相关系数大于0.97)。因此，本研究只选取上述的四种指数作为相对冻结因子算法的输入状态变量，指数的计算方法请见3-6式。

$NP{R}_{NPR}=\frac{T_b^V-T_b^H}{T_b^V+T_b^H}$ (3)

$NP{R}_{Vpol}=300-T_b^V$ (4)

$NP{R}_{STI}=600-\left(T_b^V+T_b^H\right)$ (5)

$NP{R}_{CombV}=\left(T_b^V-T_b^H\right)\cdot \left(300-T_b^V\right)$ (6)

3.3. 判别结果的精度统计

为了全面客观评价各指数应用在冻融方案的判别效果，本研究选取的精度(precision)方案分别为冻结判断精度，融化判别精度和总精度。其中，式7中state状态代表融化或冻结状态， $N_{state}^T$ 代表正确的判别样本数量， $N_{state}^F$ 代表错误判别的样本数量。

$precisio{n}^{state}=\frac{N_{state}^T}{N_{state}^T+N_{state}^F}$ (7)

4. 结果分析

4.1. 算法的阈值评估

应用FF_rel方法时，阈值的选取对冻融判别精度的影响非常显著，因此，有必要对算法中冻结阈值开展影响评估。本研究选取50˚为亮温观测角度，评估不同阈值条件下的融化、冻结及总精度，并选取最优阈值应用于各NPR指数地表冻融的判别。如图2所示，NPR、STI和CombV三指数随阈值的增加，冻结精度迅速下降，仅Vpol指数对阈值下降相较其他三种指数并不明显，最低判别精度也高于0.70。四种指数最优阈值范围也有明显差异，其中STI和CombV的最优阈值范围介于0.2~0.4之间，而Vpol对应的最优阈值为0.19。此外，四种指数的融化精度都随阈值的增大迅速饱和(接近1.0)，其中，NPR指数输入在阈值等于0.71时融化精度接近1.0，而此时对应的冻结精度仅为0.64，说明NPR有较强的低估冻结状态的趋向。

Figure 2. Frozen, thaw and total precision of four polarization indexes as changes of threshold at the angle of 50˚

图2. 四种指数在50˚观测角度下，冻结、融化和总判别精度随阈值变化(0.0~1.0)的特征

4.2. 算法精度评估

采用不同指数的判别模型获得了显著的判别精度差异。其中，Vpol获得了四种指数的最高总精度0.876和最高的冻结判定精度0.860，而STI指数获得了最高的融化判定精度，达到0.915。表1表明，四种指数均能够获得较高的融化判定精度，但对冻结判定的精度均低于融化状态。其中，广泛采用的NPR指数对冻结状态的判别精度仅为0.75，这充分表明NPR指数并不适合作为青藏高原地表冻融的判别。

Table 1. System resulting data of standard experiment (Precision)

表1. 标准试验系统结果数据(Precision)

4.3. 观测角度对阈值和判别精度的影响

不同入射角度的陆基微波辐射计观测作为地表微波遥感的重要组成部分，因此，有必要评估不同入射角度条件下的精度情况，选取最佳的冻融判别观测角度。如图3(c)所示，随着观测角度增大，不同指数的判别精度均有显著提升，特别是针对CombV指数，精度从40˚时的0.79迅速提升60˚时的0.87，另三种指数也有不同程度的明显提升。但除STI指数在65˚时还有还存在微弱的增加外，其他三种指数均在60˚时达到最高，尽管针对L波段两文观测存在明显鲁斯特角度效应，即大角度对表层土壤的冻融更加敏感，但总精度的统计结果显示，应用FF_rel方法时55˚或60˚入射角时精度已达到饱和，大于60˚的亮温观测并不能提高冻融精度的反演。

针对冻结精度的分析表明，STI和CombV指数也同样的60˚时达到最大，而NPR和Vpol指数则在55˚时已达到最高精度，且随着角度增大NPR指数的判别精度出现了明显的下降趋势，而Vpol指数并未出现明确的变化，Vpol指数的精度在大于55˚时对角度变化不敏感。而对融化状态的精度统计则表明，随着角度的增加，融化精度经历了上升–下降–再上升的过程，即45˚精度高于50˚，而55˚又大于45˚，这说明FF_rel在判定融化时存在角度陷阱效应，因此，不适于应用50˚时作为精度判断的观测角度。

Figure 3. Frozen, thaw and total precision of four polarization indexes corresponding to 40˚ - 65˚ incident angle of l-band brightness temperature measurements

图3. L波段亮温观测入射角度变化从40˚至65˚ (以5˚为间隔)对应的冻结、融化和总判别精度

4.4. 基于SMAP观测的算法验证

为了验证算法的有效性，研究利用SMAP卫星的双极化亮温观测产品，依据4.1~4.3结获取的最优冻融指数Vpol以及最优化阈值，开展了玛曲地区季节性土壤冻融过程的反演。反演结果如图4所示，基于Vpol指数的相对冻结因子算法，能够非常清晰的反映青藏高原季节性冻融过程，特别是针对降轨时刻的反演，非常明确的捕捉了5 cm深度土壤冻融的冻结，融化以及短期的反复冻结过程，并且获得了较高的精度。其中，降轨时刻的反演的精度为0.90，升轨时段的精度为0.92。但需要指出的是，尽管统计结果显示升轨时段的精度略高于降轨时刻，但考虑到土壤垂直方向冻融状态不均匀性要明显弱于升轨时刻。因此，针对星载微波辐射计观测尺度，相对冻结因子算法在降轨时刻的反演结果更加可信。

Figure 4. F/T retrieval results based on FF_rel algorithm using SMAP measurements of descending and ascending over pass times (blue/Yellow Bardenotes freeze/thaw), and concurrent soil temperature at 5 cm depth and air temperature

图4. 基于SMAP升轨和降轨观测的冻融反演结果(蓝色柱状标记冻结状态，黄色则表明融化状态)，以及同期的土壤观测网5 cm深处土壤温度平均值及气温

5. 结论与讨论

本研究针对青藏高原表层土壤冻融过程的复杂性，以青藏高原典型季节性冻土下垫面的冻融特征为研究切入点，以陆基L波段被动微波遥感观测和土壤温度的野外实测为研究基础，基于相对冻结因子算法，以四中代表性的极化指数作为输入变量，开展不同观测角度下的表层土壤冻融精度统计和分析并针对星载微波辐射计反演结果进行有效性验证。其中，基于对相对冻结因子算法的精度分析得出主要结论如下，(1) 明确了在青藏高原季节性冻土区，最适采用Vpol和STI作为相对冻结因子算法的冻融判别指数，该两种算法最大程度地减小了算法本身对冻融状态判断的混淆程度，而且受冻结/融化样本数量比例的影响较小；(2) 论证了相对冻结因子算法对入射观测角度的敏感性，得出55˚和60˚的入射角是L波段观测冻融的最优观测角度。

此外，本研究基于陆基微波辐射计获取的F/T参考值和最优化阈值，利用星载微波辐射计SMAP的双极化亮温观测，以区域尺度的土壤温度观测网数据为验证来源，展示了SMAP升轨/降轨亮温观测反演结果对季节性土壤冻融转换的效果，并统计了相对冻结因子算法的精度。结果表明，(1) 星载微波辐射计同样可以在陆基微波辐射计观测试验的基础上，获取明确的节性土壤冻融转换信号。(2) Vpol指数在星载微波辐射计像元尺度上，也能够获取较高的冻融反演精度。但是，需要指出的是，本研究同时具有一定的局限性，是因为青藏高原表层土壤冻融过程的复杂性远超现有认知，特别是冬季强冻结期间的表层冻融日循环和反复融化冻结过程，以及土壤冻结的垂直不均匀性对相对冻结因子算法的影响，都并未在本文中进行深入的探讨，因此，这将是我们未来工作重点关注的方向。