基于隨機(jī)森林算法與高分觀測的太湖葉綠素a濃度估算模型*

杭 鑫 曹 云 杭蓉蓉 徐 萌 霍 焱 孫良宵 朱士華

1 江蘇省氣候中心，南京 210019

2 中國人民解放軍61175部隊(duì)，南京 210049

3 上海翔竑信息技術(shù)有限公司，上海 202172

4 江蘇省常州市氣象局，常州 213000

提 要： 基于2018年和2019年有效觀測的高分1號(GF-1)衛(wèi)星影像和湖面原位觀測的葉綠素a濃度數(shù)據(jù)，利用隨機(jī)森林算法定量評價(jià)特征變量重要性的功能，選擇有效特征波段組合，建立了基于原位自動(dòng)監(jiān)測數(shù)據(jù)的太湖葉綠素a濃度的遙感反演模型。結(jié)果表明：綠光波段(0.52～0.59 μm)和紅光波段(0.63～0.69 μm)是遙感估算葉綠素a濃度的關(guān)鍵波段，與其他波段組合可以定量估算葉綠素a濃度；分季節(jié)構(gòu)建太湖葉綠素a濃度估算模型效果更好，春、夏、秋、冬各季模型的決定系數(shù)R2分別達(dá)0.84、0.85、0.96和0.82；太湖葉綠素a濃度夏季最高，秋、春季次之，冬季最低；春、秋和夏季葉綠素a濃度空間變化較明顯，而冬季變化不明顯，葉綠素a濃度高值區(qū)主要集中在西部沿岸區(qū)、竺山湖、梅梁湖和部分湖心區(qū)。研究表明：隨機(jī)森林模型可以客觀確定遙感反演葉綠素a濃度的有效波段，實(shí)現(xiàn)大面積內(nèi)陸水體葉綠素a濃度的高精度估算。

引 言

葉綠素a是浮游植物或藻類植物中最豐富的色素，準(zhǔn)確定量估算葉綠素a濃度對于客觀評價(jià)水體富營養(yǎng)化程度、實(shí)施水環(huán)境治理和生態(tài)保護(hù)具有重要意義(朱廣偉等，2018；尹艷娥等，2014)。太湖是一個(gè)典型的內(nèi)陸Ⅱ類水體，水質(zhì)狀況存在明顯的空間異質(zhì)性(水利部太湖流域管理局，2019)。傳統(tǒng)的葉綠素a濃度測量方法主要是人工采樣監(jiān)測，在實(shí)驗(yàn)室采用分光光度法進(jìn)行分析，最快也要近一天時(shí)間才能獲取結(jié)果，無法滿足政府防控需求。2008年以來采用了浮標(biāo)站自動(dòng)監(jiān)測水質(zhì)參數(shù)，這是一種基于熒光法的半定量分析方法，可在半小時(shí)內(nèi)獲得結(jié)果，其精密性和時(shí)效性都較高(劉苑等，2010)。但傳統(tǒng)方法缺乏時(shí)空連續(xù)性，觀測站點(diǎn)有限，觀測數(shù)據(jù)可能有缺失，難以準(zhǔn)確描述復(fù)雜水體的葉綠素a濃度與光譜特征之間的關(guān)系(王橋等，2008)。衛(wèi)星遙感技術(shù)具有監(jiān)測范圍廣、時(shí)效快和連續(xù)動(dòng)態(tài)等優(yōu)勢(張鵬等，2018)，越來越多地應(yīng)用于大面積水體葉綠素a濃度估算(He et al, 2020；趙少華等，2019)。遙感反演葉綠素a濃度的方法大體可分為經(jīng)驗(yàn)統(tǒng)計(jì)方法(Zhang et al, 2014)、半經(jīng)驗(yàn)半分析方法和分析模型(張明慧等，2018)。經(jīng)驗(yàn)方法主要通過建立水體反射率與實(shí)測葉綠素a濃度之間的定量關(guān)系來估算葉綠素a濃度，半經(jīng)驗(yàn)半分析方法則是在經(jīng)驗(yàn)方程的基礎(chǔ)上結(jié)合輻射傳輸模型，有一定的物理基礎(chǔ)。由于葉綠素a光學(xué)特征復(fù)雜和圖像的大氣校正不精確，光譜特征與組分濃度之間的關(guān)系較為復(fù)雜，經(jīng)驗(yàn)、半經(jīng)驗(yàn)?zāi)Ｐ洼^難解決復(fù)雜的非線性問題，反演精度可能較差(張玉超等, 2009b)，而分析模型雖具有較高的反演精度和較好的通用性(韓留生等，2014)，但機(jī)理復(fù)雜，建模難度大，實(shí)用性較差(潘應(yīng)陽等，2017)。總體而言，精確遙感反演渾濁水體葉綠素a濃度仍是當(dāng)前較為困難的任務(wù)(Nazeer and Nichol, 2016)。

隨著人工智能技術(shù)的發(fā)展，機(jī)器學(xué)習(xí)算法開始在水質(zhì)參數(shù)遙感反演中得到應(yīng)用。由于機(jī)器學(xué)習(xí)方法不依賴于固定的模型框架，而是通過不斷地“學(xué)習(xí)”模型校正過程中反饋的誤差，來完善自變量與因變量之間的復(fù)雜關(guān)系，因而是解決非線性回歸問題的有效方法(Lary et al, 2016；孫全德等，2019)。已有研究證實(shí)BP神經(jīng)網(wǎng)絡(luò)(BP neural network)模型在遙感反演太湖葉綠素a濃度時(shí)是可行的(朱云芳等，2017)，但神經(jīng)網(wǎng)絡(luò)模型需要確定網(wǎng)絡(luò)結(jié)構(gòu)，參數(shù)較多，且隨著訓(xùn)練樣本的增加，所構(gòu)建的最優(yōu)網(wǎng)絡(luò)可能是局部最優(yōu)，會(huì)出現(xiàn)“過學(xué)習(xí)”的現(xiàn)象。支持向量機(jī)(support vector machine,SVM)模型在反演葉綠素a 濃度時(shí)也獲取了較高的精度(Kong et al, 2017)，且與BP模型相比，SVM模型的反演精度可能更高，穩(wěn)定性、魯棒性和空間泛化能力也可能更好，但對于大規(guī)模訓(xùn)練樣本的SVM模型可能會(huì)耗費(fèi)大量的機(jī)器內(nèi)存和運(yùn)行時(shí)間，導(dǎo)致學(xué)習(xí)效率降低(張玉超等，2009a)。Breiman提出的隨機(jī)森林模型是一種組合分類智能算法，具有極高的準(zhǔn)確率、極強(qiáng)的數(shù)據(jù)挖掘能力及分析復(fù)雜相互作用分類特征的能力等很多優(yōu)點(diǎn)，并且可以給出變量重要性估計(jì)，甚至被譽(yù)為當(dāng)前最好的機(jī)器學(xué)習(xí)算法之一(Breiman, 2001；李文娟等，2018；劉揚(yáng)和王維國，2020；王可心等，2021)。近年來已有一些學(xué)者嘗試?yán)秒S機(jī)森林機(jī)器學(xué)習(xí)算法開展近岸水色遙感監(jiān)測，但在渾濁水體定量估算葉綠素a濃度的研究仍相對較少，其中葉綠素a 濃度實(shí)測數(shù)據(jù)也大多采用人工定點(diǎn)采樣實(shí)驗(yàn)室分析而來，自動(dòng)監(jiān)測數(shù)據(jù)的應(yīng)用較為少見，而水質(zhì)參數(shù)自動(dòng)監(jiān)測應(yīng)是今后的必然趨勢，因此，有必要充分利用好水質(zhì)自動(dòng)觀測數(shù)據(jù)，提高葉綠素a 濃度定量估算精度。

本文利用2018年和2019年有效觀測的GF-1影像和同步自動(dòng)觀測的葉綠素a濃度數(shù)據(jù)，采用隨機(jī)森林機(jī)器學(xué)習(xí)算法，選擇有效特征波段組合，建立太湖葉綠素a濃度的遙感估算模型，了解太湖水體葉綠素a濃度時(shí)空分布特征，試圖為藍(lán)藻水華發(fā)生發(fā)展預(yù)測預(yù)警提供重要參數(shù)，為藍(lán)藻水華防控提供科學(xué)依據(jù)。

1 數(shù)據(jù)和方法

1.1 研究區(qū)域

選取太湖作為研究區(qū)，其地理范圍在30°55′～31°30′N、119°55′～120°40′E；湖面面積約為2 338 km2，平均水深約為2 m。

1.2 遙感數(shù)據(jù)

選取GF-1衛(wèi)星搭載的Wide-Field-of-View(WFV)傳感器觀測得到的衛(wèi)星影像作為遙感數(shù)據(jù)源，所有數(shù)據(jù)均來自于高分辨率對地觀測系統(tǒng)江蘇數(shù)據(jù)與應(yīng)用中心，時(shí)間范圍為2018年1月至2019年5月。選擇有藍(lán)藻水華發(fā)生、晴朗或少云天氣時(shí)質(zhì)量較高的影像，最后共獲得18天27景影像，日期分別為，2018年：1月12日，2月5日，2月6日，2月13日，2月23日，4月8日，4月28日，5月15日，5月23日，6月25日，7月19日，7月20日，10月6日，10月27日，12月18日；2019年：1月17日，1月24日，5月3日。

由于原始影像為1級(相對輻射校正產(chǎn)品)，為準(zhǔn)確獲得水體表面反射率，需對數(shù)據(jù)進(jìn)行正射校正、輻射定標(biāo)、大氣校正和影像鑲嵌等預(yù)處理。

1.3 葉綠素a濃度實(shí)測數(shù)據(jù)

葉綠素a濃度數(shù)據(jù)取自江蘇省生態(tài)環(huán)境廳布設(shè)在太湖的19個(gè)水質(zhì)浮標(biāo)站自動(dòng)監(jiān)測數(shù)據(jù)，站點(diǎn)分布如圖1。由于衛(wèi)星影像過境時(shí)間為每日11:30(北京時(shí)，下同)左右，因此實(shí)測數(shù)據(jù)選取與衛(wèi)星觀測相同日期每日11:30的瞬時(shí)觀測值。

圖1 水質(zhì)浮標(biāo)站點(diǎn)分布Fig.1 Distribution of water quality buoy sites

1.4 隨機(jī)森林基本原理

隨機(jī)森林算法基本思路是：(1)首先確定原始總樣本集D和變量個(gè)數(shù)M；(2)基于原始訓(xùn)練樣本集D，采用重采樣技術(shù)從中抽取N個(gè)決策樹數(shù)目Ntree(number of trees)與D中樣本數(shù)量相同的子訓(xùn)練樣本集D1，D2，D3，…，Dn，分別建立Ntree個(gè)回歸樹模型，未抽取的部分組成袋外數(shù)據(jù)(OOB)作為測試樣本；(3)確定樹節(jié)點(diǎn)預(yù)選變量個(gè)數(shù)Mtry(number of variable per level)的數(shù)值，Mtry代表在決策樹節(jié)點(diǎn)做決定時(shí)所用變量數(shù)，一般Mtry須小于等于原始數(shù)據(jù)變量個(gè)數(shù)M；(4)針對每個(gè)訓(xùn)練集生長一棵分類回歸樹，按照節(jié)點(diǎn)不純度最小原則在樹的每個(gè)節(jié)點(diǎn)處，依據(jù)法則在Mtry個(gè)特征變量中選擇高分類能力的特征進(jìn)行分支生長，并且每棵樹都不做任何裁剪，任其最大限度的生長；(5)重復(fù)步驟(4)n次，得到由n棵決策樹組成的隨機(jī)森林；(6)隨機(jī)森林的最終回歸結(jié)果為每棵樹結(jié)果的平均值，預(yù)測精度則以每棵回歸樹的平均OOB誤差來確定。

隨機(jī)森林模型精度取決于Ntree和Mtry。Ntree決定了隨機(jī)森林的總體規(guī)模，Mtry決定了單棵決策樹的生長狀況，兩者分別從宏觀和微觀層面決定隨機(jī)森林的精度。在回歸模型中，Mtry值通常為變量數(shù)的三分之一，Ntree值根據(jù)模型誤差隨決策樹數(shù)目變化情況來決定。

隨機(jī)森林模型支持定量化比較各自變量之間對于模型的重要程度，在建模時(shí)，可以從大量特征變量中選取對最終結(jié)果影響較大的數(shù)目較少的特征變量，刪除一些和任務(wù)無關(guān)或者冗余的特征變量，簡化的特征數(shù)據(jù)集也常常會(huì)得到更精確的模型。變量重要性度量的主要評價(jià)指標(biāo)為精度平均減少值IncMSE (increasing the mean square error)和節(jié)點(diǎn)不純度減少值IncNodePurity (increasing the node impurity)，值越大，表明該變量越重要，反之則相對不重要。

2 結(jié)果與分析

2.1 樣本處理與潛在變量分析

從所有18天GF-1衛(wèi)星過境時(shí)刻的實(shí)測葉綠素a濃度數(shù)據(jù)中，剔除缺測、異常值及受云影響的152個(gè)數(shù)據(jù)，最終得到190個(gè)數(shù)據(jù)組成了實(shí)測葉綠素a濃度樣本數(shù)據(jù)集。再根據(jù)不同季節(jié)分組，分別得到春季(3—5月)57個(gè)，夏季(6—8月)30個(gè)，秋季(9—11月)20個(gè)和冬季(12—2月)83個(gè)實(shí)測葉綠素a濃度的樣本子集。從各樣本子集中隨機(jī)選取四分之三的數(shù)據(jù)作為訓(xùn)練樣本集，剩下四分之一數(shù)據(jù)作為測試樣本集。

基于隨機(jī)森林模型的葉綠素a濃度遙感估算模型是以GF-1 WFV影像為主要數(shù)據(jù)源。采用GF-1 WFV 4個(gè)波段的反射率進(jìn)行水質(zhì)參數(shù)反演有其合理性，但單一波段往往不能很好地反映影響因子與水質(zhì)參數(shù)之間的關(guān)系(朱云芳等，2017)。參考方馨蕊等(2019)，將GF-1 WFV 4個(gè)波段及不同組合的39個(gè)變量因子作為潛在變量進(jìn)行篩選(表1)。

表1 參與隨機(jī)森林建模的39個(gè)潛在變量因子Table 1 39 latent variable factors involved in random forest modeling

2.2 關(guān)鍵波段變量因子確定

首先，確定隨機(jī)森林回歸模型中最重要的兩個(gè)輸入?yún)?shù)Mtry和Ntree。Mtry取潛在變量總數(shù)39的三分之一，即Mtry=13，而Ntree值根據(jù)誤差分析來確定。由于參數(shù)Mtry值固定不變，那么參數(shù)Ntree值越大，則誤差Error越小或趨于穩(wěn)定，代表模型精度越高。為考察不同季節(jié)葉綠素a濃度反演的情況，分別建立了全部樣本模型(MODAll)、春季模型(MODSpr)、夏季模型(MODSum)、秋季模型(MODAut)和冬季模型(MODWin)共5個(gè)模型。5個(gè)模型在不同的Ntree和Mtry組合下所對應(yīng)的分類精度如圖2所示。由圖看出，MODAll和MODSpr的誤差在Ntree為600時(shí)趨于穩(wěn)定，MODSum的誤差在Ntree為500時(shí)趨于穩(wěn)定，MODAut的誤差在Ntree為200時(shí)趨于穩(wěn)定，MODWin的誤差在Ntree為400時(shí)趨于穩(wěn)定。綜上所述，我們確定在MODAll和MODSpr中，Ntree=600；在MODSum中，Ntree=500；在MODAut中，Ntree=200；在MODWin中，Ntree=400。

圖2 決策樹數(shù)目Ntree與模型誤差的關(guān)系(a)MODAll,(b)MODSpr,(c)MODSum,(d)MODAut,(e)MODWinFig.2 The relationship between decision tree Ntree and model error(a) MODAll, (b) MODSpr, (c) MODSum, (d) MODAut, (e) MODWin

然后，分別根據(jù)單一指標(biāo)IncMSE和IncNodePurity篩選波段因子變量。隨機(jī)森林回歸模型通過評估每個(gè)變量對總體模型預(yù)測精度提高的能力，對這些變量進(jìn)行排序，評估各變量在模型中的相對重要性。利用上述方法確定的參數(shù)Mtry和Ntree分別進(jìn)行建模和模型優(yōu)化，MODAll、MODSpr、MODSum、MODAut和MODWin模型各單獨(dú)訓(xùn)練5 000次，從中各選擇精度相對較高(相關(guān)系數(shù)R>0.85)的20個(gè)模型，對應(yīng)每個(gè)變量都可以得到20個(gè)IncMSE和IncNodePurity值。分別計(jì)算每個(gè)變量的IncMSE和IncNodePurity指標(biāo)的平均值，然后根據(jù)這兩個(gè)指標(biāo)平均值的大小分別進(jìn)行排序。將指標(biāo)平均值排名相對靠前且重要性度量曲線出現(xiàn)較明顯拐點(diǎn)前的特征變量認(rèn)為是相對重要的變量，由此得到全部樣本模型MODAll、春季模型MODSpr、夏季模型MODSum、秋季模型MODAut和冬季模型MODWin共5個(gè)模型的單指標(biāo)重要變量篩選結(jié)果：

(1) 基于IncMSE指標(biāo)的篩選結(jié)果

MODAll：DVI(2, 3)，VI(2,1,3)，RVI(2, 3)，NDVI(2, 3)，VI(2,1,3,4)，VI(2,3,4)，DVI(2, 4)，VI(2,1,4)，B1，VI(3,1,2)；

MODSpr：VI(2,1,4)，DVI(2, 3)，VI(2,1,3)，VI(2,1,3,4)，DVI(1, 2)；

MODSum：VI(3,1,2)，RVI(1, 3)，DVI(2, 3)，RVI(2, 3)，DVI(1, 3)；

MODAut：VI(2,3,4)，RVI(1, 2)，VI(1,2,3,4)，VI(1,2,3)，VI(4,1,2)，NDVI(2, 4)，VI(1,2,4)，VI(1,3,4)，DVI(1, 2)；

MODWin：RVI(2, 3)，RVI(2, 4)，VI(2,1,3,4)，VI(2,3,4)。

(2) 基于IncNodePurity指標(biāo)的篩選結(jié)果

MODAll：VI(2,1,3,4)，VI(2,3,4)，VI(2,1,3)，B4,VI(2,1,4)，DVI(2, 3)；

MODSpr：VI(4,1,3)，NDVI(1, 4)，VI(4,1,2,3)，EVI；

MODSum：RVI(1, 3)，VI(3,1,2)，VI(2,3,4)，VI(1,3,4)；

MODAut：VI(2,1,3)，RVI(2, 3)，DVI(1, 2)，VI(4,1,2,3)，VI(4,1,2)，DVI(2, 3)，VI(4,2,3)，B4,NDVI(2, 3)，EVI；

MODWin：RVI(2, 3)，VI(2,3,4)，VI(3,1,2)，RVI(2, 4)。

從上述兩類篩選結(jié)果中，發(fā)現(xiàn)同時(shí)包含B2和B3波段變量因子出現(xiàn)次數(shù)較多，在IncMSE指標(biāo)的篩選結(jié)果中出現(xiàn)19次，占比為58%，在IncNodePurity指標(biāo)的篩選結(jié)果中出現(xiàn)16次，占比為57%；而包含B2或B3波段變量因子出現(xiàn)次數(shù)更多，在IncMSE指標(biāo)篩選結(jié)果中共出現(xiàn)32次，占比為97%，僅1個(gè)變量因子不包含B2或B3波段，在IncNodePurity指標(biāo)篩選結(jié)果中累計(jì)出現(xiàn)25次，占比為89%，僅3個(gè)變量因子不包含B2或B3波段，說明B2和B3波段在所有變量因子中占主導(dǎo)地位。

最后，根據(jù)結(jié)合IncMSE和IncNodePurity的綜合因子篩選關(guān)鍵波段因子。已有文獻(xiàn)在篩選特征變量因子時(shí)，大多選擇IncMSE和IncNodePurity中的一種指標(biāo)(羅曉春等，2019)。但從以上兩種指標(biāo)的篩選結(jié)果中可以發(fā)現(xiàn)，同一種模型不同指標(biāo)篩選出的變量因子并不完全相同。相對應(yīng)于這兩個(gè)指標(biāo)，所有模型經(jīng)過篩選得到兩組不同的變量因子。對于MODAll，單用指標(biāo)IncMSE篩選得到10個(gè)變量，用指標(biāo)IncNodePurity得到6個(gè)變量，兩組變量中出現(xiàn)5個(gè)相同因子。相類似的，MODWin有3個(gè)相同因子，MODSum中有2個(gè)相同因子，MODAut中僅有1個(gè)相同因子。我們注意到模型MODSpr的兩組篩選結(jié)果中甚至沒有出現(xiàn)相同的變量因子。表明僅用IncMSE或IncNodePurity一個(gè)指標(biāo)并不能完全反映變量的重要性，存在一定的局限性。為此，考慮綜合IncMSE和IncNodePurity兩種指標(biāo)構(gòu)建一組新的變量相對重要性評價(jià)指標(biāo)(relative importance evaluation index，RIEI)，具體計(jì)算方法是：對于MODAll、MODSpr、MODSum、MODAut和MODWin這5個(gè)模型中的每個(gè)模型，首先分別將以上20個(gè)模型的IncMSE值和IncNodePurity值進(jìn)行歸一化處理，然后再求平均值，得到每一個(gè)變量的RIEI值：

(1)

式中:IncMSEi為第i個(gè)IncMSE值，IncMSEmin為20個(gè)IncMSE中的最小值，IncMSEmax為20個(gè)IncMSE中的最大值；IncNodePurityi為第i個(gè)IncNodePurity值，IncNodePuritymin為20個(gè)IncNodePurity中的最小值，IncNodePuritymax為20個(gè)IncNodePurity中的最大值。

將包含全部樣本和春、夏、秋和冬季共5個(gè)模型的RIEI值繪制成變量重要性評估曲線(圖3)。

圖3 基于RIEI值的5個(gè)模型的變量重要性評估曲線(a)MODAll,(b)MODSpr,(c)MODSum,(d)MODAut,(e)MODWinFig.3 Variable importance evaluation curve of five models based on RIEI value(a) MODAll, (b) MODSpr, (c) MODSum, (d) MODAut, (e) MODWin

將變量重要性排名相對靠前且曲線出現(xiàn)較明顯拐點(diǎn)前的特征變量認(rèn)為是相對重要的變量，分別篩選出5個(gè)模型的重要特征變量：

MODAll：VI(2,1,3,4)，VI(2,3,4)，VI(2,1,3)，B4，VI(2,1,4)，DVI(2, 3)，NDVI(2, 3);

MODSpr：DVI(1, 2)，VI(2,1,4)，VI(2,1,3,4)，RVI(1, 2)，VI(2,1,3)，DVI(2, 3)，VI(1,2,3);

MODSum：RVI(1, 3)，VI(3,1,2)，RVI(2, 3)，VI(2,3,4)，VI(1,3,4)、DVI(2, 3)，DVI(1, 3)，VI(1,2,3,4)，VI(2,1,3,4)，VI(1,2,3)，VI(2,1,3);

MODAut：VI(4,1,2)，DVI(1, 2)，VI(2,3,4)，VI(4,1,2,3)，VI(1,2,3,4);

MODWin：RVI(2, 3)，VI(2,3,4)，RVI(2, 4)，VI(2,1,3,4)，DVI(2, 3)，VI(3,1,2)。

在以上模型的變量中，同時(shí)包含B2和B3波段的變量因子出現(xiàn)25次，占比為69%，明顯高于IncMSE指標(biāo)的58%和IncNodePurity指標(biāo)的57%，表明使用綜合指標(biāo)篩選的結(jié)果明顯優(yōu)于單一指標(biāo)。同時(shí)，所有的變量因子都包含了B2或B3波段，再次證明綠光波段B2(0.52～0.59 μm)和紅光波段B3(0.63～0.69 μm)是遙感反演葉綠素a濃度的關(guān)鍵波段因子，對準(zhǔn)確估算渾濁水體葉綠素a濃度具有重要的意義。

2.3 模型建立

根據(jù)上述5個(gè)模型篩選的重要特征變量，分別重新構(gòu)建隨機(jī)森林模型，其中參數(shù)Mtry為特征變量個(gè)數(shù)的三分之一，分別取1、2、3和4四個(gè)數(shù)值，Ntree則根據(jù)前面誤差分析結(jié)果分別選取400、500和600三個(gè)數(shù)值。對應(yīng)每組參數(shù)組合(Mtry，Ntree)重復(fù)建模5 000次，從中選出各模型精度最高的參數(shù)組合(表2)。由表2可知，在MODAll、MODSpr、MODSum、MODAut和MODWin共5個(gè)模型中，秋季模型MODAut精度(R)最高達(dá)0.99，對應(yīng)的參數(shù)組合為(2，400)，包含全部樣本的模型MODAll和冬季模型MODWin精度(R)最小，均為0.84，對應(yīng)的參數(shù)組合分別為(3，400)和(2，600)。

表2 各模型精度Table 2 Each model accuracy

2.4 模型驗(yàn)證

進(jìn)一步驗(yàn)證模型的反演精度，將上述建立的5個(gè)隨機(jī)森林模型MODAll、MODSpr、MODSum、MODAut和MODWin反演的葉綠素a濃度值，分別與實(shí)測葉綠素a濃度值進(jìn)行比較，各模型葉綠素a濃度估算值和實(shí)測值擬合關(guān)系見圖4。由圖可知，各模型估算值與實(shí)測值之間均呈現(xiàn)較高的相關(guān)性，均方根誤差(RMSE)均較低。其中MODAut模型的擬合精度最高，決定系數(shù)(R2)為0.96，RMSE為2.1 mg·m-3，MODSum模型次之，R2為0.85，RMSE為2.0 mg·m-3，MODAll模型的擬合精度最低，R2為0.77，RMSE為2.2 mg·m-3。這一結(jié)果表明用所有樣本構(gòu)建的模型的擬合效果不如分季節(jié)構(gòu)建的模型，說明分季節(jié)模型估算的葉綠素a濃度值更加接近實(shí)測值，其中秋季模型擬合效果又明顯好于其他3個(gè)季節(jié)模型。

圖4 各模型葉綠素a濃度估算值和實(shí)測值之間的散點(diǎn)關(guān)系(a)MODAll,(b)MODSpr,(c)MODSum,(d)MODAut,(e)MODWinFig.4 Scatter plots between estimated and measured chlorophyll a concentration of each model(a) MODAll, (b) MODSpr, (c) MODSum, (d) MODAut, (e) MODWin

3 討 論

有效波段的選擇是高精度估算葉綠素a濃度的關(guān)鍵(姜廣甲等，2013)。利用隨機(jī)森林模型識別、量化特征變量重要性的功能，客觀篩選出相對重要的變量因子，本文確定了綠光波段(0.52～0.59 μm)和紅光波段(0.63～0.69 μm)及其組合為定量反演太湖葉綠素a濃度的關(guān)鍵波段。盡管目前對于渾濁Ⅱ類水體葉綠素a的光譜特征和敏感波段的研究還較少(潘應(yīng)陽等，2017)，但已有一些學(xué)者針對不同內(nèi)陸水體的實(shí)測光譜數(shù)據(jù)進(jìn)行了分析，得到了類似的結(jié)果。有統(tǒng)計(jì)表明水體葉綠素a濃度與0.54 μm和0.701 μm反射峰的相關(guān)系數(shù)接近1，相近波段的特征光譜能較好反演葉綠素a濃度(吳傳慶等，2009)；水體中藻類最顯著的光譜特征是0.56 μm附近的反射峰，該峰值的存在與否通常被認(rèn)為是判斷水體是否含有藻類的依據(jù)，而0.682 μm波段與葉綠素a濃度相關(guān)性最好(楊婷等，2011)；此外，葉綠素a濃度的變化也會(huì)影響浮游植物吸收峰的數(shù)值和位置，從而影響到最佳波段的選擇(潘應(yīng)陽等，2017)。這些大多是針對特定水體特定區(qū)域、選用較少的測量數(shù)據(jù)的研究結(jié)果可能并不具有普適性(馮馳等，2015)。對于較為渾濁的內(nèi)陸水體，由于存在浮游植物、懸浮物、溶解有機(jī)物等許多影響葉綠素a吸收的物質(zhì)，各組分之間彼此混合、交互作用，水體的光譜特征更加復(fù)雜，實(shí)際測量的葉綠素a反射吸收峰也會(huì)有明顯不同(羅建美等，2017)。因此，利用隨機(jī)森林模型確定葉綠素a光譜特征的關(guān)鍵波段，可以避免特定水域、特定葉綠素a濃度測量的局限性，對于定量遙感反演大面積渾濁Ⅱ類水體葉綠素a濃度不失為一次有益的嘗試。

太湖葉綠素a濃度具有明顯的時(shí)空分布特征(樂成峰等，2008)。本文利用春、夏、秋、冬四個(gè)季節(jié)的隨機(jī)森林模型，分析了太湖葉綠素a濃度的時(shí)空分布特點(diǎn)。圖5為2018年全湖葉綠素a濃度均值隨時(shí)間的變化情況，可見夏季平均葉綠素a濃度最高，冬季最低，二者平均濃度分別為9.6 mg·m-3和7.1 mg·m-3，秋季由于受到夏季高濃度的影響，葉綠素a濃度高于春季，分別為8.6 mg·m-3和7.7 mg·m-3，這與樂成峰(2008)研究較為一致。以一景影像代表各季節(jié)的空間分布情況(圖6)，可以看到冬季葉綠素a濃度較低，空間變化不明顯，這與冬季溫度降低有關(guān)(賈春燕，2008)；春季葉綠素a濃度開始出現(xiàn)較明顯的空間變化，西北部湖區(qū)，特別是梅梁湖和竺山湖附近葉綠素a濃度較高，向湖心區(qū)逐漸減小，這主要是由于存在眾多的入湖河流和高密度城市排污口，造成水體富營養(yǎng)化嚴(yán)重(劉聚濤等，2011)；而夏、秋季葉綠素a濃度空間變化最為顯著，西部沿岸區(qū)、竺山湖、梅梁湖和部分湖心區(qū)葉綠素a濃度明顯偏高，除了與富營養(yǎng)化程度有關(guān)外，還與夏、秋季盛行東南風(fēng)引起的湖流有關(guān)(秦伯強(qiáng)等，2004)。東太湖部分水域葉綠素a濃度始終呈現(xiàn)相對較高的水平，則可能是受該區(qū)域豐富的水生植物影響?？紤]到太湖各區(qū)域均存在不同程度的水生植物，且隨著季節(jié)交替，水生植物面積變化幅度較大(Zhao et al, 2013)，因此，本研究參考了相關(guān)文獻(xiàn)(楊婷等，2011；朱云芳等，2017)，沒有將水生植物剔除，可能會(huì)對葉綠素a濃度反演產(chǎn)生干擾，后續(xù)將會(huì)進(jìn)一步研究水生植物的影響。另外，本文僅采用了2018年1月至2019年5月的資料，且由于高分衛(wèi)星觀測條件限制，數(shù)據(jù)時(shí)序不連續(xù)，因此需要進(jìn)一步搜集更多數(shù)據(jù)，增加樣本數(shù)，以期構(gòu)建更為精確的模型。

圖5 2018年太湖葉綠素a濃度時(shí)間變化Fig.5 Temporal variation of chlorophyll a concentration in Taihu Lake in 2018

圖6 各季節(jié)太湖葉綠素a濃度空間分布(a)冬季，(b)春季，(c)夏季，(d)秋季Fig.6 Spatial distribution of chlorophyll-a concentration in Taihu Lake in each seasons(a) winter, (b) spring, (c) summer, (d) autumn

由于人工采樣實(shí)驗(yàn)室分析法和水質(zhì)浮標(biāo)站自動(dòng)監(jiān)測法在原理、方法和步驟等方面不同，加之太湖湖體各處的藻種、水質(zhì)存在較大差異，人工采樣實(shí)驗(yàn)室分析與自動(dòng)監(jiān)測的葉綠素a濃度之間的差異較大。據(jù)江蘇省環(huán)境監(jiān)測中心利用全湖人工觀測數(shù)據(jù)與原位觀測數(shù)據(jù)對比，表明人工采樣分析與自動(dòng)監(jiān)測數(shù)據(jù)之間的相關(guān)性較差，整體而言，人工采樣分析得到的結(jié)果約為自動(dòng)監(jiān)測的2.4倍。因此，本文利用自動(dòng)監(jiān)測數(shù)據(jù)來估算的葉綠素a濃度，比目前大多數(shù)采用人工采樣實(shí)驗(yàn)室分析數(shù)據(jù)的估算結(jié)果要明顯偏小(曹紅業(yè)等，2016；馮馳等，2015；宋挺等，2017)，但可能更客觀地反映了太湖湖體葉綠素a濃度及其分布狀況。水質(zhì)參數(shù)的自動(dòng)監(jiān)測分析應(yīng)當(dāng)是今后的必然趨勢，事實(shí)上，近年來江蘇省太湖水污染防治工作中使用的葉綠素a濃度數(shù)據(jù)就主要來自于自動(dòng)監(jiān)測結(jié)果，因此有必要在后續(xù)的研究中搞清楚自動(dòng)觀測數(shù)據(jù)和人工采樣數(shù)據(jù)之間的關(guān)系。

4 結(jié) 論

(1)用隨機(jī)森林模型可以確定綠光波段和紅光波段為遙感反演太湖葉綠素a濃度的關(guān)鍵波段，但仍需結(jié)合其他波段，這樣可以避免特定水域、特定葉綠素a濃度測量的局限性，結(jié)果更客觀。

(2)分季節(jié)構(gòu)建的太湖葉綠素a濃度隨機(jī)森林估算模型能夠得到更加接近實(shí)測值的結(jié)果，春、夏、秋、冬各季模型的決定系數(shù)R2分別為0.84、0.85、0.96和0.82，RMSE分別為1.9、2.0、2.1和1.5 mg·m-3。

(3)太湖葉綠素a濃度呈明顯時(shí)空變化特征，其中夏季葉綠素a濃度最高，秋、春季次之、冬季最低；春、秋和夏季全湖葉綠素a濃度空間變化較明顯，冬季空間變化不明顯，葉綠素a濃度高值區(qū)主要集中在西部沿岸區(qū)、竺山湖、梅梁湖和部分湖心區(qū)。

(4)模型使用實(shí)測數(shù)據(jù)來源于水質(zhì)浮標(biāo)站自動(dòng)監(jiān)測，跟人工采樣實(shí)驗(yàn)室分析法獲得結(jié)果相比明顯偏低，這是因?yàn)閮煞N監(jiān)測方式在原理、方法和步驟等方面有所不同，而且太湖湖體不同區(qū)域的藻種、水質(zhì)也存在較大差異。但水質(zhì)自動(dòng)化監(jiān)測是未來的發(fā)展方向，充分利用好自動(dòng)觀測數(shù)據(jù)具有重要的現(xiàn)實(shí)意義。