基于GOCI的渤海海域赤潮信息遙感監(jiān)測與分析

程 玉，張圣佳，李金宇，朱紅春

(山東科技大學 測繪與空間信息學院，山東 青島 266590)

赤潮又名紅潮，是由海水中的某些浮游植物、原生動物或細菌在適宜的海域條件下，短時間突發(fā)性的劇增或聚集引起水體變色的異?，F(xiàn)象[1]。赤潮不僅對漁業(yè)資源尤其是海水養(yǎng)殖業(yè)造成影響，還會給沿海旅游業(yè)造成經(jīng)濟損失。此外，赤潮還對海洋生物造成巨大威脅，人類可能因意外飲食受赤潮污染的海洋產(chǎn)品造成中毒或死亡[2-3]。2009—2019年間，渤海共發(fā)生赤潮92次，主要發(fā)生在每年的5—9月份，其中，6—8月份赤潮發(fā)生的頻率最高。赤潮已經(jīng)成為我國三大近海污染之一，防止和減少赤潮造成的危害已成為十分緊迫的任務(wù)。

赤潮的形成機理十分復(fù)雜，雖然赤潮發(fā)生的具體原因還未詳細查明，但普遍認為海域的富營養(yǎng)化是赤潮發(fā)生的物質(zhì)基礎(chǔ)，溫度、鹽度、氮磷含量是發(fā)生赤潮的主要條件，氣溫、降水、氣壓等是赤潮發(fā)生的誘導因素[4]。目前，尚無有效的手段能阻止赤潮的發(fā)生，只能通過監(jiān)測和預(yù)報赤潮盡可能減少赤潮的危害。

傳統(tǒng)的赤潮監(jiān)測主要依靠船只調(diào)查和岸濱人工觀測。常規(guī)船只監(jiān)測赤潮需要對采樣進行生物、化學分析。飛機等飛行物監(jiān)測的費用高且受天氣限制，數(shù)據(jù)覆蓋面小，同步性差。衛(wèi)星遙感技術(shù)為赤潮的快速、大范圍連續(xù)監(jiān)測提供了重要手段[5]。常用的赤潮監(jiān)測方法有葉綠素濃度法、波段比值法、赤潮指數(shù)法等[6-9]。姜德娟等[10]基于中分辨率成像光譜儀(moderate-resolution imaging spectroradiometer,MODIS)數(shù)據(jù)，比較分析了葉綠素濃度異常法、波段比值法、紅波段差值法、karenia brevis指數(shù)法等對渤海赤潮遙感監(jiān)測的有效性，結(jié)果表明幾種方法都能較好地監(jiān)測到赤潮，葉綠素濃度異常法的赤潮監(jiān)測效果最好。薛存金等[11]針對單一海洋參數(shù)對赤潮提取不足的問題，提出一種多海洋參數(shù)赤潮監(jiān)測方法，利用海表面溫度、葉綠素a濃度、懸浮泥沙濃度等海洋參數(shù)，根據(jù)專家知識設(shè)置參數(shù)閾值，設(shè)計出一套渤海赤潮災(zāi)害提取規(guī)則，成功提取渤海灣的赤潮。李陽東等[12]利用懸浮物濃度和赤潮指數(shù)，通過海洋災(zāi)害公報面積設(shè)置不同的閾值組合提取赤潮，以2011年5月29日溫州海域和2014年5月23日舟山海域的赤潮事件為例進行赤潮日變化規(guī)律研究。利用海洋參數(shù)提取赤潮，需要設(shè)置合適的閾值進行圖像分割，從而達到赤潮提取的目的，因此，分割閾值的確定是一個重點問題。但是，目前所用的閾值分割方法大部分基于經(jīng)驗統(tǒng)計閾值，或根據(jù)海洋公報面積，或假彩色圖像通過目視匹配的方法設(shè)置閾值，這些方法不僅繁瑣而且還具有較強的主觀性，因此，需要采取自適應(yīng)閾值分割方法確定最優(yōu)閾值，減少人工干預(yù)。此外，對于赤潮的研究大都集中在赤潮的提取和變化監(jiān)測方面，對于赤潮演變導致水色參數(shù)變化的研究很少。

針對這一問題，本研究探討了赤潮發(fā)生前后海洋水色參數(shù)的變化特征。通過對葉綠素a和赤潮指數(shù)圖像進行多種方法的自適應(yīng)閾值分割提取，根據(jù)效果從中選取合適的赤潮提取方法，分析赤潮的晝變化特征，并對赤潮發(fā)展演變造成的影響進行研究。

1 研究區(qū)及數(shù)據(jù)

1.1 研究區(qū)概況

渤海位于37°07′～41°00′N，117°35′～122°16′E，是半封閉內(nèi)海，僅東面與黃海相連，也是中國最大的內(nèi)海，面積約77 000 km2，海岸線全長約3 800 km，平均水深18 m。渤海處于北溫帶，受北方大陸性氣候影響，夏無酷暑，冬無嚴寒，全年區(qū)域溫差不顯著。沿岸周邊河流眾多，其中主要有黃河、海河、遼河和灤河等。河流攜帶周邊三省一市的工農(nóng)業(yè)廢水和生活污水直接流入渤海，陸域活動的影響直接導致渤海水質(zhì)下降和富營養(yǎng)化，呈現(xiàn)二類水體特征，赤潮災(zāi)害頻繁發(fā)生。

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

地球靜止海洋水色成像儀(geostationary ocean color imager, GOCI)搭載于韓國通信、海洋和氣象衛(wèi)星(communication ocean and meteorological satellite, COMS)，空間分辨率為500 m(圖1)；時間分辨率極高，可以做到每天8次，每次間隔1小時的觀測頻率；波段范圍412～865 nm，包括6個可見光波段和2個近紅外波段。在本研究中，根據(jù)海洋公報中公布的赤潮災(zāi)害信息，分別下載了2014年4月28日、5月5日、5月8日、5月12日、5月15日、5月26日、5月28日、6月9日、6月12日、6月13日和6月23日的L1B數(shù)據(jù)，用于赤潮晝變化和海洋參數(shù)變化分析。幾何校正采用自帶的經(jīng)緯度信息生成地理位置查找表文件(geographic lookup
Table, GLT)，利用GOCI數(shù)據(jù)處理軟件(GOCI data progressing software, GDPS)內(nèi)置的二類水體大氣校正算法實現(xiàn)大氣校正，獲得反射率數(shù)據(jù)，后續(xù)進行懸浮物、有色溶解有機物(chromophoric dissolved organic matter, CDOM)等海洋水色參數(shù)反演。本研究中葉綠素a濃度、懸浮物、漫衰減系數(shù)(Kd(490))和CDOM等水色參數(shù)均基于GOCI影像反演計算得到，歸一化熒光高度和海表面溫度數(shù)據(jù)來自美國國家航空航天局海洋水色網(wǎng)站[13]。

圖1 渤海區(qū)域(綠色邊界)GOCI真彩色合成影像圖

2 研究方法

2.1 赤潮指數(shù)

Ahn等[14]基于寬視場水色掃描儀(sea-viewing wide field-of-view sensor, SeaWiFS)數(shù)據(jù)，提出利用赤潮指數(shù)(red tide index, RI)來表征光學特性復(fù)雜的東北亞海域爆發(fā)赤潮的可能性。宋德彬[15]基于GOCI數(shù)據(jù)分析了赤潮和潔凈水體的反射率曲線，其中赤潮水體在550 nm處的反射率明顯較高；與潔凈水體相比較，在443和490 nm處，赤潮水體反射率與550 nm處的差異顯著增大。基于此差異，采用張賀等[16]改進的赤潮指數(shù)來進行渤海赤潮信息的提取，表達式如下：

(1)

其中：RI表示赤潮指數(shù)；Rrs(443)、Rrs(490)和Rrs(555)分別表示443、490和555 nm波段的遙感反射率。使用藍綠波段反射率來構(gòu)建赤潮指數(shù)，并通過減去443 nm波段的反射率來減少懸浮泥沙對赤潮信息提取的影響。

2.2 閾值確定方法

2.2.1 最大類間方差法

最大類間方差法(OTSU)要求選擇最佳閾值，使圖像前景和背景間的灰度值差距最大，衡量灰度差異的標準定義為類間方差，即當前景和背景間的類間方差取最大值時，獲得最佳閾值。假設(shè)圖像的像素總數(shù)為N，L個灰度級[1,2,…,L]，灰度值為i的像素總數(shù)為ni，像素i的出現(xiàn)概率

(2)

設(shè)圖像的初始閾值為t，若像素的灰度值小于t，則屬于C0類，否則屬于C1類，C0和C1分別代表圖像的背景類和前景類。圖像背景類和前景類像素出現(xiàn)的概率和灰度均值分別用ω0、ω1、μ0、μ1表示，有：

(3)

(4)

(5)

其中，σ2表示類間方差，當σ2取最大值時，獲得最佳閾值。

2.2.2 K-means聚類

K-means聚類算法依據(jù)數(shù)據(jù)中點的歐式距離將給定的數(shù)據(jù)歸入需要劃分的K類中。假設(shè)將數(shù)據(jù)集X(X中的數(shù)據(jù)個數(shù)為n)分為K個類別，將歐式距離作為相似度指標，起始階段在數(shù)據(jù)集中隨機選取k個值分別作為每一類的中心，之后計算數(shù)據(jù)集中每一個點x與各個中心u之間的歐式距離，并且依據(jù)歐式距離的大小將每個點歸入與之最近的中心所屬類別中，使各類別的聚類平方和I最小，即

(6)

之后更新每個聚類的中心位置，新的中心位置即為每一類所有點的平均值，不斷重復(fù)這一過程，當聚類平方和I達到最小時重復(fù)過程終止，此時的分類為K-means聚類的最終結(jié)果。

2.2.3 最大熵法

最大熵法通過引入信息論中最大熵原則，統(tǒng)計每個灰度值的出現(xiàn)概率，計算灰度級的熵，保證各類灰度的總熵最大，獲得使前景或背景內(nèi)部灰度盡可能均勻的最優(yōu)閾值，或是尋找使分割前后圖像的信息量差異最小的最優(yōu)閾值。假設(shè)某個灰度級出現(xiàn)概率為pi，則灰度熵的計算公式為：

(7)

設(shè)灰度值q為初始閾值，該閾值分割的背景類C0和前景類C1的累計概率用P0(q)和P1(q)表示，

(8)

背景類和前景類的概率密度函數(shù)可分別表示為：

(9)

(10)

背景類和前景類對應(yīng)的熵分別表示如下：

(11)

(12)

在該閾值下，圖像總熵為：

H(q)=H0(q)+H1(q)。

(13)

計算所有分割閾值下的圖像總熵，將最大熵對應(yīng)的分割閾值作為最佳分割閾值。

3 赤潮信息提取及晝變化分析

通過比較渤海區(qū)域云量，選取2014年5月28日云量最少的1景GOCI影像，即編號為00(8:00)的影像，分別使用OTSU、K-means聚類以及最大熵方法對葉綠素a和RI指數(shù)計算結(jié)果確定最佳的赤潮分割閾值，利用該閾值進行遙感影像二分類，得到赤潮類型的提取結(jié)果并統(tǒng)計其面積。以目視解譯結(jié)果為參考結(jié)果，采用不同方法進行閾值分割后提取的面積及誤差見表1。

表1 3種閾值分割方法面積提取及誤差分析

3種閾值分割方法對葉綠素a影像確定的閾值不斷增大，OTSU方法確定的閾值為10.5 mg·m-3，分割后得到的面積最接近目視解譯面積，為8 089.25 km2，相對誤差為9.1%；K-means聚類方法確定的閾值為11.7 mg·m-3，赤潮面積為6 619.5 km2，相對誤差25.6%；最大熵法確定的閾值為12.3 mg·m-3，提取的赤潮面積為6 100 km2，相對誤差為31.4%，是3種閾值分割算法中誤差最大的算法。對于RI圖像閾值分割結(jié)果，OTSU算法的赤潮提取結(jié)果最好，相對誤差為5.8%；K-means聚類算法提取的赤潮面積相對誤差為8.7%；最大熵法的相對誤差結(jié)果最大，為9.4%。從3種方法提取面積的誤差分析，OTSU算法對于葉綠素a影像和RI影像中赤潮提取結(jié)果最好，K-means聚類方法次之，最大熵法的結(jié)果最差。對于葉綠素a影像，3種分割方法提取赤潮面積的平均相對誤差為22%，而RI影像3種方法的平均相對誤差為7.9%，因此，RI影像更適合赤潮提取與分析。綜合來看，使用RI影像并采用OTSU閾值分割算法提取赤潮更有效，其赤潮提取結(jié)果見圖2。根據(jù)海洋公報的記錄，2014年5月28日， 河北秦皇島灤河口至遼寧綏中海域發(fā)生大規(guī)模赤潮，實驗中提取的赤潮范圍與赤潮記錄基本一致，可以證明RI影像采用OTSU閾值分割算法提取赤潮準確率較高。赤潮的發(fā)生是多因素共同作用的結(jié)果，通過查詢歷史氣象資料，2014年5月24日，河北秦皇島附近地區(qū)有中雨，最高氣溫18 ℃，降雨導致河流流量增加，攜帶更多營養(yǎng)鹽進入海洋，在5月25日之后，氣溫上升5～8 ℃，最高氣溫達30 ℃，且天氣狀況良好，充足的光照有利于浮游植物的生長，促進赤潮的大范圍暴發(fā)。本次赤潮的主要藻類是抑食金球藻，其最適宜的生長溫度為20 ℃[17]，秦皇島附近海域5月26日的海表面溫度在19 ℃左右，基本滿足抑食金球藻的溫度條件，促進了藻類的快速增殖。而未發(fā)生赤潮的遼寧錦州附近海域雖然也有降雨，但氣溫上升幅度較小，海表面溫度在14 ℃左右，不利于藻類的生長，因此錦州附近海域并未有赤潮暴發(fā)，同時也間接證明提取的赤潮結(jié)果準確性較高。

圖2 赤潮提取結(jié)果和兩個感興趣區(qū)域(黃色方框1和2)

赤潮指數(shù)圖像能反映赤潮在一天內(nèi)的變化，通過閾值分割提取赤潮面積，可大致了解赤潮加劇和消亡的過程。由于太陽高度角對于地表的太陽輻照度影響較大[18]，從而影響地表尤其是水體的表觀輻亮度[19]，進而影響分割閾值的確定，因此對2014年5月28日8景GOCI影像計算得到的RI結(jié)果采用OTSU方法分別計算了各自的分割閾值，如表2所示。分別提取每幅影像中的赤潮，選取兩個無云覆蓋的區(qū)域，分別統(tǒng)計區(qū)域內(nèi)8景影像中的赤潮面積，其面積變化結(jié)果見圖3。

圖3 2014年5月28日感興趣區(qū)域的赤潮面積變化

表2 2014年5月28日8景影像的分割閾值

區(qū)域1和區(qū)域2的赤潮面積均呈先增大后減小的變化趨勢。區(qū)域1的赤潮面積較大，一天之內(nèi)的面積變化比較劇烈。赤潮面積從8:00左右的2 411 km2開始增加，10:00左右達到一天之中的最大值，面積為2 565.25 km2，之后赤潮面積開始減少，15:00達到最小值，面積為1 258.25 km2，減少了51%。區(qū)域2中赤潮面積從8:00的715.25 km2增加到10:00的772 km2，15:00左右赤潮面積降為384.25 km2，降低了約50.2%。從面積變化上可以看出，赤潮在一天內(nèi)的變化十分明顯，10:00左右赤潮最為嚴重，分布面積最大，這可能與浮游植物的垂直遷移有很大關(guān)系。當光照充足時，浮游植物會上浮到海水表面吸收光照，進行光合作用，產(chǎn)生營養(yǎng)物質(zhì)，當光照減弱后，浮游植物會下潛至海底吸收營養(yǎng)。但10:00之后光照仍然充足，赤潮面積反而開始下降。根據(jù)氣泡調(diào)節(jié)浮力機制[20]，日出之后，浮游植物上浮至海面進行光合作用，由于小氣泡的作用細胞浮力增加，浮游植物量在10:00的時候形成上層高峰，中午強日照下的光合產(chǎn)物使得細胞內(nèi)部壓強增大導致氣泡破裂，細胞失去浮力造成浮游植物下沉，因此赤潮面積在10:00之后開始下降。根據(jù)項斯端等[21]藍藻型富營養(yǎng)化湖泊藻量的晝夜變化規(guī)律，藻量水平在一天中會出現(xiàn)兩個峰值：第一次峰值出現(xiàn)在上午10:00，中午12:00至14:00藻量大幅下降并出現(xiàn)低谷；下午藻量水平上升，18:00至20:00之間出現(xiàn)第二次高峰，之后藻量又開始下降。但由于GOCI衛(wèi)星在15:00之后不再采集影像，所以此規(guī)律還需進一步研究證明。

4 赤潮發(fā)生過程的水色參數(shù)變化分析

赤潮往往是由某種藻類生物異常繁殖導致，只有當各種理化條件基本滿足某種赤潮藻生理、生態(tài)需求時，該種藻類才可能形成赤潮。赤潮的發(fā)生是多因素共同作用的結(jié)果，與赤潮發(fā)生區(qū)域的環(huán)境密切相關(guān)。赤潮生物的生長發(fā)育和繁殖都要從環(huán)境中吸收營養(yǎng)物質(zhì)和能量，同時各個階段又受周圍環(huán)境因素的制約。赤潮會對環(huán)境因素造成影響，使得水色參數(shù)在赤潮前后發(fā)生變化，因此，赤潮與環(huán)境之間相互影響。

根據(jù)海洋災(zāi)害公報，2014年5月15日—8月7日，河北秦皇島灤河口至遼寧綏中海域發(fā)生大規(guī)模赤潮，但從2014年6月12日的葉綠素a影像和RI影像中并未提取到大范圍赤潮，因此，可以認為此次大范圍赤潮在6月12日左右就已基本消亡。本研究選取一些赤潮發(fā)生區(qū)域云量較少的影像對赤潮導致的海洋水色參數(shù)變化進行分析，按照公報顯示的赤潮記錄，4月28日、5月5日、5月8日和5月12日赤潮還未爆發(fā)，5月15日、5月26日和5月28日為赤潮快速生長繁殖時期，6月9日、6月12日、6月13日、6月23日為大規(guī)模赤潮消亡時期。由于云層的覆蓋，選取這些時段也只能提取到部分赤潮區(qū)域的連續(xù)信息。

4.1 漫衰減系數(shù)Kd(490)

水體漫衰減系數(shù)是描述輻照度隨海水深度增加而衰減的重要海洋光學參數(shù)，可以反映水體的環(huán)境狀況。如圖4(a)所示，赤潮發(fā)生前Kd(490)在0.2 m-1左右，赤潮發(fā)生后Kd(490)開始增大，赤潮發(fā)生第一天達到0.29 m-1，之后數(shù)值迅速增大，5月26日達到最大值0.68 m-1，直至6月12日赤潮基本消亡，Kd(490)值減小到0.09 m-1，此后Kd(490)值基本穩(wěn)定在0.1 m-1左右。在清澈水體中，水體下行輻射隨海水深度增加而衰減的趨勢減小且出現(xiàn)穩(wěn)定的狀態(tài)，因此漫衰減系數(shù)較??；當赤潮發(fā)生后，下行輻照度衰減程度較大，導致漫衰減系數(shù)增大。

4.2 懸浮物濃度

懸浮物濃度對浮游植物種群生長影響是多種機制疊加的結(jié)果，較低濃度的懸浮物有利于赤潮藻類的生長，但高濃度的懸浮物具有抑制作用[22]。由圖4(b)可見，5月上旬，懸浮物濃度達2 g·m-3,伴隨著赤潮的生長繁殖，濃度不斷升高，5月28日上升到3.6 g·m-3，之后濃度下降至1～1.5 g·m-3波動。從赤潮發(fā)生前后懸浮物的濃度來看，懸浮物濃度總體較低，赤潮藻類的生長并沒有受到懸浮物的影響，因此，隨著赤潮的生長，懸浮物濃度升高，赤潮消亡后，濃度隨之下降。

圖4 各水色參數(shù)時間序列圖

4.3 CDOM

CDOM為水體內(nèi)一類重要的光吸收物質(zhì)，其濃度和組成能改變水下光場，對水體生態(tài)系統(tǒng)有明顯影響。從圖4(c)分析，CDOM變化也遵循先增大后減小的趨勢，且在赤潮最嚴重時達到最大，赤潮消亡，CDOM也隨之減小。但整體變化較小，基本在0.12～0.24 m-1變化。浮游植物的大量生長聚集會伴隨一些植物自身降解，為海水中CDOM提供了充足來源，促使CDOM數(shù)值上升。

4.4 海表面溫度

赤潮爆發(fā)需要適宜溫度，同時，大量的浮游植物聚集也會對海表面溫度造成影響。但圖4(d)的海表面溫度變化情況并沒有顯示出海表面溫度隨著赤潮的生長消亡發(fā)生變化。從4月28日開始，海表面溫度一直上升，原因是隨著季節(jié)變化和氣溫的不斷升高海表面溫度也不斷上升，沒有體現(xiàn)出受赤潮影響導致的變化。不過，此次赤潮的主要種類是抑食金球藻，適宜的生長溫度為15～25 ℃，而5月26日的海表面溫度達15 ℃，促進了赤潮的大規(guī)模爆發(fā)。

4.5 歸一化熒光高度

赤潮爆發(fā)時，浮游植物大量聚集，葉綠素a濃度異常增高，其水體光譜存在一個熒光峰，可以用歸一化熒光高度來定量表示(圖4(e))。5月8日前歸一化熒光高度值約為0.35，赤潮爆發(fā)后浮游生物量增加，歸一化熒光高度值上升，6月9日歸一化熒光高度值和赤潮爆發(fā)前基本持平，之后一直維持在0.37左右。歸一化熒光高度情況基本反映葉綠素a濃度隨赤潮發(fā)生的變化情況，盡管葉綠素a濃度與歸一化熒光高度不成正比，但赤潮的生長消亡也會導致葉綠素a濃度發(fā)生變化。

綜合各種水色參數(shù)變化，赤潮的爆發(fā)會導致Kd(490)、懸浮物濃度、CDOM和歸一化熒光高度出現(xiàn)不同程度的升高，盡管本研究中的海表面溫度并沒有隨之變化，但多項研究表明，赤潮的發(fā)生會導致海表面溫度的輕微上升[23-24]。Kd(490)、懸浮物和CDOM這3種水色要素在赤潮發(fā)生后的數(shù)值低于赤潮發(fā)生前，這可能是由于赤潮的生長發(fā)育消耗了水體中的各種物質(zhì)，從而達到了凈化水體的作用，導致數(shù)值較赤潮發(fā)生前低。

5 結(jié)論

赤潮是由多種因素作用的結(jié)果，同時，赤潮也會引起多種海洋參數(shù)發(fā)生變化，為赤潮的預(yù)測預(yù)警提供了依據(jù)。

1) 在閾值分割方法上，OTSU方法在葉綠素a影像和赤潮指數(shù)影像上提取的赤潮面積精度最高，K-means聚類方法次之，最大熵法的精度最差。在赤潮監(jiān)測方法上，葉綠素a濃度圖像3種分割方法和赤潮指數(shù)圖像的平均相對誤差分別為22%和7.9%，赤潮指數(shù)圖像的提取結(jié)果更準確。因此，選擇OTSU方法結(jié)合赤潮指數(shù)圖像進行赤潮提取。

2) 以區(qū)域1和區(qū)域2的赤潮面積為例研究了赤潮的晝變化情況，赤潮面積存在明顯的晝變化特征，兩個區(qū)域的赤潮面積在一天的變化均先增大后減少，在10:00面積達到最大值，這可能與浮游植物的垂直遷移有關(guān)。

3)Kd(490)、懸浮物濃度、CDOM、歸一化熒光高度等水色參數(shù)均隨赤潮的發(fā)生而增高，赤潮消亡后數(shù)值也隨之下降。海表面溫度在赤潮發(fā)生周期一直上升，季節(jié)變化導致赤潮引起的溫度上升不明顯，海洋水色參數(shù)的變化可以為赤潮的預(yù)警預(yù)報提供依據(jù)。

由于云層遮擋等原因?qū)е鲁喑北O(jiān)測的時間序列不完整，而且GOCI數(shù)據(jù)的空間分辨率不夠精細，所以赤潮提取面積會有較大誤差，下一步可以考察更高分辨率的數(shù)據(jù)結(jié)合更多的海洋環(huán)境參數(shù)建立赤潮預(yù)警預(yù)報模型。同時對于赤潮發(fā)生后海洋水色參數(shù)的數(shù)值低于赤潮發(fā)生前的數(shù)值的原因還需進一步研究。