基于混合像元分解的香格里拉市高山松空間分布變化研究

劉蓉姣，張加龍，陳培高

(西南林業大學 林學院，云南 昆明 650224)

香格里拉市地處世界自然遺產“三江并流”的核心地帶，境內地形地貌復雜,森林資源豐富，但生態系統較為脆弱[1]。目前對其森林資源、生態環境和多樹種分布等方面已有研究，但對該區域某單一樹種分布變化的研究較少[2]。本研究圍繞香格里拉市優勢樹種高山松空間分布變化展開。高山松(Pinusdensata)是我國的特有樹種，源于第三紀云南松(P.yunnanensis)與油松(P.tabuliformis)的天然雜交種，耐寒、旱、貧瘠，廣泛分布于西藏東南部，川西部，云南西北部和青海南部2 700～4 200 m的高山地帶，多成純林分布[3-4]有較高的經濟價值和生態價值。

由于受傳感器分辨率的限制和地面地形氣候多樣性等因素的影響，遙感影像中，像元通常都是幾種地物的混合體，很少由單個均勻的地表覆蓋類組成[5]。使用高分辨率影像進行目視解譯所花費的代價太大，無法實現大范圍的分類[6]。而對于中低分辨率影像來說，混合像元的存在很大程度上降低了分類精度。如果每個混合像元都可以被分解并且其覆蓋類型組分(通常稱為端元組分)占像元的百分比(豐度)也能夠被求出，分類將會更精確。多年來，國內外學者利用混合像元分解技術提高遙感分類的準確性，并取得了有效進展。楊超[7]利用混合像元分解技術進行土地利用/覆被類型細分類。王銘[8]采用GOCI遙感影像提取時令水體及變化檢測。H.Luy[9]利用混合像元分解對衛星時間序列數據準確地估算了泛太平洋25個城市的植被。M.Dawelbaitetal[10]利用混合像元分解技術，監測North Kurdufan State (Sudan)1987-2008年的沙漠化過程，并在時間序列下分析了研究區內沙漠與植被覆蓋的空間特征及時間變化特征。陳虹兵等[11]利用像元二分模型法和線性光譜混合模型法對提取植被覆蓋度，對比分析2種模型的適用性，最后利用線性光譜混合模型生成了延慶區植被覆蓋分布圖,為縣域尺度應用Landsat 8 OLI影像估測植被覆蓋度提供參考。博文杰等[12]利用線性光譜混合模型提取莆田市常太鎮枇杷樹信息，繪出枇杷樹的(Eriobotryajaponica)豐度圖，據此估計出枇杷的種植面積。

本研究運用混合像元分解結合決策樹分類，提取香格里拉市高山松空間分布情況，分析出高山松面積的動態變化和高山松林地的發展趨勢，對香格里拉市森林資源的及時保護和進行宏觀決策提供參考依據。同時，建立的高山松分類指數對高山松分類研究具有一定參考價值。

1 材料與方法

1.1 數據源

香格里拉市是迪慶藏族自治州下轄市之一，位于云南省西北部的滇、川、藏大三角區域，地理位置99°20′-100°19′E、26°52′-28°52′N[13]。香格里拉市地形呈現西北高、東南低，最高點海拔5 545 m，最低點海拔1 503 m，海拔高差4 042 m，平均海拔3 459 m，地貌形態集山地、高原、盆地、河谷為一體[14]，市境內的土地總面積為1 141 739 hm2，林業用地面積950 911.7 hm2，占83.3%，森林覆蓋率為74.99%，據2016年全國森林資源二類調查統計，云冷杉、高山松、高山櫟、云南松為主要優勢樹種，占全市森林面積的90.8%[15](圖1)。

圖1 香格里拉市位置Fig.1 The location of Shangri-La

收集的數據有遙感數據和外業數據。遙感數據分別為香格里拉市1997、2002、2007、2012年的Landsat 5 TM數據和2017年Landsat 8 OLI影像。其中，2012年Landsat TM退役，采用時間相近的2011年的影像代替，形成每5 a間隔的時間序列影像集(表1)，研究區所在影像的軌道號為132/040、132/041、131/041，三景影像覆蓋了整個研究區，按照研究區邊界進行了裁剪，不存在拼圖問題[16]。同時收集了2016年的全國森林資源二類調查數據，香格里拉市DEM數據。外業數據收集于2018年4月，對香格里拉市42塊30 m×30 m的樣地(168個樣地點)進行坐標定位、周圍地物情況記錄、4個方位拍照和樣地內樹種進行樹高胸徑測量。

表1 研究區Landsat影像Table 1 Landsat images in the research area

1.2 研究方法

以Landsat遙感影像為數據源，首先，對5期影像進行輻射定標、大氣校正、幾何校正和地形校正的數據預處理；其次，通過最小噪聲分離、純凈像元指數計算和N維可視化，獲取植被、裸地和陰影的純凈像元；然后利用線性波譜分離和最小能量變化進行豐度提取；接下來，利用分量構建指數，獲取植被指數及DEM數據；最后，利用獲取的分類因子構建多元數據進行決策樹分類并驗證分類結果(圖2)。

圖2 技術路線Fig.2 Flowchart of methodology

1.2.1 數據預處理 利用Envi對Landsat遙感影像進行了輻射定標，大氣校正，幾何校正，地形校和鑲嵌裁剪的預處理，最終獲取研究所需影像。其中，幾何校正利用一幅精糾正過的香格里拉市SPOT 5影像為基準影像，校正中每一景影像選擇了30個控制點，平均RMS誤差控制在30 m內。香格里拉市地形和地勢落差較大，采用地形校正的方法來補償由于不規則的地形起伏而造成的地物亮度的變化[17]，校正過程采用J.Nichol et al[18]的坡度匹配模型。一副研究區影像涉及3景數據，在完成各景的地形校正后，利用研究區的矢量邊界進行切割，將裁剪后的3景圖像進行鑲嵌后形成研究區整體圖像。外業數據處理，根據外業調查數據所記錄的東西南北情況及所拍攝照片，運用ArcGIS獲取實地高山松光譜曲線。

1.2.2 基于混合像元分解的決策樹分類

1.2.2.1 高山松、土壤、陰影分量提取 樣地調查中高山松多成純林分布(圖3)，森林中郁閉度不是很高的地方存在一些裸露的土壤。在日光斜照射下，具有一定高度的地物都會形成陰影,樹木也不例外[19]，而且香格里拉市地勢落差較大，陰影在影像中的存在就更明顯了。所以在只包含某一個樹種的影像中一個像元內就可能同時包含植被、土壤、陰影多個地物的光譜信息。

圖3 外業調查Fig.3 Field survey

首先，對遙感影像進行最小噪聲分離(minimum noise fraction)判斷圖像數據維數，分離數據中的噪聲[20]，降噪后前5個波段已經集中了數據的85%以上的信息，對前5個波段進行PPI計算獲取純凈像元。其次，通過N維可視化結合圖像端元和參考端元選取高山松、裸地、陰影的端元。對照影像通過不斷調整獲取裸地和陰影的純凈端元，利用外業數據繪制的實地高山松波譜曲線和Landsat影像不斷調整，獲取高山松純凈端元。最后，在N維可視化中勾選的高山松繪制波譜曲線并與外業調查數據獲取的高山松進行對比，得到高山松與實地外業數據獲取的高山松波譜曲線對比情況(圖4)。利用IDL編寫的完全約束最小二乘算法[21]混合像元分解擴展工具對地物端元進行混合像元分解，避免了提取的豐度中出現負值的現象。

圖4 高山松光譜曲線對比Fig.4 Comparison of spectral curves of Pinus densata from 1997 to 2017

1.2.2.2 決策樹提取高山松空間分布 基于該特征運用豐度中3個分量構建植被-土壤指數、植被-陰影指數、植被-土壤-陰影指數和歸一化多分量指數4種新的指數[22-23]，新指數用于決策樹中來區別植被里的高山松、云南松、云冷杉和其他植被的關鍵數據。利用ENVI 5.3進行波段計算構建指數(表2)。

表2 新指數計算公式Table 2 New index calculation formulae

運用4種指數分離高山松與其他樹種，決策最優指數進行分類。4種指數總體精度都不差，但相比較歸一化多分量指數精度最高，歸一化多分量指數總體精度都高于87%，其他3種指數精度在84%～87%，結合歸一化多分量指數經過多次試驗通過制定規則見圖5，運用該規則進行分類，將影像分為高山松、云南松、云冷杉、其他林、非林地5類。香格里拉市云南松和云冷杉是占比較大的優勢樹種，占香格里拉市林地總面積的31.85%[24]，該規則通過海拔和其他樹種的歸一化指數先分類出了香格里拉市兩大優勢樹種，這能更準確得到高山松空間分布情況。

圖5 決策樹分類示意Fig.5 Knowledge expert decision tree

1.2.2.3 精度評價 研究采用混淆矩陣[25]進行分類結果精度驗證。通過統計樹種的面積與香格里拉市面積并計算其比例[26]，運用公式(1)計算出驗證點個數為507個。隨機篩選出507個驗證點，其中高山松的驗證點76個，云冷杉105個，云南松42個，其他樹種208個,非植被76個。對本研究分類結果進行了驗證。

(1)

式中，n為驗證點個數，B為卡方分布臨界值，由(1-α)/k決定，其中α為置信度，k為分類個數。πi為第i類地物的比例，bi為每個像元的分類絕對誤差(一般取5%)。

(2)

(3)

2 結果與分析

對分類結果進行了驗證，精度評價結果見表3。各年分類結果總體精度分別為69.42%、76.73%、81.07%、78.90%和76.53%。其中高山松的分類精度較好，其次是云冷杉和云南松，非林地與其他林錯分的情況相對多一些。

表3 精度評價結果Table 3 Accuracy evaluation results

2.1 高山松空間分布變化

運用基于混合像元分解的決策樹分類，得到分類結果(圖6)，紅色區域為其他樹種空間分布情況，其他樹種包括了如車桑(Dodonaeaviscosa)、杜鵑(Rhododendronsimsii)、華山松(Pinusarmandii)、箭竹(Fargesiaspathacea)其他等多種樹種；綠色區域為高山松空間分布情況；藍色區域為云南松空間分布情況；電氣石綠色區域為云冷杉空間分布情況。其他樹種與云冷杉空間分布情況有重疊部分，而且兩者的空間分布面積最大。

圖6 分類結果Fig.6 Classification results

對高山松總面積進行了統計(圖7)。1997年高山松的總體覆蓋面積為1 608.67 km2，占香格里拉市植被面積的17.08%。2002年高山松總面積為1 393.07km2，占香格里拉市植被面積的14.79%。2007年高山松總面積為1 358.66 km2，占香格里拉市植被面積的14.42%。2012年高山松總面積為1 480.46 km2，占香格里拉市植被面積的15.72%。2017年高山松總面積為1 541 km2，占香格里拉市植被面積的16.36%。2002年高山松覆蓋面積比1997年減少了13.40%，2002年比2007年減少了2.47%，2012年高山松覆蓋面積比2007年增加了8.96%，2017年比2012年增加了4.06%。不同年份香格里拉市高山松分布狀況不同，不同年份受降水、氣溫、人為行為、自然更替等不同因素的影響都可能影響香格里拉市高山松覆蓋面積。

圖7 各年香格里拉市高山松面積統計Fig.7 P.densata areas of the study years in Shangri-La city

應用ArcGIS中疊加分析的交集取反分析(symmetrical difference)[27]，將相鄰年份高山松空間分布進行交集取反分析，獲得相鄰年份高山松增量或減量的空間分布情況，通過分析圖8可知，1997-2002年5 a間高山松急劇減少，2002-2007年高山松總面積減少量低于1997-2002年。2012-2017年高山松面積增量小于2007-2012年的增量。

圖8 高山松分布變化Fig.8 Distribution change maps of P.densata

2.2 香格里拉市高山松隨高程的分布變化

結合DEM數據得知，每年高山松空間分布的海拔范圍情況分別為：1997年高山松分布的高程范圍是1 827.3～4 516.5 m，2002年高山松分布海拔范圍1 995.2～4 483.8 m，2007年高山松分布海拔范圍2 094.3～4 371.1 m，2012年高山松分布海拔范圍2 100.4～4 360.5 m，2017年高山松分布海拔范圍2 092.6～4 250 m。對同一海拔區間的不同年份高山松覆蓋面積占當年高山松總面積的百分比進行了簡單統計對比。海拔在1 800～2 300 m的不同年份高山松覆蓋率分別為：0.27%、0.11%、0.02%、0.04% 和0.02%。2 300～2 800區間每年高山松森林覆蓋情況分別是：4.28%、3.31%、2.38%、2.55%和2.59%。2 800～3 300 m區間高山松面積占當年高山松總面積的比率分別為：30.93%、29.95%、30.02%、35.61%和36.81%。3 300～3 800 m海拔區間內的高山松覆蓋率是最高的，該區間高山松面積占當年高山松總面積的比率分別為：56.14%、58.00%、57.60%、58.45%和58.36%。3 800～4 300 m海拔區間的高山松面積與當年高山松總體面積之間的比值分別為8.30%、8.54%、3.23%、3.36%和2.22%。4 300～4 800 m海拔區間的高山松面積與當年高山松總體面積之間的比值分別為0.08%、0.09%、0.004%、0.002%和 0.003%。該區間是每年中高山松分布最少的區間。

每年高山松分布區間最低海拔為1 827.3 m，最高海拔為4 483.8 m。利用ArcGIS軟件的空間分析結合重分類將高山松分布區域海拔高度劃分6個高程區，分別為1 800～2 300、2 300～2 800、2 800～3 300、3 300～3 800、3 800～4 300、4 300～4 800 m。由圖9可知，不同海拔區間對比高山松分布情況為：3 300～3 800 m>2 800～3 300 m>3 800～4 300 m>1 800～2 300 m>4 300～4 800 m。

圖9 各年各海拔區間高山松空間分布Fig.9 Spatial distribution of P.densata at different altitudes in each year

3 結論與討論

以香格里拉市高山松為研究對象，利用1997年至2017年的Landsat遙感影像，通過混合像元分解結合決策樹分類，獲取了1997年至2017年近20 a高山松空間分布，分析其長時間序列變化情況。1997、2002、2007、2012年和2017年分類總體精度分別為：69.42%、76.73%、81.07%、78.90%和76.53%。

通過多次試驗將分類的閾值從一年推廣到多年進行分類。其中，NDVI植被指數區分林地非林地閾值為0.45，高山松歸一化多分量指數閾值為0.333，云南松歸一化多分量指數閾值為0.208，區分云冷杉的歸一化多分量指數值為0.362，高程為3 700～4 350 m。

香格里拉市高山松覆蓋率在1997-2002年急劇減少，2002年高山松覆蓋面積比1997年減少了13.40%，2007年高山松覆蓋面積比2002年減少了2.47%，從2007-2012年和2017年香格里拉市高山松覆蓋率都在逐年增加，增幅沒有很大，2012年高山松覆蓋面積比2007年增加了8.96%，2017年比2012年高山松覆蓋面積增加了4.06%。其中1997-2002年變化最顯著。

結合高程分析了高山松的空間分布情況。1997、2002、2007、2012年和2017年香格里拉市高山松主要分布在海拔高度為2 800～3 800 m的地區，高山松生長分布最低海拔1 800 m，最高海拔為4 400 m。不同海拔區間對比高山松分布情況為：3 300～3 800 m>2 800～3 300 m>3 800～4 300>1 800～2 300 m>4 300～4 800 m。其中高山松最適宜生存海拔為3 300～3 800 m。總之，不同年份香格里拉市高山松分布狀況不同，不同年份的降水、氣溫、人為行為、自然更替等不同因素都會影響香格里拉市高山松覆蓋面積。

隨著遙感技術的不斷發展，對樹種分類的探究也越來越被需要。研究樹種分布空間變化首先要能準確進行樹種分布的提取，這就取決于分類結果，然而影響分類精度的關鍵因素之一是分類方法的選取。由于本研究區面積過大、地形差異大、植被覆蓋率高、生物多樣性豐富以及遙感影像分辨率低等因素的影響，運用光譜信息、坡度、坡向、紋理等特征難以區分出不同的樹種，因此采用混合像元分解結合決策樹進行分類，在解決多光譜數據由于大氣散射和空間分辨率低等問題造成的光譜混合現象后，結合前人經驗運用決策樹進行分類，獲取了香格里拉市高山松空間分布情況，為后序分布變化提高可靠數據。今后希望能用分辨率較高的數據對林分區別較大的區域進行混合像元結合決策樹準確分類更多樹種。