盧小鵬，葉慶衛(wèi)，呂翠蘭

（1.寧波海洋環(huán)境監(jiān)測中心站，浙江 寧波315012；2.寧波大學(xué) 信息科學(xué)與工程學(xué)院，浙江 寧波315211）

0 引言

當(dāng)前常用的潮流分析及預(yù)報(bào)方法［1－3］以調(diào)和分析為主，但是在潮流預(yù)報(bào)的實(shí)踐檢驗(yàn)中發(fā)現(xiàn)調(diào)和分析方法用于潮流復(fù)雜海域的潮流預(yù)報(bào)誤差較大。仔細(xì)分析認(rèn)為造成這樣的原因主要是該種海域的潮流除了受潮汐力的影響外，還受到其他因素的影響，如周邊地形的曲折岸線、島嶼的阻攔、海底的復(fù)雜水深、以及其他水道潮水的補(bǔ)充等，這都會(huì)造成這些水域的潮流沒有外海潮流表現(xiàn)的周期性完美。在這些復(fù)雜水域進(jìn)行調(diào)和分析計(jì)算時(shí)，雖然在一些高頻分潮上會(huì)有所體現(xiàn)，但從預(yù)報(bào)結(jié)果來看并不能很好地表現(xiàn)水域的實(shí)際潮流情況，對非潮汐力干擾較大的潮流不能準(zhǔn)確描述，這樣就表現(xiàn)出誤差較大的情況，特別是一些較大的極值不能很好把握。

針對這種情況，本文引入信號稀疏表示［4－8］來克服部分非線性效應(yīng)、噪聲因素等對預(yù)報(bào)的干擾?？紤]從實(shí)測數(shù)據(jù)本身出發(fā)，將這些數(shù)據(jù)看成一系列的時(shí)間序列信號，選用AR模型即自回歸模型，分析這些數(shù)據(jù)的歷史關(guān)聯(lián)性。從表面來看，這只是對實(shí)測信號進(jìn)行分析，好像完全顛覆了調(diào)和分析的原理，與傳統(tǒng)方法毫無聯(lián)系。實(shí)際分析會(huì)發(fā)現(xiàn)該方法對傳統(tǒng)調(diào)和分析原理也存在一定的繼承，因?yàn)樵汲绷鲾?shù)據(jù)本身就包含這個(gè)調(diào)和原理，再通過AR模型對這些潛在信息進(jìn)行自回歸研究，尋找信號自身的自適應(yīng)性，尋求信號間的變化規(guī)律。所以不管被預(yù)測的系統(tǒng)受到多少因素的影響而產(chǎn)生復(fù)雜信號，都能夠通過建立一種較為通用的模型來自適應(yīng)獲取這些因素的規(guī)律，并進(jìn)行預(yù)測。

AR模型方法在經(jīng)濟(jì)學(xué)、信息學(xué)以及自然現(xiàn)象的預(yù)測上已廣泛運(yùn)用，也是氣象要素序列常用的分析方法［9］，在水文預(yù)報(bào)中主要用于洪峰預(yù)測，但涉及潮流方面的運(yùn)用分析研究還較少。本文在常規(guī)AR模型基礎(chǔ)上引入稀疏優(yōu)化進(jìn)行了改進(jìn)，提出一種稀疏AR模型。把常規(guī)的AR模型看成是一種過完備稀疏表示，采用重復(fù)多次抽樣的方法構(gòu)建多個(gè)欠定方程組，然后利用稀疏優(yōu)化求解方法（如 OMP、BP等［4－8］）求解欠定方程組，獲得稀疏的AR系數(shù)解。多次優(yōu)化求解不同的稀疏抽樣欠定方程組后進(jìn)行平均來消除噪聲因素的干擾，獲得穩(wěn)定精確的AR系數(shù)解。

1 基于稀疏AR建模的基本原理

1.1 基于 AR 模型的過完備稀疏基構(gòu)造［10－11］

假設(shè)信號符合AR模型：

上式中xk為時(shí)間序列參數(shù)，al為相關(guān)系數(shù)。由（1）式可以看出AR模型具有較強(qiáng)的時(shí)間相關(guān)性，即第k時(shí)刻的信號xk只與它前p個(gè)時(shí)刻的信號xk－l（l＝1，2，…，p）有關(guān)，與其他時(shí)刻的信號和激勵(lì)無關(guān)，則第k時(shí)刻的信號xk可以由其前p個(gè)時(shí)刻的信號線性表示。

對某一信號采樣了n點(diǎn)，因?yàn)樾盘柗螦R模型，則可以得到方程組：

將方程組（2）用矩陣形式表示：

令矩陣

其中我們把矩陣Z定義為過完備稀疏基，并且要求np＞p。因?yàn)樾盘柗螦R模型，方程（4）得出的解→a最多也只有p個(gè)非零元素，即→a最壞的情況是p為稀疏的。而n－p＞p，則向量→x′∈Rn－p可以由矩陣Z的列向量稀疏表示。所以我們可定義矩陣Z為過完備稀疏基。

1.2 稀疏表示

因?yàn)榉螦R模型［10－11］的信號具有較強(qiáng)的時(shí)間相關(guān)性，則第k時(shí)刻信號xk可以由其前p時(shí)刻的信號稀疏表示。根據(jù)這一特點(diǎn)，可以構(gòu)建一個(gè)矩陣Φ，用OMP算法求得信號的稀疏系數(shù)→a。

首先從過完備稀疏基中隨機(jī)抽取m行，其中m＜p，得到欠定方程組：

其中i1，i2，…，im是從過完備稀疏基Z中任意選取的m行的行號，并令：

因?yàn)榉匠探M（5）是欠定方程組，有無窮多個(gè)解，我們需要的解是最稀疏解，則方程組可以轉(zhuǎn)換為求：

其中：ε為實(shí)數(shù)。因?yàn)槭剑?）的求解是NP－h(huán)ard問題，若Φ滿足約束等距性，并且δ＜－1（RIP，即約束等距性，就是對任意K－稀疏向量s→，觀測矩陣Φ滿足式1－，其中δ∈ （0，1）），則可以將式（6）轉(zhuǎn)換成求?1最小范數(shù)，即式：

求?1最小范數(shù)下最優(yōu)化問題，本文采用OMP算法求解式（7）得到稀疏系數(shù)→a。因?yàn)棣凳请S機(jī)抽取過完備稀疏基Z的m行構(gòu)造的m行p列矩陣，所以每次選取后構(gòu)造的結(jié)果和求的稀疏系數(shù)結(jié)果都不穩(wěn)定，為了穩(wěn)定實(shí)驗(yàn)結(jié)果，可以重復(fù)上述步驟N次，對最后求得的N個(gè)稀疏系數(shù)（j＝1，2，…，N）進(jìn)行求和平均得到→a＾。

用過完備稀疏基Z乘以最終求得的稀疏系數(shù)→a＾可以重構(gòu)第p＋1到n的采樣信號。即：

而信號的向前預(yù)測一點(diǎn)的計(jì)算表達(dá)式由下式表示：

向前預(yù)測多點(diǎn)可以多次調(diào)用式（9）來獲得。其計(jì)算流程如圖1所示。

圖1 基于稀疏AR建模的信號模擬算法流程圖Fig.1 The flow chart of signal imitation based on sparse AR model

2 實(shí)驗(yàn)驗(yàn)證

這里的主要研究對象是潮流信號，為了更好地驗(yàn)證AR模型在潮流預(yù)報(bào)中的優(yōu)越性，建議選擇在流況復(fù)雜區(qū)域采集潮流信號。另外為了保證該模型的可靠性和穩(wěn)定性，以及實(shí)驗(yàn)和驗(yàn)證的需要，要求數(shù)據(jù)的時(shí)間序列長度不可太短，要大于1 000以上，實(shí)驗(yàn)選用的采樣時(shí)間為10min。因?yàn)槌绷鲾?shù)據(jù)是矢量，首先需將其分為北分量和東分量：

2.1 建立模型

這里選用MATLAB進(jìn)行實(shí)驗(yàn)仿真，建立AR自回歸模型和計(jì)算處理算法公式。以時(shí)間序列潮流的北分量或東分量信號｛χt｝（t＝1，2，…，n）為實(shí)驗(yàn)數(shù)據(jù)，建立多個(gè)自回歸AR模型?？紤]信號本身可能存在干擾，這里以小波去噪對原始信號進(jìn)行預(yù)處理，采用的是啟發(fā)式閾值函數(shù)形式（heursure），并設(shè)定為軟門限閾值（s）處理，閾值門限是由小波去噪函數(shù)自動(dòng)確定的。只根據(jù)第一層小波分解（sln）稀疏估計(jì)噪聲水平，用sym8小波對信號作5層分解完成對信號的去噪。

在AR模型基礎(chǔ)上，運(yùn)用randn函數(shù)建立一個(gè)標(biāo)準(zhǔn)正態(tài)分布的矩陣作為測量矩陣，并進(jìn)行歸一化處理，目的是使該測量矩陣限定于一定范圍，保證矩陣運(yùn)行的收斂。然后根據(jù)壓縮感知理論建立這個(gè)稀疏求解模型，以上述歸一化處理矩陣為測量矩陣，將AR模型中建立的矩陣為基矩陣，與向量信號一起建立等式。

然后求解這個(gè)基矩陣的系數(shù)，這里選用OMP正交匹配追蹤算法求解。通過該算法獲取稀疏矩陣中與向量的相關(guān)度最大的列數(shù)，再通過多次循環(huán)，逐次求解與殘差向量相關(guān)度最大的列數(shù)，最終求得這個(gè)稀疏矩陣的系數(shù)。

2.2 實(shí)驗(yàn)流程

實(shí)驗(yàn)以寧波光明碼頭某點(diǎn)的5m層的潮流為例，該碼頭位于寧波市北侖區(qū)穿山半島北岸，該水域潮流漲落主要受螺頭水道影響，因受局地凹槽地形因素影響較大，故流況很是復(fù)雜。

下面運(yùn)用前文所建立的模型對該站數(shù)據(jù)進(jìn)行實(shí)驗(yàn)驗(yàn)證（圖2），實(shí)驗(yàn)過程中分別對模型階數(shù)M選用了500、400、350、300、250和200，稀疏度K 分別為200、150、130、100、90、80、70、60和50，重復(fù)求解次數(shù)Dn從2到10分別進(jìn)行調(diào)整，并通過對預(yù)報(bào)結(jié)果的方差進(jìn)行比對，尋找最佳參數(shù)。

圖2 光明站點(diǎn)AR模型驗(yàn)證Fig.2 AR model verification at Guangming station

從上述實(shí)驗(yàn)結(jié)果來看，在不考慮階數(shù)、稀疏度和求解次數(shù)變化的條件下，運(yùn)用該計(jì)算方法的北分量的預(yù)報(bào)均方差多數(shù)控制在9～11cm／s之間，出現(xiàn)大于11cm／s的值很少，大于11cm／s的情況出現(xiàn)過11次。因?yàn)槭芨咚闺S機(jī)測量矩陣的影響，該方差會(huì)受隨機(jī)矩陣的變動(dòng)存在變化，但總體來看，預(yù)報(bào)均方差值較為穩(wěn)定。

通過計(jì)算會(huì)發(fā)現(xiàn)，稀疏度越小，計(jì)算速度就越快，當(dāng)然求解次數(shù)越少，計(jì)算速度也越快。從階數(shù)分析，在階數(shù)變化中，方差結(jié)果存在變化，但變化不顯著，基本也在9～11cm／s之間變動(dòng)，也就是說回報(bào)結(jié)果基本穩(wěn)定可靠。

3 AR模型與調(diào)和分析方法對比

AR模型預(yù)報(bào)方法是從“時(shí)域”角度進(jìn)行分析，而調(diào)和分析方法是從“頻域”角度進(jìn)行分析，與前者分析的角度截然不同。調(diào)和分析方法將潮流信號看成是由多個(gè)不同振幅、不同頻率、不同相位的正弦函數(shù)組成，該方法中潮流信號的長度及潮流信號的時(shí)間間隔將影響到函數(shù)頻率的選擇。

運(yùn)用調(diào)和分析方法對光明站點(diǎn)進(jìn)行回報(bào)（圖3），回報(bào)均方差結(jié)果為15.207 3cm／s，比上述結(jié)果大4～5cm／s。僅從這兩者結(jié)果來看，前者優(yōu)于后者，可以適用。再對其他站點(diǎn)資料運(yùn)行實(shí)驗(yàn)，回報(bào)均方差結(jié)果也顯示兩者相差不明顯，AR模型優(yōu)于調(diào)和分析。

圖3 光明站點(diǎn)調(diào)和模型預(yù)報(bào)驗(yàn)證Fig.3 Harmonic forecast model verification at Guangming station

從細(xì)節(jié)圖可知（圖4），兩種方法對極值點(diǎn)的預(yù)報(bào)，基于AR模型的預(yù)報(bào)方法明顯優(yōu)于調(diào)和分析方法，前者更遵從于實(shí)測數(shù)據(jù)，對極值潮位特別是多峰值的把握更為準(zhǔn)確，而調(diào)和分析方法對多峰值則是進(jìn)行光滑處理，從圖形上看是更為漂亮，但對實(shí)際情況的描述則是略顯弱勢。對其它站點(diǎn)的實(shí)驗(yàn)也驗(yàn)證了這一特點(diǎn)。

圖4 光明站點(diǎn)AR模型（a）和調(diào)和模型（b）驗(yàn)證細(xì)節(jié)圖Fig.4 Detailed drawing of AR model（a）and harmonic model（b）verification at Guangming station

4 小結(jié)

通過多個(gè)站點(diǎn)的實(shí)驗(yàn)說明，基于AR模型的潮流預(yù)報(bào)方法具有一定的可行性，它與調(diào)和分析方法是對潮流分析的不同角度的解釋，但它對潮流的描述在保留天體因素信號的基礎(chǔ)上，更遵從于客觀事實(shí)，特別是對潮流出現(xiàn)多峰值的把握更為準(zhǔn)確。因此對于流況復(fù)雜水域，如回流明顯的區(qū)域，運(yùn)用基于AR模型的潮流預(yù)報(bào)方法有明顯的優(yōu)越性。

（References）：

［1］FANG Guo－h(huán)ong，WEI Ze－xun，WANG Yong－gang.Development of tide and tidal current regional prediction in China［J］.Advances in Earth Science，2008，23（4）：331－336.

方國洪，魏澤勛，王永剛.我國潮汐潮流區(qū)域預(yù)報(bào)的發(fā)展［J］.地球科學(xué)進(jìn)展，2008，23（4）：331－336.

［2］CHEN Zong－yong.Tidology［M］.Beijing：Science Press，1980.

陳宗鏞.潮汐學(xué)［M］.北京：科學(xué)出版社，1980.

［3］FANG Guo－h(huán)ong，ZHENG Wen－zhen，CHEN Zong－yong，et al.A－nalysis and prediction of tides and tidal currents［M］.Beijing：China Ocean Press，1986.

方國洪，鄭文振，陳宗鏞，等.潮汐和潮流的分析與預(yù)報(bào)［M］.北京：海洋出版社，1986.

［4］ZHAO Yi－jiu.Study on compressive sampling and recovery algorithm of sparse analog singnal［D］.Chengdu：University of Electronic Science and Technology of China，2012.

趙貽玖.稀疏模擬信號壓縮采樣與重構(gòu)算法研究［D］.成都：電子科技大學(xué)，2012.

［5］PATI Y C，REZAIIFAR R，KRISHNAPRASAD P S.Orthogonal matching pursuit：Recursive function approximation with applications to wavelet decomposition［R］.Piscataway，USA：IEEE，1993：40－44.

［6］DONOHO D.Compressed sensing［J］.IEEE Transactions on Information Theory，2006，52（4）：1 289－1 306.

［7］TSAIG Y，DONOHO D L.Extensions of compressed sensing［J］.Signal Processing，2006，86（3）：549－571.

［8］YE Qing－wei，SUN Yang，WANG Xiao－dong，et al.An improve LLE algorithm with sparse constraint［J］.Journal of Computational and Theoretical Nanoscience，2013，10（12）：2 872－2 876.

［9］HUANG Jia－you.Meteorological statistical analysis and forecasting method［M］.Beijing：Meteorological Press，2004：3.

黃嘉佑.氣象統(tǒng)計(jì)分析與預(yù)報(bào)方法［M］.北京：氣象出版社，2004：3.

［10］NIE Shu－yuan.The historical process of analysis on Yule’s established AR（P）model［J］.Statistics and Decision，2011，237（3）：4－7.

聶淑媛.尤爾建立時(shí)間序列線性自回歸AR（P）模型的歷史過程探析［J］.統(tǒng)計(jì)與決策，2011，237（3）：4－7.

［11］YU Chun－mei.Review on sparse optimization algorithms［J］.Computer Engineering and Applications，2014，50（11）：210－217.

于春梅.稀疏優(yōu)化算法綜述［J］.計(jì)算機(jī)工程與應(yīng)用，2014，50（11）：210－217.