王曉紅,辛守英,陽麗虹,馬明浩,焦琳琳
(1.華北理工大學(xué) 礦業(yè)工程學(xué)院,河北 唐山 063210;2.河北省礦區(qū)生態(tài)修復(fù)產(chǎn)業(yè)技術(shù)研究院,河北 唐山 063210;3.河北省礦山綠色智能開采技術(shù)創(chuàng)新中心,河北 唐山 063210;4.國家知識產(chǎn)權(quán)局專利局專利審查協(xié)作廣東中心,廣東 廣州 510700)
植被是陸地生態(tài)系統(tǒng)的重要組成部分,在氣候調(diào)節(jié)、水土保持和生態(tài)系統(tǒng)穩(wěn)定方面發(fā)揮著關(guān)鍵作用[1]。近些年,氣候變化和人類活動(dòng)對植被的影響不斷加劇,而植被對氣候變化和人類活動(dòng)的響應(yīng)也變得更加敏感[2]。因此,監(jiān)測植被覆蓋的時(shí)空分布變化以及了解氣候因素和人類活動(dòng)的驅(qū)動(dòng)機(jī)制,對改善生態(tài)環(huán)境和提升生態(tài)服務(wù)功能具有重要的意義。
遙感技術(shù)以其覆蓋范圍廣、時(shí)空連續(xù)性強(qiáng)的優(yōu)勢,成為常用于監(jiān)測植被覆蓋分布變化的重要工具[3]。歸一化差異植被指數(shù)(normalized difference vegetation index,NDVI)作為一種高分辨率且具有良好時(shí)效性的遙感數(shù)據(jù),被廣泛用于區(qū)域植被覆蓋變化的動(dòng)態(tài)監(jiān)測[4]。如賈云飛等[5]利用2000—2019年的NDVI數(shù)據(jù)對黃土高原延河流域及其溝壑區(qū)植被覆蓋變化進(jìn)行分析。吳晶晶等[6]利用長時(shí)序NDVI數(shù)據(jù),對祁連山地區(qū)的生態(tài)修復(fù)前后植被覆蓋時(shí)空變化進(jìn)行分析。同時(shí),針對NDVI對植被覆蓋變化的驅(qū)動(dòng)機(jī)制的研究,如于璐等[7]運(yùn)用回歸分析和殘差分析法,分離出了不同時(shí)期京津風(fēng)沙源氣候變化和人類活動(dòng)對植被變化的影響。石淞等[8]運(yùn)用趨勢分析、Hurst指數(shù)、偏相關(guān)分析、殘差分析及相對作用分析法,深入解析黑龍江省植被時(shí)空變化特征及其對氣候變化和人類活動(dòng)的響應(yīng)機(jī)制。由此,植被覆蓋時(shí)空變化主要受氣候和人類活動(dòng)兩大因素的影響[9]。其中,氣溫和降水是能夠量化影響植被覆蓋變化的主要?dú)夂蛞蛩豙10]。而人類活動(dòng)因素對植被覆蓋變化的影響較難直接定量,主要通過獲取人類活動(dòng)的相對貢獻(xiàn)率來評估,目前定量區(qū)分人類活動(dòng)相對貢獻(xiàn)率應(yīng)用最廣泛的是殘差趨勢法[11],其能夠從氣候因素中有效地分離出人類活動(dòng)對植被覆蓋變化的影響。此外,對植被覆蓋的NDVI進(jìn)行植被覆蓋等級劃分,這有助于揭示在不同干擾情況下區(qū)域植被的動(dòng)態(tài)變化特征,并更好地探索植被覆蓋的動(dòng)態(tài)變化的驅(qū)動(dòng)機(jī)制[12]。
塞罕壩機(jī)械林場作為京津冀乃至華北地區(qū)防風(fēng)固沙和涵養(yǎng)水源的重要生態(tài)屏障,是集森林、草原、草甸、水體、濕地等為一體的復(fù)雜生態(tài)系統(tǒng),其植被覆蓋變化多受自然氣候環(huán)境和人類活動(dòng)等的綜合影響[13]。然而,以往對塞罕壩的研究主要關(guān)注植被景觀格局時(shí)空尺度效應(yīng),對氣候變化和人類活動(dòng)與不同等級植被覆蓋區(qū)相互作用的研究相對較少[14]。基于此,本文采用氣溫、降水和MODIS-NDVI數(shù)據(jù),通過趨勢分析、偏相關(guān)性分析、殘差分析和相對作用分析等方法,逐像元地分析2000—2020年塞罕壩林場的生長季植被NDVI的空間分布,并進(jìn)行植被覆蓋等級劃分,在不同植被覆蓋等級的基礎(chǔ)上分析氣候因素和人類活動(dòng)的相對貢獻(xiàn)率,揭示氣候變化和人類活動(dòng)與植被覆蓋變化的相互作用機(jī)制,闡釋塞罕壩植被覆蓋現(xiàn)狀特征和趨勢變化驅(qū)動(dòng)機(jī)制,為塞罕壩植被保護(hù)與森林經(jīng)營管理提供科學(xué)的理論指導(dǎo)。
塞罕壩機(jī)械林場(42°02′~42°36′ N、116°51′~117°39′ E)位于河北省承德市圍場滿族蒙古族自治縣的最北部(以下稱為塞罕壩林場),如圖1所示。該林場地處內(nèi)蒙古高原和冀北山地的接壤地帶,是省屬的大型國有人工林林場[15],總面積926.347 km2,森林覆蓋率為75.5 %[16];總場下有三道河口、千層板、北曼甸、陰河、第三鄉(xiāng)和大喚起6個(gè)分場。林場屬于寒溫帶半干旱半濕潤大陸性季風(fēng)氣候[17],年均氣溫-1.4℃,年降水量490 mm;冬季受內(nèi)蒙古高壓控制,盛行西北風(fēng)且風(fēng)力大,氣候寒冷;夏季受太平洋氣壓影響,多東南風(fēng),風(fēng)力小,降雨頻繁。塞罕壩林場地形復(fù)雜,分為壩上和壩下兩部分。其中,北部地區(qū)屬于壩上內(nèi)蒙古高原的一部分,土壤類型主要為風(fēng)沙土,伴隨有草甸土和沼澤土;中部和南部地區(qū)屬于壩下山地,土壤類型以棕壤和灰色森林土為主[14]。
圖1 研究區(qū)位置與地形
1.2.1 MODIS-NDVI數(shù)據(jù)
本文采用的是MOD13Q1 NDVI(https://ladsweb.modaps.eosdis.nasa.gov)遙感數(shù)據(jù),時(shí)間跨度為2000—2020年;時(shí)間分辨率是16 d,空間分辨率為250 m,行列號為h26v04。原始數(shù)據(jù)通過MRT(MODIS Reprojection Tools)工具批量進(jìn)行格式轉(zhuǎn)換、投影轉(zhuǎn)換、裁剪等操作,即將HDF格式轉(zhuǎn)換為GeoTIFF,坐標(biāo)投影為WGS84坐標(biāo)系,裁剪為塞罕壩林場區(qū)域;并使用ArcGIS軟件剔除NDVI數(shù)據(jù)中的無效值,留下NDVI真值。為了消除云層、大氣、衛(wèi)星軌道偏移與太陽高度角等的影響,利用最大值合成(MVC)方法合成月尺度數(shù)據(jù)。根據(jù)以往研究,將每年的4—10月作為生長季[18]。為了更好地反映研究區(qū)植被覆蓋生長狀況,計(jì)算生長季的NDVI均值進(jìn)行研究[19]。
1.2.2 氣象數(shù)據(jù)
本文使用的氣溫和降水氣象數(shù)據(jù)來源于國家地球系統(tǒng)科學(xué)數(shù)據(jù)中心——黃土高原分中心(http://loess.geodata.cn),在ArcGIS軟件中根據(jù)塞罕壩林場矢量邊界裁剪獲得2000—2020年分辨率為1 km逐月平均氣溫和降水?dāng)?shù)據(jù),并計(jì)算每年的生長季平均氣溫和降水?dāng)?shù)據(jù);利用ArcGIS軟件對氣氣溫和降水?dāng)?shù)據(jù)重采樣為與MODIS-NDVI數(shù)據(jù)像元大小和位置相同的柵格數(shù)據(jù)。
1.3.1 趨勢分析
本研究使用SPSS軟件對NDVI柵格數(shù)據(jù)逐像元進(jìn)行一元線性回歸分析2000—2020年塞罕壩林場的時(shí)空變化趨勢。其具體計(jì)算公式如下:
①
式中:Slope為NDVI變化趨勢斜率即NDVI年際變化速率;i為研究時(shí)段內(nèi)第i個(gè)年份;n為研究時(shí)段的年數(shù)(n=21);NDVIi為第i年的生長季NDVI均值。根據(jù)回歸特性,當(dāng)Slope<0,表示NDVI隨時(shí)間減少,即區(qū)域植被覆蓋呈減少趨勢;當(dāng)Slope>0,表示NDVI隨時(shí)間增加,即區(qū)域植被覆蓋呈增加趨勢;總之Slope的絕對值越大,NDVI隨時(shí)間變化越快,即區(qū)域植被覆蓋變化趨勢越快。此外,通過F顯著性檢驗(yàn)方法檢驗(yàn)NDVI的變化趨勢的可置信程度,檢驗(yàn)的具體計(jì)算公式如下:
②
1.3.2 偏相關(guān)性分析
本研究利用偏相關(guān)分析是逐像元計(jì)算2000—2020年塞罕壩林場年均NDVI與氣溫和降水的氣候因子的偏相關(guān)系數(shù)。其具體計(jì)算公式如下:
③
式中:rxy,z為在控制變量z的情況下x和y之間的偏相關(guān)系數(shù),即xy相關(guān)分析中剔除z的影響;rxy、rxz和ryz為不同變量間的相關(guān)系數(shù);x、y、z為NDVI與氣溫與和降水的氣候因子。其不同變量間的相關(guān)系數(shù)具體計(jì)算公式如下:
④
⑤
式中:r為NDVI與氣溫和降水的氣候因子的偏相關(guān)系數(shù);n為研究時(shí)段的年數(shù);k為自變量的個(gè)數(shù)。
1.3.3 殘差分析
殘差分析是區(qū)分氣候變化和人類活動(dòng)對NDVI的影響[21]。該方法通過對NDVI與氣候因素的氣溫和降水逐像元進(jìn)行多元線性回歸分析,擬合得到NDVI的預(yù)測值,將其視為氣候因素對NDVI的影響,進(jìn)而利用NDVI的預(yù)測值與從遙感影像獲取的NDVI觀測值做差,差值即為殘差,是表示人類活動(dòng)對NDVI的影響。其具體計(jì)算公式如下:
NDVICC=a×T+b×P+c
⑥
NDVIHA=NDVIobs-NDVICC
⑦
式中:NDVICC為多元線性回歸模型中NDVI預(yù)測值;NDVIobs為從遙感影像獲取的NDVI觀測值;NDVIHA為殘差;T和P分別為平均氣溫和降水;a和b分別為回歸系數(shù),c為回歸方程的常數(shù)項(xiàng)。
1.3.4 相對作用分析
Sun等應(yīng)用相對作用分析方法確定了氣候因素和人類活動(dòng)在華北地區(qū)植被NDVI變化過程中的相對貢獻(xiàn)率[22]。該方法是基于一元線性回歸趨勢分析和殘差分析得到NDVICC和NDVIHA的變化趨勢,在不同情景下,氣候因素和人類活動(dòng)對生長季NDVI相對貢獻(xiàn)率的計(jì)算方式不同[23],如表1所示。情景1和4表示氣候因素和人類活動(dòng)的共同作用導(dǎo)致的NDVI變化,情景2和5表示僅氣候因素導(dǎo)致的NDVI變化,情景3和6表示僅人類活動(dòng)導(dǎo)致的NDVI變化。
表1 不同情景下植被NDVI變化的驅(qū)動(dòng)因素相對貢獻(xiàn)率計(jì)算
根據(jù)塞罕壩林場植被覆蓋特征,本文將2000—2020年生長季NDVI均值劃分為3個(gè)等級。具體劃分標(biāo)準(zhǔn)為:低于0.55的區(qū)域被歸為低等植被覆蓋區(qū),高于0.65的區(qū)域歸為高等植被覆蓋區(qū),而0.55~0.65的區(qū)域則劃分為中等植被覆蓋區(qū);如圖2展示了這些不同等級植被覆蓋區(qū)在塞罕壩林場內(nèi)的空間分布。
圖2 不同植被覆蓋等級空間分布
根據(jù)圖2可知,低等植被覆蓋區(qū)主要分布在三道河口林場的邊緣地帶以及少量分布在千層板、北曼甸、陰河和大喚起分林場;中等植被覆蓋區(qū)主要分布在西部的三道河口和千層板、東部的北曼甸和陰河林場、南部的第三鄉(xiāng)和大喚起分林場;高等植被覆蓋區(qū)主要分布在南部的第三鄉(xiāng)和大喚起分林場以及少量分布在其他分林場。因此,塞罕壩林場植被NDVI覆蓋具有明顯的東高西低、南高北低的空間分布特征。這種空間分布受多方面原因的綜合影響:一方面受到土壤類型的影響,研究區(qū)東部主要為黑土,為冷濕草原草甸植被下形成的土壤,土壤較濕潤,有機(jī)質(zhì)含量多,土層深厚,適宜植被生長;西部土壤類型主要是風(fēng)沙土,土壤養(yǎng)分低,植被生長條件嚴(yán)苛[24]。另一方面水熱條件也在植被覆蓋變化中發(fā)揮重要作用,東部降水明顯多于西部,東南部氣溫高于東北部[25],導(dǎo)致了研究區(qū)植被覆蓋東高西低、南高北低的分布特點(diǎn)。此外,根據(jù)植被覆蓋分級的像元空間分布占比可知,低等植被覆蓋區(qū)面積遠(yuǎn)遠(yuǎn)小于中高等植被覆蓋區(qū),中高等植被覆蓋率高達(dá)94.62%,而低等植被覆蓋率僅為5.38%,如表2所示。
表2 不同植被覆蓋等級像元占比
本文對2000—2020年的塞罕壩林場植被NDVI覆蓋等級化數(shù)據(jù)逐像元進(jìn)行趨勢分析,并通過F檢驗(yàn)的分析結(jié)果,如圖3所示。
圖3 低(a)、中(b)、高(c)植被覆蓋區(qū)NDVI變化趨勢
根據(jù)圖3可知,從不同植被覆蓋等級的NDVI變化速率來看,低等植被覆蓋區(qū)NDVI增長速率較快(圖3a),平均增長速率為0.005 6/a,中等植被覆蓋區(qū)NDVI增長速率存在著空間差異(圖3b),西部NDVI增長速率較快,而中部和東部增長速率較緩,平均增長速率為0.004 3/a。這可能是由于西部早期植被覆蓋度不高,地勢平坦,造林工作主要致力于西部沙地[26]。高等植被覆蓋區(qū)NDVI增長速率較緩(圖3c),平均增長速率為0.003 1/a。總之,NDVI變化速率與植被覆蓋情況大致相反,高植被覆蓋區(qū)NDVI增長速率慢,低植被覆蓋區(qū)NDVI增長較快。此外,通過2000—2020年逐年NDVI植被覆蓋等級化的數(shù)據(jù)逐像元進(jìn)行趨勢分析,并通過F顯著性檢驗(yàn)的結(jié)果如圖4所示。
圖4 低(a)、中(b)、高(c)植被覆蓋區(qū)歸一化差異植被指數(shù)年平均變化
根據(jù)圖4可知,從塞罕壩林場逐年平均植被NDVI的變化情況來看,生長季年均NDVI呈顯著波動(dòng)上升趨勢(Slope>0,P<0.05)。近20年以來,塞罕壩林場植被NDVI增長趨勢可分為2個(gè)階段,2000—2007年為顯著上升階段(Slope>0,0.01
0,P<0.01),即2007年之后,塞罕壩林場全域植被NDVI平均增長速率明顯高于2000—2007的植被NDVI平均增長速率。從不同植被覆蓋等級的植被的變化速率來看,低植被覆蓋區(qū)NDVI增長速率快(圖4a),平均增長速率為0.005 1/a,中植被覆蓋區(qū)NDVI平均增長速率為0.003 3/a(圖4b),高植被覆蓋區(qū)NDVI增長較緩(圖4c),平均增長速率為0.002 7/a。
2.3.1 氣候因素對植被NDVI的影響
通過逐像元計(jì)算2000—2020年的生長季平均NDVI與平均氣溫和降水的氣候因素的偏相關(guān)系數(shù)并得到通過F顯著性檢驗(yàn)的結(jié)果,如圖5所示。
圖5 植被NDVI與氣溫(a)和降水(b)的偏高相關(guān)系數(shù)空間分布
根據(jù)圖5可知,植被NDVI與氣溫和降水的相關(guān)性主要呈現(xiàn)正相關(guān),且降水與NDVI的相關(guān)性更為顯著,其偏相關(guān)系數(shù)大都大于0.5。此外,在不同植被覆蓋區(qū)統(tǒng)計(jì)了NDVI與氣溫和降水的平均偏相關(guān)系數(shù),如表3所示。
表3 不同植被覆蓋區(qū)NDVI與氣溫和降水的平均偏相關(guān)系數(shù)
根據(jù)表3可知,隨著植被覆蓋程度的改善,NDVI與氣溫的相關(guān)性越高,與降水的相關(guān)性越低。此外,從不同植被覆蓋等級方面,評估氣候因素對塞罕壩林場植被NDVI覆蓋變化的相對貢獻(xiàn)率,如圖6所示。
圖6 低(a)、中(b)、高(c)植被覆蓋區(qū)氣候因素的相對貢獻(xiàn)率
根據(jù)圖6可知,低植被覆蓋區(qū)氣候因素對NDVI的相對貢獻(xiàn)率多低于50%(圖6a),其平均相對貢獻(xiàn)率約為39.22%。中高等植被覆蓋區(qū)東北部氣候因素的相對貢獻(xiàn)率相對較高(圖6b和6c),其相對貢獻(xiàn)率均在50%以上。這主要是因?yàn)槿眽瘟謭鰱|北部屬壩上地區(qū),海拔較高,人類活動(dòng)較少。
2.3.2 人類活動(dòng)對植被NDVI的影響
從不同植被覆蓋等級方面,評估人類活動(dòng)對塞罕壩林場植被NDVI覆蓋變化的相對貢獻(xiàn)率,如圖7所示。
圖7 低(a)、中(b)、高(c)植被覆蓋區(qū)人類活動(dòng)的相對貢獻(xiàn)率
由圖7可知,低植被覆蓋區(qū)的人類活動(dòng)對NDVI相對貢獻(xiàn)率多在50%~75%之間(圖7a),其平均相對貢獻(xiàn)率約為60.78%。中高等植被覆蓋區(qū)的西部和東南部人類活動(dòng)的相對貢獻(xiàn)率相對較高(圖7b和7c),其相對貢獻(xiàn)率均在50%以上。
自2000年以來,塞罕壩林場的不同植被覆蓋區(qū)的NDVI呈上升趨勢,其中低等植被覆蓋區(qū)的NDVI逐年平均增長速率為0.005 1/a,中等植被覆蓋區(qū)的NDVI逐年平均增長速率為0.003 3/a以及高等植被覆蓋區(qū)的NDVI逐年平均增長速率為0.002 7/a。與賀軍亮等對承德市全域的2000—2018年植被的NDVI整體呈上升趨勢研究結(jié)果[27]和王倩等對張家口和承德區(qū)域的2001—2020年植被的NDVI呈增加趨勢研究結(jié)果[28]一致。此外,2000—2020年,塞罕壩林場植被NDVI增長趨勢由2007年轉(zhuǎn)折分成2個(gè)階段,2000—2007年為顯著上升階段;2008—2020年為極顯著快速上升階段,這很大程度上是與2007年塞罕壩林場晉升成了國家級自然保護(hù)區(qū)、進(jìn)入了塞罕壩林場國家級自然保護(hù)區(qū)一期總體規(guī)劃項(xiàng)目實(shí)施有關(guān)[29]。
根據(jù)偏相關(guān)性分析可知,塞罕壩林場植被的NDVI時(shí)空變化是受氣候因素和人類活動(dòng)的共同影響。在氣候因素中,降水和溫度是與植被覆蓋增長直接相關(guān)的氣候因素,低等植被覆蓋區(qū)的氣候因素對NDVI的平均相對貢獻(xiàn)率約為39.22%,中高等植被覆蓋區(qū)東北部氣候因素的相對貢獻(xiàn)率相對較高,這主要因?yàn)槿眽瘟謭鰱|北部屬壩上地區(qū),海拔較高,人類活動(dòng)相對較少。趙明偉等[30]發(fā)現(xiàn),降水是驅(qū)動(dòng)華北平原北部和內(nèi)蒙古大部分區(qū)域植被覆蓋時(shí)空動(dòng)態(tài)變化的重要因素。邵雅琪等[31]研究結(jié)果表明,NDVI與降水呈顯著正相關(guān),與氣溫呈顯著負(fù)相關(guān)。而本研究中塞罕壩林場植被的NDVI與降水和氣溫的相關(guān)性均呈正相關(guān),且降水與NDVI的相關(guān)性更為顯著,即降水對植被的NDVI的影響大于氣溫,這可能是因?yàn)槿眽瘟謭鰵鉁乜臻g分布差異小、降水空間分布差異大所致。本研究在氣候因子中僅考慮了降水和氣溫,并未考慮相對濕度和極端氣候等其他氣候因素的影響,因此,在今后的研究中應(yīng)考慮該類因素影響使研究內(nèi)容更加準(zhǔn)確充實(shí)。研究期間,人類活動(dòng)在低植被覆蓋區(qū)的NDVI平均相對貢獻(xiàn)率約為60.78%。中高等植被覆蓋區(qū)的西部和東南部人類活動(dòng)相對貢獻(xiàn)率相對較高,這是因?yàn)槲鞑康貏萜教?是塞罕壩林場造林工程的核心區(qū)域,東南部海拔較低,屬于壩下山地區(qū),易受到人類活動(dòng)的影響。
本文基于2000—2020年長時(shí)序NDVI數(shù)據(jù)劃分植被覆蓋等級進(jìn)行塞罕壩林場植被覆蓋變化監(jiān)測,同時(shí)以氣溫和降水?dāng)?shù)據(jù)為基礎(chǔ)分析不同植被覆蓋等級的植被覆蓋變化受氣候和人類活動(dòng)的主要驅(qū)動(dòng)機(jī)制,得如下結(jié)論:
根據(jù)塞罕壩林場植被覆蓋特征,將林場內(nèi)NDVI低于0.55的區(qū)域歸類為低植被覆蓋區(qū),NDVI高于0.65的區(qū)域歸類為高植被覆蓋區(qū),NDVI為0.55~0.65的區(qū)域?yàn)橹械戎脖桓采w區(qū)。在不同等級植被覆蓋等級區(qū),植被NDVI的增長速率有所差異,低等植被覆蓋區(qū)NDVI增長速率較快,其次是中等植被覆蓋區(qū)和高等植被覆蓋區(qū)。尤其是2007年之后,塞罕壩林場全域植被NDVI平均增長速率明顯高于2000—2007的植被NDVI平均增長速率。此外,在不同等級植被覆蓋等級區(qū),NDVI的主要驅(qū)動(dòng)因素也有所不同,低等植被覆蓋區(qū)NDVI主要受到人類活動(dòng)的影響,中高等植被覆蓋區(qū)主要驅(qū)動(dòng)因子存在著空間異質(zhì)性,西部和東南部低海拔區(qū)的人類活動(dòng)是主要的驅(qū)動(dòng)因素,東北部高海拔區(qū)植被覆蓋變化則主要受到氣候因素的影響。