基于InVEST 模型的珲春天保工程区水源涵养功能时空动态变化研究

郭玉亭，王英伟，杨成江

（1.东北林业大学林学院，黑龙江 哈尔滨 150040；2.吉林省重点国有林技术服务中心，吉林 长春 130022）

【研究意义】森林生态系统是维持社会经济发展和健康生态系统的基础［1］。21 世纪以来，由于经济、社会的迅速发展和人口的急剧增长，人类活动剧烈地改变了地表自然环境，产生了资源约束趋紧、土地退化、环境污染、生态系统破坏等一系列问题。目前，全球约有60%的生态系统已经或正在退化，导致生态系统服务供应量锐减，人类福祉遭受严重损害［2-3］。为减少这些问题的影响，国务院2000 年开始实施世界上最大和最全面的生态保护工程——天然林保护工程（简称“天保工程”）。天保工程森林面积占我国天然林面积的53%，主要目标是在一定程度上提升森林生态系统调节气候、固碳释氧、涵养水源、保水固土、保护生物多样性等服务功能，保护环境与促进社会和谐发展［4］。【前人研究进展】自天保工程实施后，众多学者对国内天保工程的生态服务功能及其对经济社会的影响等方面进行评析。例如，王慧等［5］通过森林生态连清体系对吉林省森林局天保工程进行效益评价，结果发现生物多样性保护和涵养水源两项服务的价值量比例超过50%，为生态系统主导功能；黄龙生等［6］通过应用分布式测算方法，对东北和内蒙古重点国有林区天保工程保育土壤物质量及价值量进行评估，结果表明天保工程的实施提升了研究区森林生态系统保育土壤的功能及效益；李娜娜等［7］采用森林生态系统定位数据及森林二类调查数据，选取涵养水源、保育土壤、固碳释氧3 项指标对山西省天保工程生态效益进行评价，发现山西省天然林保护工程自实施以来，有效地保护和恢复了天然林资源，发挥了天然林的巨大生态效益。从研究方法来看，生态系统服务功能的评估方法有很多，如InVEST 模型、ARIES 模型、SolVES模型、TESSA 模型等。其中，InVEST 模型具有驱动数据简便容易掌握、定量化评价精确、评估流程与结果空间表述明确、情景仿真预测功能强大等优点，已被国内外学者广泛应用于生态系统研究中［8］。例如，王茜等［9］利用InVEST 模型研究秦岭产水量变化对土地利用变化的响应，结果表明秦岭的产水量由秦岭南坡至北坡依次递减，不同的土地利用类型产水功能差异较大，与各地类植被覆盖率、潜在蒸散发及面积等因素有关；Bai 等［10］选取水质净化、水源涵养、土壤保持3 个指标通过InVEST 模型对美国肯塔基州进行评估，发现气候与土地利用变化为水质净化、水源涵养、土壤保持的主要影响因素，其中气候变化对水源涵养的影响大于土地利用变化对水源涵养的影响。【本研究切入点】珲春天保工程区是国家“一带一路”重要实施区，区内生态系统类型多样、自然资源丰富，被称为“生态宝库”［11］。然而，关于该区域的水源涵养服务功能在长时间尺度的动态变化还鲜有针对性研究。因此，基于InVEST 模型科学评价天保工程实施以来珲春天保工程区水源涵养服务功能的时空变化趋势，并分析其影响因素非常有必要。【拟解决的关键问题】为探究近20 年天保工程给该区带来的生态效益，本研究基于InVEST 模型，在多年林业资源档案数据、卫星遥感数据的基础上，从涵养水源、水质净化两个方面对近20 年天保工程水源涵养服务时空演变情况进行评价与分析，并探讨水源涵养服务时空动态变化的影响因素，以期为水资源利用、水环境污染治理的研究奠定基础，并为今后天然林保护修复制度的实施提供决策依据。

1 数据来源与研究方法

1.1 研究区概况

珲春天保工程区位于吉林省珲春市（包含少部分图们市），地处中、俄、朝三国交界的边境地区（地理位置129°52'～131°18' E、42°25'～43°30' N）。在地形地貌方面，研究区属于长白山系东部老爷岭山脉中低山区，地势起伏大，高差悬殊，整体呈东南、东北、西部较高，中部、西南部较低的簸箕形盆地。在气候方面，研究区属近海中纬度温带海洋性季风气候，年均气温5.65 ℃，年均降水量为618.1 mm，年均蒸发量894.2 mm。在自然资源方面，研究区内水资源丰富，年均水资源总量20.58×108m3，地下水总量3.86×108m3［12］。

1.2 数据来源

天保工程2000 年开始实施、2020 年结束，且在2015 年采取全面停止天然林商业性采伐政策，故选取2000、2010、2015、2020 年为时间节点进行研究。本研究数据包括土地利用数据、气象数据、土壤属性数据和数字高程模型（DEM）数据和潜在蒸散发数据（表1）。其中，将土地利用数据重分类为耕地、林地、草地、水体、建设用地及湿地6 类；2000—2020 年降水量数据源自中国气象数据网，并在ArcGIS 中通过反距离权重法插值得到；各栅格数据的精度统一为30 m，投影坐标系统一采用WGS_1984_UTM_Zone_51N。

表1 研究数据来源情况Table 1 Research data sources

1.3 研究方法

1.3.1 转移矩阵及动态度分析方法 土地利用转移矩阵法是MARKOV 模型运用在土地利用变化上的方法［13］，是分析研究系统的状态与转移量的定量过程，包括转出面积与转入面积，计算公式为：

式中，T1、T2为研究时段内土地利用转移矩阵的两种表示方式，T2称为MARKOV 转移概率矩阵；Sij是研究时段内由地类i转移为地类j的面积；Pij是土地利用类型由研究期初地类i转移为研究期末地类j的概率；n是土地利用类型的数量。

动态度分析法用于测算研究区整体或单一土地利用变化速度，其中综合土地利用动态度用P表示、单一土地利用变化速度用K表示，计算公式［14］为：

式中，ΔUi-j为某一研究期内地类i向其他地类转出的面积之和，Ui为研究期初地类i利用类型的面积，t为研究期间隔年数。

式中，Ua、Ub分别为期初和期末某土地利用类型的面积，t为研究期间隔年数。

1.3.2 水源涵养分析方法 InVEST 模型的水源涵养模块，是在Budyko 水热耦合平衡与年均降水量的基础上，同时考虑实际蒸散发。产水量计算公式［15］如下：

式中，Y（x）为单元x年产水量（mm），AET（x）为单元格x实际蒸散量（mm），P（x）为单元x年降水量（mm）。

基于InVEST 模型中的Water Yield 模块结果，运用ArcGIS 计算得到珲春天保工程区的水源涵养量，计算公式［10］如下：

式中，WRij为LUCCj上像元i的年水源涵养量（mm），Y（x）为单元格x年产水量（mm），Runoffij为LUCCj上像元i的年地表径流量（mm），Pij为LUCCj上像元i的年均降水量（mm），Cj为地类j的地表径流系数。

1.3.3 水质净化分析方法 水质净化服务通过InVEST Nutrient Delivery Ratio 模块进行计算，模块基于土地利用/覆被类型在年尺度上计算栅格单元N、P 营养盐的保持量和输出量，利用N、P 含量表征水质状况，TN、TP 输出量越大表明栅格单元水质净化能力越弱，计算公式［16］如下：

式中，ALVx为栅格x各污染物年输出值（kg）［14］，HSSx为在栅格x上水文敏感性得分，POLx为栅格x污染物输出系数。

1.3.4 灰色关联度分析方法 灰色系统理论提出对各系统进行灰色关联度分析的概念，灰色关联度是以因素之间发展趋势的相似或相异程度作为衡量各系统因素关联程度的一种方法［17］。灰色关联度分析法多用于模糊的影响因素与生态服务功能之间，由于影响因素与生态服务功能之间关系预先并不清楚，故本研究运用灰色关联度法，选取自然、经济、社会以及政策方面的9 个因子，对水源涵养服务功能的影响因素进行量化分析。灰色关联度的计算公式为：

式中，γ为灰色关联系数，i=1，2，3，...m；k=1，2，3，...n。

2 结果与分析

2.1 土地利用类型变化特征

2000—2020 年林地和耕地是珲春天保工程区的主要土地类型（图1、图2），且林地整体呈波动上升趋势，耕地呈减少趋势。此外，土地利用动态变化度可一定程度上展现土地利用类型变化的剧烈程度［18-19］。草地面积呈现先增加后减少的趋势（3.62%～-6.84%），且在2015—2020 年期间下降速率较快；耕地面积呈现逐渐减少态势（-1.00%～-1.41%），减少速率较平稳；建设用地呈较剧烈的增长（-0.43%～11.30%），增长速率呈倍速上升；林地面积增长较为剧烈（-0.05%～0.47%）；湿地面积在工程实施期间大面积减少（-17.26%～0.20%）；水体面积在20 年间稳步增长，增长速率前快后慢（2.47%～1.96%）。2010—2020 年综合土地利用变化度（0.48%）远大于2000—2010 年（0.19%），前者是后者的2.23 倍。

图1 2000—2020 年土地利用类型Fig.1 Land use status from 2000 to 2020

图2 2000—2020 年不同土地利用类型面积及单一土地利用动态变化度Fig.2 Area of different land use types and dynamic change degree of single land use from 2000 to 2020

2000—2020 年不同时段中耕地向草地或林地转移率最大，特别是2000—2010 年耕地转草地的转移比率达62.6%。2000—2010 年，研究区土地利用类型具有实质性转移的共有25 种（171.37 km2），按面积大小对此类转移进行排序，前7 种累计转移比率达89.15%（表2）；面积占比较大的是耕地转草地（转移面积107.29 km2，转移比率62.60%），其次是林地转耕地（转移面积20.04 km2，转移比率11.69%）。2010—2015 年，研究区土地利用类型具有实质性转移的共有30 种（745.79 km2），比上一阶段增加574.42 km2，按面积大小进行排序后发现前7 种累计转移比率达86.01%；本阶段土地利用变化以耕地和草地间转换为主（转移面积150.54 km2，转移比率20.19%），其次为林地转耕地（转移面积127.82 km2，转移比率17.14%）。2015—2020年，研究区土地利用类型具有实质性转移的共有30 种（78 056.71 km2），按面积大小进行排序后发现前7 种累计转移比率达83.24%；转移面积占比最大的是耕地与林地间的转换（转移面积150.33 km2，转移比率19.26%），其次为草地向耕地的转化（转化面积140.74 km2，转移比率18.03%）。

表2 2000—2020 年不同土地利用类型主要变化面积和转移比率Table 2 Major changes in area and transfer rates of different land use types from 2000 to 2020

2.2 水源涵养服务功能时空变化特征

2.2.1 水源涵养量时空变化特征 研究区水源涵养量呈东北高、西南低的分布格局（图3），2000—2020 年研究区水源涵养量总体呈先减少后增加趋势，各期水源涵养总量分别为32.08×108、6.41×108、25.39×108、30.36×108m3。2000—2010年水源涵养总量减少5.67×108m3，降幅17.67%；2010—2015 年水源涵养总量减少1.02×108m3，减少趋势较一期缓慢；2015—2020 年水源涵养总量大幅增加4.97×108m3，涨幅19.57%。

图3 水源涵养量空间分布Fig.3 Spatial distribution of water conservation

2000—2020 年，各土地利用类型中，以林地的水源涵养贡献率最大（表3）。2000 年，林地水源涵养量最大，占总量82.51%；2010 年，林地水源涵养量占总量的79.86%，较2000 年减少20.32%，其余地类水源涵养量无过大波动；2015年，林地水源涵养量达21.38×108m3，相比2010年增长较少，草地水源涵养量有所减少；2020 年，林地水源涵养量占总量82.18%，涨幅16.70%，其余地类水源涵养量除草地增加外无过大变动。总体来看，水源涵养量大小呈林地＞耕地＞草地＞水体＞湿地。

表3 不同土地利用类型的水源涵养量（×108 m3）Table 3 Water conservation of different land use types(×108 m3)

2.2.2 水质净化时空变化特征 水质净化能力与N、P输出量相反。研究区N、P输出呈西南部高、中部及东北部低的分布格局（图4、图5），2000—2020年，研究区N、P输出量年均值总体呈波动中减少趋势（表4），表明水质净化能力呈增加趋势。2000—2020年，研究区各期N输出量年均值分别为0.1806、0.1834、0.1678、0.1703 kg/hm2，P输出量年均值分别为0.0206、0.0209、0.0175、0.0200 kg/hm2。2000—2010年，N输出量年增加0.0028 kg/hm2、涨幅1.55%，P输出量年增加0.0003 kg/hm2、涨幅1.46%；2010—2015年，N、P输出量分别年均减少0.0156、0.0034 kg/hm2，减幅分别为8.51%、1.62%。2015—2020年，N、P 输出量稍有增多，N 输出量年均增加0.0025 kg/hm2、涨幅1.49%，P输出量年均增多0.0025 kg/hm2、涨幅1.43%。

表4 不同土地利用类型的单位面积N、P 输出量Table 4 Average N and P output of different land use types

图4 N 输出量空间分布Fig.4 Spatial distribution of N output

图5 P 输出量空间分布Fig.5 Spatial distribution of P output

从各地类N、P 输出情况来看，耕地是N、P输出量的最大贡献者。2000、2000、2015、2020年各地类对水质净化能力大小表现为水体＞湿地＞建设用地＞草地＞林地＞耕地。

2.3 水源涵养服务功能影响因素

本研究采用灰色关联分析法探究天保工程区水源涵养与水质净化功能间的影响因素及它们之间的关系，结果发现在各影响因素中，自然因子影响作用最大。自然因子中，对水源涵养、水质净化功能影响最大的是降水量、蒸散发、森林蓄积3 个因素；经济因子中，对水源涵养、水质净化功能影响最大的因素是第一产业值（表5）。2000 年天保工程在吉林省珲春市正式展开实施，珲春辖区范围内的森林质量不断提升，使得森林蓄积在20 年间持续增长。根据当地林业资源档案可知，2000—2020 年研究区森林蓄积量由388.68 m3增长到426.74 m3，增幅达9.79%，珲春林区职工收入与林区管护人员的增加使森林质量与水源涵养服务功能均有提升，实现了天保工程“看住林，养活人”的初衷，同时也提升了人们对于森林资源的保护意识，破坏森林的现象少有发生。

表5 水源涵养服务功能影响因子及关联度Table 5 Influencing factors and correlation degree of water conservation function

3 讨论

2000 年天保工程实施以来，珲春天保工程区的林地面积呈先减少后增加趋势，这是由于其他土地利用类型向林地大面积转移导致的，特别是耕地向林地转移。其中，2010—2020 年综合土地利用变化动态度较高，土地利用变化较大，这主要是天保工程二期快速推进引起的。与此同时，珲春天保工程区的水源涵养服务功能也发生一系列变化。从时间上看，其水源涵养总量在一定程度上有所提升，且与林地变化趋势一致，这主要与林地面积变化及降水量变化有关；且该结果与东北内蒙古等重点国有林区水源涵养功能的研究结果［20］一致。从空间上看，研究区水源涵养主要表现为林地占据高值、草地和耕地等占据低值，这是由于林冠密度高，大气降水接触叶面积截留的水分较其他地类截留水量多，并且天然林林地中有较厚腐殖质落叶层，降水流动速度也比其他地类慢，以上两个因素均可增强截留效率，起到涵养水源的作用［21］。对于水质净化功能，N、P输出量呈逐年减少趋势，研究区内水质净化功能逐渐加强，与西南地区天保工程区内水质净化净化规律变化［22］一致。上述结果主要是由于天保工程实施后，森林面积显著增加，森林植被通过生化反应、水循环作用，在大气、土壤和降水中吸收N、P 等营养物质并贮存在内，减少了营养物质向低谷地区输移，因而森林植被量的变化直接影响水质净化的能力。此外，耕地输出N、P最多，这是由于化肥和农药大量使用的结果；建设用地N、P 输出量较耕地、林地小，这可能是由于研究区养殖畜牧业发展没有种植业发达所导致；而水作为人类生存发展必不可少的资源，则不具备N、P 营养盐输出的能力。

本研究分析了珲春天保工程区在2000—2020年间水源涵养、水质净化的时空变化特征，但各生态服务功能综合分析及各生态服务功能之间的权衡与协同关系在区域可持续发展、环境保护与管理等政策制定与实施方面也至关重要，因此需要进一步详细研究。此外，生物多样性维护、固碳释氧、水土保持等也是工程区重要的生态系统服务，但由于数据资料的限制，本研究并未涉及这些方面的生态效应，研究结果不够全面。另外，我国近年来气候呈现暖湿变化趋势，有利于植被恢复与功能提升，但如何量化评估生态工程与气候变化的相对贡献亦是未来的研究重点。

4 结论

天保工程实施后的20 年间，吉林省珲春市研究区林地总面积整体呈波动上升趋势。其中，2000—2020 年耕地向草地或林地转移率最大；水源涵养量呈现出先减小后增加的趋势，高值区集中分布在研究区外围的林地，低值区集中分布在以建设用地和耕地为主的沟谷区；水质净化功能呈现波动增长的趋势，耕地为N、P 输出高值区，林地和草地为N、P 输出低值区。天保工程通过森林蓄积方式在水源涵养、水质净化等森林生态系统服务功能的提升中作出了重要贡献。