楊傳訓(xùn),李 勇,楊 驥,舒思京

(廣東省科學(xué)院廣州地理研究所,廣東 廣州 510070)

水體懸浮泥沙濃度大小決定水體透度、渾濁度和水色等光學(xué)性質(zhì),對河口自然環(huán)境、生態(tài)環(huán)境及港航通道等具有重要影響[1-4],河口海岸作為連接淡水水生系統(tǒng)和海洋水生生態(tài)系統(tǒng),是最重要的地球生態(tài)系統(tǒng)之一,對全球生物元素循環(huán)和氣候具有公認(rèn)的重要性[5-6]。

當(dāng)前監(jiān)測懸浮泥沙濃度監(jiān)測主要以現(xiàn)場實測和遙感反演兩種方式為主[7-8]。站點式監(jiān)測能夠在時間尺度上連續(xù)、準(zhǔn)確獲取懸浮泥沙含量,但站點建設(shè)和維護成本高,站點數(shù)量少、分布稀疏,無法進行大面積的空間連續(xù)監(jiān)測[8-9];遙感反演手段能夠?qū)崿F(xiàn)大范圍的空間連續(xù)監(jiān)測,能夠較好地體現(xiàn)懸浮泥沙空間連續(xù)變化[10-11]。在珠江口懸浮泥沙遙感反演模式研究上,文獻[12]在總結(jié)和對比前人線性關(guān)系、對數(shù)關(guān)系、Gordon、負(fù)指數(shù)關(guān)系基礎(chǔ)上,提出并推導(dǎo)了懸浮泥沙遙感定量的統(tǒng)一模式,并利用 TM 影像數(shù)據(jù)反演珠江口懸浮泥沙,證明了統(tǒng)一模式與其他模式相比,精度較高;文獻[13]建立了TM 輻射亮度值與懸浮泥沙濃度的多波段對數(shù)模式;文獻[14]推導(dǎo)出懸浮泥沙定量遙感綜合模式,該模式在條件簡化下包含前人研究的半經(jīng)驗?zāi)Ｐ?并在珠江口取得成功應(yīng)用。

當(dāng)前,機器學(xué)習(xí)方法在水質(zhì)遙感監(jiān)測中已被廣泛使用[15-17],本文以珠江口懸浮泥沙水質(zhì)采樣數(shù)據(jù)及高光譜數(shù)據(jù)為基礎(chǔ),構(gòu)建珠江口人工神經(jīng)網(wǎng)絡(luò)、支持向量回歸機、隨機森林懸浮泥沙反演模型,實現(xiàn)對2006年12月19日TM影像、2013年9月10日Landsat OLI影像及2020年4月29日Sentinel-2影像數(shù)的反演與分析。

1 研究區(qū)域與試驗數(shù)據(jù)

1.1 研究區(qū)域

珠江入海口為研究區(qū)域,該區(qū)域多年平均含沙量約為0.126～0.334 kg/m3,斷面最大含沙量為1.23～4.08 kg/m3。流域內(nèi)年均降水量約2000 mm,具有典型的汛期和枯水期,汛期(4—9月)降水量占全年的67.6%～83.5%。珠江口水域?qū)儆诘湫偷牟灰?guī)則半日混合潮型,大約12 h為一個平均潮汐周期,為弱潮河口,高潮位與低潮位之間的潮位差較小,屬于典型的低潮差入?？凇?/p>

1.2 試驗數(shù)據(jù)

水質(zhì)及光譜數(shù)據(jù):為研究珠江口潮汐效應(yīng)改正的珠江口懸浮泥沙遙感反演與預(yù)測分析,共在珠江口的多處地方采集了229個點的水質(zhì)參數(shù)及地面同步高光譜數(shù)據(jù)。

影像數(shù)據(jù):采用的遙感影像數(shù)據(jù)為2006年12月19日TM影像數(shù)據(jù)、2013年9月10日Landsat OLI影像數(shù)據(jù)及2020年4月29日Sentinel-2影像數(shù)據(jù)。

1.3 泥沙光譜特征分析

從229個采集點位中選擇出18條典型點位采集的光譜曲線(如圖1所示),懸浮泥沙濃度值在1.6～456.2 mg/L中逐漸增大,結(jié)果表明在400～900 nm波段范圍,水體的反射率值隨著懸浮泥沙濃度增加而增加。從實際采樣的珠江口高光譜實測數(shù)據(jù)發(fā)現(xiàn)珠江口懸浮泥沙存在兩個反射峰值,第一反射峰位置在綠光、紅光波段(570～690 nm),第二反射峰位置在近紅外波段(780～830 nm),在含沙量較低時,第一反射峰值集中在綠光波段(560 nm)且顯著高于第二反射峰值。隨含沙量濃度的增加,第一反射峰值波段集中綠光、紅光波段(550～670 nm),同時第二反射峰的反射率逐漸升高,第一反射峰值與第二反射峰值差逐步縮小。

圖1 ASD地物光譜儀實測的水體光譜數(shù)據(jù)(N=18)

2 研究方法與模型構(gòu)建

2.1 研究方法

技術(shù)路線如圖2所示。

圖2 技術(shù)路線

2.1.1 人工神經(jīng)網(wǎng)絡(luò)

使用多層人工神經(jīng)網(wǎng)絡(luò),該模型包括1個輸入層,1～4個隱藏層和1個輸出層。輸入層是診斷光譜變量,系水質(zhì)光學(xué)特征參數(shù)(OAWQ)與原始光譜波段和波段比值之間的高相關(guān)系數(shù)確定。第一隱藏層神經(jīng)元的值等于輸入?yún)?shù)乘以其連接權(quán)重加上偏差。第二層除了使用第一隱藏層神經(jīng)元的值作為輸入?yún)?shù)外,計算過程類似于第一層。輸出層為OAWQ。假設(shè)輸入與第一隱藏層的轉(zhuǎn)換關(guān)系為Ai,兩個隱藏層的關(guān)系為AAi。Ai和AAi可表示為

(1)

式中,B1和B2為激活函數(shù);b1和b2為擬合過程中兩層之間的偏差(Biasi);n為節(jié)點數(shù);j為前一層的第j個節(jié)點;i為當(dāng)前層中的第i個節(jié)點;k為兩層之間的第k次迭代計算。在第k次迭代中,是連接輸入層第j個節(jié)點和第一隱藏層第i個節(jié)點之間的權(quán)重;是兩個隱藏層之間的連接權(quán)重。在OAWQi估算中,第二隱藏層和輸出層之間的轉(zhuǎn)換關(guān)系可表示為

(2)

2.1.2 支持向量回歸機

支持向量機是一種相對簡單的監(jiān)督機器學(xué)習(xí)算法,用于解決分類或回歸問題。支持向量回歸機在長江口等其他流域區(qū)域在懸浮泥沙研究中也都表現(xiàn)出不錯的結(jié)果[18]。常用的核函數(shù)有高斯徑向基核、多項式核、sigmoid核、RBF核函數(shù)。本文采用RBF核函數(shù),公式為

(3)

式中,γ為核參數(shù)。

2.1.3 隨機森林模型

隨機森林模型是一種基于分類樹的機器學(xué)習(xí)方法,具有可輸入大量的變數(shù),分類規(guī)則清晰,結(jié)果容易理解,準(zhǔn)確度高、計算量較小,實現(xiàn)速度快等優(yōu)點,能夠用于回歸、分類、預(yù)測等分析。

隨機森林模型原理是采用bootsrap重抽樣方法從原始樣本中抽取N個樣本,對N個樣本分別建立N個決策樹模型,通過組成多棵決策樹最終通過投票得出模擬結(jié)果[19]。

2.2 模型構(gòu)建

人工神經(jīng)網(wǎng)絡(luò)(artificial neural network, ANN)回歸算法在懸浮泥沙濃度(total suspended solids, TSS)遙感反演估算結(jié)果如圖3(a)所示。結(jié)果顯示訓(xùn)練集精度和測試集精度R2均大于0.85(訓(xùn)練集R2=0.86,測試集R2=0.88),訓(xùn)練集RMSE=7.6 mg/L,MAE=4.76 mg/L,測試集RMSE=8.1 mg/L,MAE=5.54 mg/L,PRD=2.31。結(jié)合散點圖發(fā)現(xiàn),散點能較均勻分布1∶1線兩側(cè),擁有接近于1的Slope(0.88 vs 0.76),但是訓(xùn)練、測試采樣點存在部分飽和現(xiàn)象。

圖3 基于ANN、SVR、RF的TSS遙感反演估算

支持向量回歸機(support vector regression,SVR)算法在TSS遙感反演估算結(jié)果如圖3(b)所示。結(jié)果顯示測試集R2略優(yōu)于訓(xùn)練集(訓(xùn)練集R2=0.84,測試集R2=0.86)。訓(xùn)練集RMSE=5.96 mg/L,MAE=5.01 mg/L,測試集RMSE=6.53 mg/L,MAE=5.22 mg/L,PRD=2.3。

隨機森林(random forest, RF)算法在TSS遙感反演估算如圖3(c)所示。結(jié)果顯示測試集R2均大于0.90(訓(xùn)練集R2=0.91,測試集R2=0.95)。訓(xùn)練集RMSE=4.68 mg/L,MAE=3.89 mg/L,測試集RMSE=5.41 mg/L,MAE=4.18 mg/L,PRD=2.12。結(jié)合散點圖發(fā)現(xiàn)散點能較均勻分布1∶1線兩側(cè),擁有接近于1的Slope(0.90 vs 0.80)。

相對人工神經(jīng)網(wǎng)絡(luò)和支持向量回歸機方法,隨機森林模型預(yù)測方法預(yù)測精度高。3種機器學(xué)習(xí)模型精度對比見表1,與ANN相比,R2由0.88提升至0.95,RMSE降低2.69 mg/L,MAE降低1.36 mg/L。與SVR相比,R2由0.86提升到0.95,RMSE降低1.12 mg/L,MAE降低1.04 mg/L,主要原因為ANN、SVR黑箱模型內(nèi)部的運行機理復(fù)雜,對影響因子權(quán)重難以確定,加上這兩種模型對原始數(shù)據(jù)量要求高、計算量大等因素。

表1 各模型精度對比

隨機森林模型分類規(guī)則清晰,結(jié)果容易理解,準(zhǔn)確度高、計算量較小,實現(xiàn)速度快,且相關(guān)系數(shù)和預(yù)測精度明顯優(yōu)于ANN、SVR等其他模型。與ANN和SVR模型相比,RF模型在TSS水質(zhì)參數(shù)中的模擬精度均有提升。因此,隨機森林模型更適用于水體懸浮物監(jiān)測遙感估算。

結(jié)合地面采樣的229個光譜數(shù)據(jù)和水質(zhì)數(shù)據(jù),找到紅光、近紅外為懸浮泥沙敏感波段,而常用的傳感器TM5、ETM+、OLI、Sentinel-2都存在綠光與近紅外通道,可用于懸浮泥沙濃度的反演。但是多光譜傳感器波段較寬(幾十至幾百nm),直接將構(gòu)建的高光譜波段(1～20 nm)模型應(yīng)用于多光譜傳感器可能存在較大誤差,因此需要分析寬波段下懸浮泥沙反演模型的適用性及構(gòu)建。TM5、ETM+、OLI、Sentinel-2光譜響應(yīng)函數(shù)如圖4所示。綠光波段整體范圍為0.50～0.65 μm,中心波長為0.55 μm。近紅波段整體范圍為0.71～0.93 μm。TM與ETM+光譜響應(yīng)函數(shù)差異較小,Sentinel-2與OLI差異明顯。

圖4 綠光與近紅外多傳感器光譜響應(yīng)函數(shù)(TM5、ETM+、OLI、Sentinel-2)

由于ASD光譜儀采集的數(shù)據(jù)光譜分辨率較高、通道窄,而不同遙感影像每個通道的光譜分辨率不一致,為了研究多傳感器模型效果,在模型校正與驗證前,將野外實測的229條反射率采集到TM5、ETM+、OLI、Sentinel-2波段傳感器。具體計算公式為

(4)

式中,λ為波長;λmin為通道的起始波長;λmax為通道的終止波長;r(λ)為對應(yīng)λ波長的地表反射率;f(λ)為光譜響應(yīng)函數(shù)。

隨機篩選152條光譜重采樣數(shù)據(jù),基于經(jīng)驗方法分析不同傳感器通道響應(yīng)對模型的適用性,經(jīng)驗?zāi)Ｐ托ＵY(jié)果如圖5所示。上述模型擬合效果均較好(R2>0.8),說明將高光譜采集到傳感器之后近紅與綠光波段的比值能夠用于懸浮泥沙的反演,模型都滿足二次函數(shù)關(guān)系,模型形式一致,發(fā)現(xiàn)近紅與綠波段的比值能削弱外界影響,構(gòu)建模型精度較高,適用多源傳感器較穩(wěn)定。但是4種傳感器模型系數(shù)差異顯著,因此在應(yīng)用模型時,不同傳感器須獲取對應(yīng)的光譜響應(yīng)函數(shù),將實測高光譜數(shù)據(jù)重采樣到對應(yīng)傳感器通道,對模型進行重新校正。

圖5 多傳感器最優(yōu)懸浮泥沙反演模型

將重采樣剩余的77個樣本數(shù)據(jù)應(yīng)用于校正的經(jīng)驗?zāi)Ｐ偷木仍u價,評價結(jié)果如圖6所示,實測懸浮泥沙與模型計算懸浮泥沙R2均大于0.85,RMSE低于27.46,模型都具有較高的精度,其中Sentinel-2與ETM+傳感器優(yōu)于TM5與OLI傳感器。驗證散點圖斜率為0.80～0.85之間,斜率小于1。由于不同傳感的響應(yīng)函數(shù)存在差異,將高光譜重采樣到對應(yīng)傳感器波段并對模型校正到相應(yīng)的傳感器通道后4種傳感器均具有較高的精度,模型計算結(jié)果略低于實測結(jié)果。

圖6 多光譜傳感器的經(jīng)驗?zāi)Ｐ途仍u價

在機器學(xué)習(xí)中選擇模型精度較好的隨機森林模型中隨機抽選40個樣本數(shù)據(jù)進行驗證,評價結(jié)果如圖7所示,隨機森林模型的精度與經(jīng)驗半經(jīng)驗?zāi)Ｐ拖啾扔幸欢ㄌ岣?在TM5、TM7、OLI、Sentinel-2 不同傳感器上模型精度R2均大于0.9。TM5傳感器RMSE=17.86 mg/L,TM7傳感器RMSE=7.00 mg/L,傳感器RMSE=21.22 mg/L,Sentinel-2傳感器RMSE=11.38 mg/L,整體誤差較小,可用于懸浮泥沙濃度反演監(jiān)測。

圖7 隨機森林模型精度評價

3 結(jié)果分析

應(yīng)用模型評價精度最高的隨機森林學(xué)習(xí)算法,基于2006年12月19日TM影像數(shù)據(jù)、2013年9月10日Landsat OLI影像數(shù)據(jù)及2020年4月29日Sentinel-2影像數(shù)據(jù)實現(xiàn)對珠江口懸浮泥沙反演、精度評價及空間制圖,反演結(jié)果如圖8所示。

圖8 珠江口懸浮物濃度反演

由圖8可看出,2006年12月20日冬季珠江口懸浮泥沙濃度整體平均濃度較低,處于0～152 mg/L范圍,平均值76.34 mg/L,最大值152 mg/L。珠江口懸浮泥沙整體分布西高東低,珠江東側(cè)支流較少且東江來水量小,導(dǎo)致虎門到入?？谔帒腋∧嗌碀舛戎递^小,處于0～40 mg/L范圍。珠江口西部支流眾多,橫門水道、洪奇門水道因受西江來水量大,加上兩岸社會經(jīng)濟發(fā)達人口眾多,最終攜帶泥沙含量較大堆積至入海口導(dǎo)致入?？诟浇鼞腋∧嗌碀舛却?珠江口西部橫門水道、洪奇門水道、東部交椅灣懸浮泥沙濃度值較高,處于0～160 mg/L范圍,最大值為152.5 mg/L。東部交椅灣海岸由于工程施工,植被覆蓋遭到破壞,懸浮泥沙受地表徑流攜帶入海,導(dǎo)致交椅灣懸浮泥沙濃度較大。

2013年9月珠江徑流量較冬季增多,珠江干流懸浮泥沙濃度整體平均濃度較高,處于0～242 mg/L范圍,平均值125.8 mg/L,最大值242 mg/L。珠江干流較2006年12月中間段懸浮泥沙濃度值較高平均為131.3 mg/L,入海口懸浮泥沙濃度平均為48.1 mg/L。珠江口懸浮泥沙濃度的分布在東莞區(qū)域呈現(xiàn)出從近岸到離岸逐漸遞減趨勢。珠江口西部橫門水道、洪奇門水道、東部交椅灣懸浮泥沙濃度值仍為范圍內(nèi)高值區(qū)域。秋季,珠江口呈現(xiàn)入?？诤臀鞑績善瑧腋∧嗌持递^高區(qū)域,主要原因是秋季珠江口沿岸上升流減弱,西北季風(fēng)還未開始,西南側(cè)遠(yuǎn)離珠江口頂,致使懸浮泥沙濃度未能出現(xiàn)冬季僅西南側(cè)高現(xiàn)象。而珠江口入?？谝虼思竟?jié)收到臺風(fēng)、潮汐涌動影響,在珠江口形成高濃度值區(qū)域。

2020年4月處于汛期初,珠江口懸浮泥沙濃度整體平均濃度較高,處于0～248.2 mg/L范圍,平均值139.4 mg/L,最大值248.2 mg/L。但此季節(jié),珠江干流懸浮泥沙濃度值較高,濃度值均值為131.6 mg/L,深圳西側(cè)海域懸浮泥沙呈現(xiàn)大面積濃度值較高現(xiàn)象,均值高達168.2 mg/L,主要受深中通道開通工程施工影響,懸浮泥沙較大。東部交椅灣因工程竣工植被恢復(fù)懸浮泥沙明顯較2006年、2013年降低,項目工程施工對懸浮泥沙濃度影響非常大。珠江口西部橫門水道、洪奇門水道受自然徑流大影響懸浮泥沙濃度值仍較高,橫門水道南支受上游輸沙率減少懸浮泥沙濃度值也發(fā)生較大減少。

4 結(jié) 語

本文對比了機器學(xué)習(xí)人工神經(jīng)網(wǎng)絡(luò)、支持向量回歸機及隨機森林模型,發(fā)現(xiàn)隨機森林模型在珠江口懸浮泥沙遙感監(jiān)測模型精度最高。通過高光譜數(shù)據(jù)重采樣到多光譜波段上,不同傳感器構(gòu)建對應(yīng)的光譜響應(yīng)函數(shù),將實測高光譜數(shù)據(jù)重采樣到對應(yīng)傳感器通道,能實現(xiàn)多傳感器的懸浮泥沙遙感監(jiān)測。采用模型精度最高的隨機森林模型對珠江口2006年Landsat TM5、2013年Landsat ETM、2020年Sentinel-2進行了懸浮泥沙濃度反演,反演結(jié)果發(fā)現(xiàn)珠江口岸懸浮泥沙濃度呈現(xiàn)西高東低,從近岸到離岸逐漸遞減的趨勢。主要原因是珠江口地形為喇叭形,在季風(fēng)和潮汐共同作用下,珠江口頂部區(qū)域受潮汐和風(fēng)向混合作用強烈;在遠(yuǎn)離珠江口頂部區(qū)域地形較為開闊,各支流徑流的懸浮泥沙匯聚堆積在珠江口西南區(qū)域,導(dǎo)致珠江口西南區(qū)域懸浮泥沙濃度高。基于地面實測站點數(shù)據(jù)對影像反演結(jié)果進行了評價,影像反演結(jié)果與實測數(shù)值相關(guān)系數(shù)較高,模型具有較高精度;同時發(fā)現(xiàn)港珠澳大橋的修建對珠江口懸浮泥沙的波動起著較大影響。