1、第39 卷 第13期2602023年农业工程学报7月Transactions of the Chinese Society of Agricultural EngineeringVol.39No.13July 2023基于减沙效益和经济效益的流域退耕还林方案优化喻丹1.2,董晓华1.2*,彭涛12,刘冀1.2,刘超3,刘芳4(1三峡大学水利与环境学院,宜昌4430 0 0;2 三峡库区生态环境教育部工程研究中心,宜昌4430 0 0;3宜昌市夷陵区水资源利用服务中心,宜昌4430 0 0;4宜昌市水文局,宜昌4430 0 0)摘要:实施退耕还林,是控制中国水土流失、改善生态环境的有效途径。如何
2、制定最具成本-效益的退耕还林方案,以平衡生态、经济和粮食安全之间的矛盾,是保证退耕还林工程可持续发展的关键。该研究以淮河上游的息县流域为研究区,在土地利用现状的基础上,一方面通过建立分布式水文模型SWAT(s o i l a n d w a t e r a s s e s s me n t t o o l)依次模拟各子流域的退耕还林操作得到泥沙削减系数,另一方面通过GDP(g r o s s d o me s t i c p r o d u c t)与土地利用现状的空间叠置分析得到各子流域进行退耕还林的GDP损失系数,据此分别构成退耕还林的减沙效益目标和经济效益目标,采用多目标遗传算法NSGA
3、-II优化求解子流域尺度的退耕还林方案。研究结果表明:1)建立的SWAT模型对研究区径流和泥沙的模拟精度较高,Nash系数分别在0.9 0 和0.7 0 以上,决定系数均大于0.8 0,且百分比偏差均控制在-2 0%2 0%以内,可认为SWAT模型能够用于评估退耕还林的泥沙削减效果;2)子流域泥沙削减系数范围为2 6.7 0 2 6 7 5.8 5t/km,并表现出从上游到下游逐渐减小的趋势,说明在流域上游的河源区实施单位面积的退耕还林能够取得更好的泥沙控制效果;3)子流域GDP损失系数在空间上呈现出较大的差异性,既有子流域出现了GDP的增加也有子流域出现了GDP的减小,对比发现在行政市或县主
4、要居民点所在的子流域进行退耕还林需要付出更大的经济代价;4)多目标优化求解得到的退耕还林方案集将人均耕地面积维持在1.0 410-31.5410-km,明显高于粮食安全的警戒水平,同时该方案集能够在仅损失30.13%37.6 7%经济产值的同时将泥沙产量削减53.54%6 9.8 6%,并达到区域可持续发展的土壤侵蚀水平。该研究提出的基于生态减沙效益和经济效益的子流域尺度退耕还林优化方法可为流域水土保持、退耕还林工程的科学规划提供借鉴和指导。关键词:土地利用;泥沙;退耕还林;SWAT模型;空间叠置分析;多目标遗传算法;淮河息县流域doi:10.11975/j.issn.1002-6819.20
5、2304044中图分类号:S19喻丹,董晓华,彭涛,等.基于减沙效益和经济效益的流域退耕还林方案优化J.农业工程学报,2 0 2 3,39(13):2 6 0-270.doi:10.11975/j.issn.1002-6819.202304044YU Dan,DONG Xiaohua,PENG Tao,et al.Optimization of Grain for Green Program based on sediment reduction and economicbenefitsJ.Transactions of the Chinese Society of Agricultural
6、Engineering(Transactions of the CSAE),2023,39(13):260-270.(in Chinese with English abstract)0引 言土壤是人类赖以生存和发展的基础自然资源,由于特殊的地理和历史条件,中国水土流失量大且面积广,是突出的生态环境问题之一。在众多导致水土流失的因素之中,耕地又被认为是水土流失的主要策源地。为了遏制水土流失、改善生态环境,中国于1999年启动实施了退耕还林工程,在长江、黄河和淮河等流域的长期实践看,其控制水土流失的效果十分显著。2 0 14年国家林业局启动了新一轮的退耕还林工程,旨在进一步推进中国的生态文明建设
7、。退耕还林是一项涉及社会、经济和生态多个方面的复杂系统性工程,工程的实施既要重视生态效益,又要兼顾社会经济效益,因此客观评价退耕还林工程的成本收稿日期:2 0 2 3-0 4-0 7 修订日期:2 0 2 3-0 6-11基金项目:国家自然科学基金项目(52 10 90 58)作者简介:喻丹,博士,讲师,研究方向为水文模拟及变化环境下的水文响应。Email:y u d a n 8 8 0 418 16 3.c o m通信作者:董晓华,博士,教授,博士生导师,研究方向为水文学及水资源。Email:文献标志码:Adoi:10.11975/j.issn.1002-6819.202304044文章编号
8、:10 0 2-6 8 19(2 0 2 3)-13-0 2 6 0-11http:/www.tcsae.org有效性已成为衡量其是否成功的关键指标。目前已有学者围绕退耕还林实施后的生态效益评估2-5、经济效益评估6 以及综合效益评估7 进行了大量的研究。然而由于受到气候、下垫面特征以及区域经济产业结构的影响,退耕还林工程的各项效益存在着一定的空间差异8-10 。仅根据退耕还林工程实施后所进行的“事后评价”难以直接给与管理决策者用于规划退耕还林工程空间布置所需要的“事先评价”的信息。如何预先制定最具成本-效益的退耕还林方案,以平衡生态和经济之间的矛盾,已成为保证退耕还林工程可持续发展的关键。多
9、数研究采用基于情景分析(scenario-based)的方法对退耕还林工程进行优化配置,即依据一定的准则条件形成几个退耕还林空间配置方案,再根据评价结果进行方案比选、择优。黄康等-13 均根据不同的坡度设置了多个退耕还林情景,李越强等14-15 则通过设置极端的情景,即假设全流域的耕地全都进行退耕还林。于航等16 以陡坡和河岸带的耕地为目标区域,设置了4种不同的退耕还林情景。显然这种基于情景分析的方法对于情景的假设是有限的,无法全局覆盖退耕还林的所有情景,从而丢失潜http:/www.tcsae.org第13期在的更优的配置方案。因此,需要探索基于全局搜索寻优的方法实现退耕还林空间配置方案的优
10、化。淮河流域是中国中东部地区重要的生态安全屏障,水土流失一直都是流域生态保护所面对的主要问题。2021年度全国水土流失动态监测成果显示淮河流域的水土流失表现为水力侵蚀为主,且主要分布在沂沭泗河、沙颍河以及淮河干流的上游区域17。在淮河典型流域开展退耕还林的空间配置研究对淮河流域水土流失综合治理以及生态环境保护具有重要的意义。本文选择淮河上游的息县流域为研究区,先运用生态水文模型模拟和GIS空间叠置方法分别构造退耕还林的生态效益目标和经济效益目标,再以子流域为水土流失治理单元构成多目标优化问题,采用具有全局搜索能力的智能算法求解子流域尺度的退耕还林优化方案,为淮河流域水土保持、退耕还林工程的科学
11、规划提供借鉴和指导。1研究区概况与数据来源1.1石研究区概况淮河流域(东经11155 12 12 5;北纬30 55 3636)位于中国东部、长江和黄河两大流域之间。淮河干流发源于河南省南阳市桐柏县和湖北省随州市随县交界处的桐柏山太白顶,流经河南、湖北、安徽和江苏,全长10 0 0 km,流域总面积为2 7 万km。本文选取淮河干流上游的息县水文站以上流域作为研究区(下文简称息县流域),研究区概况如图1所示。高程Elevation/m高 High:1 110低Low:33流域边界Basin boundary子域边界Subbasinboundary水系Rivernetwork雨量站Rainfal
12、lstation0510203040km图1研究区息县流域概况图Fig.1 Overview of Xixian River Basin as the study area息县流域面积为10 2 19.6 7 km,流域内以山区和丘陵区为主,地势由西南和北部的山区逐渐向东降低。根据中科院资源环境科学数据中心提供的中国多期土地利用遥感监测数据集(CNLUCC)以及沈玲怡等18 对息县流域2 0 0 0 一2 0 2 0 年的土地利用变化分析可知,2000年之后息县流域并未发生明显的土地利用变化,且耕地始终是流域主要的土地利用方式。以淮河干流为分界线,干流以北的耕地种植着玉麦和玉米为主的旱地作物,
13、干流以南的耕地多为水田,主要种植水稻。其次是喻丹等:基于减沙效益和经济效益的流域退耕还林方案优化N分析系统CFSR计算得到。2研究方法2.1退耕还林的减沙效益在降雨的驱动作用下,土壤侵蚀剥离的泥沙伴随着流域水文循环过程而迁移。流域生态水文模型不仅能对流域的水沙过程进行有效模拟,还能定量评价水土保持措施对流域水沙的影响。SWAT(s o i l a n d w a t e r a s s e s s m e n ttool)模型被广泛应用于预测不同土壤类别、土地利用类型以及管理措施下流域的水量、泥沙及水质变化,具有较强的适用性。根据研究区DEM数据,SWAT模型将全流域共划分为49个子流域,如图
14、1所示。本研究选择的土地利用数据时期为2 0 10 年,土地利用分类依据的是中国土地资源分类系统,为了与SWAT模型所含的土地利用数据库对应,如图2 a所示将其重分类为耕地、林地、灌木、草地和水域,其中耕地和林地为研究区的主要土地利用类型,分别占全流域的55.41%和2 2.0 8%。基于HWSD土壤提取本研究区的土壤类型空间分布如图2 b所示,土壤类型的划分依据FAO-90土壤分类系统,其中堆积土和不饱和始成土为研究区的两种主要的土壤类型,所占比例分别为2 8.6 6%和2 1.9 5%。由于261林地,主要分布在流域河源区域的桐柏山一带。林地垂直分带明显,可分为中山落叶阔叶林,低山常绿落叶
15、阔叶混交林,丘陵常绿针叶林及经济林19。然而,人口的增长和社会经济的发展导致了息县流域土地资源的不合理利用,流域内许多地区的自然植被遭到了破坏,加之耕地降低了土壤的抗侵蚀能力,因此引发了水土流失这一生态环境问题2 0 。息县流域已被水利部划分为国家级水土流失重点预防区,此外还被纳入第四期和第五期的全国水土保持重点工程建设项目规划2 1。1.2数据来源本研究使用的数据包括:研究区地形数据DEM(d ig it a le le v a t io n m o d e l),来源于美国太空总署和国防部国家测绘局联合测量的SRTM-DEM数据(http:/srtm.csi.cgiar.org/),其空间
16、分辨率为90 m;土地利用空间分布情况从中国科学院资源环境科学数据中心(https:/ m;土壤数据来源于联合国粮农组织和维也纳国际应用系统研究所所构建的世界土壤数据库Harmonized World Soil Database-HWSD(https:/www.fao.org/soils-portal/soil-survey/soil-maps-and-databases/harmonized-world-soil-database-v12/en/),空间分辨率为1km;本研究采用中国GDP(g r o s s d o m e s t ic p r o d u c t)空间分布公里网格数据来表
17、征研究区的社会经济状况,该数据来源于中国科学院资源环境科学数据中心(https:/ 个雨量站的日降雨数据以及流域出口息县水文站的流量观测数据均由淮河水利委员会提供,时间范围为19 9 5一2 0 10 年;其他气象参数通过美国国家气象预报中心提供的全球气候再262SWAT模型自带的土壤数据库没有覆盖到本文的研究区,故根据魏怀斌等2 3 提出的方法计算土壤参数。根据相关研究2 4 推荐的HRU(h y d r o l o g i c r e s p o n s e u n i t)划分阅值(坡度10%、土地利用0、土壤10%),最终在研究流域划分了52 8 个HRU。N耕地Farmland林地F
18、orest灌木Shrub草地Grassland水域Water城镇Urban land堆积土Cumulic anthrosols不饱和始成土Dystriccambisols饱和始成土Eutriccambisols石灰性冲积土Calcaricfluvisols饱和冲积土Eutricfluvisols饱和潜育土Eutricgleysols潜育士Calcicgleysols饱和浅层土Eutricleptosols弱发育淋溶土Haplicluvisols不饱和粘磐土Dystricplanosols饱和粘磐土Eutricplanosols石灰性粗骨土Calcaricregosols不饱和粗骨土Dystri
19、cregosols饱和粗骨土Eutricregosols不饱和变性土Dystricvertisols饱和变性土EutricVertisols水体Water图2 研究区土地利用空间分布和土壤类型空间分布Fig.2 Land use and soil types throughout the study area采用本研究收集到的研究区雨量站点的降雨数据和天气发生器的气象参数驱动SWAT模型,以研究流域出口息县水文站实测的径流和泥沙数据作为参照,以19 9 5年为预热期,19 9 6 2 0 0 3年为率定期,2 0 0 42 0 10 年为验证期,采用SWAT-CUP的SUFI-2算法分别进行月
20、尺度的径流参数率定以及月尺度的泥沙参数率定,选取统计指标Nash系数(Ens)、决定系数(R)和百分比偏差(PBIAs)对SWAT模型率定期和验证期的模拟精度进行评定。流域管理措施产生的效果会随着配置地点地形、水文气候、土壤等特征的不同而不同2 5。本研究利用SWAT模型来量化空间不同位置进行退耕还林起到的水土保持效果。将建立SWAT模型采用的2 0 10 年土地利用为基准,以子流域为退耕还林的基本单元,每次仅选择一个子流域(如子流域i)进行退耕还林操作,再利用率定好的SWAT模型进行模拟,得到该子流域退耕还林后流域的产沙值S,依次遍历所有子流域,根据式(1)计算子流域退耕还林的泥沙削减系数S
21、c。So-S;Sci=A耕i式中S为基准期土地利用条件下SWAT模型模拟的流域农业工程学报(http:/www.tcsae.org)0.1020406080kma.土地利用a.Land use0.1020406080kmb.土壤类型b.Soil type2023年产沙值,t;S,为子流域退耕还林后流域的产沙值,t;A耕i为基准期土地利用条件下第i个子流域中退耕还林的面积,即是基准期该子流域耕地的面积,km。2.2退耕还林的经济效益土地是区域经济发展的最基本物质载体,土地利用和经济发展有着密不可分的关系,社会经济结构的变化直接决定了土地利用方式的变化。淮河流域作为中国重要的粮食生产基地,其农业产
22、值和人均农业产值均高于全国同期水平,流域内的耕地资源有着极其重大的经济价值。因此可认为在息县流域实施退耕还林会给区域经济造成一定的影响。本研究采用资源环境科学与数据中心提供的中国GDP公里网格分布数据来表征研究流域的经济发展水平,该数据是在综合分析了人类活动形成的土地利用格局与GDP的空间互动规律的基础上,建立一二三产业GDP数量与土地利用类型的空间相关性模型来实现GDP的空间化,且认为耕地、林地、草地、水体这些土地利用类型只贡献第一产业GDP,第二三产业GDP值仅与城乡工矿居民地有关2 6 。因此本研究针对的是退耕还林后第一产业的GDP损失。将2 0 10 年的中国GDP公里网格分布数据与基
23、准期2010年的土地利用数据进行空间叠置分析,分别得到基准期土地利用条件下耕地和林地覆盖的GDP产值空间分布,再结合子流域空间分布可统计第i个子流域内耕地和林地分别对应的第一产业GDP,分别记为G耕;和G林i并假设在同一个子流域内,单位耕地或者单位林地产生的GDP是均匀的。那么在基准期土地利用条件下,第i个子流域的第一产业GDP值Go可表示为Goi=G耕i+G林i选择第i个子流域进行退耕还林,退耕面积为A耕,同时保持其他子流域的土地利用不变。退耕还林后该子流域原有耕地产生的GDP产值就全部损失了,林地面积在原有的基础上增加了。因此第i个子流域的第一产业GDP值 G,可表示为G;=GM(A#i+
24、Aw:)A林i式中A林;为基准期土地利用条件下第i个子流域内林地的面积,km。由式(2)和式(3)即可推求退耕还林过程中第i个子流域第一产业GDP产值的变化,记为G:AG;=G;-Goi=为了消除退耕还林面积的影响,定义GDP损失系数Gc为单位退耕还林面积的GDP变化,具体表达式为G;_G林i_G耕iGci=A耕iA林iA耕i重复上述操作依次遍历所有的子流域,即可得到空间不同子流域的 GDP损失系数。Gc0表示GDP在基准期的基础上增大,Gc29.2282.7782.77136.74136.74189.95189.95331.87149子流域边界及子流域编号Sub-basin boundary
25、 andsub-basin No.图4泥沙削减系数Sci在子流域空间的分布Fig.4 Spatial distribution of sediment reduction coefficient Seiat sub-basin scale3.3子流域GDP损失系数结合现状土地利用和GDP网格数据叠置分析得到基准期条件下研究流域总GDP值为515.6 1万元,且各子流域内林地和耕地覆盖的GDP值如表3 所示,对比可知对大部分子流域来说,耕地比林地产生了更高的GDP,这也与淮河流域农业经济占主导地位这一事实相符。仅有部分子流域,如3 4、41和44,其林地产生的GDP略大于耕地。此外,子流域GDP
26、损失系数计算结果如表3所示,可以看出实施退耕还林后子流域GDP的变化有着较大的空间差异,不仅出现了GDP的降低,同时还出现了GDP的增加。根据式(5)可知,从数值上看计算的Gc越大,说明在该子流域实施退耕还林操作后GDP的损失越小,从经济角度上看就越适合进行退耕还林。表3中Gci值大于0 的子流域有3 4和41,说明这些区域原有的耕地产生的经济效益较低,在这些子流域进行退耕还林能够在不降低区域GDP水平的前提下降低流域出口的泥沙含量,那么这些子流域内的耕地应尽可能多地参与退耕还林。退耕还林导致的GDP损失最严重的为子流域3 0,其Gc值为-57 56.8 3 元/km,表明该子流域内每1km的
27、退耕还林会使GDP产值降低57 56.8 3 元。从地理位置上来看子流域3 0 位于信阳市市辖区附近,说明该区域内的耕地为当地居民的主要经济来源,具有极其重要的经济价值,盲目地实施退耕还林会带来较大的经济损失。除此之外,子流域4、11、15和2 4的Gc值均小于-10 0 0 元/km,从地理位置上看这些子流域内聚集着桐柏县(子流域11)、信阳市(子流域4和2 4)以及息县(子流域15)的主要居民点,因此在这些区域进行退耕还林时需要权衡经济损失和泥沙控制之间的平衡关系。喻丹等:基于减沙效益和经济效益的流域退耕还林方案优化N耕地面积林地GDP耕地GDP12子流域编号Area of10242627
28、3339364341331.87537.70537.70782.82782.821085.001085.001429.781429.782675.850.510 20 30 40,265表3 子流域GDP损失系数Gc计算结果Table 3Results of GDP loss coefficient Gc at the sub-basin scale林地面积Area of16192230313538374049GDPofSubbasinNo.forest land/km10.0020.003141.2640.0050.00632.14427894710114812131415161718192
29、0km2122232425262728303134353637383940414244454647493.4子流域退耕还林方案优化结果遗传算法参数的选择对优化的最终结果有着一定的影响,为了确定合适的算法参数,对种群大小、迭代次数、交叉概率和变异概率这4个参数分别进行敏感性分析。如表4所示设定不同的参数变化情景,采用OAT(o n e-a t-a-t i m e)方法,即每次更改一个参数的值同时将其他所有参数值固定在初始值。各参数情景的优化结果如图5所示。由图5a可以看出当迭代次数从10 0 增加到6 0 0 时,优化出的Pareto前沿逐渐向着两个目标函数的最小值逼近,但是当迭代次数增加到10
30、 0 0 时,Pareto前沿并没有进一步被优化。图5b中显示的种群大小大于以10 0 为分界线,种群大小超过10 0 得到的Pareto前沿明显优于种群大小为2 0 和50,当种群大小从2 0 0 增加到50 0 时Pareto前沿在两个目标函数的寻优表现上并没有显著的提高。当交叉概率从0.1增加到0.7,图5c所示的Pareto前沿也越来越优,但是当交叉概率进一步增加到0.9时,Pareto前沿出现了反向移动而变得更差。GDP offarmland/km92.48110.81107.58322.58186.5132.900.0077.100.00234.7040.93138.3530.48
31、180.06203.86139.350.00491.300.0028.670.00124.460.0018.540.00133.780.0012.120.00115.180.00323.930.0030.970.00538.630.00337.410.005.320.00313.0616.2029.92204.01234.5389.2026.190.005.770.0038.100.00429.6346.4214.120.00569.7852.9830.0252.8459.530.00472.6114.1121.4768.69164.0185.2428.7414.8887.7669.7432.
32、620.00196.520.008.5060.1823.320.00194.33Ge/forest land/104元0.000.000.690.000.000.060.000.000.280.060.500.000.000.000.000.000.000.000.000.000.000.000.000.000.070.850.390.000.000.000.220.000.260.300.000.090.432.070.130.530.000.000.100.006.426.364.4734.506.430.097.3512.837.104.4319.5122.272.244.021.958
33、.460.356.7520.781.5723.2913.750.0331.341.405.380.600.1321.9342.620.0527.890.383.2920.770.254.170.314.100.284.190.140.352.78-693.93-573.97-366.29-1069.57-344.97-9.04-953.10-546.46-444.93-225.07-1375.68453.37-780.47-322.67-1051.29-632.68-287.90-585.95-641.43-507.38-432.48-407.38-57.91-1 001.13-423.24-
34、187.83-183.97-222.51-5 756.83-992.009.30-489.50-76.75494.53-439.55-55.81-191.21136.26-380.41-7.98-213.19-160.13-134.32-142.92266图5d中所示的变异概率设置为0.0 0 0 5得到的Pareto前沿与初始值0.0 0 0 1较为接近,且都劣于其余设定值得到的Pareto前沿。综合判断最终确定遗传算法的参数设定为:迭代次数6 0 0,种群大小40 0,交叉概率0.7,变异概率0.0 1。表4遗传算法参数变化情景Table 4Change scenarios of the
35、genetic algorithm parameters参数情景送代次数种群大小ParameterIterationscenariosnumber初始值100120024003600480051000106选代次数Iteration number2.0:8 8:48 0.800+6001.5:1000牛+1.06789101112泥沙产量Sedimentyield/10ta.代次数的敏感性分析a.Sensitivity analysis of iterationnumber106交叉概率Crossoverprobability2.0口8.0K1.51.067.89101112泥沙产量Sedim
36、ent yield/10stc.交叉概率的敏感性分析c.Sensitivity analysis of crossoverprobability图5遗传算法参数敏感性分析结果Fig.5 Sensitivity analysis results of the genetic algorithmparameters根据设定的遗传算法参数进行多目标优化,得到的Pareto前沿图6 所示,构成了子流域尺度退耕还林多目标优化问题的最优方案集,该方案集一方面将GDP损失控制在10 5万元至196 万元,损失率为3 0.13%至3 7.6 7%;另一方面将流域的产沙量控制在6 3.5万t至97.9万t,相较
37、于基准期的泥沙产量分别削减了6 9.8 6%和53.54%。PIMENTEL34 将土壤的形成速率0.51.0 t/hm定义为区域可持续发展的土壤侵蚀水平,认为土壤侵蚀的速度若大于这一临界水平时,土壤有机体的多样性和丰富度会大幅度降低从而不利于该区域农业的可持续发展。据此标准换算出研究区息县流域可持续发展水平下的土壤侵蚀量为51.2 10 2.5万t,不难看出对于优化出的Pareto前沿上的所有解均没有该水平的超过上限值,说明即使是在经济效益目标的牵制下,优化得到的退耕还林方案集的生态效益目标仍然能够达到较为满意的程度。农业工程学报(http:/www.tcsae.org)2.01.91.81
38、7221.21.1交叉概率变异概率PopulationCrossoversizeprobability500.9200.11000.32000.44000.55000.7106种群大小Population size2.0383881.556789101112泥沙产量Sedimentyield/105tb.种群大小的敏感性分析b.Sensitivity analysis of populationsize152.5106Mutation probability2.08.000 18.00050.0018.0051.5.8.018.100-0-8678910111213泥沙产量Sediment y
39、ield/10td.变异概率的敏感性分析d.Sensitivity analysis of mutationprobability2023年106Pareto前沿ParetofrontK n e e 解Knee solution1.0Mutation6.06.57.07.58.08.59.09.510.0probability泥沙产量Sedimentyield/10st0.000 1图6 多目标优化结果的Pareto前沿及Knee解0.0005Fig.6 Pareto frontier and Knee solution of multi-objective0.0010.0050.010.1变
40、异概率optimization results退耕还林直接减少了粮食生产的基本物质基础,对区域粮食安全产生的影响不容忽视。人均耕地面积这一指标常被用来分析地区或区域的粮食安全问题。经统计得到Pareto前沿解集中全流域的退耕面积为2 50 8.53889.0km,根据研究流域所辖行政区统计年鉴统计研究流域的人口数量(2 7 53 9 6 0 人),由此计算退耕还林后剩余的人均耕地面积为1.0 410-3 1.5410-km,均高于陈百明等3 5 预测的近期人均耕地面积阅值0.7 410 km和远期人均耕地阈值0.6 410-3 km,说明本研究所求的退耕还林方案解集即使是减少了研究区耕地的面积
41、也仍然能够满足人均粮食供给。另外,基于本研究构造的优化问题,即流域产沙量最小和GDP损失最小,在优化算法求解的过程中会引导着退耕还林面积更多地偏向于土壤易侵蚀且耕地产值较低的区域配置。时亚坤等3 6-3 7 的研究均强调,一方面退耕还林实施区的耕地质量差且粮食产量本就较低,另一方面耕地面积的减少促使农户从粗放式的生产模式转变为集约式,进而有效改善了区域环境和土壤肥力,因此可认为退耕还林的实施有助于提高粮食产量。综合分析来看,本研究求解的退耕还林方案集能够发挥水土保持和粮食供给的协同作用。针对Pareto前沿上的Knee解(如图6 所示),在土地利用现状上确定全流域的退耕还林方案,如表5所示按照
42、子流域退耕还林面积百分比从大到小排列。可以看出上游子流域退耕还林面积比例一般较大,基本上都超过了50%。退耕还林面积占比最小的几个子流域(子流域19、3 0、4和2 4)主要集中在信阳市辖区及周边。根据表5统计全流域退耕还林总面积为2 98 8.8 7 km,占现状土地利用耕地面积的44.19%。该退耕还林方案下对应的两个目标函数分别为泥沙产量7 5.1万t和GDP损失13 0 万元,经计算可知该方案以牺牲约2 5%的经济损失换取了约6 5%的泥沙削减量。盛莉等3 8 从土壤保持价值和涵养水源价值两个方面估算了中国水土保持的生态功能服务价值,在全国平均水平得出单位重量(t)的土壤保持量可提供5
43、8.9 6 元的生态服务价值。据此可以推断本研究中削减的泥沙量能够提供约7 9 9 4万元的生态服务价值,远大于估算的GDP损失值,因此可认为优化出的退耕还林方案对保护研究区生态环境有着巨大的潜力,具有较高的成本效益。第13 期表5Knee解对应的子流域尺度的退耕还林面积Table 5 GGP area at sub-basin scale corresponding tothe Knee solution退耕面积退耕百分比AreaofPercentage排序子流域returnedof returnedRank Subbasinfarmlandfarmland/km/%149233454265
44、1165145.1572181094110471144126134014391514169172718361925203721462213考虑到退耕还林工程需要面对退耕耕地种类和还林林地种类的选择,本研究遵循研究区自然规律,坚持因地制宜的原则:子流域内面积不占绝对主导的耕地类型优先参与退耕;若子流域内水田旱地面积相当,则根据淮河干流南北的气候差异,干流以北的耕地优先考虑退水田,干流以南的区域优先考虑退旱地;根据退耕区高程差异,高程8 0 0 m以上种植落叶阔叶林,高程50 0 800m种植常绿阔叶林及落叶阔叶林,高程50 0 m以下种植常绿针叶林及经济林。据此原则各子流域具体的退耕还林方式可分
45、为五类:河源子流域11内的水田零散分布且海拔较高,因此该子流域内优先退水田还落叶阔叶林;相邻的河源子流域2 6 内面积较少的水田优先退耕,退常绿阔叶林及落叶阔叶林;子流域10 和子流域3 6 分别位于桐柏县和信阳狮河区,面积不占主导的旱地优先退耕,还常绿阔叶林及落叶阔叶林;淮河干流以北的大部分子流域(19、12、13、16、18、19、2 1和2 4)以及干流以南的部分子流域(2 5、2 7、3 5和3 8)退水田还针叶林及经济林;位于淮河干流以南的剩余子流域则退旱地还常绿针叶林及经济林。4结 论本研究结合分布式水文模型模拟、空间叠置分析以及多目标优化方法,以流域泥沙产量和经济损失最小构造多目
46、标优化问题,求解得到子流域退耕还林最优解集,并从中选取了最具成本效益的退耕还林方案。本文的主要研究结论有:1)SW A T 模型径流模拟率定期(1996 一2 0 0 3 年)和验证期(2 0 0 4一2 0 10 年)的Nash系数Ens均大于0.90,决定系数R均大于0.90,百分比偏差PBIAs分别为-4.2 6%喻丹等:基于减沙效益和经济效益的流域退耕还林方案优化退耕面积退耕百分比AreaofPercentage排序子流域returnedRank Subbasinof returnedfarmlandfarmland/km/%167.1486.0189.2482.94162.3482.
47、61186.1979.39109.8578.8377.82412.1276.51133.9874.4120.7672.2416.1069.0221.9767.3521.1764.32104.7463.8613.5863.2878.6863.2286.3662.4216.3462.3617.8859.5817.6859.0935.1759.084.8957.5316.3056.85267和-5.7 3%;对于泥沙模拟率定期和验证期的Nash系数Ens均大于0.7 0,决定系数R均大于0.8 0,百分比偏差PBIAs 均控制在2 0%之内。由于本文在数据收集过程中受到一定的限制,只针对历史时期的径
48、流和泥沙过程进行了模拟。但是考虑到研究区的气象和下垫面条件没有发生较大的突变,因此认为基于历史时期模型率定和验证232824152522642271728729343016312032223318341352336383712388393540314124424433044193.2510.2659.4346.476.4239.816.9865.8914.97114.8137.0629.331.44127.91131.1661.83123.6785.5953.3244.735.0041.9456.3355.3453.6352.9552.9351.6449.4349.2548.3234.0332
49、.1731.7227.1127.0626.7026.3421.7019.9217.0313.8713.1412.95得到的参数仍能够沿用到近期。因此,精度评定结果能够说明SWAT模型对研究流域径流和泥沙的模拟具有较好的适用性。2)泥沙削减系数在空间上呈现出上游向下游逐渐减小的过渡趋势,说明在上游进行退耕还林操作能取得更好的泥沙控制效果。3)GDP损失系数其在空间上呈现出较大的差异性,在居民点聚集的区域进行退耕还林操作会带来较大的经济损失。4)子流域尺度退耕还林最优方案集能够在经济损失了3 0.13%3 7.6 7%的同时将流域泥沙产量削减53.54%6 9.8 6%,且能够平衡水土保持、经济损
50、失和粮食供给三者之间的矛盾。同时,最具成本效益的退耕还林方案对应的退耕面积占总耕地面积的44.19%且主要集中于上游子流域,由此带来的泥沙削减率约为6 5%,经济损失约为2 5%。参考文献1王晓燕,宁佳,史文娇.南方红壤丘陵区耕地生态问题识别与修复分区J.农业工程学报,2 0 2 2,3 8(2 4):197-2 0 6.WANG Xiaoyan,NING Jia,SHI Wenjiao.Identification ofecological problems and restoration zoning in the hilly red soilregion of southern Chin