牛 婷,文 方*,常夢迪,高 潔,李蘊輝,程 帆,杜姿影,陳曉斐
(1.新疆維吾爾自治區環境保護科學研究院,新疆 烏魯木齊 830011;2.新疆環境污染監控與風險預警重點實驗室,新疆 烏魯木齊 830011;3.新疆清潔生產工程技術研究中心,新疆 烏魯木齊 830011;4.國家環境保護準噶爾荒漠綠洲交錯區科學觀測研究站,新疆 烏魯木齊 830011;5.成都信息工程大學 資源環境學院,四川 成都 610225)
富含大量有機污染物的工業、農業以及居民生活廢水排入湖水中,造成湖泊水質嚴重富營養化,從而導致藻類以及其他浮游生物的大量爆發,形成水華,而葉綠素作為藻類的主要組成部分,且葉綠素濃度便于人工測定。因此葉綠素濃度常常作為反映水體富營養化程度的重要評價指標[1]。常規的水質監測是通過采集水樣進行實驗室分析[2],比較費時費力,對大面積水體的水質狀況難以全面反應,遙感技術可以實現湖泊水質的高頻、大范圍、準實時監測[3]。葉綠素主要分為葉綠素a和葉綠素b,其中葉綠素a又是藻類植物中最主要的色素,所以葉綠素a濃度是水色遙感監測中最重要的參數之一。利用遙感數據反演葉綠素的濃度主要有三種模型,即經驗模型、半分析模型和分析模型,不同的反演模型因為水體組成成分的不同而使得反演精度也不同。經驗模型基于統計分析,算法比較成熟、過程簡單,但缺乏物理依據,反演出的模型精度較低;半分析模型通過水體組成成分的光譜特征與統計分析相結合,有一定的物理依據,模型精度較高,外推性適宜性較好;分析模型基于大量的地面實時采集數據,物理機理復雜,模型精度高,但模型建立難度大。
Gordon在1975年提出了一種新型模型稱為分析模型,它主要依賴于水體組分、固有光學量、表觀光學量之間的關系,通過代數方程直接求解葉綠素濃度,具有極高的物理性與相關性。但是相較于經驗/半經驗模型法,分析模型需要大量水體各種組分的固有光學特性數據,具體實施較為困難,建立算法的難度比較大[4]。劉忠華[5]、RUNDQUIST等[6]通過大量的實測數據,表明葉綠素a濃度在特定波長處于地表反射率具有較高的相關性。祝令亞[7]和溫新龍等[8]以太湖為例采用波段組合算法,分別建立了葉綠素a與MODIS數據、環境一號衛星CCD數據的反演模型。李旭文等[9]基于Landsat5 TM數據和地表實測數據建立了經驗模型,對梅梁湖區藍藻生物量進行了估算,證明葉綠素a濃度和差異植被指數的相關性較高,同樣也有諸多研究表明TM不同波段組合與葉綠素濃度具有較好的相關性[10-16]。
博斯騰湖位于內陸河口,是焉耆盆地各種河流的聚集地,受人類活動影響,屬于典型的二類水體。本文通過對博斯騰湖水體葉綠素a濃度的反演,研究博斯騰湖中水體葉綠素的濃度分布,從而更好地應對水體富營養化對博斯騰湖水資源的影響,解決博斯騰湖水體富營養化問題。基于以上分析,本課題選用經驗模型和半分析模型對博斯騰湖水體葉綠素a濃度進行反演建模研究。
博斯騰湖(Bosten Lake),維吾爾語意為“綠洲”,位于中國新疆維吾爾族自治區焉耆盆地東南面博湖縣境內,是中國最大的內陸淡水吞吐湖(見圖1)。湖區介于東經86°19′-87°28′,北緯41°46′-42°08′之間,湖面的平均海拔高度約為1048米,總面積1646平方公里,水域面積約為800多平方公里。其既是位于其上游區域開都河、黃水溝等水系的尾閭,又是位于其下游區域孔雀河水系的源頭,相當于一個巨大的“調節水庫”,是開孔河流域的“心臟”。湖區南邊和北邊海拔較高,受湖泊水位變化影響較小,西側和東側地勢平坦更容易受到湖泊水位變化的影響,周圍生長著廣茂的蘆葦,是中國四大集中產葦區之一。

圖1 博斯騰湖地理位置圖
本研究采樣點的布設沿用巴州環境監測站采樣點分布(見圖2),采樣點位坐標見表1。對照影像數據,選用2010年6月9日、2010年9月1日采樣點實測數據進行反演研究。

表1 博斯騰湖水質監測點位坐標

圖2 博斯騰湖水環境監測采樣點分布圖
TM成像儀是Landsat4和Landsat5攜帶的傳感器,它是一種改進的多光譜掃描儀。為保證反演的精度,選擇天氣晴朗、研究區上方天空少云或者無云及衛星過境時間最靠近實測數據時間的TM影像。根據地面水質監測數據,選取2010年6月10日、2010年8月29日的影像進行建模。
遙感器觀測目標物的輻射或者反射的電磁能量時,不可避免地受到大氣分子、氣溶膠和云粒子等大氣成分的吸收與散射的影響,再加之遙感器本身的光電系統特征、地形以及太陽高度等因素,所得到的測量值與目標物真實的光譜轄射亮度或者反射率等是有一定差異的。地表參數的遙感定量反演需要糾正目標輻射的不確定性信息。
對遙感影像數據的預處理主要包括幾何校正、輻射定標、大氣校正。本研究針對這三個方面的要求對獲取的遙感影像進行相應的預處理。
2.2.1 TM影像的幾何精校正
利用已知坐標的遙感圖像數據,選擇適當的控制點,然后基于多項式糾正模型的校正。本研究以谷歌影像數據為底圖,在ENVI軟件中對遙感影像進行幾何精糾正,糾正誤差控制在0.5個像元以內。
2.2.2 TM影像的輻射定標
本研究基于ENVI軟件,利用ENVI的Landsat Calibration模塊,讀入TM數據各個波段的定標系數,將原始數據DN值轉換為大氣外層表面反射率(輻射亮度值)。
2.2.3 TM影像的大氣校正
大氣校正就是將轄射亮度或者表觀反射率轉換為地表實際反射率,目的是消除大氣散射、吸收、反射引起的誤差。本文采用ENVI軟件的FLAASH大氣校正模塊對TM輻射定標后的數據進行大氣校正,計算影像的真實反射率。需要輸入的參數中幾何參數可以從影像數據的元數據獲得,包括影像中心經緯度、太陽高度角、傳感器類型和遙感影像獲取具體時間;大氣組分參數包括能見度(氣溶膠光學厚度,由FLAASH根據暗目標法,根2200nm暗目標的反射率與660nm暗目標的反射率之間的比值關系,反演出TM影像的一個平均的能見度數據)、氣溶膠模式和大氣模式,大氣模式和氣溶膠類型參數則根據博斯騰湖的地理位置(經緯度)和遙感影像數據成像時間進行設置。
2.2.4 采樣點光譜信息的提取
首先對需要建模的影像數據采取3*3窗口的滑動平均處理,基于IDL-ENVI軟件,通過編程將采樣點的經緯度信息轉化為TM反射率數據的文件坐標,從滑動平均后的影像中提取出對應坐標的柵格像元值,得到TM數據各波段對應的采樣點處的光譜信息。
2.3.1 APPEL模型法
APPEL模型是El-Alem等針對Modis波段通道范圍與三波段模型不完全匹配的問題而提出[17]。葉綠素a在綠光波段的光譜特征表現為高反射,根據此特征可以用來獲取最大的葉綠素信息量;有色可溶性有機物(CDOM)在藍光波段表現為強反射,因此可以使用藍光波段來消除CDOM的影響;同時根據水體懸浮物在紅光波段表現為反射的光譜特征,選用紅光波段來最小化懸浮物的影響;近紅外波段是葉綠素的敏感波段,可以通過近紅外波段來去除藍光波段的葉綠素信息、藍光波段的后向散射影響和紅光波段的葉綠素信息;因為在紅、近紅外波段,水體的光譜特征表現為強吸收,因此后向散射的影響可以忽略不計。
APPEL模型光譜指數:

以此光譜指數來構建APPEL模型的一般回歸模型:

公式(1)中bN、bB和bR分別為近紅外、藍光以及紅光波段;公式(2)中A、B為回歸模型的相關系數,Ca為反演的葉綠素濃度。
2.3.2 經驗模型法
經驗模型主要是以葉綠素a濃度和遙感參數之間的統計關系為基礎來實現對水體葉綠素a濃度的遙感反演,是一種較為廣泛的葉綠素a濃度反演模型[18]。
利用葉綠素a濃度在TM數據波段范圍內的光譜特征差異,采用不同的波段組合方式,以擴大葉綠素a吸收峰與反射峰之間的差異,從而最大化地獲取葉綠素a濃度信息量,用于提取葉綠素a濃度。
本研究通過對TM數據多個波段組合的反射率光譜值與博斯騰湖葉綠a濃度實測數據進行Pearson相關性分析,選擇相關系數較高的波段組合方式,將最佳波段組合方式計算出的反射率光譜值與博斯騰湖葉綠素a濃度實測數據進行一元線性回歸分析,最終得出反演葉綠素a濃度的經驗模型。
2.3.3 對比分析
根據經驗模型和APPEL模型反演出的葉綠素濃度模擬值分別與實測數據計算,得出模型反演精度的評價指標均方根誤差(RMSE),綜合兩種模型的各方面信息,選擇最佳模型作為博斯騰湖水體葉綠素a濃度含量反演模型。
3.1.1 經驗模型反演葉綠素a濃度結果
為提升經驗模型反演精度,全面研究TM數據波段光譜值與葉綠素a濃度的關系,本研究參考國內外學者對基于TM影像的葉綠素a濃度經驗模型反演的研究成果,選取TM影像中的TM1-TM4單波段及其各波段的線性組合,選擇采樣點數據中的23個采樣點數據(2010年6月10日的P1-P11點、2010年8月29日的P1-P12點)作為反演數據建立回歸模型,利用SPSS軟件分析其光譜值與博斯騰湖采樣點實測數據之間的PEARSON相關系數,統計結果見表2。

表2 波段光譜值與葉綠素a濃度之間的相關系數
由表2可以看出,單波段中相關系數最高的兩個波段分別為b3和b4,達到了0.717、0.715,波段組合中相關系數最高的三個組合分別為b3*b4、b3+b4和b2+b3,其相關系數分別為0.81、0.76和0.718。可以看出多波段組合的相關系數高于單波段,選擇b3*b4、b3+b4和b2+b3這三個波段組合方式作為因子分別與實測數據通過SPSS軟件做線性回歸分析,建立3個葉綠素a濃度的反演模型如下:
模型一:

模型二:


模型三:三個公式中Y代表葉綠素a濃度,波段組合代表相應波段組合的真實地表反射率值。
3.1.2 APPEL模型反演葉綠素a濃度結果
水體中的葉綠素a濃度、懸浮物、可溶性有機物(CDOM)以及后向散射對水體的光譜特征有著較大影響。APPEL模型是通過分析上述因子對水體光譜特征反射、吸收的相互關系得出。APPEL模型最早是由ELAlem等提出,原模型是根據葉綠素在近紅外波段表現為強反射來提取最大的葉綠素信息量,但是作者在實際實驗過程中發現,使用近紅外波段提取葉綠素最大信息量的反演模型精度并沒有使用綠光波段的反演模型精度高,其相關系數(表3)相差較大。

表3 APPEL模型光譜值與葉綠素a濃度之間的相關系數
研究表明,當水體中的葉綠素a濃度低于100ug/L時水體反射率紅光波段高于近紅外波段,這種現象隨著葉綠素a濃度的降低而更加明顯。博斯騰湖中的葉綠素a濃度較低,綠光波段的反射率大于紅光波段和近紅外波段,所以用綠光波段能提取到更大的葉綠素信息,因此本研究根據葉綠素在綠光波段強反射的特征來提取葉綠素的最大信息量。
綠光波段APPEL模型光譜值:

公式(6)中,模型中bG、bN、bB和bR分別為綠光、近紅外、藍光以及紅光波段的反射率。根據綠光波段APPEL模型的光譜值與實測數據建立線性回歸模型如下:
APPEL模型:

公式(7)中,X為APPEL模型的光譜指數Fa,Yappel為博斯騰湖的葉綠素a濃度。
3.1.3 最佳反演模型選擇
本研究因為2010年6月10日P14點位處于陸地上,故剔除了2010年6月10日P14這一個異常點,所以用剩下的10個采樣點數據進行模型的精度評價。將10個點的光譜值代入上述四個模型中得到葉綠素a濃度的模擬值,將實測值與模擬值進行比較(如圖3)。通過圖3可以看出,APPEL模型驗證的擬合效果較好。

圖3 葉綠素a濃度反演模型與檢驗
四個模型中,三個經驗模型的方程線性擬合度較高,從表4中可以看出b3*b4波段組合模型的R2最大為0.6638,可是其均方根誤差較其他模型最大為3.1;估測值與實測值之間相關系數最大的兩個分別為APPEL模型和(b2+b3)模型,分別為0.8887和0.8889,幾乎相同;但是APPEL模型的均方根誤差最低為0.8883,(b2+b3)模型的均方根誤差為1.4514;且APPEL模型原理是基于水體中葉綠素a、有色可溶性有機物(CDOM)、懸浮物等光譜特征的半經驗半分析模型,模型外推適宜性較強,因此綜合考慮選擇APPEL模型作為反演葉綠素a濃度的最佳模型。

表4 模型擬合效果對比
在可見光到近紅外波段范圍內,水體的吸收強度逐漸增大,在近紅外波段的吸收強度非常大,反射率極低,而植被在近紅外和綠光波段的反射強度非常大,所以本研究采用綠光波段與近紅外波段的比值來最大化提取水體信息。
通過ENIV軟件,導入大氣校正后的反射率數據,反演歸一化水體指數。歸一化水體指數(NDWI)的反演公式為(8):
公式(8)中B2、B4分別為TM影像的綠光波段和近紅外波段,根據水體條件指數設立閾值,其中閾值大于0則為水體。利用Arcgis的按屬性提取模塊提取出NDWI值大于0的像元作為博斯騰湖水體的范圍,提取結果如圖4所示。

圖4 博斯騰湖水域范圍示意圖
利用IDL編程建立APPEL數學模型,并根據APPEL模型計算葉綠素a濃度。在2010年6月-2010年10月選取博斯騰湖天空無云的TM遙感影像,反演其葉綠素a濃度,分析葉綠素a濃度的時空變化情況。反演2010年6月10日、2010年8月29日、2010年9月5日、2010年9月30日和2010年10月16日的博斯騰湖水體葉綠素濃度分布圖(圖5)。

圖5 2010年6月、8月、9月、10月葉綠素a濃度分布圖
葉綠素是藍藻、水華的主要組成成分,從葉綠素a濃度反演結果分布圖(如圖6)上可以看出,葉綠素a濃度的空間分布受到風向、水流方向的影響,而呈現出條紋狀,所有時間段內葉綠素a在湖區靠近岸邊的地方分布較大,湖中心的葉綠素a分布較少。同時表現出6月份至8月份隨著時間增長,北部靠近岸邊的區域葉綠素a濃度分布逐漸減少,而南部和東部靠近湖岸區域葉綠素a濃度分布逐漸增加,這種現象的可能原因是博斯騰湖位于焉耆盆地,地形總體趨勢是北高南低,由于水體的流動,葉綠素a濃度分布自北向南隨著時間逐漸增加。9月初博斯騰湖水體中葉綠素a濃度達到2010年6月至10月的頂點。到了10月初湖岸區域葉綠素分布逐漸消失,湖區葉綠素含量逐漸減少。該反演結果與博斯騰湖對應時間的葉綠素a濃度實測數據濃度分布趨勢基本相同,能夠反映出博斯騰湖葉綠素a在同一時間下的空間分布狀況以及不同時間下的時空分布狀況。

圖6 博斯騰湖不同時期葉綠素a濃度平均值
本研究基于葉綠素a濃度實測數據與同時期多期TM影像數據,分別通過經驗模型和APPEL模型建立反演值與實測值關系模型,經過精度評價以及對比分析,選擇最優模型用于博斯騰湖TM影像葉綠素a濃度的反演,得出如下結論:
(1)TM影像單波段與實測葉綠素a濃度相關性較高的為b3、b4波段,為0.72左右,但b2、b3、b4的多波段組合與實測葉綠素a濃度相關性高于單波段,尤其是b3、b4的乘法組合,達到0.82左右。
(2)在本研究中,采用APPEL模型在綠光波段反射率大于紅光波段和近紅外波段,因此采用綠光波段APPEL模型的光譜值建立與實測數據線性回歸模型。
(3)通過對比三個多波段組合經驗模型與一個APPEL模型發現,四個模型估測值與實測值之間相關系數均較高,在0.85-0.89之間,多波段組合經驗模型(b3*b4)R2(0.6638)最高,但均方根誤差(3.100)也最大,APPEL模型均方根誤差最小為0.8883,同時估測值與實測值之間相關系數也較高(0.8887),同時APPEL模型外推適宜性較強,因此綜合考慮選擇APPEL模型作為反演葉綠素a濃度的最佳模型。
(4)湖岸區域的葉綠素濃度較高,湖心區域較低;在6月至10月葉綠素由湖區北部向湖區東南部轉移,并伴隨著時間葉綠素濃度逐漸減小。