散亂數(shù)據(jù)擬合的自適應分片逆尺度空間算法

鐘軼君， 李崇君

(1. 浙江理工大學理學院數(shù)學科學系，浙江 杭州 310018；2. 大連理工大學數(shù)學與科學學院，遼寧 大連 116024)

散亂數(shù)據(jù)擬合(函數(shù)擬合/逼近)是一個在信號處理、計算機圖形學中被廣泛研究的問題。其通過給定的散亂的、帶噪音的點集來擬合一個較好的曲線或曲面，需要解決2個基本問題：①如何選擇一個簡單的函數(shù)空間作為描述對象函數(shù)的逼近空間；②其空間要有很好的逼近度。為此，隨著散亂點曲面擬合工具多樣性的不斷探索，有學者開始利用探索多層結構性質來改進逼近方法，如LEE等[1]提出的基于B樣條的多層逼近方法；CASTANO等[2]提出了基于小波的多層正則化方法等。JOHNSON等[3]進一步探討了在平移不變空間(principal shift invariant，PSI)中基于 B樣條和小波的散亂點曲面擬合問題。以上方法均力求在所構建的函數(shù)類中與研究對象函數(shù)本身誤差最小。對于第②個問題，學者們希望在其空間中得到散亂數(shù)據(jù)的稀疏表示解，因此推動了稀疏逼近和優(yōu)化方法在函數(shù)擬合領域的廣泛應用。隨著深層挖掘曲面中的稀疏性，幾何建模中很多的啟發(fā)式方法得到了進一步發(fā)展和推廣。近些年，信號的稀疏性研究也取得了豐碩的成果，其中細化研究信號的稀疏分布(非零元素的分布)成為了熱點，如HUANG[4]提出了結構化稀疏的概念，即一類非零元具有一定結構性的稀疏向量，如廣泛研究的非零元是成塊出現(xiàn)的塊稀疏[5-9]，在統(tǒng)計中也稱為Group Lasso[10-17]，除此之外還有樹稀疏[18]、圖稀疏[19]。劉翼鵬[19]討論了向量服從全局稀疏、局部分布不同的凸優(yōu)化模型，討論了全局稀疏局部稠密，即非零元素呈簇狀分布的稀疏向量。提取出了可以描述整個信號的結構性信息，有利于使用更少的樣本和運算量來恢復信號，且能夠決定信號是否可以有效地得到恢復。對于不同的實際場景，信號的稀疏表示向量具有不同的結構，更好地研究這些結構有利于理解和解決實際問題。隨著對稀疏性的細化探索，學者們嘗試在曲面重構，尤其是在PSI空間中的曲面重構問題所具有的結構型稀疏表示進行研究。HAO 等[20]提出了利用凸優(yōu)化模型來考慮PSI空間中基于B樣條的散亂點曲面擬合問題，利用多塊LASSO (MLASSO)得到散亂點曲面擬合問題的稀疏解；文獻[21]得到了曲面在PSI空間中每一層 B樣條上的稀疏表示解?；贛LASSO的散亂點曲面擬合問題不僅可以得到稀疏解(稀疏表示系數(shù))而且擬合效果較好，但其多個正則化參數(shù)的選取是一個待解決的問題。OSHER等[22]提出 Bregman逆尺度空間算法(inverse scale space，ISS)，即

此方法是貪婪算法和凸優(yōu)化算法的一種結合，不僅規(guī)避了凸優(yōu)化算法中參數(shù)的選取問題而且運行速率較快。之后LI和ZHONG[23]提出了一種帶有刪除機制的分片逆尺度空間算法，核心是將分片稀疏性引入到ISS算法中。在稀疏恢復問題上，相較于ISS算法，該算法能夠更好地保護信號的小尺度非零元素，在應用于散亂數(shù)據(jù)擬合問題中得到了較好的擬合效果。分片稀疏性是一類描述向量中非零元素散亂分布現(xiàn)象的稀疏性。DONOHO和HUO[24]針對片稀疏信號對應觀測矩陣具有的分塊結構，引入分塊矩陣互相干(mutual coherence)條件，可將分塊正交矩陣的精確恢復信號的稀疏度理論上界推廣到一般分塊矩陣之中, 改進了稀疏信號可精確恢復的充分性條件。其結果不僅給出了分片稀疏恢復的理論保證，也提升了L0模型和壓縮感知恢復整體稀疏信號的可信度。向量的分片稀疏結構所對應的觀測矩陣A的分塊結構恰好符合PSI空間基于 B樣條所生成基函數(shù)的層數(shù)，因此引入分片稀疏性到稀疏優(yōu)化算法(模型)中，并分析散亂點曲面擬合問題是有意義的。

然而分片ISS算法(P_ISS)[23]也有其缺陷，算法中的刪除機制增加了新的參數(shù)選取問題。目前考慮基于散亂點的曲面擬合的稀疏表示的文獻較少?；诖耍疚奶岢鲆环N自適應的分片逆尺度空間算法(adaptive piecewise ISS，aP_ISS)，由此得到散亂點曲面擬合問題的稀疏表示。此算法不僅規(guī)避了帶有刪除機制的分片逆尺度算法中需要人工選取參數(shù)的問題，而且通過對“分片符號一致性”的分析，使得算法具有自適應性，可更好地應用到散亂點曲面逼近問題中。

本文主要貢獻如下：

(1) 通過將分片稀疏的理念引入到散亂數(shù)據(jù)擬合問題中，拓廣了稀疏技術在函數(shù)逼近中的應用。

(2) 改進了文獻[22]中的分片逆尺度空間算法，規(guī)避了算法中的刪除機制。利用“分片符號一致性”使得算法更具自適應性。

(3) 稀疏算法作為高效算法已經成功地應用到了曲面擬合問題中，本文建立了一個稀疏表示系數(shù)的先驗信息，從而可以更好擬合曲面的細節(jié)信息，達到更好地逼近效果。與已有的稀疏重建算法相比，不僅有更松弛的理論保證并且擬合效果更好。

1 預備知識

1.1 散亂數(shù)據(jù)曲面擬合問題

典型的散亂點曲面擬合問題是指：存在一個散亂點集合X={x,x,… ,x} ?Ω?Rd和相對應的函數(shù)值f|X= {f1,f2, …,fn}，且需要找到一個屬于空間H的函數(shù)g能夠擬合給定的數(shù)據(jù) { (xk,fk) }nk=1。經典的光滑樣條散亂點重構模型為

由于g屬于Beppo-Levi空間，所以當數(shù)據(jù)點規(guī)模變大時計算量也會變大。為了解決該問題，文獻[3]考慮在一個由緊致支撐函數(shù)移位擴張生成的空間(稱為 PSI空間)中求解以上正則化問題。利用PSI空間不僅避免了計算量增大的問題，也由于其和小波、B樣條的聯(lián)系使得稀疏逼近數(shù)據(jù)函數(shù)成為可能。

利用分塊矩陣A= [A1, … ,AN]表示由定義的觀測矩陣

向量bσ表示被噪音污染的數(shù)據(jù){fi}，則曲面擬合問題g(xi) ≈fi可以記為bσ=Ax+e，其中x為曲面表示系數(shù)；e為噪音。

1.2 符號說明

定義1[24].觀測矩陣A的列互相干為

其中，ai為矩陣A的第i列。

定義 2[23].設矩陣A= [A1, … ,AN]是N個矩陣的拼接，定義第i塊的子陣互相干為

命題1[23].第i塊的子陣互相干μi,i滿足

其中，αi∈[0,1]為衡量第i個子矩陣塊內互相干與矩陣整體互相干比例關系的參數(shù)。

Oracle估計[21]是最小二乘解在真實支撐集S上的部分，oracle的非零元素可表達為

分片oracle估計[22]是oracle估計限制在子空間Si上的部分

1.3 逆尺度空間算法

文獻[25]中提出了一個自適應的逆尺度空間算法(adaptive inverse scale space，aISS)，實現(xiàn)算法步驟如下：

輸入：A，bσ,threshold ≥ 0 。

輸出：x(tk)。

步驟1.初始化

步驟2. 計算x(tk)為問題

在Sk=Ik上的解。

步驟3.更新時間變量為

步驟4.更新次梯度為

其中，t=tk+1，ejn∈R為第j個分量為1其余分量為0的向量。

步驟5.更新指標集為

文獻[21]中給出了如下ISS系統(tǒng)(1)的理論保證條件。

假設 1[21].限制強凸性(restricted strong convexity)：存在γ∈ ( 0,1]，使得

假設 2[21].不可表示條件或精確恢復條件(irre-presentable condition or exact recovery condition)：存 在η∈ ( 0,1)， 使 得其中

定理1[21]. 設假設1和2均成立，那么存在一個時刻τ使得Bregman ISS路徑(1)在時刻τ之前沒有錯誤選擇，即?t≤τ，supp(x(t)) ? supp(x)以大于的概率成立。若

的概率有sign(x(τ))=sign(x)，即Bregman ISS求解的最優(yōu)解達到了符號一致性。

1.4 分片稀疏

對于給定的向量

(1) 全局稀疏。x稱為s-稀疏的若x至多含有個非零元素；

(2) 塊稀疏[4-9]。x稱為塊s-稀疏的若x至多含有s個含有非零元素的塊，即塊l0范數(shù)2,0||x||=不超過s；

(3) 分片稀疏[23]。向量可以分成N個部分，且每部分xiT∈Rni都是一個稀疏向量，若已知si=||xi||0(i=1,…,N)，則稱x是分片(s1,… ,sN)-稀疏向量。

1.5 分片逆尺度空間算法(P_ISS)

分片逆尺度空間系統(tǒng)(piecewise inverse scale space，P_ISS)[22]旨在恢復分片稀疏向量中每一個子向量∈Rni的支集，即P_ISS方法分別且同時求解xi，再由所有的xi共同構成x，即

算法1[22].分片逆尺度空間算法

輸入：觀測矩陣A= [A1, … ,AN]，觀測向量bσ，噪音程度參數(shù)σ。

輸出：分片稀疏向量x=(x1T, … ,xTN)T∈Rn。

步驟1.初始化。

步驟2.最小二乘。更新為

的解。其中I為指標集；IP為在指標集I上的投影算子。

步驟3.并行時間路徑更新。

步驟4.并行的次梯度向量更新。

步驟5.并行的指標集更新。

步驟6.刪除機制。設置向量

滿足不等式

的分量為0。引入分片稀疏性后提出的帶刪除機制的分片Bregman逆尺度空間(P_ISS)算法。因該算法結構類似于正交匹配算法，所以其運行速率較快，又因求解的是分片凸優(yōu)化問題，故其收斂到l1稀疏解。數(shù)值實驗表明相比于 ISS算法，分片Bregman逆尺度空間算法能夠更好地保護小尺度元素不被噪音混淆。

2 自適應分片逆尺度空間算法(aP_ISS)

2.1 P_ISS算法中刪除機制的缺陷

在改良的貪婪算法中，比如壓縮采樣匹配追蹤算法 (compressive sampling MP，CoSaMP)[26]，利用刪除機制對算法當前恢復的解進行自適應的修正是較為常用且有效的方法。CoSaMP在算法最后一步保留最大的s(s=||x||0)個非零元素，其余的元素均設置為0以克服類似StOMP算法[27]在支集添加冗余元素的缺陷。在算法P_ISS迭代過程中，由于分片稀疏向量x整體到達符號一致性[22]的時刻其中iτ為第i片到達符號一致性的時刻。即所有的“片”到達了符號一致性算法才會停止更新支集。在所有的N個子向量都到達符號一致性之前，早到達符號一致性的“片”會由于不斷的迭代產生錯誤或冗余的非零元素，所以在P_ISS算法中設置刪除機制可以有效地避免冗余元素的產生。但是，刪除機制中引入的參數(shù)C使得算法增加了選取問題，即參數(shù)C與數(shù)據(jù)密切相關，需要經驗性選取，消耗大量時間。

2.2 分片符號一致性

本文考慮對算法進行分片符號一致性的改進，即當某一“片”到達符號一致性時, 停止此“片”的支集更新。從而保證在最小二乘計算中，該“片”不會隨著晚到達符號一致性的“片”的更新而增加冗余的元素。該算法自行地對每一“片”實行分片符號一致性的“監(jiān)督”，使得算法自適應性更強。

2.3 自適應分片逆尺度空間算法(aP_ISS)

算法2.自適應分片逆尺度空間算法

輸入：觀測矩陣A= [A1, … ,AN],觀測向量bσ，噪音程度參數(shù)σ。

輸出：

步驟1.初始化。

步驟2.最小二乘。更新為

的解。其中，I為指標集；PI為在指標集I上的投影算子。

步驟3.并行時間路徑更新。

步驟4.并行的次梯度向量更新。

步驟5.并行的自適應的指標集更新。若第i片未達到符號一致性，則

實際算法操作中，通過算法當前迭代向量和oracle向量的符號一致性，即考慮第i片的向量xi和分片oracle估計向量是否達到了符號一致性

并衡量是否更新指標集。下面給出 aP_ISS的性能保證定理。首先給出一些相關的引理。

引理1.矩陣AS= [AS1,…,ASN]的奇異值Λ滿足

其中，

證明：在每個子矩陣Ai中，選擇指標集為Si的列向量組成的子矩陣ASi。考慮Grassmannian矩陣

其中，Φij為ΦS中的元素，且ΦS中對角線元素滿足非對角線元素滿足|Φij| ≤μ。因此

證畢

引理 2.設矩陣A= [A1, … ,AN]的每一個子矩陣Ai的互相干參數(shù)是αi，假設1成立，則

其中，

定理2.設假設條件1和2均成立，若分片稀疏信號滿足

其中，smax為si=|Si|的最大值；

3 數(shù)值實驗

3.1 實驗設計

本文的散亂數(shù)據(jù)曲面擬合問題，對應于分片稀疏模型，矩陣A= [A1, … ,AN]是由定義的觀測矩陣帶噪音的觀測數(shù)據(jù){fi}是用向量bσ表示。此部分通過給定的散亂點集{xk,yk}和其相對應的觀測數(shù)據(jù){f(xk,yk)}擬合曲面并得到其稀疏表示g(x,y)。其逼近函數(shù)是由張量積 B樣條基函數(shù)(x,y)(i= 1 ,…,N;j= 1,…,ni)表示，其中N為多層B樣條的層數(shù)，基函數(shù)的總數(shù)為n=n1+…+nN，于是待擬合曲面為

考慮到N層B樣條基函數(shù)之間是相關的，于是此問題轉變?yōu)榱艘粋€求解分片稀疏向量(曲面的系數(shù)向量)，即

的分片稀疏逼近問題。

本文通過對一個函數(shù)f添加高斯噪音來人工生成一個帶噪音的數(shù)據(jù)集{((xk,yk) ,f(xk,yk)):k=1,… ,9 00}，即

對比以下算法求解散亂數(shù)據(jù)曲面擬合問題的稀疏表示的數(shù)值效果：

(1) 逆尺度空間算法(aISS)。

(2) 帶有刪除機制的分片逆尺度空間算法(P_ISS)，其中刪除機制中參數(shù)C=0.3。

(3) 多塊LASSO (MLASSO)利用Jacobi型的逼近交替方向乘子算法(alternating direction method of multipilers，ADMM)[28]求解MLASSO，其中正則化參數(shù)的值取為γ1=0.3，γ2= 0 .15，γ3=0.2，γ4=0.3。

(4) 本文提出的自適應分片逆尺度空間算法(aP_ISS)。通過均方根誤差(root-mean-square，RMS)、分片稀疏度及算法運行時間(CPU)進行4種算法比較，RMS為

其中，

3個噪音分別取值為σ1= 0 .25,σ2= 0 .1,σ3= 0 .01。

3.2 實驗結果分析

從表1～3的數(shù)值結果和圖1～3的擬合效果可以發(fā)現(xiàn)，若不考慮分片稀疏性的a_ISS算法，本文算法(aP_ISS)較其他 2種算法在近似散亂點集上優(yōu)勢較明顯，同時在保護分片稀疏度上效果更好，即每一層支集均不為空集且稀疏。相較于分片稀疏性的MLASSO算法，aP_ISS算法在運行效率上有明顯的優(yōu)勢，MLASSO模型存在事先人工輸入和調節(jié)多個正則化參數(shù)的缺陷，本文算法的自適應性保證了算法可以更快地運行，且 MLASSO算法近似的曲面效果遠不如本文算法。P_ISS算法的主要缺陷在于刪除機制中參數(shù)C的選取依賴數(shù)據(jù)、噪音等因素，需要花費大量時間，而本文算法不僅避免了參數(shù)的選取而且能夠得到和 P_ISS算法同等的近似效果。

另從表1～3的結果還可以觀察到，分片稀疏逼近的目的是為了在每一層上都有非零元素表示，且是稀疏表示，本文算法除了近一步改進已有的分片算法 P_ISS，更重要的是為了通過“分片符號一致性”來體現(xiàn)分片稀疏的思想。實驗結果表明本文的算法不僅在擬合效果上更好，并且達到了分片稀疏表示的效果。

表1 觀測函數(shù)f1對應曲面擬合的數(shù)值結果

表2 觀測函數(shù)f2對應曲面擬合的數(shù)值結果

表3 觀測函數(shù)f3對應曲面擬合的數(shù)值結果

圖1 觀測函數(shù)f1對應4種算法擬合出的近似曲面的對比圖

圖2 觀測函數(shù)f2對應4種算法擬合出的近似曲面對比圖

圖3 觀測函數(shù)f3對應4種算法擬合出的近似曲面對比圖

4 結 束 語

本文利用 aP_ISS處理散亂點曲面擬合問題,通過引入分片稀疏性，可以有效地保護在PSI空間中多層支集的分片稀疏性；進一步通過對分片逆尺度空間系統(tǒng)的分片符號一致性的考慮，避免了分片算法中參數(shù)的選取問題。實驗結果表明，該算法近似散亂點曲面效果較好，驗證了利用稀疏優(yōu)化算法處理散亂點曲面擬合問題的優(yōu)勢。

分片稀疏思想應用到散亂數(shù)據(jù)擬合等問題仍然欠缺較為完善的研究，通過對稀疏優(yōu)化算法的自適應分片處理來求解曲面擬合問題，將是未來研究的重點。