基于GOCI 數據的黃河口高懸浮物濃度海域的葉綠素a 提取指數構建與濃度反演

張曦元，萬劍華，劉善偉，宋冬梅

（中國石油大學（華東）海洋與空間信息學院，山東 青島 266580）

海水的固有光學特征是水體中所含各類物質對光的吸收和散射的結果，光譜特性復雜且具有區域性差異[1]。GORDON H R 等[2]依據水質成分將海水分為Ⅰ類水體和Ⅱ類水體，Ⅰ類水體的光學特性主要受葉綠素a（Chlorophyll-a，Chl-a） 及其分解產物的影響；Ⅱ類水體的光譜特征受Chl-a、總懸浮物（Total Suspended Matter，TSM）及黃色物質等多重影響，其光學特性也更為復雜。渤海西南部黃河口附近及萊州灣、渤海灣近岸海域為典型的Ⅱ類水體，黃河攜帶大量泥沙匯入萊州灣，在萊州灣、渤海灣海域隨洋流擴散，高濃度TSM 的不透光性干擾了Chl-a 濃度的光譜量測，使得該海域水體光譜特征復雜，難以準確估算Chl-a 濃度[3]。采用GOCI（Geostationary Ocean Color Imager） 衛星遙感技術構建Chl-a 濃度反演模型，能夠動態監測高TSM 濃度海域的Chl-a 濃度，打破傳統方法在時空上的局限。

目前，國內外學者針對復雜渾濁水體Chl-a 濃度的反演，提出了多種方法來構建其反演模型。三波段因子（Three-band factor，TBA）模型可用來估算渾濁水體的Chl-a 濃度，該方法具有明確的物理意義且精度較高[4]，因而被研究者廣泛應用[5-7]；四波段半解析模型能夠提高太湖高濁度水體Chl-a 濃度的估算精度[8]；改進三波段半分析模型可用于黃河口高渾濁水體Chl-a 濃度的反演，該算法需滿足λ1和λ2的波長在660～690 nm 之間，且λ3的波長要大于750 nm 的條件[9]。修正型葉綠素指數是針對實測水體光譜的近紅外與紅光波段的反射率進行比值修正，能夠削減懸浮物對葉綠素光譜信息的影響，有效提高了太湖懸浮物主導型水體的Chl-a 濃度的反演精度[10]；綜合葉綠素指數能夠提高長江口渾濁水域Chl-a 濃度的估算精度，但其只能對MERIS（Medium Resolution Imaging Spectrometer Instrument）等高光譜數據使用[11]。上述方法針對高濃度懸浮物水體Chl-a 濃度反演模型的研究取得了較高的精度，但其普適性較差，對衛星數據波段區間的要求較高，對于時間分辨率較高的GOCI 數據等并不成立。

因此，本文針對黃河口近岸、萊州灣及渤海灣高TSM 濃度的渾濁海域，分析了實測水體光譜曲線的特征，提出一種新的Chl-a 提取指數（Chlorophyll Extraction Index，CEI），該指數旨在高TSM濃度的水體中，通過因子組合來削弱TSM 對Schl的影響。

1 研究區與數據

1.1 研究區

黃河是世界上泥沙含量最高的河流，每年向海洋輸入巨量的懸浮泥沙以及各類營養鹽等物質，因黃河下游流速的降低易造成入海口物質沉積[12]。隨著經濟的快速發展和排水區人為活動的強烈影響，近二十年來，污水排放的增加、農業化肥的過量使用以及工業廢水的排放，使得黃河口近岸以及萊州灣、渤海灣海域的海水富營養化十分嚴重，藻類、浮游生物的大量繁殖導致海洋生態環境嚴重惡化[13]。因此，本文以高TSM 濃度的黃河口近岸及萊州灣、渤海灣海域為研究區，利用遙感技術監測該區域Chl-a 濃度的變化有利于保護渤海的生態環境。

1.2 實測數據

本文實測樣本數據包括光譜反射率數據和水質參數濃度數據，其觀測站點位置如圖1 所示。水體光譜反射率數據是2011 年7 月4 日于黃河口近岸測得的，共18 組；水質參數濃度數據包括Chl-a 和TSM 濃度，2011—2018 年共收集有效樣本數據619組，采集的樣本數據如表1 所示。取樣時間為當地時間9∶00 到15∶00，天氣晴朗、海上無大風。

表1 取樣數據說明

由實測數據可知，黃河口近岸海域葉綠素濃度范圍為0.73～15.4 μg/L，均值為3.79 μg/L；懸浮物濃度范圍為225.36～2 291.84 mg/L，均值為536.78 mg/L。渤海灣南部海域葉綠素濃度范圍為0.51～18.1 μg/L，均值為3.93μg/L；懸浮物濃度范圍為3.08～1 086.45mg/L，均值為116.41 mg/L。萊州灣海域葉綠素濃度范圍為0.45～40.7 μg/L，均值為4.46 μg/L；懸浮物濃度范圍為4.17～1 760.86 mg/L，均值為141.92 mg/L。由此可見，黃河口近岸懸浮物濃度遠遠高于其他海域，與順布日等[14]研究結果較為一致，嚴重干擾了Chl-a濃度對光譜波段敏感性的響應。

1.3 衛星數據及預處理

GOCI 是世界上首顆專門監測海洋水色產品的靜止軌道衛星[15]，其參數信息如表2 所示。GOCI衛星監測范圍能夠完全覆蓋我國渤海海域，每日可獲取8 景影像，在動態變化監測方面具有較大優勢。

表2 GOCI 衛星波段參數信息表

為獲取與現場取樣時間匹配的影像數據，本文根據具體取樣時刻選擇同步的GOCI 影像。影像預處理包括大氣校正、幾何校正、云掩膜等操作，其中大氣校正精度是影響反演模型質量的主要因素。由于渤海近岸水體的渾濁度高，其復雜的光學特性嚴重影響到大氣校正的精度。GDPS 軟件內置的算法難以降低渤海渾濁水體對大氣校正帶來的干擾[15]，何賢強基于紫外波長提出了一種適用于海洋水色圖像的大氣校正算法（UV-AC），能夠解決光譜匹配算法在渾濁水域失效的問題[16]。因此，本實驗采用UV-AC 法對GOCI 數據進行大氣校正。

2 研究方法

本文針對黃河口及其附近海域高TSM 濃度的渾濁水體，通過分析水體實測光譜數據變換曲線與Chl-a 濃度之間的相關性，確定GOCI 衛星數據所對應的敏感波段的反射率值。根據最大葉綠素指數（Maximum Chlorophyll Index，MCI） 的波段基線內差值原理構建葉綠素敏感指數（Chlorophyll-a Sensitivity Index，Schl）、懸浮物敏感指數（Total Suspended Matter Sensitivity Index，STSM），進一步構建葉綠素提取指數，該方法能夠削弱TSM 對Schl的影響。根據采樣點同步的GOCI 衛星數據得到CEI，以CEI 和實測Chl-a 濃度為樣本集，基于隨機森林算法，采用訓練樣本構建Chl-a 濃度的反演模型，而后利用測試樣本進行模型精度檢驗。

2.1 水體光譜數據的特征分析

水體的光譜反射特征受到純水、Chl-a、TSM 等有色物質的共同影響。本研究采集的海水光譜信息位于黃河口近岸，大量泥沙的匯入使得TSM 濃度相對較高；內陸污水的排放加劇了海水的富營養化程度，浮游植物的大量繁衍使得Chl-a 濃度相對較高；此外，儀器、光照等外界因素也會對數據采集產生一定的影響，因此該區域水體的光譜特征較為復雜。現場測量的海表光譜曲線如圖1 所示。

圖1 黃河口近岸海域實測光譜反射率曲線

黃河口近岸海域水體的光譜曲線基本符合Ⅱ類水體的特征，但又有其獨特性。400～500 nm 區間由于TSM 和Chl-a 的共同強吸收作用，光譜值較低；550 nm 附近由于TSM 散射作用形成較大的反射峰；660 nm 附近出現弱反射峰，此處受TSM 吸收作用和Chl-a 反射作用的共同影響；680 nm 附近出現吸收谷，該谷值通常被認為是判定水體是否含有Chl-a的依據；750～850 nm 范圍內因TSM 對水體光譜的后向散射作用而形成微小的反射峰，當水中TSM含量增大時，其反射峰也會向長波方向移動，也稱之為紅移現象；900 nm 后的光譜反射率存在較大噪聲，對此不做具體分析。

為了減少外界因素對采集光譜數據所造成的影響，本文采用二階微分變換的方式對原始反射率進行處理。原始光譜反射率及二階微分處理后的光譜反射率與Chl-a 和TSM 濃度的相關關系如圖2 所示。圖2（a）為原始光譜反射率與Chl-a 和TSM 濃度的相關性分析，圖2（b）為二階微分處理后光譜反射率與Chl-a 和TSM 濃度的相關性分析。

由圖2（a）可見，原始光譜的反射率與Chl-a 和TSM 濃度最大相關系數分別為-0.32 和0.98，分別出現在737 nm 和388 nm 附近。二階微分處理后的光譜反射率與Chl-a 和TSM 濃度的相關關系如圖2（b）所示，經二階微分處理后的光譜反射率與Chl-a 濃度的相關性較原始光譜反射率有明顯提升，但其與TSM 濃度的相關性略有下降，Chl-a 和TSM 濃度最大相關系數分別為0.88 和-0.85，分別出現在745 nm 和550 nm 附近，該相關系數能夠滿足實驗要求。

圖2 光譜反射率與Chl-a 和TSM 濃度的相關性分析

此外，經二階微分處理后的光譜反射率與Chl-a和TSM 濃度的相關性分析曲線大致呈對稱分布，555 nm、660 nm、680 nm 及745 nm 處呈明顯對稱分布。Chl-a 濃度與二階微分光譜反射率的相關性在這四處均呈現正相關性，555 nm 附近光譜與Chl-a濃度由負相關迅速轉為強正相關，且與TSM 濃度呈現明顯的負相關；660 nm 附近光譜與Chl-a 濃度正相關系數處于較低水平，與TSM 濃度的負相關系數處于較高水平，其原因可能是此處TSM 吸收作用大于Chl-a 反射效果；680 nm 附近光譜與TSM濃度的負相關系數由-0.78 驟減為-0.25，其與Chl-a濃度的正相關系數也有所下降，可能是TSM 的強反射與Chl-a 的強吸收作用的結果；745 nm 附近光譜與Chl-a 濃度的正相關系數達到最高值，與TSM 濃度的負相關系數超過-0.6。

2.2 葉綠素a 提取指數構建

MCI 利用2 個端點波段的遙感反射率，構建一條跨2 個端點波段區間的光譜基線。水體的Chl-a在信號波段之間具有光譜反射峰，其遙感反射率值與波段基線處的內差值之差為MCI[17]。

根據MCI 的波段基線內差值的原理，首先確定Chl-a 與TSM 的光譜吸收波段及多個敏感波段。通過分析海表光譜曲線特征以及Chl-a、TSM 與二階微分處理后的光譜反射率之間的相關性可知：555 nm處的光譜反射率與TSM 濃度由強正相關突然轉為強負相關；660 nm 附近TSM 對光譜的吸收效果較Chla 更為明顯；680 nm 附近Chl-a 的吸收作用較強，且該波段是判斷水中Chl-a 濃度高低的重要依據；745 nm 處因TSM 對水體光譜的后向散射作用而形成微小的反射峰，與Chl-a 濃度顯示強正相關。在此基礎上，分別計算Chl-a 與TSM 吸收波段的光譜基線內差值，本文將該值定義為Schl和STSM，其具體表達式如下。

式中，R 表示GOCI 數據在該波段的反射率；λ表示該波段對應的中心波長。

如圖3 所示，Schl是從680 nm 處的反射率值到660 nm 和745 nm 之間基線的反射率值的距離；STSM是從660 nm 處的反射率值到555 nm 和745 nm 之間基線的反射率值的距離。假設Chl-a 濃度是固定的，因為660 nm 處TSM 對光譜的吸收作用明顯，660 nm 的反射率值將會隨著TSM 濃度的增加而減少，則Schl將增大、STSM將減小。如果僅使用Schl對來反演Chl-a 濃度，Schl的增大將導致Chl-a 濃度的反演值偏高。

圖3 基于GOCI 遙感反射率和Chl-a 濃度構建的葉綠素提取指數示意圖

因此，當水體TSM 濃度較高時，Schl難免受到STSM的干擾。通過因子組合Schl-STSM來削弱TSM 對Schl的影響，能夠更好地突顯Chl-a 對光譜反射率特征變化的敏感性，該因子組合即為CEI。

2.3 基于隨機森林算法構建Chl-a 反演模型

隨機森林（Random Forests，RF） 是一種聯合多個回歸決策樹的算法，利用bootstrap 重采樣技術，對訓練集采取隨機有放回的重復抽取n 個樣本，并在bagging 的基礎上對每棵決策樹進行隨機特征的選擇，最后，采用簡單平均法對多個分類回歸樹得到的回歸結果進行算術平均，得到最終的模型輸出[17]。隨機森林因其隨機性地引入，使得結果不容易過擬合，具有良好的非線性數據集處理能力，可以達到反演模型的精度要求。

本實驗利用2011—2018 年多次采集的黃河口近岸、萊州灣及渤海灣水域的實測Chl-a 濃度及同期的GOCI 衛星遙感影像，以Schl和Chl-a 濃度為樣本集，根據空間均勻分布的原則，采用交叉驗證法對不同時期所測樣本進行篩選，分別得到591 個訓練樣本和200 個測試樣本。利用591 個訓練樣本基于隨機森林算法構建黃河口近岸高TSM 濃度海域的Chl-a 濃度反演模型，而后采用決定系數（R2）及平均相對誤差（MRE）對200 個測試樣本集進行模型精度檢驗。精度評價指標計算公式如下。

式中，yi和y 分別為Chl-a 濃度的實測值和反演值；n 為樣本數量；決定系數R2在0～1 之間，值越大，模型的反演精度越高；MRE 值越小，模型的反演精度越高。

3 結果與分析

3.1 基于不同指數法構建Chl-a 反演模型

本文將619 組樣本的實測Chl-a 濃度分別與三波段因子（TBA）、歸一化葉綠素指數（NDCI）、Schl以及CEI 進行相關性分析，結果如圖4 所示，決定系數I（R2）為Schl＞CEI ＞NDCI ＞TBA。

圖4 不同指數與實測Chl-a 濃度的相關性分析

利用419 個訓練樣本及同期的GOCI 衛星遙感影像，基于隨機森林算法分別建立TBA、NDCI、Schl及CEI 與Chl-a 濃度的反演模型，而后采用200個測試樣本依次對其模型進行精度檢驗并對比分析，結果如表3 所示。

目前，葉綠素濃度反演精度，對于Ⅰ類水體誤差不大于30%，而對于Ⅱ類水體誤差40%～50%。如表3 所示，無論是訓練樣本還是測試樣本，基于CEI 和Schl指數反演得到的Chl-a 濃度值與實測Chl-a濃度值之間的誤差均在40%～50%之間，其精度評價參數均優于TBA 與NDCI。因此，認為CEI 和Schl指數能夠滿足Ⅱ類水體反演精度的要求，可在渤海南部海域開展Chl-a 濃度的反演研究。

3.2 針對不同類別水體構建Chl-a 反演模型

由表3 可知，基于CEI 構建的反演模型精度優于Schl，其原因可能是本實驗在黃河口及其附近海域采集的樣本多位于高TSM 濃度的海域中，在高TSM 濃度的海域中，采用CEI 構建的Chl-a 反演模型的精度較高，而Schl更適合TSM 濃度低的海域。因此，將研究區分為TSM 主導型水體和Chl-a 主導型水體兩類，分別采用CEI 和Schl構建Chl-a 反演模型。

表3 訓練樣本及測試樣本的評價參數結果

3.2.1 確定分區標準

GOCI 衛星每天可獲取8 景影像，云霧的干擾會使得Chl-a 濃度數據出現缺失或異常，下午（14時）海霧發生頻數最低[18]，因此本實驗選取的GOCI影像均為同步于實測樣本數據當天的14 時的影像數據。根據GOCI 衛星光譜數據采用況潤元提出的光學分類法[19]，利用實測樣本數據的采集時間及同步的GOCI 衛星數據，將黃河口近岸海域、萊州灣及渤海灣海域劃分為Chl-a 主導型水體和TSM 主導型水體兩部分。此時公式如下。

式中，R6、R4分別為GOCI 的第6 波段和第4波段的反射率值；λ6、λ4分別對應其中心波長680nm和555nm。將（R6-R4）/λ6-λ4＞0 的海域歸為TSM主導型水體，反之則為Chl-a 主導型水體。

3.2.2 構建不同類型水體的Chl-a 反演模型

根據空間均勻分布的原則，采用交叉驗證法對不同時期所測樣本進行篩選，符合Chl-a 主導型水體的采樣點共有123 個，符合TSM 主導型水體的采樣點共有496 個。

基于隨機森林算法構建Chl-a 濃度反演模型時，Chl-a 主導型水體采用Schl，為進一步驗證CEI 在高TSM 濃度海域的精度，TSM 主導型水體采用Schl與CEI 進行對照，而后利用R2及MRE 分別對其測試樣本進行模型精度檢驗，結果如圖5、表4 所示。

圖5 不同類別水體的Chl-a 反演模型精度對比

對比表3 和表4 可以看出，基于不同類別水體構建的反演模型的MRE 明顯減小，R2均得到明顯提高，因此，在TSM 濃度較高的海域開展Chl-a 反演模型研究時，將樣本分為Chl-a 和TSM 主導型水體后，再分別采用Schl與CEI 構建Chl-a 反演模型能夠提高Chl-a 的反演精度。

表4 訓練樣本及測試樣本的評價參數結果

3.2.3 Chl-a 反演模型時空分析

本文基于隨機森林訓練模型對渤海海域開展了連續8 年（2011—2018 年）的Chl-a 濃度反演，實驗每月選取一幅晴朗天氣14 時的影像數據進行反演，按空間像素點進行平均，獲得Chl-a 濃度的年平均值。年平均反演結果如圖6 和圖7 所示，2011—2018 年渤海海域Chl-a 濃度的年際變化無顯著規律性，但其空間分布情況具有明顯的區域性特征，渤海灣沿岸及黃河口近岸的萊州灣海域的濃度值最高。

圖6 渤海灣海域的Chl-a 濃度反演圖

圖7 萊州灣海域的Chl-a 濃度反演圖

通過分析渤海海域地理位置及查閱相關資料發現，渤海灣潮差小、流速低、水體交換能力較弱，內陸廢水注入的營養鹽有助于浮游生物的生長繁殖，因而導致Chl-a 含量增加；萊州灣Chl-a 濃度偏高的原因可能是黃河徑流為近岸海域帶來豐富的營養物質，近岸水深淺且光照充足，為浮游植物提供適宜的生長環境，從而導致Chl-a 濃度偏高。

4 結 論

（1）通過分析水體光譜變換曲線與Chl-a、TSM濃度的相關性，根據MCI 原理構建了Schl和STSM，進一步構建了CEI，CEI 通過Schl-STSM的差值組合大大降低了TSM 對Schl的影響。

（2）基于隨機森林算法，根據GOCI 衛星的反射率與Chl-a 濃度值，分別采用TBA、NDC、Schl及CEI 構建Chl-a 濃度的反演模型，結果表明：采用TBA 和NDCI 所構建的Chl-a 反演模型的濃度值與實測值之間的MRE 均大于50%，不能滿足Ⅱ類水體反演精度的要求；基于Schl構建的反演模型精度優于CEI，在高TSM 濃度的海域中，采用CEI 能夠提高Chl-a 濃度反演模型的精度。

（3）根據波譜斜率將渤海南部的萊州灣、渤海灣海域劃分為Chl-a 主導型水體和TSM 主導型水體兩類，分別采用Schl與CEI 構建Chl-a 主導型水體和TSM 主導型水體的Chl-a 反演模型，結果顯示：針對不同類別水體構建的反演模型的精度較分類前有了明顯提高，且TSM 主導型水體構建的反演模型的MRE 與Chl-a 主導型水體的結果較為接近，證實了在高TSM 濃度的海域，采用CEI 指數法能夠有效削弱TSM 對Schl的影響。