一種改進的局部線性回歸估計器及其在雷達高度計海況偏差估計中的應用

蔣茂飛 許 可 劉亞龍 王 磊

?

一種改進的局部線性回歸估計器及其在雷達高度計海況偏差估計中的應用

蔣茂飛①②③許 可*①②劉亞龍④王 磊①②

①(中國科學院微波遙感技術重點實驗室 北京 100190)②(中國科學院國家空間科學中心 北京 100190)③(中國科學院大學 北京 100049)④(國家海洋局煙臺海洋環境監測中心站 煙臺 264000)

在建立雷達高度計海況偏差(Sea State Bias, SSB)非參數模型時，通常會用到局部線性回歸(Local Linear Regression, LLR)估計器，而傳統的局部線性回歸估計器涉及高維矩陣運算，當建模的數據量較大時，估計海況偏差需要大量的時間，從而使得非參數估計方法很難用于高維海況偏差模型。該文提出一種改進的局部線性回歸(Improved Local Linear Regression, ILLR)估計器，可以避免傳統的LLR估計器所需的高維矩陣運算，在不影響海況偏差估計結果的條件下，將局部線性回歸估計器獲取加權函數的時間復雜度由降低為，從而大幅地降低估計海況偏差所需的時間，為實現高維非參數海況偏差模型的實時運算奠定了基礎。

雷達高度計；海況偏差；非參數估計；LLR估計器

1 引言

全球海平面上升已越來越引起人們的關注，雷達高度計的一個重要應用就是測量平均海面高度[1, 2]。但是由于海況偏差的存在，高度計測量得到的平均海面會低于實際的平均海面，必須對它進行校正。海況偏差主要是由于海面的非高斯分布特性引起，包括電磁偏差和偏斜度偏差[3]。利用現有方法得到的海況偏差的不確定度能夠達到2 cm[4]，已經成為Jason系列測高衛星最大的誤差源[5, 6]。

目前，業務化運行的雷達高度計的海況偏差校正主要采用以高度計測得的風速和有效波高作為輸入的非參數SSB模型[6,7]。非參數SSB模型的構建通常需要用到局部線性回歸估計器，而傳統的局部線性回歸估計器需要大量的矩陣運算，從而消耗大量的運算時間，極大地限制了非參數估計在高維SSB模型中的應用。

本文提出了一種改進的局部線性回歸估計器，該估計器避免了高維矩陣運算，將從局部線性回歸估計器獲取加權函數的時間復雜度由降低為。數值實驗結果表明，在其它條件都相同的條件下，和傳統的局部線性回歸估計器相比，改進的局部線性回歸估計器不會影響SSB估計的結果，但可以大幅地減低估計SSB所需的時間。

2 局部線性回歸估計器及其在海況偏差非參數估計中的應用

GASPAR等人[8]在1998年首次將非參數估計的方法用于SSB模型的構建，和參數模型[4, 9]相比，非參數模型的精度更高，從而被更廣泛地應用于高度計測量海面高度的海況偏差校正。用于非參數SSB模型的估計器主要有Nadaraya-Watson (NW)估計器[8]和局部線性回歸估計器[10]，NW估計器在接近邊界處會出現較大偏差，而局部線性回歸估計器則可以很好解決這一問題[11,12]，所以局部線性回歸估計器得到了更廣泛的使用。

2.1 局部線性回歸估計器

(2)

(4)

最小化式(4)，可以得到

(6)

(8)

2.2局部線性回歸估計器在海況偏差非參數估計中的應用

通常使用交叉點處的海面高度不符值來建立SSB模型，衛星繞地球運行一周的軌跡包括兩個弧段(pass)，分別為上升軌道和下降軌道，上升軌道和下降軌道的交點就稱為交叉點，而交叉點海面高度不符值是指上升軌道和下降軌道在交叉點處測得的海面高度之差[13]。

(9)

(11)

將式(12)轉化為矩陣的形式，可以得到

(13)

(15)

在估計SSB時，可制作一個SSB估計值關于風速、有效波高SWH的2維查找表，再通過插值的方法得到任意高度計測量點的SSB估計值。為了得到SSB估計值的查找表，通常是先利用一個cycle的數據計算查找表中每個網格點的輸入變量所對應的SSB估計值，再對多個cycle得到的結果進行平均。每個cycle的交叉點的有效個數一般可達5000~7000個，如果將所有交叉點都用于建模，那么式(15)中矩陣求逆需要大量的運算，消耗大量計算時間。在文獻[8]的模型中，從每個cycle隨機抽取500個有效的交叉點來建立SSB估計模型。為了避免了式(15)中矩陣求逆運算，GASPAR等人[10]在2002年利用LSQR算法來求線性方程式(14)的解，并將每個cycle的所有交叉點數據都用于SSB的建模。LSQR算法是PAIGE等人[14]在1982年提出的一種解稀疏線性方程組的迭代方法，它被用于求線性方程組的解或者使最小的解。為了進一步提高求解線性方程式(14)的速度，本文使用了LSMR算法來解方程式(14)，LSMR算法是FONG等人[15]在2010年提出的一種解線性方程組的迭代算法。和LSQR算法一樣，LSMR算法也是用于求線性方程組的解或者使最小的解，但和 LSQR算法相比，使用LSMR算法可以更快和更安全的達到收斂條件。

2.3改進的局部線性回歸估計器

非參數SSB模型最終得到的是一個關于風速和有效波高SWH的2維查找表，需要計算查找表中每個網格點的輸入所對應的SSB估計值。雖然利用LSQR和LSMR算法來解線性方程式(14)可以避免式(15)式矩陣求逆運算，但式(11)中的加權函數和式(12)中加權矩陣中的元素都是通過局部線性回歸估計器得到，從式(5)和式(8)也可以看出通過局部線性回歸估計器得到加權函數時也需要用到高維矩陣運算，本文希望避免高維矩陣運算，在對結果影響不大的條件下，能大大降低運算時間。

對于傳統的利用風速()和有效波高(SWH)來建立SSB估計模型時，,在的鄰域，回歸函數的局部線性模型為

(17)

其中，,,和可以表示為

(19)

(21)

所以，權值函數可以表示為

(23)

將式(23)替換到式(8)，就可以將ILLR估計器用于非參數SSB模型的構建，和LLR估計器相比，ILLR估計器避免了高維矩陣運算。

由式(8)和式(23)可知，兩種估計器得到的權函數有所區別，ILLR估計器和LLR估計器在本質上是相同的，都是通過對式(3)求偏導，并使其偏導為零得到，如果將式(5)重寫為，其展開后就得到了式(18)，所以只要不是奇異矩陣，ILLR估計器和LLR估計器得到的結果就完全相同。在LLR估計器中，由于是對角矩陣，和在運算過程中會產生很多結果為零值,它們會繼續參與和的矩陣運算，而在ILLR估計器中，LLR估計器中和在計算過程中產生的那些中間結果為零的值不會參與式(18)，式(9)的運算，這就大大節約了計算時間。從算法時間復雜度方面分析，在LLR估計器中，權函數由得到，是的矩陣，是的矩陣，是的矩陣，是的矩陣，是的矩陣，所以和相乘需要次乘法運算，和相乘需要次乘法運算，矩陣求逆需要次乘法運算，和相乘需要次乘法運算，和相乘需要次乘法運算，所以利用LLR估計器獲得權函數的時間復雜度為。由式(18)可知，計算,,和分別需要進行, 2, 2和次乘法運算，是的矩陣，所以利用ILLR估計器獲得權函數的時間復雜度為。

3 數值實驗

3.1數據來源和處理

本文使用Jason-2高度計的地球物理數據集(GDR)產品，Jason-2衛星于2008年6月發射，它是TOPEX/POSEIDON, Jason-1測高衛星的后繼星，被公認為目前精度最高的測高衛星。本文在建立非參數SSB模型時使用了cycle 51~100 共50個cycle的數據，而在評估ILLR估計器對SSB估計結果的影響時使用了2009(cycle19~cycle55), 2010 (cycle56~cycle92)和2011(cycle93~cycle128)共3年的數據。

建立SSB模型需要用到風速、有效波高、交叉點處未經SSB校正的海面高度不符值。計算海面高度所需的各種校正項都需要插值到交叉點處，本文采用3次樣條插值的方法，為了得到更好的插值效果，只有在交叉點兩側各有4個連續的測量點時得到的交叉點才有效，這樣可以保證較好的高度計數據質量。

3.2核函數和帶寬的選擇

核函數和帶寬的選擇對于非參數SSB模型非常重要。目前，最常用的核函數有高斯核函數和Epanechnikov核函數，而用于非參數SSB模型的帶寬包括全局帶寬和局部帶寬兩種，全局帶寬是指對于所有的輸入，風速的帶寬和有效波高的帶寬都為固定值，而局部帶寬是指和會根據數據密度自動進行調節。

和高斯核函數相比，Epanechnikov核函數可以將線性系統式(14)中矩陣轉換稀疏矩陣，從而大大提高了計算效率，但在數據稀疏區域估計的效果比較差。和全局帶寬相比，局部帶寬可以根據數據的密集程度調節帶寬，從而在數據稀疏區域優勢較為明顯。本文的數值實驗只是為了評估改進的局部線性估計器的計算效率以及對估計結果的影響，所以對核函數和帶寬的選擇沒有太多的要求。本文采用了兩種核函數和帶寬的組合，第1種是采用高斯核函數和全局帶寬，對應于文獻[10]的模型，第2種是Epanechnikov核函數和局部帶寬，對應于文獻[12]的模型。

3.3實驗結果

3.3.1 LLR估計器和ILLR估計器估計SSB效率比較

首先比較利用LLR估計器和ILLR估計器估計SSB所用的時間，表1給出了利用高斯函數和全局帶寬得到的結果，隨著每次用于SSB估計的交叉點數的增大，LLR估計器里的矩陣的維度就相應地增大，所需的計算時間也會明顯增大。ILLR估計器中矩陣的維度不會隨著的增大而改變，雖然其所需的時間也會隨著的增大而增大，但其增大的幅度相對較小。表2給出了利用Epanechnikov核函數和局部帶寬估計SSB所用的時間，正常情況下，Epanechnikov核函數可以將線性系統式(14)中的矩陣轉換為稀疏矩陣，從而提高計算效率，降低估計SSB所需的時間，但由于采用了局部帶寬，在估計SSB時，對于每一輸入，都需要先確定其所在網格所對應的帶寬，這樣就會增加計算量，從而使得計算時間增大。和表1一樣，隨著的增大，利用LLR估計器估計SSB所需的計算時間大幅增加，而利用ILLR估計器估計SSB所需的計算時間增幅相對較小。從表1和表2可以看出，和LLR估計器相比，利用ILLR估計器大可以大幅度地降低估計SSB所需的時間，越大，改進的效果越明顯。

表1利用高斯核函數和全局帶寬估計SSB所用的時間

表2利用Epanechnikov核函數和局部帶寬估計SSB所用的時間

非參數SSB模型最終得到的是一個關于和SWH的查找表，需要計算每個查找表中每個網格點所對應的SSB估計值，通過和SWH只能得到有限的SSB的信息，將新的參量(如波周期參數)引入SSB模型的構建已經成為SSB估計的一個新的研究方向[6, 16, 17]。在建立高維非參數SSB模型時，最終得到的結果將是一個高維查找表，其所對應的網格點數將遠遠多于2維查找表的網格點數，使用LLR估計器所需消耗的時間將非常巨大，并且隨著維度的增加，必須要有更多的數據參與模型的構建才能保證模型的有效性，這對于LLR估計器是巨大的挑戰，而ILLR估計器則可以很好地解決這一問題。

3.3.2 LLR估計器和ILLR估計器估計SSB結果比較

圖1是利用高斯核函數和全局帶寬得到的結果，根據數據的密集程度將數據劃分為兩部分，圖中陰影部分表示高度計測量值個數小于100的區域，代表數據稀疏區域，而另外一部分則代表數據密集區域。圖1(a)和圖1(b)的等值線分別表示利用LLR估計器和ILLR估計器得到的SSB查找表在(U, SWH)平面的分布，等值線的單位為cm。很明顯，利用兩種估計器得到的SSB隨和SWH的變化趨勢完全相同，在相同的條件下，SSB的幅值都隨著SWH的增大而增大；而在SWH相同的條件下，隨著的增大，SSB的幅值都是先增大，而當增大到一定值時，SSB的幅值都逐漸減小。圖1(c)的等值線表示利用LLR估計器和ILLR估計器得到的SSB之差在(,SWH)平面的分布，可以看出，在數據密集區域，兩種估計器得到的SSB之差大約在量級，基本上可以忽略。

圖1 利用高斯核函數和全局帶寬得到的結果在(U, SWH)平面的分布

圖2是利用Epanechnikov核函數和局部帶寬得到的結果，圖2(a)和圖2(b)的等值線分別表示利用LLR估計器和ILLR估計器得到的SSB查找表在(,SWH)平面的分布，等值線的單位為cm，同圖1一樣，利用兩種估計器得到的SSB隨和SWH的變化趨勢也是相同的，由于使用了Epanechnikov核函數，在低風速高海況和高風速低海況的數據稀疏區域，由于使用到的數據比較少，得到的結果不如圖1中利用高斯核函數得到的結果平滑。圖2(c)的等值線表示利用LLR估計器和ILLR估計器得到的SSB之差在(,SWH)平面的分布，在數據密集區域，利用兩種估計器得到的SSB之差大約在cm量級，可以認為利用兩種估計器得到的SSB基本相同。

圖2 利用Epanechnikov核函數和局部帶寬得到的結果在(U, SWH)平面的分布

圖3給出了利用LLR估計器和ILLR估計器得到的SSB估計值的散點圖，圖中的結果是通過將兩種估計器得到的SSB查找表分別用于2011年Jason-2高度計SSB的估計得到，可以看出，兩種估計器得到的SSB估計值一致性非常好。再將兩種估計器得到的SSB的查找表分別用于2009~2011年Jason-2高度計SSB的估計，分別利用每年的數據進行統計，表3給出了統計結果，當利用高斯核函數和全局帶寬時，兩種估計器得到的SSB估計值之差的最大值在量級，平均值在量級；而當利用Epanechnikov核函數和局部帶寬時，兩種估計器得到的SSB估計值之差的最大值在cm量級，平均值在量級。所以不管是哪種核函數和帶寬的組合，ILLR估計器對SSB的估計結果的影響都可以忽略。

圖3 利用LLR估計器和ILLR估計器得到的SSB估計值的散點圖

表3利用兩種估計器得到的SSB之差(cm)統計

目前，國際上主要采用解釋方差來對SSB模型進行評估。解釋方差是指未經SSB校正的海面高度不符值的方差減去經過SSB校正后的海面高度不符值的方差，解釋方差可以理解為海面高度不符值的方差中SSB能夠解釋的部分，解釋方差越大，說明SSB模型越有效。對兩種不同的核函數和帶寬組合，分別使用2009~2011 3年的數據計算利用兩種估計器得到的SSB的解釋方差，表4給出了利用兩種估計器得到的SSB的解釋方差之差。當使用高斯核函數和全局帶寬時，兩種估計器得到的SSB的解釋方差之差在量級，而當使用Epanechnikov核函數和局部帶寬時，兩種估計器得到的SSB的解釋方差之差在量級，所以不管是哪種核函數和帶寬的組合，和LLR估計器相比，ILLR估計器對SSB模型有效性的影響都可以忽略，從而也不會影響高度計最終的測高精度。

表4 利用兩種估計器得到的SSB的解釋方差之差(cm2)

4 結論

非參數模型由于其較高的精度已被廣泛用于雷達高度計的SSB校正，在利用非參數估計方法估計SSB時，通常會采用LLR估計器來求加權矩陣里的權函數，而LLR估計器涉及高維矩陣運算，當數據量較大時需要大量的計算時間。

本文提出了一種ILLR估計器，該估計器可以避免高維矩陣運算，從而極大地提高計算效率。本文分別使用高斯核函數加全局帶寬以及Epanechnikov核函數加局部帶寬這兩種組合來估計SSB，在兩種組合中，利用ILLR估計器得到的SSB和利用LLR估計器得到的SSB之差基本上可以忽略；ILLR估計器可以大幅度地降低估計SSB所需要的時間，當數據量很大時，這種改進的效果更為明顯。目前業務化運行的高度計都是采用基于風速和有效波高的2維SSB模型，將新的參量引入SSB模型已經成為SSB經驗模型研究的一個新方向，本文提出的ILLR估計器可以在不影響SSB估計結果條件下大幅度降低估計SSB所需要的時間，這為實現高維非參數SSB模型的實時運算奠定了基礎。

[1] 王磊. 高精度衛星雷達高度計數據處理技術研究[D]. [博士論文], 中國科學院(國家空間科學中心), 2015.

WANG L. Study on the data processing for high precision satellite radar altimeter[D]. [Ph.D. dissertation], National Space Science Center, Chinese Academy of Sciences, 2015.

[2] 王磊, 許可, 史靈衛, 等. 一種消除合成孔徑雷達高度計延遲校正中殘余誤差的新算法及仿真驗證[J]. 電子與信息學報, 2015, 37(11): 2713-2718. doi: 10.11999/JEIT150282.

WANG L, XU K, SHI L W,A new range migration correction algorithm and its simulation for SAR altimeter[J].&, 2015, 37(11): 2713-2718. doi: 10.11999/JEIT150282.

[3] 劉亞龍. HY-2雷達高度計海面高度定標技術研究[D]. [博士論文], 中國海洋大學, 2014.

LIU Y L. Calibration technology for HY-2 radar altimeter sea surface height[D]. [Ph.D. dissertation], Ocean University of China, 2014.

[4] GASPAR P, OGOR F, LETRAON P Y,. Estimating the sea-state bias of the topex and poseidon altimeters from crossover differences[J].-, 1994, 99(C12): 24981-24994. doi: 10.1029/ 94JC01430.

[5] PRANDI P, PHILIPPS S, PIGNOT V,. SARAL/AltiKa global statistical assessment and cross-calibration with Jason-2[J]., 2015, 38(10):297-312. doi: 10.1080/01490419.2014.995840.

[6] TRAN N, VANDEMARK D, CHAPRON B,. New models for satellite altimeter sea state bias correction developed using global wave model data[J]., 2006, 111(C9): 141-152. doi: 10.1029/2005JC003406.

[7] BONNEFOND P, EXERTIER P, LAURAIN O,. SARAL/AltiKa absolute calibration from the multi-mission corsica facilities[J]., 2015, 38(10):171-192. doi: 10.1080/01490419.2015.1029656.

[8] GASPAR P and FLOREANS J P. Estimation of the sea state bias in radar altimeter measurements of sea level: Results from a new nonparametric method[J]., 1998, 103(C08): 15803-15814. doi: 10.1029/ 98JC01194.

[9] CHELTON D B. The sea-state bias in altimeter estimates of sea-level from collinear analysis of topex Data[J]., 1994, 99(C12): 24995-25008. doi: 10.1029/94JC02113.

[10] GASPAR P, LABROUE S, OGOR F,. Improving nonparametric estimates of the sea state bias in radar altimeter measurements of sea level[J]., 2002, 19(10): 1690-1707. doi: 10.1175/ 1520-0426(2002)019<1690:INEOTS>2.0.CO;2.

[11] WANG J X and XIAO Q X. Local polynomial estimation of time-dependent diffusion parameter for discretely observed SDE models[J]., 2014, 28(4): 871-878. doi: 10.2298/ FIL1404871W.

[12] SU L and ULLAH A. Local polynomial estimation of nonparametric simultaneous equations models[J]., 2008, 144(1): 193-218. doi: 10.1016/j.jeconom. 2008.01.002.

[13] DETTMERING D, SCHWATKE C, and BOSCH W. Global calibration of SARAL/AltiKa using multi-mission sea surface height crossovers[J]., 2015, 38(10):206-218. doi:10.1080/01490419.2014.988832.

[14] PAIGE C C and SAUNDERS M A. LSQR: an algorithm for sparse linear-equations and sparse least-squares[J]., 1982, 8(1): 43-71. doi: 10.1145/355984.355989.

[15] FONG D C L and SAUNDERS M. LSMR: An iterative algorithm for sparse least-squares problems[J]., 2010, 33(5): 2950-2971. doi: 10.1137 /10079687X.

[16] TRAN N, VANDEMARK D, LABROUE S,. Sea state bias in altimeter sea level estimates determined by combining wave model and satellite data[J]., 2010, 115(C03020): 1-7. doi: 10.1029/ 2009JC005534.

[17] FENG H, YAO S, LI L Y,. Spline-based nonparametric estimation of the altimeter sea-state bias correction[J]., 2010, 7(3):577-581. doi: 10.1109/LGRS.2010.2041894.

Improved Local Linear Regression Estimator and Its Application to Estimation for Radar Altimeter Sea State Bias

JIANG Maofei①②③XU Ke①②LIU Yalong④WANG Lei①②

①(Key Laboratory of Microwave Remote Sensing, Chinese Academy of Sciences, Beijing 100190, China)②(National Space Science Center, Chinese Academy of Sciences, Beijing 100190, China)③(University of Chinese Academy of Sciences, Beijing 100049, China)④(Yantai Marine Environmental Monitoring Center Station, State Oceanic Administration, Yantai 264000, China)

The Local Linear Regression (LLR) estimator is usually used when developing a nonparametric model for radar altimeter Sea State Bias (SSB). However, the conventional LLR estimator contains matrices with high dimension. When a large number of data are used to develop the SSB model, the SSB estimation costs too much time. Therefore, the nonparametric estimation method can hardly be used to develop high-dimensional SSB models. This paper presents an Improved LLR (ILLR) estimator, complexity fromtowhich can avoid high-dimensional matrix operations. The improved LLR estimator can greatly reduce the time for SSB estimation without affecting the estimated accuracy. So the improved LLR estimator can laid the foundation for high-dimensional SSB models.

Radar altimeter; Sea state bias; Nonparametric estimation; Local Linear Regression (LLR) estimator

TN953

A

1009-5896(2016)09-2314-07

10.11999/JEIT151280

2015-11-17；

2016-05-05；

2016-06-12

許可 xuke@mirslab.cn

蔣茂飛： 男，1989年生，博士生，研究方向為雷達高度計數據處理.

許 可： 男，1967年生，博士，研究員，博士生導師，主要研究方向為星載雷達高度計系統技術、合成孔徑雷達高度計系統技術和信號處理技術.

劉亞龍： 男，1985年生，博士，研究方向為雷達高度計數據處理.

王 磊： 男，1986年生，博士，研究方向為雷達高度計信號處理