一種利用低秩矩陣填充技術恢復氣象數據的方法

臧 芳

(湖南大學信息科學與工程學院 湖南 長沙 410012) (湖南電子科技職業學院機械與電子信息工程分院 湖南 長沙 410205)

一種利用低秩矩陣填充技術恢復氣象數據的方法

臧 芳

(湖南大學信息科學與工程學院 湖南 長沙 410012) (湖南電子科技職業學院機械與電子信息工程分院 湖南 長沙 410205)

矩陣填充技術是繼研究壓縮感知技術之后又一備受關注的新技術，并已被應用于科學與工程等諸多領域。矩陣填充技術主要研究如何通過已知的部分矩陣數據來達到恢復整個矩陣的目的。氣象數據為天氣預報和開展各種氣象服務提供依據，具有空間相關性和時間穩定性，加上氣象數據矩陣的低秩性，為研究如何利用低秩矩陣填充技術對氣象數據的連續獲取提供了條件。最后通過真實氣象數據進行了仿真實驗，證明了該方法的有效性。

低秩性 氣象數據 矩陣填充 氣象預測

0 引 言

壓縮采樣[1]理論一經提出，就引起信息論[2-5]、圖像處理[6-9]、無線通信、大氣等領域的廣泛關注。其中，矩陣填充[10-22]是壓縮采樣領域的重要理論，主要作用是能利用已知矩陣的部分數據來恢復缺失數據。

氣象數據采集技術的發展與更新為當前氣象數據的采集、處理、預測提供了強有力的支持。氣象數據的采集通過無線傳感器完成，涉及地區廣，問題規模大。整個采集過程受采集器的精確度、存儲容量、處理速度、能耗等方面的影響。同時氣象數據又有時間與空間的相關性優勢，即地理位置相鄰的同一時段區域其天氣往往相差不多。因此，采集完整的氣象數據既不經濟也不現實。如果氣象數據所對應的矩陣低秩合理，那么氣象數據的完整性可以用矩陣的恢復技術來實現。

文中主要針對氣象數據的完整性研究提出了利用矩陣填充技術來實現的方法。最后采用了真實氣象數據進行實驗驗證，結果證明：文中提出的矩陣填充算法能較好地恢復氣象數據，并能有效解決傳感器存儲容量過大、能耗過大等問題。

1 相關研究

工作實際中碰到的矩陣大多具有低秩性或類似低秩性，而低秩性是該矩陣具有唯一解的必備條件[12]。某一缺失元素的矩陣，如果能在原有數據的基礎上根據矩陣的低秩結構來補全矩陣元素，稱此過程為矩陣恢復。

矩陣填充理論研究側重研究矩陣的可填充性。而僅當研究矩陣元素數量大于一定數才有可填充性。Candès等研究者[12,28]用大量實例證明了當采樣矩陣的奇異值及已知矩陣元素數量達到一定的標準時，絕大多數低秩陣可求凸優化問題而完好地還原矩陣。但經驗證該方法并不十分可行。

近年來，國內外許多學者針對該問題設計出了精確低秩矩陣填充重構算法[23-27]，取得了一系列可喜的成果。涉及精確低秩矩陣填充的重構算法除了奇異值閾值SVT(Singular Value Thresholding)算法[29-31]，還有SET算法[13]和OptSpace 算法[22]。

1.1 SVT算法

Cai[14]等提出了一種簡單的適用于較大規模矩陣填充的方法：SVT算法。利用SVT算法解矩陣恢復問題，可描述為：

min‖X‖*subject toPΩ(X)=PΩ(M)

算法1矩陣填充的SVT算法

1:初始化Y0=0，

2：while not converged do

3: Xk=Dτ(Yk-1),

4: Yk=Yk-1+δkPΩ(M-Xk).

5: end while

1.2 SET算法

SET(Subspace Evolution And Transfer)算法[13]是由Dai等提出的一種新算法。SET算法在已知秩的情況下采用查找列空間來匹配矩陣。主要有子空間演化(使用線性查找程序細化列空間)和子空間轉換構成。SET算法針對此問題在子空間轉換部分，設置用于檢測障礙的機制，并將列空間轉換位置。實踐證明，該算法對超低秩矩陣非常有效，缺點是收斂慢，在若干的迭代次數后，因檢測機制遠離列空間而易失效。

1.3 OptSpace 算法

針對矩陣恢復問題，Keshavan等提出了OptSpace算法[22]。其主要使用譜方法在格拉斯曼流形[32-33]上進行目標函數優化，可從隨機觀測中恢復低秩矩陣的算法。重點對目標矩陣進行盡可能恢復，在對修復后的矩陣進行投影時借助縮放因子來補充矩陣中檢測到的項數相對于原始矩陣而言的較小均值，最后運用格拉斯曼流形空間上的梯度下降法來求解獲得最小目標函數，即獲得矩陣恢復最小誤差值。文獻[34]給出了該算法可行的理論支撐。與原始矩陣填充區別在于：必須事先獲得秩。

與SET算法有著相同之處：兩者均作用于格拉斯曼流形上，且通過流形上的線性查詢得到目標值。不同之處在于：SET算法要求秩必須已知，利用列空間來獲得目標矩陣；而OptSpace算法行空間和列空間均需查找，增大了數據修復的難度，而且線性搜索不可確保收斂性，也是由于搜索路徑影響該算法得到最佳答案。

2 氣象數據特征分析

為證明氣象數據恢復能利用矩陣填充技術進行研究，首先對氣象數據的特性進行詳細分析。

2.1 平穩性

由實際氣象數據分析可知，同一站點在較短時間內氣象變化一般波動不大。在此，利用相鄰時間段的氣象數據差值分布來探究氣象數據在一段時間內的相對平穩性。用s表示氣象觀測站點號，t表示某個時間點，則相鄰間隔差值函數可表示為：

(1)

由于氣象數據量龐大，實驗只選取每年各季節的一段數據集。圖1描述了春夏秋冬各季度氣象數據相鄰間隔差值分布，橫軸表示歸一化的差值大小，縱軸表示累計分布百分比。從下面四個組圖中得知，接近1的相鄰間隔差值都保持在0.1之內，證明了相鄰時間段的氣象數據數值有較好的時間平穩性。

圖1 春夏秋冬氣象數據的相鄰間隔差值分析圖

2.2 低秩性

相鄰地區如株洲與長沙，經常能同時有晴天或同時下雨。這說明氣象數據除了具有時間上的平穩性，還具有空間上的相關性。為此，假設氣象數據構成的矩陣具有低秩性，如果能證明其具有低秩性這證明氣象數據的恢復適合用低秩矩陣填充技術實現。以下通過分析真實氣象數據矩陣的最大L個奇異值累計貢獻率來驗證其低秩性。

矩陣A∈Cm×n的奇異值分解成A=U∑VT。σi為矩陣A的奇異值，n為奇異值總個數，最大L個奇異值累積貢獻率：

(2)

由于數據集較大，實驗分別選取一年四季數據中各自一段數據集進行驗證。圖2描述了12天，288個測量時間點四季氣象數據矩陣的最大L個奇異值累積貢獻率。圖2的四個組圖說明最大的50個奇異值，基本蓄積了全部奇異值接近99%的能量。充分證明氣象數據具有很好的低秩性，具有能利用矩陣填充技術的條件。

圖2 四季氣象數據最大L個奇異值累積貢獻分布

3 問題的解決

3.1 數據采集

為了解決傳感器能量消耗過大的問題，就不能在每個時間點采集所有站點的氣象數據。氣象數據的時間與空間相關性使采集的部分數據對剩余部分利用矩陣填充技術進行推斷成為可能。Candes等提出伯努利模型[12]證明了矩陣元素可獨立采樣，但要求獲取的元素要滿足在矩陣的每一行和每一列都有分布。該實驗利用伯努利模型對氣象數據進行采集。假設氣象數據的矩陣為A∈Cm×n，m表示站點數目，n為總測量時間點。對應站點對應時間點的氣象數據測量值aij表示矩陣元素的大小。每一列數據表示某一時間點采集到全部站點的氣象數據。通過[0,1]投影矩陣決定在該時刻是否采集。則相應的采樣矩陣表示成：

Bm×n=P[A]Ω

(3)

采樣描述如圖3所示。

圖3 氣象數據的采集模型

3.2 數據預測

圖3中純深灰色區域為未采集氣象數據，利用填充技術得到如圖4(b)中深灰底白杠部分數據。

圖4 氣象數據的預測模型

計算數學基礎文獻[12]中指出：矩陣具有低秩性的情況下，根據矩陣恢復技術定義，可將矩陣恢復問題可轉化為如下：

minimizerank(X)subject toXij=Aij(i,j)∈Ω

(4)

式中的Aij為采集到的值，Ω為采樣矩陣元素下標的集合。由上述公式可以推斷：矩陣填充問題就是在保證矩陣低秩性的良好結構的前提下進行數據補全。填充數據時秩(rank)用核范數來替代。核范數定義為：

(5)

式中σl代表降序的第l個奇異值。由此，解決的優化問題由式(4)導出式(6)：

minimize ‖X‖*subject toXij=Aij(i,j)∈Ω

(6)

文獻[12]中指出，對于滿足不相關假設的秩為r并且大小為n1×n2的A低秩矩陣，測量矩陣和稀疏矩陣間不包含相關元素，當采樣數目達到一定時就可精確還原。但這類優化問題無法處理大型矩陣，如氣象數據。研究者針對此種凸優化問題給出了若干應對算法，如文中提到的SVT算法[25]、IALM算法[35-36]等、貪婪追蹤方法和基于格里斯曼流形空間的方法等。

3.3 連續預測

本實驗對于氣象數據的連續獲取，采用基于移動窗口形式來保持算法中被填充矩陣的大小不變來實現。首先對一段連續時間點均勻采集部分數據，獲得初始窗口。然后利用矩陣填充對未采集的數據進行恢復，得到完整數據集。在得到前一階段完整數據之后，同樣只采集下一時刻的部分數據，保持矩陣窗口大小不變，向前移動一個時刻點，重新組成采樣矩陣，然后利用矩陣填充恢復其他數據。如圖5所示。

圖5 氣象數據連續預測窗口模型

圖5中橫向1～10表示首次采集時間段，2～11表示下一個采集時間段。該過程重要步驟為在矩陣右側增添一列數據，然后去掉原有矩陣左側首列數據。

4 連續獲取數據算法

先確定初始窗口的大小，得到均勻采集部分數據后，應用矩陣填充技術獲得全部數據。基于移動窗口形式保持運算矩陣大小不變，然后更新采樣矩陣來獲取后續的氣象數據。描述如算法2。

算法2基于矩陣填充的氣象數據獲取

1.初始化：t=0，窗口列大小=n，采樣比例=v。

2.應用均勻采樣，選取n長的時間段，按照采樣比例，對數據進行連續采樣，獲取采樣矩陣A(t)。

3.矩陣填充算法應用于A(t)，獲得恢復矩陣X(t)。

4.t=t+1，進入下一窗口，按照采樣比例，對下一時刻的數據進行采樣，獲取新增一列的采樣位置集合C，去掉前一窗口第一列，重新組成新的采樣矩陣A(t)。

5.矩陣填充算法應用于A(t)，獲得恢復矩陣X(t)。

6.未結束轉到4，否則停止。

5 實 驗

5.1 實驗過程

本測試利用均方值誤差計算法來衡量預測氣象數據值與實際氣象數據值之間的誤差。如果誤差越小，說明該算法精度越高。

實驗數據采用2014年到2015年，每年中四個季度的氣象降雨量數據，通過不同數據采樣率(60%、50%)，對文中提到的SVT、SET、OptSpace三種算法進行性能權衡。通過柏努利模型對若干個站點均勻采集氣象數據，柏努利采樣模型，可以保證獲取的元素滿足在矩陣的每一行和每一列都存在分布，有助于矩陣填充算法對氣象數據的重構。實驗中還借助MATLAB軟件進行構圖，利用低秩矩陣填充算法求解。模型中設置n=24×12, 12代表12天的時間測量點。

5.2 實驗結果

圖6為采樣率為50%連續預測恢復誤差對比圖，恢復誤差都小于0.04，算法精度較高。圖7描述了60%的數據連續預測恢復誤差對比實驗結果。與圖6相比，圖7隨著采集數目增多，填充精度也有了較大的提高。因此可以推斷：如果充分隨機的采集一定量的氣象數據，是完全能夠填充其他氣象數據的。柏努利采樣模型，可以保證獲取的元素滿足在矩陣的每一行和每一列都存在分布，有助于矩陣填充算法對氣象數據的重構。測試結果充分證明了利用低秩矩陣填充技術恢復氣象數據的方法的實效性。

圖6 采樣率為50%的氣象數據預測誤差對比圖

圖7 采樣率為60%的氣象數據預測誤差對比圖

6 結 語

本文首先介紹了幾種精確低秩矩陣填充的重構算法，如適合于較大規模矩陣恢復的SVT算法、使用譜方法在格拉斯曼流形上進行目標函數優化的OptSpace算法和采用尋找列空間來匹配采樣矩陣，秩必須已知的SET算法。然后重點對幾種矩陣填充方法用于低秩氣象數據矩陣填充問題進行數據測試，最后選取不同參數實施矩陣重建，并驗證了矩陣重建效果。

從以上實驗得知：如果充分獲取一定比例的隨機氣象數據，可通過柏努利采樣模型恢復得到絕大部分其他氣象數據。柏努利采樣模型，可以保證獲取的氣數據在矩陣的每一行和每一列都有分布。根據氣象數據具有的時間與空間的相關性，氣象數據構成的矩陣還具有低秩特點，基于矩陣填充算法重構出氣象數據矩陣。

文中提出的矩陣填充算法能較好地恢復氣象數據，并能有效解決傳感器存儲容量過大、能耗過大等問題。該方法實用有效。

[1] Donoho D L.Compressed Sensing[J].Information Theory,Ieee Transactions on,2006,52(4):1289-1306.

[2] 葉蕾,郭海燕,楊震.基于壓縮感知重構信號的說話人識別系統抗噪方法研究[J].信號處理,2010,26(3):321-326.

[3] 季云云,楊震.基于自相關觀測的語音信號壓縮感知[J].信號處理,2011,27(2):207-214.

[4] 孫林慧,楊震.基于壓縮感知的分布式語音壓縮與重構[J].信號處理,2010,26(6):824-829.

[5] 高暢,李海峰,馬琳.面向內容的語音信號壓縮感知研究[J].信號處理,2012,28(6):851-858.

[6] 朱明,高文,郭立強.壓縮感知理論在圖像處理領域的應用[J].中國光學,2011,4(5):441-447.

[7] 王芳.基于壓縮感知的圖像處理與重構方法[J].光電子·激光,2012(1):196-202.

[8] 焦鑄.壓縮感知在圖像處理中的應用[D].西安理工大學,2011.

[9] 丁偉.基于壓縮感知的水下圖像處理[D].中國海洋大學,2013.

[10] Balzano L,Nowak R,Recht B.Online Identification and Tracking of Subspaces From Highly Incomplete Information[C]//Communication,Control,and Computing (allerton),2010 48th Annual Allerton Conference on,2010:704-711.

[11] Keshavan R H,Montanari A,Oh S.Matrix Completion From a Few Entries[J].Information Theory,Ieee Transactions on,2010,56(6):2980-2998.

[12] Candès E J,Recht B.Exact Matrix Completion via Convex Optimization[J].Foundations of Computational Mathematics,2009,9(6):717-772.

[13] Dai W,Milenkovic O.Set:an Algorithm for Consistent Matrix Completion[C]//Acoustics Speech and Signal Processing (icassp),2010 Ieee International Conference on,2010:3646-3649.

[14] Cai J,Candès E J,Shen Z.A Singular Value Thresholding Algorithm for Matrix Completion[J].Siam Journal on Optimization,2010,20(4):1956-1982.

[15] Candes E J,Plan Y.Matrix Completion with Noise[J].Proceedings of the Ieee,2010,98(6):925-936.

[16] Chen Y.Incoherence-optimal Matrix Completion[J].IEEE Transactions on Information Theory,2013,61(5):2909-2923.

[17] Geelen J F.Maximum Rank Matrix Completion[J].Linear Algebra and Its Applications,1999,288:211-217.

[18] Candès E J,Tao T.The Power of Convex Relaxation:Near-optimal Matrix Completion[J].Information Theory,Ieee Transactions on,2010,56(5):2053-2080.

[19] Ma S,Goldfarb D,Chen L.Fixed Point and Bregman Iterative Methods for Matrix Rank Minimization[J].Mathematical Programming,2011,128(1):321-353.

[20] Chen Y,Bhojanapalli S,Sanghavi S,et al.Coherent Matrix Completion[J].Proceedings of the 31st International Conference on Machine Learning,Beijing,China,2014:674-682.

[21] Recht B.A Simpler Approach to Matrix Completion[J].The Journal of Machine Learning Research,2011,12(4):3413-3430.

[22] Keshavan R,Montanari A,Oh S.Matrix Completion From Noisy Entries[C]//Advances in Neural Information Processing Systems,2009:952-960.

[23] 趙玉娟,鄭寶玉,陳守寧.矩陣填充及其在信號處理中的應用[J].信號處理,2015(4):423-436.

[24] 張瑋奇,張宏志,左旺孟,等.基于加權核范數最小化的矩陣填充模型[J].計算機科學,2015,42(7):254-257.

[25] 方茂中.關于矩陣填充和非負矩陣的研究[D].華東師范大學,2008.

[26] 韋仙.基于矩陣填充技術重構低秩密度矩陣[J].武漢工程大學學報,2015,37(2):72-76.

[27] 劉旭東,陳德人,王惠敏.一種改進的協同過濾推薦算法[J].武漢理工大學學報(信息與管理工程版),2010,32(4):550-553.

[28] Candès E J,Plan Y.Matrix Completion with Noise[J].Proceedings of the Ieee,2009,98(6):925-936.

[29] 陳峰峰.奇異值閾值算法在Netflix問題中的應用研究[D].清華大學,2011.

[30] 孟樊,楊曉梅,周成虎.遙感影像低秩信息的矩陣填充復原方法[J].測繪學報,2014,43(12):1245-1251.

[31] 王平.基于凸優化的矩陣重建問題算法的研究[D].海南師范大學,2014.

[32] Ngo T,Saad Y.Scaled Gradients on Grassmann Manifolds for Matrix Completion[C]//Advances in Neural Information Processing Systems,2012:1412-1420.

[33] Boumal N,Absil P A.Low-rank Matrix Completion Via Preconditioned Optimization on the Grassmann Manifold[J].Linear Algebra & Its Applications,2015,475:200-239.

[34] Wong R K W,Lee T C M.Matrix Completion with Noisy Entries and Outliers[D].University of California,2016.

[35] Lin Z,Chen M,Ma Y.The Augmented Lagrange Multiplier Method for Exact Recovery of Corrupted Low-rank Matrices[J].Arxiv Preprint Arxiv:1009.5055,2010.

[36] Goldfarb D,Ma S,Scheinberg K.Fast Alternating Linearization Methods for Minimizing the Sum of Two Convex Functions[J].Mathematical Programming,2013,141(1):349-382.

METHODOFRESTORINGMETEOROLOGICALDATAUSINGLOW-RANKMATRIXFILLINGTECHNIQUE

Zang Fang

(CollegeofComputerScienceandElectronicEngineering,HunanUniversity,Changsha410012,Hunan,China) (InstituteofMechanicalandElectronicInformationEngineeringBranch,HunanVocationalCollegeofElectronicandTechnology,Changsha410205,Hunan,China)

Matrix filling is another attractive new research field after compressive sensing. It has been applied in the fields of science and engineering. Matrix filling principally studies how to restore the whole matrix elements by the known partial matrix elements. The meteorological data have the spatial correlation，time stability, and low rank of meteorological data matrix. All these advantages provide the continuous acquisition of meteorological data for researching how to use the low rank matrix filling technique. Finally, the simulation experiments are performed by real weather data, and the effectiveness of the method is proved.

Low-rank Meteorological data Matrix filling Meteorological forecast

TP391

A

10.3969/j.issn.1000-386x.2017.09.063

2016-11-11。臧芳，講師，主研領域：軟件工程。