大柳塔煤矿区植被恢复适度人工干预与调控研究

2018-08-29于昊辰牟守国王小予雷科玉

生态与农村环境学报 2018年8期

于昊辰, 牟守国, 王小予, 雷科玉

(中国矿业大学环境与测绘学院, 江苏徐州 221116)

煤炭资源高强度开采会损毁大面积土地资源,而当前煤矿区土地复垦率与植被覆盖率均较低,易加剧区域水土流失和土地沙化[1-2]。西方发达国家早在20世纪初便开始开展矿区土地恢复工作[3],包括分类种植、废弃矿山生态景观、环境监测与监督管理[4-6]等。国外不同学者也从植被演替方向、肥料使用、土壤改造等方面着手,促进了矿业废弃地植被自然恢复、植被恢复技术等领域的研究。而我国矿区土地复垦与生态环境治理起步较晚,在矿区植被恢复方面的研究集中在植被恢复与重建技术[7-8]、植被物种筛选[9]、植被恢复对生态系统中其他生物因子与非生物因子的影响[10-11]等方面。

干旱半干旱地区因水资源匮乏、土壤条件贫瘠以及生态环境脆弱等原因,受损生态系统难以通过自我修复完成植被恢复。故自然条件下该地煤矿区土壤恢复周期一般达到数十年以上,而进一步完成生态恢复更需数十年到数百年的时间。采煤后通过必要的补种措施等进行适当人工干预调控,能促进煤矿区人工植被恢复。而人为干预煤矿区植被恢复,无论是在空间、时间上,还是在干预类型、频度、强度等方面,都需要综合全面考虑。当前植被恢复人工干预的主要研究包括必要性、模式、条件、影响等方面[12-13],但人工如何干预植被恢复方面的深入研究较少。中度干扰假说提出,适度人工干扰有利于形成与保护物种多样性,可促进生态系统的演化更新与持续发展[14]。基于此,笔者以大柳塔煤矿区为例,根据1988—2016年研究区植被状况划分人工干预等级,通过比较不同人工干预等级下植被覆盖变化趋势及恢复状况,构建人工干预等级与植被恢复的关系模型并提出相应措施,试图为煤矿区植被恢复适度人工干预与调控以及相似区域植被恢复提供借鉴。

1 数据与方法

1.1 研究区域

大柳塔煤矿区(图1)位于陕西省榆林市神木县城西北约52 km处,地理坐标为39°13′53″～39°21′32″ N,110°12′23″～110°22′54″ E,面积126.55 km2。研究区为黄土高原与沙漠过渡区域,气候干燥,降水量少,蒸发量大,地势西北高东南低,海拔1 000～1 250 m,水系发育受限。研究区煤层埋深较浅,长期的采矿活动造成地表塌陷、台阶、地裂缝、植被退化等问题,导致矿区局部土地沙漠化现象严重,生态环境维稳性、抗干扰能力和自我恢复能力较低。

图1 大柳塔煤矿与不同干预等级样点图位置示意Fig.1 Location of Daliuta coal mine and its samples of different intervention levels (Ⅰ-Ⅳ)

研究区半固定沙地、固定沙丘和沙地沙丘间低地面积较大,其主要植被类型为沙生植被,以沙柳为主;另外在洼地、河流滩地和湖泊周围有少量湿地植物,大多浅根耐旱。根据其自然地理特点及对照实验设定,查询研究区卫星云图与相关自然地理资料,综合考虑影响植被恢复过程各因素,解译遥感影像图,分析在植被恢复过程中人工干预痕迹及程度差异,确定干预等级Ⅰ～Ⅵ样点的具体位置(图1)。基于前人对人工干预相关研究[15],结合实地调研与遥感影像解译,划分各样点对应的干预等级(表1)。依据大柳塔煤矿人工植被特点[16],其主要人工植被为成排的沙柳,间有少量沙蒿,故主要参考人工种植沙柳的面积比例与成排程度划分样点干预等级。

表1人工干预等级Ⅰ～Ⅵ划分标准

Table1Classificationdescriptionofhumaninterventionlevels(Ⅰ-Ⅵ)

干预等级样点描述干预程度 Ⅵ人工栽植沙柳面积近60%,且成排栽植强度干预 Ⅴ人工栽植沙柳面积约40%,且成排栽植较强干预 Ⅳ人工栽植沙柳面积约25%,基本成排栽植中等干预 Ⅲ人工栽植沙柳面积约15%,基本成排栽植较轻干预 Ⅱ人工栽植沙柳面积低于10%,未成排栽植轻度干预 Ⅰ基本无人工栽植沙柳,基本不受人类干预微度干预

1.2 数据来源与处理

选用遥感影像主要为美国陆地卫星影像Landsat TM/OLI,其中Landsat TM影像是包含7个波段的Landsat-5遥感影像数据,Landsat OLI为包括9个波段以及高分辨率全色波段的Landsat-8遥感影像数据。处理分析1988—2016年间7期无云、雾、条带等影响，植被生长状况良好的遥感影像(表2),并作图像预处理,包括遥感图像波段合成、辐射校正和裁剪等。

表2研究数据

Table2Researchdata

获取时间数据类型数据标识 1988-09-24Landsat 5 TMLT51270331988268BJC00 1992-07-17Landsat 5 TMLT51270331992199BJC01 1995-09-12Landsat 5 TMLT51270331995255CLT00 1998-07-02Landsat 5 TMLT51270331998183ULM00 2005-09-07Landsat 5 TMLT51270332005250BJC00 2010-08-20Landsat 5 TMLT51270332010232IKR00 2016-09-21Landsat 8 OLILC81270332016265LGN00

1.3 NDVI提取与植被覆盖度计算

提取植物光谱中的近红外(对树叶健康状况最敏感)与红外波段(被植被叶绿素强吸收)2个波段值,运用ENVI 5.3软件运算,基于光谱信息提取研究时点NDVI值(INDV)并分级,计算公式为

(1)

式(1)中,ρNIR为近红外波段地表反射率;ρR为可见光红波段地表反射率。

采用卫亚星等[17]提出的分级法处理分析图像数据,可较直观准确地反映研究区土地覆盖类型变化规律。由于研究样点区域范围较小,根据NDVI值将样点研究区域分为3类:NDVI值在>0.5～1之间的区域为a类(植被生长状况较好);NDVI值在0.3～0.5之间的为b类(植被生长状况一般);NDVI值<0.3为c类(植被生长状况较差)。其中高等级土地单元比低等级的植被生长状况好(a>b>c)。

提取NDVI值后,运用像元二分模型法估算植被覆盖度,计算公式为

(2)

式(2)中,FC为植被覆盖度,取值范围为0～1;INDV,soil为像元内完全裸土或无植被覆盖的NDVI值;INDV,veg为像元内完全被植被所覆盖的NDVI值。

从理论上讲,裸土INDV,soil取值应接近0,INDV,veg取值应接近1。但因大气、地表湿度和土壤类型等因素影响,INDV,veg值会随着空间不同而发生变化;另外,因Landsat TM/OLI影像空间分辨率为30 m,单像元实地面积为900 m2,而研究区内纯植被面积大于900 m2的区域较少,因此INDV,veg取值不为1。INDV,veg和INDV,soil一般取一定置信度范围内植被覆盖度的最大值和最小值,并根据经验估算[18]。结合1988—2016年干预等级Ⅰ～Ⅵ各样点NDVI值与对比试验设定,对其取值分别确定为0.01与0.6[19]。运用ENVI 5.3软件中Band Math模块计算得出像元值介于[0,1]的植被覆盖度图。依据前人研究[20]并结合研究区实际,将植被覆盖度划分为5个等级:Fc≤0.1,无植被覆盖(裸地或水体,N);0.1

1.4 植被恢复适度人工干预与调控模型

运用SPSS 22.0软件对1988—2016年T值与干预等级分别进行线性、指数、对数、二次以及三次线性方程拟合,确定拟合方程及相关系数。考虑普遍实用性,通过修正方差改进煤矿区植被恢复适度人工干预与调控模型。

2 结果

2.1 NDVI值变化

运用ENVI 5.3与Excel 2016软件统计并整理a、b、c类区域像元数,得出不同时期不同植被等级所占面积及其比例(表3)。在1988年,c类区域NDVI值所占比例均较高,除Ⅳ级干扰为84.8%外,其余均达90%以上,说明研究初始年NDVI值均较低。而随着年份的推移,c类区域逐渐向a和b类区域转移。特别是2016年,Ⅰ～Ⅳ级干预的c类区域面积比例均低于1%。另外,与初始1988年相比,Ⅳ和Ⅲ级干预a类区域面积比例增幅最大,分别为66.0和53.5百分点;而Ⅵ、Ⅴ、Ⅱ、Ⅰ级干预b类区域面积比例增幅最大,依次为63.7、59.7、52.6和74.7百分点。

2.2 植被覆盖度变化

根据植被覆盖度划分等级,运用ENVI软件中Density Sliced功能对植被覆盖度图进行密度分割,得到1988—2016年研究区干预等级Ⅰ～Ⅵ的植被覆盖度等级图(图2)。

表31988—2016年干预等级Ⅰ～Ⅵ的各区域NDVI面积占比

Table3TheproportionofNDVIareaintheinterventionlevels(Ⅰ-Ⅵ) %

干预等级1988年1992年1995年1998年2005年2010年2016年 abcabcabcabcabcabcabc Ⅵ0.00.0100.00.00.0100.00.021.878.20.01.498.62.894.72.50.086.813.217.463.718.9 Ⅴ0.58.990.60.110.889.11.215.283.51.814.184.11.847.850.40.968.530.725.068.66.5 Ⅳ0.015.284.80.025.874.20.038.761.30.031.069.117.167.415.40.095.74.366.033.10.9 Ⅲ0.06.493.61.020.878.32.257.540.31.335.063.741.258.80.019.880.20.053.546.20.3 Ⅱ0.07.093.00.011.888.30.749.350.00.016.483.60.793.06.30.398.90.840.159.60.3 Ⅰ0.01.898.20.01.698.40.02.597.50.02.897.20.613.585.90.944.354.823.576.50.0

a、b、c分别表示NDVI值>0.5～1、0.3～0.5、<0.3的区域。

图2 1988—2016年干预等级Ⅰ～Ⅵ植被覆盖度等级图Fig.2 Vegetation coverage grade map of intervention grade Ⅰ-Ⅵ in 1988-2016

运用ENVI统计干预等级Ⅰ～Ⅵ在不同时期各植被覆盖等级所对应的像元数,整理统计结果可得植被覆盖等级所占面积及比例,并对各植被覆盖等级赋予相应权重,计算不同干预等级下植被覆盖度的加权平均值T。干预等级Ⅵ时，各植被覆盖等级面积与T值见表4。各研究时点不同干预等级下T值见表5。

3 讨论

3.1 NDVI值与植被覆盖度变化分析

3.1.1NDVI值变化分析

1988—2016年,研究区域各类别植被覆盖区域面积存在一定变化。相同气候、降水等自然条件下,不同干预等级的植被覆盖变化趋势相近;早期a类区域面积接近于0,随着时间的推移,c类区域面积逐渐减少,a、b类区域面积逐渐增多,其中绝大部分增加的植被覆盖面积均属于b类;2010—2016年间,a类区域面积增幅十分明显。而不同干预等级下NDVI变化存在差异,分析各研究时点不同干预等级下NDVI变化,可发现一定干预等级下,随干预等级增加，c类区域面积所占比例逐渐减少,a、b类区域逐渐增加,但干预等级达到某一强度值后,随干预等级增加，c类区域面积所占比例总体呈增加趋势,a、b类区域呈下降趋势;在2005与2010年,随干预等级增加,c类区域面积占比出现2次下降。

表4干预等级Ⅵ时不同植被覆盖等级面积与T值

Table4AreaandTvalueofvegetationcoverageatinterventionlevelⅥ hm2

植被覆盖等级1988年1992年1995年1998年2005年2010年2016年 F0.000.600.776.136.1312.1713.45 H3.213.745.4516.7711.759.759.45 M5.629.3610.141.116.132.091.11 L15.1810.317.650.000.000.000.00 N0.000.000.000.000.000.000.00 T1)14.9218.7321.0431.0530.0032.1032.57

F为高植被覆盖; H为较高植被覆盖; M为中植被覆盖; L为低植被覆盖; N为无植被覆盖(裸地或水体);T为植被覆盖度的加权平均值。1)单位为%。

表51988—2016年不同干预等级T值汇总

Table5TvalueofinterventionlevelsⅠ-Ⅵin1988-2016%

干预等级1988年1992年1995年1998年2005年2010年2016年 Ⅵ14.9218.7321.0431.0530.0032.1032.57 Ⅴ15.1717.5520.8428.5729.6531.7632.46 Ⅳ16.4119.7124.3530.5331.4932.3133.54 Ⅲ12.7217.4823.0130.3631.3232.5433.43 Ⅱ11.2113.6518.1425.4428.7731.0632.71 Ⅰ7.789.4611.9918.9325.3628.0330.80

T为植被覆盖度的加权平均值。

3.1.2植被覆盖度变化分析

1988—2016年间,相同程度的人工干预等级下,随着时间推移植被覆盖均有明显好转,植被覆盖度总体呈现上升趋势;且随着干预等级增加,T值大体呈现倒“U”型趋势,即随着干预等级的增加,T值呈现先增加后降低的趋势,但是干预等级超过某一值时,也出现了植被覆盖度升高的现象;2010年干预等级为Ⅲ、Ⅳ、Ⅵ级时,植被覆盖度基本已达区域最佳值。

3.2 植被恢复适度人工干预与调控模型

运用Excel 2016软件制作1988—2016年不同干预等级下T值散点图(图3)。通过散点图可以看出,在Ⅲ、Ⅳ级干预程度下,各研究时点T值散点分布最为密集;而其他干预等级下,随着干预程度与Ⅲ级和Ⅳ级差距增大,T值散点图越来越分散。表明干预等级为中等时更利于植被恢复。

T为植被覆盖度的加权平均值。

不同研究时点T值变化与干预等级变化趋势基本一致,但仍存在一定差异。因此需要对各年份T值与干预等级的相关关系分别进行研究。运用SPSS 22.0软件对1988—2016年T值与干预等级关系散点图分别进行线性、指数、对数、二次以及三次线性方程拟合,确定拟合方程及相关系数,探究煤矿区植被恢复人工干预模型。以1998年为例,T值与干预等级Ⅰ～Ⅵ方程拟合相关系数见表6。

表61998年T值与干预等级方程拟合相关系数

Table6CorrelationcoefficientofequationfittingforTvalueandinterventionlevelin1998

方程拟合模型摘要参数估计值R2FP值常量b1b2b3 线性0.657.310.0520.462.01 对数0.8420.360.0120.426.44 二次0.9012.790.0312.527.96-0.85 三次0.9720.700.053.6519.12-4.550.35 指数0.636.830.060.0620.420.08

比较5种拟合方程,三次方程的决定系数(R2)最高且均为0.95以上,且显著性均较好。故说明三次线性方程可较客观反映植被覆盖度变化。考虑普遍实用性,对各三次线性方程作一定的修正改进,构建矿区植被恢复适度人工干预与调控模型。

不管何种干预等级,2005、2010和2016年研究区T值都基本达到阈值,而实际上煤矿区植被恢复工作大部分处在治理初期，或已破坏但未开始治理。该研究旨在探究人工如何干预与调控,能耗费较低成本且快速实现煤矿区植被恢复,并为相似区域生态恢复提供理论支撑,故研究时点的重点应在植被恢复初期和中期(1988—1998年)。基于此,以干预等级为自变量,以T值为因变量,构建煤矿区植被恢复适度人工干预模型,公式为

(3)

式(3)中,Tf为结束年份T值，%;T0为初始年份T值，%;ΔYt为年期修正;yb与ya分别为结束年份与初始年份;x为人工干预等级。

为评价模型合理性,估算1988—1998年不同干预等级下的T值,并与实际T值进行比较,采用下式计算误差系数EC(CE):

(4)

式(4)中,VE和VR分别为T的估算值与实测值,%。

据表7可知,1988—1998年不同干预等级T估算值与实测值的EC均小于13.5%,平均误差系数为-0.12%,平均误差系数绝对值低于5.8%,能较为客观地反映不同人工干预等级下T值随时间推移的变化,该模型可应用于煤矿区植被恢复最佳人工干预程度的选择。实际应用时,可根据研究区域的总土地面积通过式(3)计算T值,利用模型判断何种人工干预等级可低耗高效地达到预期植被恢复效果。

表7T的估算值与实测值误差分析

Table7TheerroranalysisofestimationandmeasuredTvalue

VE和VR分别为T的估算值与实测值; EC为误差系数。**表示估算值与实测值的误差系数大于10%。

3.3 措施建议

干预等级与T值关系见图4。

T为植被覆盖度的加权平均值。

大柳塔煤矿区面积共126.55 km2,运用ArcGIS软件矢量化估算大柳塔煤矿区村庄与工矿用地面积约为35 km2,要使全区植被恢复达到较好状态(植被覆盖度达0.35),综合考虑年期修正影响,最佳人工干预等级应介于Ⅲ～Ⅳ级之间(图4),因此大柳塔煤矿区植被恢复可参照Ⅲ～Ⅳ级干预开展植被恢复。即当人工栽植沙柳面积比例为15%～25%,且基本成排栽植时,可较低成本、较高效完成煤矿区植被恢复。

从干预等级及其分布情况看,同Ⅴ级与Ⅳ级干预相比,Ⅲ级与Ⅱ级干预样点距村庄与工矿用地较远,人类活动相对较少;而同Ⅵ级和Ⅰ级干预相比,Ⅲ级和Ⅱ级干预样点距村庄与工矿用地较近,人类活动也相应较多。

距离人类活动较近区域在植被恢复时投入成本较大,干预强度较大,其植被恢复效果会相对较好,但植被恢复工作不够经济;距离人类活动区域较远的区域投入成本较少,植被恢复效果较差;而Ⅳ级与Ⅲ级干预植被恢复效果较好且成本较为适中,因此综合考虑经济性与恢复效果,Ⅲ级干预模式更适于大柳塔矿区,可在开展该地区植被恢复工作时作为参考。

另外,大柳塔煤矿区人工干预植被恢复区域大多集中于村庄附近,且经人工干预后植被恢复效果均较好,而远离村庄区域基本无人工干预痕迹。推测这是导致大柳塔煤矿区总体植被覆盖状况不佳的原因之一,建议应参考Ⅲ级与Ⅱ级干预,将植被恢复逐步向距离人类活动区(村庄、工矿用地等)较适中的区域转移。