基于SOM的海表溫度遙感數據集的EOF算法重構

王建樂 ，朱建華 ，何宇清 ，宋占杰

(1.天津大學電視與圖像信息研究所，天津 300072；2.國家海洋技術中心，天津 300112)

基于SOM的海表溫度遙感數據集的EOF算法重構

王建樂1，朱建華2，何宇清1，宋占杰1

(1.天津大學電視與圖像信息研究所，天津 300072；2.國家海洋技術中心，天津 300112)

針對海表溫度數據集的數據缺失，提出了一種基于自組織映射算法（SOM）和經驗正交函數算法（EOF）有機結合的重構缺失值的新方法。該方法應用了SOM的非線性估計，能夠很好的反映數據集的非線性結構，并把SOM估計的結果用于EOF算法的初始化，克服了EOF對數據集初始化敏感的問題。在處理過程中，對奇異值分解使用了lanczos算子分解矩陣，提高了程序運行效率。此外，該方法還引入蒙特卡羅交叉校正集，確定最佳重構的EOF模態數，最終高精度計算出重構誤差。使用AQUA遙感衛星海表溫度數據進行實驗，結果表明該方法能夠很好地重構出缺失率高達83.23%的數據集，且重構精度高。

遙感；海表溫度SST；SOM算法；EOF算法；lanczos算子

中分辨率成像光譜儀 (MODIS)是美國國家航空航天局(NASA)對地觀測系統(EOS)計劃中最有特色的傳感器之一。MODIS每兩天連續提供地球上任何地點白天反射圖像和晝夜的發光光譜圖像數據，包括對地球陸地、海洋和大氣觀測的可見光和紅外波譜數據。MODIS在紅外大氣窗的5個波段可以用來反演SST[1]。相比其他海洋表面溫度的獲取方式，衛星遙感具有大面積同步觀測、動態與長期觀測、實時或者準實時等特點，同時具備費用相對較低，可以監測船舶、浮標不易到達的海區等諸多優點[2]。但是受氣候隨機變化的影響，數據的缺失給遙感產品的使用帶來了很大的限制。多數情況下傳統檢測雖然可以容忍部分數據的缺失，但隨科技的發展，許多部門對完整的遙感數據集的需求越來越大。目前已經發展了許多種方法來解決填充缺失值的問題，諸如最優插值法（OI）、經驗正交函數分解法（EOF）、自組織映射法（SOM）、本征模態分解法（POD）、期望最大化法（EM）、奇異譜分析法（SSA）等等。這些方法大致可以分為兩類：確定性的方法和隨機性的方法[3]。本文旨在利用搭載在AUQA衛星上的MODIS傳感器測得的遙感數據對我國長江口附近海域的海表溫度缺失數據進行反演重構研究，進而分析該海域的物理現象。

自組織映射算法[3]（Self-Organizing Maps，SOM）是一種非線性、隨機、可以對相同特性的個體進行分類，強調各類別之間的相鄰結構的算法。它模擬人腦中處于不同區域的神經細胞分工不同的特點，即不同區域具有不同的響應特征，而且這一過程是自動完成的。該算法基于無監督學習法則，其訓練完全是隨機的、數據之間的。自組織映射算法允許高維數據向低維網格數據的估計。盡管該估計聚焦在數據的拓撲特性上，但是SOM算法可以對缺失數據進行非線性插值。

經驗正交函數算法[2]（Empirical Orthogonal Functions，EOF）是一種線性、確定性、無參數，并利用了遙感數據集時間和空間的相關性，能夠對高維空間數據進行線性估計的算法。EOF算法可以對缺失數據進行連續地線性插值，而且與傳統的最優插值（OI）相比，不需要遙感數據集的先驗信息。但是由于其運行時數據集不能有缺失值，所以對于數據集中缺失值的初始化比較敏感，這給EOF的應用增加了負擔。

本文提出了一種基于SOM的EOF算法來重構遙感數據集缺失值，該方法兼顧了自組織映射算法和經驗正交函數算法的優點。首先利用SOM算法的非線性特性對原始數據去噪，隨后利用其運行結果作為EOF算法處理時對數據集的初始化，再利用EOF算法對數據集進一步去噪，最后利用其連續插值特性有效地重構缺失數據。

1 海表溫度重構

1.1 SOM算法[3]

SOM算法是一種基于無監督的學習原則的算法，即訓練完全是基于數據自身的，不需要輸入數據的信息的算法。針對此方法，本文使用了一個二維的網絡，該網絡由形狀是矩形格子框架的c個單元（或編碼向量）組成。網絡的每一個單元都有數量和學習數據樣本xn(n=1，2，…，N)的長度T一樣多的權重。網絡的所有單元可以由一個權重矩陣m(t)=[m1(t)，m2(t)，…，mc(t)]組成，其中 mi(t)是在時間 t的 i單元的T維權重向量，t代表學習過程的步驟。每一個單元都通過鄰域函數λ(mi，mj，t)和其相鄰單元相連接，該鄰域函數定義了在時間t時鄰域的形狀和大小。該鄰域在整個學習過程中可以不變，也可以改變。

首先，學習過程開始于隨機的初始化網絡節點的權重。然后，再隨機的選擇樣本xt+1，進而計算權重最接近樣本的神經元——最佳匹配單元BMU（Best Matching Unit）。BMU定義為：

式中：I是網絡節點索引的集合。BMU指出了最佳匹配節點的索引，‖·‖是標準的歐式范數。

但是，如果隨機選擇的樣本包含缺失值，我們就不能完全地求解BMU。此時，使用Cottrell和Letremy[4]的調整的SOM算法來代替。含有缺失值的隨機選擇的樣本xt+1分為兩個子集，其中前者是沒有缺失值的子集，后者是含有缺失值的子集。在沒有缺失值的子集上定義范數：

式中：xt+1，k，k=[1，2，…，T]表示已選擇的數據向量的第 k個值；mi，k（t），k=[1，2，…，T]，i=[1，2，…，c]表示第i個編碼向量的第k個值。k在沒有缺失值的子集遍歷所有的索引。然后，再計算BMU。

在找到了最佳匹配單元BMU后，網絡權重調整為：

式中：ε（t）是調整后獲得的參數，值域為[0，1]，隨時間的增加而減小。在更新權重時，需要考慮的神經元數依賴于鄰域函數λ（mi，mj，t），需要權重調整的神經元數通常是隨著時間減少的。

在更新權重后，下一個樣本被隨機地從數據矩陣中選出，程序通過尋找樣本的BMU被再次執行。當SOM算法收斂的時候，就停止遞歸的學習程序。此時，得到了填充缺失值后的數據集。

1.2 EOF算法

Beckers[2]等在2003年提出了一種無參數的基于經驗正交分解（EOF）方法來重構時間序列數據中缺失值。相比經典的最優插值(OI)法，有著不需要相關函數矩陣，信噪比(SNR)和相關長度等先驗信息的優勢。Alvera-Azcarate[5]等在2005年提出了基于EOF分解的數據插值方法——DINEOF（Data Interpolating Empirical Orthogonal Functions），對1995年5月9日到10月22日的6個月的AVHRR衛星數據（剔除了其中數據缺失比例小于5%的圖像）進行了缺失值重構，并且引進了交叉校正集，允許建立EOF過程中的最優截斷和對缺省值的估計誤差（從交叉校正集中獲得）。并且與最優插值法（Optimal Interpolation，OI）進行了比較，DINEOF 法與 OI法有相似的重構精度，但是后者程序的運行時間比前者高了30多倍。

本文使用了Alvera-Azcarate[5-6]的DINEOF，其工作原理如下：

假設X0為m×n維的二維數據矩陣，其中m＞n（m是空間維，n是時間維），X0中包括一些缺失值，本文用NaN表示缺失值。

數據的去均值和初始化。計算X0中有效數據的均值0，令 X=X0-0；隨機挑選部分有效數據點集 Xcv作為獲取最佳重構模態數的交叉校正集（本文使用蒙特卡羅交叉校正集[7]）；對X中的處于交叉校正集的位置的值賦值為NaN；對X中所有的NaN用0代替，使缺失值的初始值為數據集的無偏估計值。同時，令k=1。

對X使用式（5）進行奇異值分解SVD，得到最主要的k個模態，使用式（6）計算缺失點的重構值。

式中：U，S，V分別對應SVD分解后數據集的空間模態矩陣、奇異值矩陣和時間模態矩陣；T表示矩陣的轉置。

式中：i，j為矩陣的空間與時間下標；up和vp分別是空間模態U和時間模態V的第p列；ρp為相應的奇異值，p=1，2，…，k。該過程迭代N次（N為運行程序前設定的最大的迭代次數），并且計算交叉校正集Xcv的重構值與原始值的均方根誤差RMS。

令k=2，…，kmax，重復上述步驟，計算出對應的均方根誤差，比較得出均方根誤差值最小時對應的模態數P，其中，kmax是根據時間維數n確定的（kmax≤n）。缺失點集的值用P模態時計算的重構值替換，交叉校正點集Xcv處的值用原始值替換，令 k=P，重復公式（5）、（6）對 X 的分解與重構，計算出所有點的重構值，仍記為X。再令X=X+0，此時就得到了重構的數據集。

針對EOF算法執行過程中奇異值分解SVD大矩陣效率低的問題，本文引入lanczos算子[8]加速SVD分解過程。又針對EOF算法的迭代過程效率低（一般預先設定的迭代次數都比實際使用的次數多）的問題，使用丁又專[9]提出的迭代收斂準則，實驗表明，該迭代準則能夠很大程度上提高程序運行效率。

1.3 疊加方法

本文將上述兩種方法進行有機結合。首先運行非線性估計的SOM算法重構缺失值，然后用SOM估計的結果作為運行EOF方法的初始化值。

對于SOM算法，需要選擇最優的網格大小c；對于EOF算法，需要選擇奇異值和奇異向量的最優數p。這兩個參數的選擇是為了使用蒙特卡羅交叉校正集[7]，并且對于c和p這兩個參數的所有組合都使用相同的交叉校正集。最后，使用給出最小交叉校正誤差的SOM和EOF的參數組合來最終重構數據的缺失值。

該算法中使用的蒙特卡羅交叉校正集的每一個校正集的均方根誤差的計算公式為：

式中：dataFilled表示插值后的測試數據；testdata表示插值前的測試數據；N表示插值數據的總數。

2 遙感數據來源和研究區域

本文使用MODIS衛星數據反演長江入海口附近海域，從2010年8月1日到2010年11月8日的海表溫度（SST）數據集（http://oceancolor.gsfc.nasa.gov/cgi/browse.pl?sen=am）。選擇裝載在AQUA遙感衛星上的MODIS L2海表溫度產品。為了更好地利用數據的時間空間相關性，鑒于海表溫度的保守性，數據是每隔約3 d取1次的衛星遙感海表溫度資料，總共37幅圖像。選擇長江入海口作為本文研究區域進行試驗，其坐標范圍為 26°N～34°N，120°E～128°E，觀測數據像素點的分辨率是1 km。而且，為了避免白天海洋表面熱氣等因素的影響，增加重構精度，研究中只選擇使用夜間的數據。

首先對最原始的下載的數據壓縮包進行解壓縮，然后利用seadas軟件導出每一張圖像的數據。由于系列圖像的數據值的經緯度坐標存在偏差，需要利用surfer軟件對數據進行網格化操作來校正坐標。另外，原始數據中包含有陸地，還需要剔除陸地值以減少數據的重構誤差。本文利用ArcGIS軟件做出研究區域的陸地模板，剔除每一天圖像中的陸地部分，就構成了可以進行處理的海表溫度的原始數據集。處理后的數據包含37幅含有云覆蓋的圖像，每幅圖像包含121×121 個像素，覆蓋區域是 26°N～34°N，120°E～128°E。

經過上述處理后，數據集中的云的平均覆蓋率高達83.23%，有些時刻的數據缺失率接近100%，亦即該天的遙感數據基本沒有，最高的也只不過是42.00%。與J MBeckers(2003)的數據預處理不同，本文沒有剔除云覆蓋率在95%以上的圖像。

3 實驗結果及其分析

3.1 基于SOM的EOF算法重構結果

利用該算法重構長江口附近海域的SST。從原始數據中選取五重交叉校正集作為蒙特卡羅交叉校正集，每個交叉校正集的數據占原始數據集的比例是5%，初始化方法是把每一時間維的平均值作為該列數據中的缺失值的代替值。預設定的最大EOF數是10（由于使用的數據集的缺失率很高，故選取的EOF數比較小，但是符合客觀情況）。

圖1是重構前后的對比圖[1]。其中，圖1A為2010年8月1日的SST，圖1B為2010年9月21日的SST，圖1C為2010年10月27日的SST；圖1D、1E和1F分別為對應的數據重構后的結果。9月到10月的時間段正是海表溫度逐漸降低的過程，從圖1E和圖1F可以看出，重構后的海溫分布圖很好地體現了這一海溫漸變的過程，同時海洋的流場形態能夠保持，未因重構后的插值而破壞。而更值得注意的是這3張原始圖像，每一幅圖的數據缺失率都是相當高的，分別是42.00%（本實驗中缺失率最低的）、81.39%和97.53%。尤其是在存在大面積空缺時，傳統的插值方法難以達到這樣的補缺效果，特別是圖1C，采用一般的重構算法無法進行補缺。而本文中提出的算法在首先經過SOM插值結果對EOF的初始化，再經過經驗正交函數來提取數據在空間及時間域上的物理特征，通過保留最佳的模態數可對缺失的物理觀測數據進行有效地重構，同時保持了數據的時空分布特征，而在本實驗中數據集的總體缺失率高達83.23%，相對缺失率比較低的數據集存在更多噪聲的情況下，也能夠達到這樣的效果。

3.2 模態分析[9]

圖1 重構前后對比圖

本實驗數據集的最佳模態數為2，兩個模態對原始數據中總方差的解釋比例分別為89.92%和2.93%。前2個模態已經解釋了原始數據中總方差的92.85%。最主要的2個空間模態如圖2所示。第1空間模態（圖2A）的SST在空間上大致呈現出西北方向海域至東南海域方向數值依次增加的趨勢，正好反映了該研究區域內溫度分布從西北到東南海域方向依次遞增的客觀情況，而且在30°N的水平方向上呈現出臺灣暖流與東海水團相遇的狀態。而這種情況符合該海域的水溫變化趨勢。第2空間模態（圖2B）還刻畫了江蘇沿岸流和長江入海口對海溫的影響。冬季陸地溫度比海洋溫度低，而沿岸水溫受大陸影響，溫度比鄰近海域偏低。

雖然SOM可以單獨重構數據，但本文給出的新算法重構出的交叉校正誤差更小，見圖3。

圖2 本試驗SST的2個最主要特征模態圖

圖3 SOM與SOM+EOF的對比

圖4 測試數據和重構后數據的散點圖

3.3 誤差分析

由于該區域空間變異較大，使用幾個空間模態來重構，而且該原始數據集的缺失率相當高，質量受到多方面的干擾而較差，使得重構后的數據集表現出來的特征不明顯，對空間變異大的地方解釋能力也不強。

在選擇測試數據和校正數據后，這些數據在原始數據集中都當做缺失值處理，而且原始數據中也包含噪聲，這樣勢必增加了數據處理時的數據的缺失率，又因為原本數據集的缺失率就很高（高達83.23%），所以在重構的結果中就反映出來了。但是，不可否認，本算法作用于該遙感數據集后（37幅圖像的數據集，總共有數據點909 577個，缺失數據點總共798 239個），重構的誤差為0.408 3℃，這個數值的確已經滿足需求（Alvera-Azcárate[5]的算法，其誤差為0.744 0℃）。測試數據和重構后數據的散點圖更好地說明了實驗結果，如圖4（由于實驗數據點數量很大，所以只選取了前100個）。

4 討論與展望

SOM通過尋找最優參考矢量集合來對輸入模式集合進行分類。每個參考矢量為一輸出單元對應的連接權向量。SOM能夠通過其輸入樣本學會檢測其規律性和輸入樣本相互之間的關系，并且根據這些輸入樣本的信息自適應調整網絡，使網絡以后的響應與輸入樣本相適應。SOM神經網絡通過學習能夠識別成組的相似輸入向量，使網絡層中緊鄰的神經元對相似的輸入向量產生響應。自組織映射神經網絡不但能學習輸入向量的分布情況，還可以學習輸入向量的拓撲結構，其單個神經元對模式分類不起決定性作用，要靠多個神經元的協同。

傳統意義上講，SOM和EOF都能各自的填充缺失數據，但是實驗表明聯合使用這兩種方法的效果更好。首先，SOM算法可以非線性估計。在這種意義上，即使數據集有復雜的非線性結構，SOM的編碼向量也能夠成功地捕捉數據的非線性特性。同時，該估計法是在低維的網格（本文是二維的）操作的，這很可能失去數據的內在信息。其次，EOF算法是使用奇異值分解SVD進行的線性估計。正因為如此，EOF算法將不能反映數據集的非線性結構，但是其估計的空間維數可以像輸入數據一樣高，而且該算法還是連續估計的。鑒于此，采取這兩種算法的聯合的新算法揚長避短，最終使得數據集的重構效果更好。

在下一步的研究工作中,如需要處理更大型的數據矩陣時,還需要做一些附加的優化工作,例如整個算法處理過程中，數據集的數據點的有效性的判斷以減少整個數據集的污染，SOM算法的訓練法則的改進以及EOF收斂準則的改善以提高整個程序的效率等。在將來的研究工作中，還會考慮利用該算法應用于其他一些海洋遙感要素進行重構，例如海水鹽度、懸浮物濃度、葉綠素濃度等。

[1]盛崢,石漢青,丁又專.利用DINEOF方法重構缺測的衛星遙感海溫數據[J].海洋科學進展,2009,27(2):243-249.

[2]BECKERSJ,RIXENM.EOF calculations and data fillingfromincomplete oceanographic datasets[J].Journal ofAtmospheric and Oceanic Technology,2003,20(12):1839-1856.

[3]SORJAMAAA,MERLINP,MAILLETB,et al.SOM+EOF for findingmissingvalues[C]//European Symposiumon Artificial Neural Networks.Bruges,Belgium:d-side publication.,2007:115-120.

[4]COTTRELL M,LETREMY P.Missingvalues:Processingwith the kohonen algorithm[C]//Dans ASMDA 2005 CD-ROMProceedings-ASMDA 2005,Brest,France,2005:489-496.

[5]ALVERA-AZCARATE A,BARTH A,RIXEN M,et al.Reconstruction of incomplete oceanographic data sets using Empirical Orthogonal Functions:Application tothe Adriatic Sea surface temperature[J].Ocean Modelling,2005,9(4):325-346.

[6]ALVERA-AZCARATE A,BARTH A,BECKERS J，et al.Multivariate reconstruction ofmissingdata in sea surface temperature,chlorophyll,and wind satellite fields[J].Journal ofGeophysical Research,2007,112:C03008.

[7]LENDASSE A,WERTZ V,VERLEYSEN M.Model selection with cross-validations and bootstraps-application to time series prediction with rbfn models[J].LNCSSpringer-Verlag,2003,2714:573-580.

[8]TOUMAZOU V,CRETAUX J.Using a Lanczos eigensolver in the computation of empirical orthogonal functions[J].Monthly Weather Review,2001,129(5):1243-1250.

[9]丁又專.衛星遙感海表溫度與懸浮泥沙濃度的資料重構及數據同化試驗[D].南京:南京理工大學,2009.

EOF Method to Reconstruct SST Remote Sensing Dataset Based on SOM

WANG Jian-le1,ZHU Jian-hua2,HE Yu-qing1,SONG Zhan-jie1
(1.Institute of TV and Image Information,Tianjin University,Tianjin 300072,China;2.National Ocean Technology Center,Tianjin 300112,China)

For the missing data of SST dataset,a method based on SOM and EOF algorithms to reconstruct the missing data was proposed,which could use the nonlinear estimation of SOM to reflect the nonlinear structure of dataset and could initialize the input of EOF algorithm utilizing the result of SOM to get rid of the sensitivity problem of EOF initialization to dataset.In the process of SVD,Lanczos operator was used to decompose matrix so as to enhance the efficiency of procedure.Monte-Carlo cross validation was introduced to assure the optimal mode EOF number and calculate the error of reconstruction.After using AQUA satellite remote sensing SST data,the result shows that this method can reconstruct the dataset with missing rate of 83.23%and have a high reconstruction precision.

remote sensing;SST;SOM algorithm;EOF algorithm;Lanczos operator

TP75

A

1003-2029（2012）01-0067-05

2011-05-24

國家自然科學基金資助項目(60872161)

王建樂（1986-），男，碩士生，主要研究方向為隨機信號處理和海態參數估計。

何宇清（1973-），男，博士，講師，主要研究方向為數字信號處理和數字圖像處理。