LAI和FVC植被参数对VIC模型土壤含水量模拟的影响研究

王 慧，孙亚勇 ，黄诗峰 ，马建威，杨永民

（1.中国水利水电科学研究院，北京 100038；2.水利部防洪抗旱减灾工程技术研究中心，北京 100038）

1 研究背景

土壤含水量是陆-气水循环及地表能量平衡中的一项重要参数，在水文过程、气候变化、生态平衡、农业生产等不同领域中都有着重要作用。分布式水文模型是获取大范围、高时空分辨率土壤含水量数据的重要手段。VIC水文模型（Variable Infiltration Capacity Macroscale Hydrologic Model，VIC）是一种基于水量平衡与陆-气间能量平衡模拟的分布式水文模型，在径流、蒸散、土壤水模拟等方面获得了广泛的应用。Nijssen等［1］利用VIC模型生成了1980—1993年全球范围2°分辨率的逐日土壤含水量分布图；Andreadis等［2］利用基于VIC模型模拟的0.5°分辨率的土壤水分和径流资料，重现了历史上北美大陆1920—2003年的干旱过程；Wu等［3］利用VIC模型模拟了1971—2005全国范围30 km分辨率的逐日土壤含水量；Zhang等［4］利用VIC模型建立了1952—2012年全中国0.25°分辨率的地表水文通量与状态数据集。

植被作为下垫面的重要要素，是土壤水分影响地表能量的重要因素，而植被覆盖度及植被健康状况对地表蒸散、近地表大气水分分布等都有一定影响［5-6］。植被参数是分布式水文模型中的重要参量，其中叶面积指数（Leaf Area Index，LAI）与植被覆盖指数（Fractional Vegetation Cover，FVC）两种主要的植被参数很大程度影响着植被对土壤水分的作用能力。VIC模型中，默认情况下，每种植被每个月的LAI为固定的常数，不存在年际变化；每个网格的FVC为固定的常数，不存在年内与年际变化。然而，实际的LAI及FVC具有时空差异性，会随区域和时间的改变而变化。VIC模型在默认情况下的静态植被参数无法准确反映出现实植被的时空差异性，影响了模型的模拟精度。目前，已有一些学者尝试使用精度更高的植被参数作为VIC模型的输入参数，并进行了相关应用分析。Mao等［7］利用陆面数据同化系统（Land Data Assimilation System，LDAS）的LAI产品制作了VIC模型的植被参数数据；Yuan等［8］使用甚高分辨率辐射计（Advanced Very High Resolution Radiometer，AVHRR）的土地分类数据及LDAS提供的植被数据制作了VIC模型的静态植被参数数据，并以其作为VIC模型的输入参量，能较好地模拟中国汉江流域的径流量。遥感对植被有良好的动态监测能力，能够有效提供动态LAI及FVC数据。MODIS（Moderate-resolution Imaging Spectroradiometer，中分辨率成像光谱仪）数据能够提供全球的动态植被产品，对VIC模型的应用能提供重要的数据支撑。Tang等［9］发现基于MODIS的动态LAI数据能使得VIC模型在墨西哥西部北美季风区的土壤含水量及蒸散模拟精度提高；Zhou等［10］开展了静态植被参数与基于MODIS的动态植被参数在模型应用中的对比研究，发现使用动态植被参数在中国褒河流域能获得更合理的蒸散及径流量模拟结果。国内对MODIS植被数据产品在VIC模型中模拟土壤含水量的研究还很少。

基于此，本文利用MODIS LAI及FVC，开展LAI及FVC对VIC模拟土壤含水量的影响研究，研究两者对VIC模型模拟表层土壤含水量的影响，为提高VIC模型土壤含水量的模拟精度提供参考。

2 研究区及数据

2.1 研究区 本文研究区位于西藏自治区那曲县，经度为91.68～92.44°E，纬度为31.06～31.82°N。该地区平均海拔4650 m，属高原丘陵地形，坡度较为平缓。研究区地表覆盖类型绝大部分为草地，属高原草甸，植被主要是冷生、中生、多年生草本植物。该地区属亚寒带半湿润高原季风气候，高寒缺氧，空气干燥，太阳辐射强度较大，年平均气温为-2℃，全年大风日100 d左右，年降水量约500 mm，分冷暖季。六月至九月相对温暖，受南亚季风影响，集中了75%的年降水量；十月至次年五月为冰冻期。由于低温冻土等原因，青藏高原的高原草甸生态系统极其脆弱，对气候变化十分敏感，抗干扰能力很低［11］。

2.2 数据及预处理 为开展VIC模型的土壤含水量模拟研究，本文收集了MODIS植被产品数据、气象产品数据、野外站点土壤含水量实测数据等，并对数据进行了预处理，以此进行VIC模型土壤含水量模拟结果的相关分析与精度评价。

2.2.1 MODIS植被产品数据 MODIS即中分辨率成像光谱仪，搭载于Terra、Aqua两颗卫星上，是美国地球观测系统（EOS）计划中用于观测全球生物和物理过程的重要传感器，具有光谱范围广、更新频率高的优点。它每1～2d实现一次全球覆盖，并提供陆海温度、陆地表面覆盖、云、水汽、火情等多种应用产品数据的免费服务。本文选择时间范围为2008年1月1日至2011年8月31日的MODIS LAI 8日产品MOD15A2以及由MODIS NDVI 16日产品MOD13Q1计算的FVC，公式如下：

式中：NDVIs、NDVIv分别为裸土与植被的NDVI值，根据郭鹏［12］的研究，取NDVIs=0.05，NDVIv=0.8。

为研究动态植被参数对VIC模型表层土壤含水量模拟的影响，将试验点上基于MODIS的8日LAI数据与16日FVC数据使用下式分段线性插值为单日数据：

式中：k=1，2，…n-1为植被数据根据时间顺序排列的序号；tk为第k条数据对应的时间点；Vk为tk时的植被数据值；t为插值时间；V()t为t时的植被数据值。2.2.2 土壤含水量观测数据 本文使用中国科学院青藏高原研究所的青藏高原多圈层模拟与数据同化中心（Data Assimilation and Modeling Center for Tibetan Multi-spheres，DAM）在西藏自治区那曲县周边设立的青藏高原中部地区土壤湿度与温度监测网（Central Tibetan Plateau Soil Moisture and Tempera⁃ture Monitoring Network，CTP-SMTMN）所采集的土壤含水量数据［13-16］。数据每30 min输出一次，埋设的土壤水分传感器能给出深度5、10、20和40 cm处的土壤湿度数据。

考虑到观测数据时长与数据的连续性，本文选择CTP-SMTMN的L08（31.662°E，91.795°N）、L21（31.172°E，92.197°N）、L27（31.664°E，92.342°N）、L34（31.259°E，91.799°N）站点2010年9月1日至2011年8月31日之间的数据进行研究（图1）。

图1 研究区及土壤含水量观测站点位置分布

站点观测所得5、10、20和40 cm处的土壤含水量数据，与VIC模型模拟出的0～15 cm土壤含水量数据不能进行直接对比。为方便对比，对站点观测的土壤含水量根据以下公式进行计算，折算为0～15 cm土壤层的平均体积含水量：

式中：Z0～5、Z5～10、Z10～15分别为0 ～ 5、5 ～ 10、10 ～ 15 cm土壤层厚度，均为5cm；θ0～15、θ0～5、θ5～10、θ10～15分别为0～15、0～5、5～10、10～15cm土壤体积含水量，其中：

化简后得：

式中：θ5、θ10、θ20分别为站点观测的5、10、20 cm处土壤体积含水量。

2.2.3 其他数据 气象驱动数据选择了2008年1月1日至2011年8月31日的中国气象局陆面数据同化系统第一版本［17］（CMA Land Data Assimilation System Version1.0，CLDAS-V1.0）日产品数据，包括降水量、气温、风速、气压、短波辐射、湿度等要素。

土壤参数数据部分由世界土壤数据库（Harmonized World Soil Database，HWSD）提供，部分由SWAT模型数据库辅助软件SPAW计算得到，包括土壤饱和体积含水量、土壤水下渗速率、水力传导率、土壤饱和水势、土壤厚度、土壤水扩散系数、饱和容量曲线形状参数等。

3 模型及方法

3.1 VIC模型VIC模型是基于水量平衡与陆-气间能量平衡的分布式水文模型，由华盛顿大学、普林斯顿大学及加利福尼亚大学伯克利分校共同研制开发［18-19］。它弥补了传统水文模型对能量过程描述的不足，同时考虑冠层蒸发、植被蒸腾和裸土蒸发，并可以进行土壤冻融、积雪融雪、湖泊、沼泽相关参量的计算。此外，VIC模型也考虑到了单一网格内存在的地表覆盖类型、高程、土壤蓄水容量等属性的不同［18］。以下为VIC模型模拟土壤含水量的主要物理过程原理。

VIC模型采用下式计算上两层土壤的含水量［20］：

式中：θ为土壤体积含水量，zi为第i层的土壤深度；I为下渗速率；E为蒸散量；为水力传导率；D()θ为土壤水扩散系数。

I、需通过土壤实验进行测算，本文这些参量从由制备好的土壤参数数据直接给出。

蒸散量在很大程度上受到植被数据的影响。VIC模型中的蒸散量由三部分组成：冠层蒸发量Ec、植被蒸腾量Et及裸土蒸发量即：

式中：LAI会对冠层蒸发Ec、植被蒸腾Et有影响。冠层最大截留水量为LAI的0.2倍［21］，在降水充沛时，冠层最大截留水量与Ec成正比；冠层阻力与LAI成反比，与Et也成反比；则LAI与Ec、Et成正比关系［20，22，23］。

每个网格上的总蒸散量E，是Ec、Et、El根据网格内不同地表覆盖所占面积分数计算的加权和：

式中：Cn为第n种植被类型在该网络的植被覆盖比例；CN+l为裸土覆盖比例，且

在VIC模拟表层土壤含水量的过程中，植被参数中的LAI与FVC在降水截留、冠层蒸发及植被蒸腾过程中都起到了一定的作用。其中，LAI会影响到植被冠层最大截留水量与冠层阻力，进而影响到冠层蒸发与植被蒸腾，LAI与二者成正比；由式（7），植被覆盖比例越大，相应植被的冠层蒸发与植被蒸腾越大；由式（8），网格内所有植被覆盖比例总和的大小（即该网格FVC）也会影响到裸土面积，进而影响到该网格裸土蒸发的强度。

3.2 植被数据组合方案为了分析LAI与FVC对VIC模拟表层土壤含水量的影响差异，制定了5种植被数据组合方案（表1），分别就每种方案进行土壤含水量的模拟。

其中，方案5由VIC模型默认LAI计算FVC，与默认LAI组合。通过下式利用默认LAI计算FVC：

3.3 统计指标为了全面评估VIC使用不同植被数据组合方案所模拟的表层土壤含水量的准确性，使用偏差（Bias）、均方根误差（RMSE），相关系数（R）这三种统计指标来进行评估［24］。偏差描述了VIC模型模拟的表层土壤含水量与站点观测值之间的绝对误差，由下式计算：

表1 植被数据组合方案

式中：为线性平均；t为表层土壤含水量的模拟/观测时刻；θsim()t为VIC模型模拟的t时刻的表层土壤体积含水量；θob()t为站点观测的t时刻的表层土壤体积含水量。

均方根误差描述了VIC模型模拟的表层土壤含水量相对于站点观测值的精度，RMSE由下式计算：

相关系数描述了VIC模型模拟的表层土壤含水量与站点观测值在时间上的一致性及相对准确性，由下式计算：

4 结果与分析

使用VIC模型，分别采用表1所示5种植被数据组合方案，对2008年1月1日至2011年8月31日4个试验点的土壤含水量进行模拟，采用2010年9月1日至2011年8月31日的测站数据对模拟结果进行精度验证。

4.1 动态LAI对模拟结果的影响由图2可知，相比方案4，使用方案2在4个试验点上的表层土壤含水量模拟值更高。其中，2010年10月至2011年5月的冰冻期间，使用方案2模拟的表层土壤含水量在L21、L27、L34站点相比方案4有较明显的高估；2010年9月与2011年6—9月的非冰冻期间，方案2模拟的表层土壤含水量在4个站点都比方案4略有高估。由表2可知，相比方案4，使用方案2在4个试验点的土壤含水量模拟值的平均Bias从0.101 cm3·cm-3上升至0.145 cm3·cm-3，RMSE从0.135 cm3·cm-3上升至0.166 cm3·cm-3，R从0.483上升至0.792，表明误差增大，模拟精度降低，但时间上的一致性及相对准确性提高，即变化趋势更加符合站点观测值的表现。

由图2可知，相比模型默认LAI值（2.0～3.8），MODIS LAI值（0.1～1.7）整体较低。根据VIC模型对冠层最大截留水量的定义，LAI降低，使得冠层最大截留水量减少，冠层蒸发量降低；LAI降低使冠层阻力增大，植被蒸腾也会随之降低；根据式（4）、式（5），冠层蒸发量与植被蒸腾量的减少会导致地表蒸散量的减少，从而使得表层土壤含水量升高，令土壤含水量的模拟结果高估、误差增大、精度下降。比起默认LAI，MODIS LAI的走势更符合试验区的气候特征与植被变化规律，这可能是使用MODIS LAI的方案2的模拟结果的变化趋势比使用方案4更准确的原因。

4.2 动态FVC对模拟结果的影响由图2可知，方案3的模拟效果比方案4有明显改善。尤其在冰冻期间，使用方案3的模拟结果能呈现出冰冻期间表层土壤含水量的低值现象；在10～12月、5～6月，方案3能够较准确地模拟出冻融变化过程。由表2可知，相比方案4，使用方案3在4个试验点的土壤含水量模拟值的平均Bias从0.101 cm3·cm-3下降至0.031 cm3·cm-3，RMSE从0.135 cm3·cm-3下降至0.071 cm3·cm-3，R从0.483上升至0.834，误差大幅降低，模拟精度明显提高，变化趋势也更加符合站点观测值的表现。

由图2可知，基于MODIS NDVI计算的FVC值（0.02～0.91）大部分时间远低于默认FVC值（恒定为1），且冰冻期稳定在0.04以下。在冰冻期，试验区的植被几乎全部枯萎（基于MODIS NDVI计算的FVC值始终低于0.04，MODIS LAI始终低于0.2），蒸散以土壤蒸发为主，根据式（4）～式（7），默认FVC恒定为1会导致土壤蒸发全年为0，使得方案4在冰冻期的表层土壤含水量模拟值持续偏高；基于MODIS NDVI计算的FVC使模拟中冰冻期的土壤蒸发量处于合理范围，因此能正确模拟出冰冻期表层土壤含水量的低值现象及冻融变化过程。

4.3 动态LAI、FVC对模拟结果的综合影响由图2可知，方案1的模拟效果比方案4有明显改善，能模拟出冰冻期间液态土壤含水量的低值现象，并能够较准确地模拟出冻融变化过程。由表2可知，相比方案4，使用方案1在4个试验点的土壤含水量模拟值的平均Bias从0.101 cm3·cm-3下降至0.032 cm3· cm-3，RMSE从0.135 cm3·cm-3下降至0.071 cm3· cm-3，R从0.483上升至0.836，误差大幅降低，模拟精度明显提高，变化趋势更符合站点观测值的情况。方案1的三项统计数据与方案3持平，使用MODIS LAI仅使得方案1的Bias比方案3高0.001 cm3·cm-3，R比方案3高0.002。FVC对VIC模型在试验区对表层土壤含水量的模拟精度的影响远远高于LAI。

图2 使用方案1-5模拟所得表层土壤含水量结果对比

综合以上分析结果，提出了在缺乏准确的植被数据源时的替代方案，即利用默认LAI计算FVC参数进行VIC模型模拟（即方案5）。由图2可知，方案5在冰冻期的模拟效果比起方案4有所改善，在冰冻期土壤含水量有下降的趋势。由表2可知，相比方案4，使用方案5在4个试验点的土壤含水量模拟值的平均Bias从0.101 cm3·cm-3下降至0.067 cm3·cm-3，RMSE从0.135 cm3·cm-3下降至0.100 cm3·cm-3，R从0.483上升至0.720，误差降低，模拟精度提高，变化趋势比方案4符合站点观测值的表现，但没有模拟出冰冻期液态土壤含水量的低值现象与冻融变化过程。

由图2可知，方案5使用基于模型默认LAI计算所得的FVC值（0.63～0.85）低于1，土壤蒸发量的模拟得到改善，进而冰冻期的表层土壤含水量模拟值有所降低，比仅使用默认数据的方案4所得模拟结果更加符合实际。不过由于FVC值过高，土壤蒸发量的模拟精度仍然较低，且基于默认LAI计算所得的FVC的变化趋势不符合试验区气候植被特征，故该方案尚有发展的潜力，如使用更加准确的LAI或NDVI计算FVC作为模型输入，进行表层土壤含水量的模拟。

表2 VIC使用不同植被数据组合方案在各试验点模拟表层土壤含水量结果验证统计

5 结论

本文以MODIS遥感数据估算的动态植被参数数据驱动VIC模型，模拟西藏那曲县2008年1月1日至2011年8月31日的表层土壤含水量（0～15cm），以此分析动态植被参数LAI及FVC对VIC模拟表层土壤含水量精度的影响。结果表明：（1）使用MODIS动态LAI数据驱动VIC模型，R从0.483上升至0.792，其表层土壤含水量的模拟值与站点观测值在时间上的变化趋势一致性更好，模拟值的相对准确性提高；（2）使用MODIS动态FVC数据驱动VIC模型，平均Bias从0.101 cm3·cm-3下降至0.031 cm3·cm-3，RMSE从0.135 cm3·cm-3下降至0.071 cm3·cm-3，其表层土壤含水量的模拟精度显著提升，并且R从0.483上升至0.834，模拟值与站点观测值在时间上的变化趋势一致性显著提高，较准确地模拟了表层土壤含水量的相对变化，尤其能更准确地模拟出冻融过程中表层土壤含水量的变化趋势；（3）相比静态植被参数，使用MODIS遥感数据估算的动态植被参数显著提升了VIC模拟表层土壤含水量的精度，较准确的模拟出冰冻期液态土壤含水量的低值，较好的模拟出10—12月及5—6月的土壤冻融变化过程。

致谢：感谢由中国气象局提供的陆面数据同化系统第一版本（CMA Land Data Assimilation System Version1.0，CLDASV1.0）生产的气象数据。

参考文献：

［1］NIJSSEN B，SCHNUR R，LETTENMAIER D P.Global retrospective estimation of soil moisture using the vari⁃able infiltration capacity land surface model，1980-93［J］.Journal of Climate，1999，14（8）：1790-1808.

［2］ANDREADIS K M，CLARK E A，WOOD A W，et al.20th century drought in the conterminous United States［J］.Journal of Hydrometeorology，2005，6（6）：985-1001.

［3］WU Z，LU G，WEN L，et al.Thirty-five year（1971-2005）simulation of daily soil moisture using the variable infiltration capacity model over China［J］.Atmosphere-Ocean，2007，45（1）：37-45.

［4］ZHANG X J，TANG Q，PAN M，et al.A long-term land surface hydrologic fluxes and states dataset for China［J］.Journal of Hydrometeorology，2014，15（5）：2067-2084.

［5］MCPHERSON R A.A review of vegetation--atmosphere interactions and their influences on mesoscale phenome⁃na［J］.Progress in Physical Geography，2007，31（3）：261-285.

［6］ SR R A P.Influence of the spatial distribution of vegetation and soils on the prediction of cumulus convective rain⁃fall［J］.Reviews of Geophysics，2001，39（2）：151-177.

［7］MAO D Z，CHERKAUER K A.Impacts of land-use change on hydrologic responses in the Great Lakes region［J］.Journal of Hydrology，2009，374（1）：71-82.

［8］YUAN F，XIE Z，LIU Q，et al.An application of the VIC-3L land surface model and remote sensing data in sim⁃ulating streamflow for the Hanjiang River basin［J］.Canadian Journal of Remote Sensing，2004，30（5）：680-690.

［9］TANG Q，VIVONI E R，MU OZARRIOLA F，et al.Predictability of evapotranspiration patterns using remotely sensed vegetation dynamics during the North American Monsoon［J］.Journal of Hydrometeorology，2012，13（1）：103-121.

［10］ZHOU S，LIANG X，CHEN J，et al.An assessment of the VIC-3L hydrological model for the Yangtze River ba⁃sin based on remote sensing：a case study of the Baohe River basin［J］.Canadian Journal of Remote Sensing，2004，30（5）：840-853.

［11］罗亚勇，孟庆涛，张静辉，等.青藏高原东缘高寒草甸退化过程中植物群落物种多样性、生产力与土壤特性的关系［J］.冰川冻土，2014，36（5）：1298-1305.

［12］郭鹏.基于MODIS的青藏高原地表参数时空变化研究［D］.北京：中国石油大学，2010.

［13］CHEN Y，YANG K，QIN J，et al.Evaluation of AMSR-E retrievals and GLDAS simulations against observations of a soil moisture network on the central Tibetan Plateau［J］.Journal of Geophysical Research-Atmospheres，2013，118（10）：4466-4475.

［14］QIN J，YANG K，LU N，et al.Spatial upscaling of in-situ soil moisture measurements based on MODIS-derived apparent thermal inertia［J］.Remote Sensing of Environment，2013，138（138）：1-9.

［15］YANG K.A multi-scale soil moisture and Freeze-Thaw Monitoring Network on the Tibetan Plateau and its appli⁃cations［J］.Bulletin of the American Meteorological Society，2013，94（12）：1907-1916.

［16］ZHAO L，YANG K，QIN J，et al.Spatiotemporal analysis of soil moisture observations within a Tibetan meso⁃scale area and its implication to regional soil moisture measurements［J］.Journal of Hydrology，2013，482（5）：92-104.

［17］SHI C，JIANG L，ZHANG T，et al.Status and Plans of CMA Land Data Assimilation System（CLDAS）Project［C］//Proceedings of the EGU General Assembly Conference，2014.

［18］LIANG X，LETTENMAIER D P，WOOD E F，et al.A simple hydrologically based model of land surface water and energy fluxes for general circulation models［J］.Journal of Geophysical Research Atmospheres，1994，99（D7）：14415-14428.

［19］NIJSSEN B，LETTENMAIER D P，LIANG X，et al.Streamflow simulation for continental-scale river basins［J］.Water Resources Research，1997，33（4）：711-724.

［20］GAO H，TANG Q，SHI X，et al.Water Budget Record from Variable Infiltration Capacity（VIC）Model［Z］.2010，

［21］DICKINSON R E.Modeling evapotranspiration for three-dimensional global climate models［J］.Washington Dc American Geophysical Union Geophysical Monograph，1984，29：58-72.

［22］BLONDIN C.Parameterization of Land-Surface Processes in Numerical Weather Prediction［M］.Springer New York，1991.

［23］DUCOUDR N I，LAVAL K，PERRIER A.SECHIBA，a New Set of Parameterizations of the Hydrologic Ex⁃changes at the Land-Atmosphere Interface within the LMD Atmospheric General Circulation Model［J］.Journal of Climate，1993，6（2）：248-273.

［24］SUN Y，HUANG S，MA J，et al.Preliminary evaluation of the SMAP radiometer soil moisture product over China using in situ data［J］.Remote Sensing，2017，9（3）：292.