楊明芬吳旭

（1.西藏自治區(qū)科技信息研究所，西藏 拉薩 850000；2.成都理工大學(xué)，四川 成都 610059）

生態(tài)環(huán)境是人類生存和發(fā)展的重要物質(zhì)基礎(chǔ)之一，對(duì)生態(tài)環(huán)境脆弱性變化展開(kāi)研究，是制定區(qū)域可持續(xù)發(fā)展戰(zhàn)略的重要前提。對(duì)生態(tài)環(huán)境進(jìn)行檢測(cè)和評(píng)價(jià)分析，有助于及時(shí)了解環(huán)境質(zhì)量變化的規(guī)律并采取科學(xué)的應(yīng)對(duì)措施。生態(tài)環(huán)境脆弱性評(píng)價(jià)是區(qū)域生態(tài)綜合評(píng)價(jià)的重要內(nèi)容和生態(tài)規(guī)劃的重要依據(jù)，近年來(lái)在地理學(xué)、生態(tài)學(xué)和環(huán)境科學(xué)研究中廣泛應(yīng)用。

以RS和GIS技術(shù)支持的生態(tài)環(huán)境脆弱性研究被更多的學(xué)者采用。生態(tài)環(huán)境脆弱性評(píng)價(jià)可分為定性評(píng)價(jià)和定量評(píng)價(jià)，定量方法對(duì)環(huán)境的脆弱性、穩(wěn)定性、敏感性和外部環(huán)境脅迫對(duì)系統(tǒng)造成的影響進(jìn)行定量描述。［1］利用生態(tài)數(shù)學(xué)語(yǔ)言提出了生態(tài)脆弱帶基礎(chǔ)判定的數(shù)學(xué)方法。［2］用統(tǒng)計(jì)模型對(duì)生態(tài)環(huán)境脆弱性密度分布進(jìn)行了初步測(cè)算。［3］通過(guò)多源遙感數(shù)據(jù)下植被覆蓋變化研究對(duì)脆弱性評(píng)價(jià)進(jìn)行了研究。由于生態(tài)環(huán)境脆弱性研究涉及地學(xué)、生態(tài)學(xué)等多個(gè)領(lǐng)域的研究，具有很強(qiáng)的綜合性，因此傳統(tǒng)評(píng)價(jià)方法使用的研究方法、評(píng)價(jià)指標(biāo)體系和分級(jí)標(biāo)準(zhǔn)呈現(xiàn)多樣化。各種研究方法的重點(diǎn)集中在確定生態(tài)脆弱性評(píng)價(jià)指標(biāo)、權(quán)重設(shè)定、指定各指標(biāo)穩(wěn)定度的標(biāo)準(zhǔn)方面。這些方法需要較長(zhǎng)的時(shí)間來(lái)處理前期數(shù)據(jù)資料，研究周期較長(zhǎng)，無(wú)法快速定期對(duì)生態(tài)環(huán)境脆弱性現(xiàn)狀和變化趨勢(shì)進(jìn)行評(píng)價(jià)監(jiān)測(cè)，而且評(píng)價(jià)結(jié)果不夠全面。

文章在前人建立的生態(tài)環(huán)境脆弱性定量評(píng)價(jià)方程基礎(chǔ)上，研究一種新的區(qū)域態(tài)環(huán)境脆弱性定量評(píng)價(jià)模型?；诙嘣炊鄷r(shí)相的衛(wèi)星遙感影像數(shù)據(jù)提取植被指數(shù)、土地覆蓋分類、DEM坡度與坡向數(shù)據(jù)作為評(píng)價(jià)模型輸入的原始數(shù)據(jù)。采用多元非線性回歸方法提取生態(tài)環(huán)境脆弱性定量評(píng)價(jià)所需的模型參數(shù)。通過(guò)評(píng)價(jià)模型的計(jì)算實(shí)現(xiàn)研究區(qū)生態(tài)環(huán)境脆弱性的快速、定量評(píng)價(jià)，形成脆弱性評(píng)價(jià)分級(jí)圖，并對(duì)模型評(píng)價(jià)的精度和穩(wěn)定性進(jìn)行了驗(yàn)證。

1 研究方法

1.1 基本思路

文章采用生態(tài)環(huán)境脆弱性的因子相關(guān)分析方法，參考［4］劉燕華和王經(jīng)民的經(jīng)驗(yàn)公式,構(gòu)建脆弱性定量評(píng)價(jià)模型：

式（1）中F為生態(tài)環(huán)境脆弱性值，U為土地覆蓋類型的變化梯度，f為變化區(qū)域的植被覆蓋度，其中：

式（2）中U為變化的土地類型面積；Si為該類型區(qū)域前一年土地總面積為土地覆蓋類型變化梯度；i為土地覆蓋類型，參考前人研究基礎(chǔ)并結(jié)合與植被有關(guān)的特性，劃分為林地、草地、濕地、水體、人工用地、沙地等六種類型；

式（3）中f是估算的指標(biāo)覆蓋度，DNVI是歸一化植被指數(shù)，NDVImax和NDVImin分別是研究區(qū)同一年內(nèi)的最大值和最小值。

通過(guò)遙感監(jiān)督分類提取NDVI、土地覆蓋變化范圍，并通過(guò)編寫(xiě)IDL程序計(jì)算土地變化的面積和總面積的最大值與最小值，求出變化梯度，根據(jù)NDVI求出植被覆蓋度，代入式（4）求出由于變化梯度導(dǎo)致的生態(tài)環(huán)境脆弱性值。在IDL環(huán)境中，對(duì)研究區(qū)的數(shù)字地形進(jìn)行高程、坡度與坡位的綜合分析，選取與植被相關(guān)性最大的因素代入評(píng)價(jià)方程中，計(jì)算得出脆弱性評(píng)價(jià)結(jié)果。將數(shù)字地形模型與評(píng)價(jià)結(jié)果疊加，得出區(qū)域生態(tài)環(huán)境脆弱性綜合評(píng)價(jià)等級(jí)。實(shí)驗(yàn)流程如圖1。

圖1 實(shí)驗(yàn)流程圖

1.2 實(shí)驗(yàn)數(shù)據(jù)

研究區(qū)位于四川省若爾蓋縣，四川、甘肅、青海三省交界處，地處青藏高原東部邊緣，海拔3400-3900m，東西與南北最大距離約150km，面積約10620km2，地理坐標(biāo)為 102°08′E 至 103°39′E，32°56′至34°19′N之間，是長(zhǎng)江和黃河流域上游地區(qū)重要的水源涵養(yǎng)地，也是中國(guó)殘存面積最大的高原泥炭沼澤濕地，整個(gè)研究區(qū)形成一個(gè)相對(duì)完整的地理自然區(qū)，對(duì)于研究區(qū)域生態(tài)脆弱性變化極具代表性。

1.2.1 Landsat數(shù)據(jù)。Landsat1-8系列衛(wèi)星是由美國(guó)航空航天局 NASA（National Aeronautics and Space Ad?ministration）發(fā)射的對(duì)地觀測(cè)衛(wèi)星。文章采用Land?sat5 TM數(shù)據(jù)提取土地覆蓋分類信息，并計(jì)算土地類型變化率。TM的覆蓋范圍為之間的所有陸地區(qū)域，東西掃描范圍約183km，南北掃描范圍約為170km，數(shù)據(jù)更新周期為16天，空間分辨率30m。文章獲取2001-2010年研究區(qū)內(nèi)的Landsat5 TM遙感數(shù)據(jù)（每年4景，共40景數(shù)據(jù)）?？紤]到研究每一景影像可能存在多云覆蓋的問(wèn)題，我們選擇的單景影像平均云量小于10%的影像數(shù)據(jù)，成像時(shí)間集中在8-9月。選擇出每一年內(nèi)成像質(zhì)量較好的影像數(shù)據(jù)，對(duì)各景數(shù)據(jù)進(jìn)行鑲嵌，根據(jù)研究區(qū)矢量地圖裁剪得出研究區(qū)影像數(shù)據(jù)，對(duì)裁剪的研究區(qū)數(shù)據(jù)進(jìn)行大氣校正、正射校正、去除薄云處理。其中投影為UTM、坐標(biāo)系為WGS-1984，如圖2。

圖2 研究區(qū)歷年影像

1.2.2 MODIS數(shù)據(jù)。文章采用MOD13A1 NDVI數(shù)據(jù)產(chǎn)品，共獲取2001、2003、2006、2008年4個(gè)年份，16天間隔的產(chǎn)品數(shù)據(jù)（每年23景，共計(jì)92景影像數(shù)據(jù)，空間分辨率500m）。在ENVI軟件環(huán)境中對(duì)每年數(shù)據(jù)及十年數(shù)據(jù)進(jìn)行HANTS濾波處理，生成長(zhǎng)時(shí)空序列數(shù)據(jù)集，采用波段運(yùn)算工具，對(duì)NDVI的月均、年均值進(jìn)行計(jì)算。對(duì)數(shù)據(jù)集進(jìn)行Hants濾波處理，然后利用波段運(yùn)算工具計(jì)算2001-2010年均NDVI的數(shù)據(jù)集，如圖3。

圖3 NDVI年均數(shù)據(jù)集

1.2.3 地形數(shù)據(jù)。本文選取的DEM數(shù)據(jù)為SRTMDE?MUTM（90m分辨率），提取高程數(shù)據(jù)；坡度數(shù)據(jù)為SRTMSLOPE（90m分辨率）；坡向數(shù)據(jù)為SRTMAS?PECT（90m分辨率），如圖4所示。這三種數(shù)據(jù)用于分析其與NDVI的相關(guān)性，估算未知區(qū)域NDVI的取值，為定量評(píng)價(jià)模型提供所需參數(shù)。

圖4 地形數(shù)據(jù)

1.3 實(shí)驗(yàn)步驟

1.3.1 計(jì)算土地覆蓋類型的變化梯度。研究區(qū)典型地物類型包括：河流、湖泊、沼澤、沼澤化草甸、濕草甸、灌叢、灌木林、森林、裸地、裸巖、沙地和建筑用地等［5］。其中河流、湖泊、泥炭沼澤、沼澤化草甸和濕草甸是研究區(qū)的主要濕地類型［7］。文章參考《全國(guó)遙感監(jiān)測(cè)土地利用/覆蓋分類體系》、前人研究成果，并考慮模型計(jì)算的簡(jiǎn)便性，把土地覆蓋類型分為：草地、林地、裸地（裸土、沙地、裸巖、積雪）、水體（河流、湖泊）、沼澤、人工用地（居民地、耕地、道路）六類［6］。為保證分類精度和信息提取的快捷，文章采用ENVI軟件下的支持向量機(jī)分類（Support Vector Machine，SVM）方法對(duì)遙感圖像進(jìn)行監(jiān)督分類，再通過(guò)人工操作修改錯(cuò)誤分類。分類過(guò)程中結(jié)合野外調(diào)查的先驗(yàn)知識(shí)選擇訓(xùn)練樣本。影像中難免存在不可去除的云層干擾等因素，由此產(chǎn)生錯(cuò)誤分類，文章通過(guò)與影像觀測(cè)時(shí)間最接近的數(shù)據(jù)進(jìn)行比較人工修改分類。把相同年份的遙感預(yù)處理數(shù)據(jù)與分類結(jié)果數(shù)據(jù)疊加，在遙感預(yù)處理影像中隨機(jī)選取各類型純凈像元，采用ENVI的分類后處理工具中的誤差混合矩陣進(jìn)行分類精度評(píng)價(jià)，分類結(jié)果精度如表1所示（以2001年為例）。

表1 土地覆蓋類型分類精度表

從各時(shí)相分類精度可以看出，總體精度在85%以上，可以滿足應(yīng)用需求。其中沼澤與水體、沼澤與草地、裸地與草地誤分的現(xiàn)象明顯，這主要是由于高植被覆蓋沼澤與草地的光譜信息比較接近，地覆蓋草地與裸露地表、沙地等裸地在圖像上比較相似，從而導(dǎo)致誤分；同時(shí)少量區(qū)域由于厚云根本上無(wú)法提出，對(duì)林地、草地也存在誤分。

文章在獲取2001、2003、2006、2008年土地覆蓋類型的分類數(shù)據(jù)基礎(chǔ)上，采用的分類統(tǒng)計(jì)工具計(jì)算每種覆蓋類型相鄰年份的變化值，計(jì)算出土地覆蓋類型的變化梯度，對(duì)土地覆蓋類型的動(dòng)態(tài)變化進(jìn)行分析。通過(guò)計(jì)算結(jié)果獲得變化的類型、面積、像元數(shù)量，并統(tǒng)計(jì)出土地覆蓋類型變化的百分比。

1.3.2 植被覆蓋度估算。然后基于多源遙感數(shù)據(jù)的生態(tài)環(huán)境脆弱性評(píng)價(jià)，評(píng)價(jià)模型參數(shù)組要依賴于植被覆蓋度和土地覆蓋類型的變化梯度。但一個(gè)地區(qū)的生態(tài)環(huán)境脆弱性往往是多種因素共同的結(jié)果。研究植被指數(shù)、地形高度、坡度以及坡向變化的規(guī)律有利于區(qū)域生態(tài)環(huán)境脆弱性定量評(píng)價(jià)的準(zhǔn)確性。文章對(duì)獲取的2001-2010年NDVI數(shù)據(jù)集與研究區(qū)地形數(shù)據(jù)的相關(guān)性分析，步驟如下：

將NDVI數(shù)據(jù)集與研究區(qū)高程、坡度、坡向數(shù)據(jù)疊加，并配準(zhǔn)；

編寫(xiě)程序，在研究區(qū)內(nèi)指定的矩形窗口范圍內(nèi)生成隨機(jī)采樣點(diǎn)，獲取每一年份同一經(jīng)緯度的NDVI值、海拔高程、坡度值、坡向值；

通過(guò)程序自動(dòng)提取對(duì)應(yīng)坐標(biāo)點(diǎn)的植被指數(shù)、值高程值、坡度和坡向值，生成NDVI-高程、NDVI-坡度、NDVI-坡向的數(shù)據(jù)記錄文件，以便進(jìn)一步分析它們之間的關(guān)系。

1.3.3 NDVI與地形數(shù)據(jù)的相關(guān)性分析。將采集到的10個(gè)時(shí)相NDVI數(shù)據(jù)與地形高程進(jìn)行研究：利用IDL編程繪制出對(duì)應(yīng)坐標(biāo)點(diǎn)的NDVI隨地形高度變化的散點(diǎn)圖，對(duì)每一時(shí)相2000多個(gè)采樣點(diǎn)進(jìn)行統(tǒng)計(jì)分析。從繪制出的散點(diǎn)圖可以觀察出：NDVI值隨著高度的增加，是先增加后減少，其變化規(guī)律和植被分布的垂直地帶性是相符合。坡度為0的時(shí)候，ND?VI的值在坡度0.4-0.5區(qū)間分布均勻，這表明在草地這樣的平坦地形，NDVI的值是沒(méi)有反映坡度的變化規(guī)律。因此NDVI隨坡度的變化表現(xiàn)平直，沒(méi)有明顯的升降變化。同時(shí)，NDVI與坡向的變化分布均勻，沒(méi)有明顯的規(guī)律性變化，數(shù)據(jù)與各種擬合函數(shù)相關(guān)性極低。NDVI與坡度和坡向的相關(guān)系數(shù)大都小于0.2，因此跟研究中只針對(duì)NDVI與高度值進(jìn)行回歸分析和函數(shù)擬合。根據(jù)散點(diǎn)圖顯示，NDVI與高度之間有明顯的規(guī)律性變化，但無(wú)明顯的線性或非線性函數(shù)關(guān)系，考慮對(duì)NDVI值進(jìn)行某種轉(zhuǎn)換，研究分析轉(zhuǎn)換后的值與高程之間的相關(guān)關(guān)系，如表2所示。

表2 NDVI與地形數(shù)據(jù)的相關(guān)性分析

文章嘗試對(duì)NDVI與NDVI均值差的絕對(duì)值取自然對(duì)數(shù)，然后分析轉(zhuǎn)換后的值與高程的關(guān)系、轉(zhuǎn)換值，通過(guò)編寫(xiě)程序生成與高程值的散點(diǎn)圖，如圖5所示。根據(jù)散點(diǎn)圖的分布情況進(jìn)一步分析擬合函數(shù)的形態(tài)。

圖5 NDVI與高程的非線性分析

1.3.4 基于多元非線性回歸的NDVI估值。對(duì)多年份數(shù)據(jù)的轉(zhuǎn)換值Y與高度值的散點(diǎn)分布進(jìn)行分析，發(fā)現(xiàn)散點(diǎn)的分布特征十分接近e的冪函數(shù)類型，因此設(shè)計(jì)一種多元非線性回歸方程對(duì)數(shù)據(jù)點(diǎn)進(jìn)行擬合，并通過(guò)IDL編程計(jì)算回歸參數(shù)：Y=C*eax+b，其中x為高度值，a、b、c為待定參數(shù)。通過(guò)以上計(jì)算阿對(duì)未知點(diǎn)的NDVI進(jìn)行估值，從而估算出研究區(qū)的植被覆蓋度。圖6描述了計(jì)算的待定參數(shù)形成的擬合函數(shù)與Y的關(guān)系。

擬合函數(shù)對(duì)NDVI估計(jì)值的精度檢驗(yàn)如表2所示。

表2 NDVI估值精度檢驗(yàn)

圖6 NDVI非線性函數(shù)擬合

1.4 定量評(píng)價(jià)模型

非線性回歸模型表示為：

其中f為期望函數(shù)，xn為第n個(gè)響應(yīng)的回歸向量或自變量向量。模型的期望響應(yīng)是關(guān)于參數(shù)的非線性函數(shù)。θ表示未知參數(shù)，當(dāng)分析一個(gè)特定的數(shù)據(jù)集時(shí)，假定向量xn，，n=1，…N，N是固定的，并幾種考慮期望響應(yīng)關(guān)于θ的相依關(guān)系。構(gòu)造N維向量η（θ），其第n個(gè)向量為：ηn（θ）=f（xn，θ），n=1，2…，N；記非線性回歸模型為：

對(duì)于非線性的情形，利用迭代的方法求解最小二乘估計(jì)量。利用期望函數(shù)的線性近似，通過(guò)迭代改善的初始值θ0，直到迭代無(wú)變化為止。

將期望函數(shù)f（xn，θ）在θ0處進(jìn)行一階泰勒展開(kāi)：

與η（θ0）相比，現(xiàn)在點(diǎn)y。因此取更好的參數(shù)值θ1=θ0+δ0，開(kāi)始另一次迭代，繼續(xù)計(jì)算新的殘差z1=y-η（θ），新的導(dǎo)數(shù)矩陣V1和新的增量。重復(fù)上述過(guò)程，直至達(dá)到收斂為止。

2 實(shí)驗(yàn)結(jié)果與分析

在生態(tài)環(huán)境評(píng)價(jià)中，由于評(píng)價(jià)單元不同，評(píng)價(jià)結(jié)果會(huì)存在差異,因此評(píng)價(jià)單元的選擇與研究區(qū)范圍密切相關(guān)［7］。文章以遙感影像中的柵格數(shù)據(jù)為載體，以每一格柵格單元為評(píng)價(jià)的最小單位，在模型計(jì)算的基礎(chǔ)上，按照指數(shù)的數(shù)值大小進(jìn)行等級(jí)劃分，確定其優(yōu)劣程度。把脆弱性指標(biāo)值劃分為5個(gè)等級(jí)：無(wú)脆弱性、輕微脆弱性、一般脆弱性、中度脆弱性和嚴(yán)重脆弱性。文章采用此方法對(duì)2003、2006、2008三年的評(píng)價(jià)結(jié)果進(jìn)行分級(jí)，以此分級(jí)閾值作為評(píng)價(jià)，最終得到研究區(qū)的生態(tài)環(huán)境脆弱性定量評(píng)價(jià)結(jié)果，如圖7所示。

圖7 區(qū)域生態(tài)環(huán)境脆弱性定量評(píng)價(jià)計(jì)算結(jié)果

由計(jì)算結(jié)果可知研究區(qū)評(píng)價(jià)值主要集中在0.2-0.8之間，生態(tài)環(huán)境脆弱性嚴(yán)重的區(qū)域主要集中在北部和中部部分地方。這里主要是裸露地表、沼澤邊緣和建筑用地、耕地等受人類活動(dòng)影響較為頻繁的區(qū)域。從圖7可以看出，各評(píng)價(jià)等級(jí)相間分布，各等級(jí)間沒(méi)有明顯的界線。

3 結(jié)論

文章從少量因子出發(fā)，先抽取和植被生態(tài)方面有關(guān)的因子（土地覆蓋類型變化率、植被覆蓋度），研究生態(tài)環(huán)境脆弱性和它們之間的關(guān)系。再研究植被指數(shù)和地形變化的規(guī)律，逐個(gè)因子添加不僅可以使研究的復(fù)雜問(wèn)題變的簡(jiǎn)單而且有助于揭示它們之間的聯(lián)系。主要研究成果包括：實(shí)現(xiàn)了生態(tài)環(huán)境脆弱性的遙感定量研究。解決了如何從遙感影像中提取出反映生態(tài)環(huán)境脆弱性的因子，并把這種參數(shù)帶入到生態(tài)環(huán)境脆弱性的方程中，從而實(shí)現(xiàn)生態(tài)環(huán)境脆弱性的遙感定量研究是問(wèn)題的關(guān)鍵。從植被的角度出發(fā)把造成生態(tài)環(huán)境脆弱性的因子劃分為：植被覆蓋度、土地覆蓋類型變化梯度、地形高度、坡度以及坡向。從遙感影像中提取出土地覆蓋類型變化的梯度，利用NDVI和植被覆蓋度的關(guān)系，定量表達(dá)植被覆蓋度。通過(guò)模型求出生態(tài)環(huán)境脆弱度，實(shí)現(xiàn)生態(tài)環(huán)境脆弱性的遙感定量研究，研究具有創(chuàng)新意義。揭示了研究區(qū)植被指數(shù)隨地形高度、坡度、坡向變化的規(guī)律。研究中對(duì)每一時(shí)相采集到2000多個(gè)具有典型意義的點(diǎn)，數(shù)據(jù)豐富，資料可靠。通過(guò)研究發(fā)現(xiàn)，NDVI和高程存在著明顯的相關(guān)關(guān)系，通過(guò)IDL編程求出了具體的定量表達(dá)關(guān)系式。實(shí)驗(yàn)表明植被指數(shù)隨地形坡度和坡向的變化規(guī)律不明顯，只簡(jiǎn)單地分布在均值周圍?；诙嘣蔷€性回歸方法，建立了生態(tài)環(huán)境脆弱性植被遙感定量研究的模型，通過(guò)改造傳統(tǒng)意義上的脆弱度方程，實(shí)現(xiàn)了參數(shù)的替代。該模型中脆弱度被定義為土地覆蓋類型變化的梯度和植被覆蓋度的商。通過(guò)反映生態(tài)環(huán)境脆弱性本底值的方法來(lái)評(píng)價(jià)劃分生態(tài)環(huán)境脆弱性等級(jí)。