同倫攝動(dòng)稀疏正則化方法及其在

1、同倫攝動(dòng)稀疏正則化方法及其在初始地形地貌重構(gòu)中的應(yīng)用 摘要：稀疏正則化方法在參數(shù)重構(gòu)中起到了越來越重要的作用。與傳統(tǒng)的正則化方法相比，稀疏正則化方法能較好地重構(gòu)稀疏變量。由于稀疏正則化的不可微性，需要對(duì)已有的經(jīng)典算法進(jìn)行改進(jìn)。本文構(gòu)建同倫攝動(dòng)稀疏正則化方法克服標(biāo)準(zhǔn)稀疏正則化的不可微性，并應(yīng)用該方法應(yīng)用到項(xiàng)目問題中，能夠有效地重構(gòu)山體初始表面。數(shù)值實(shí)驗(yàn)表明，所提出的方法是收斂和穩(wěn)定的。 關(guān)鍵詞：稀疏正則化；同倫攝動(dòng)稀疏正則化方法；參數(shù)重構(gòu)；重構(gòu)山體初始表面1. 引言 參數(shù)重構(gòu)在許多應(yīng)用中起到重要作用，如波動(dòng)率和政策參數(shù)重構(gòu)1, 2，在其他工程實(shí)踐領(lǐng)域也有重要的應(yīng)用，例如：圖像去噪3, 4，地震信

2、號(hào)重構(gòu)5, 6以及心電信號(hào)重構(gòu) 7, 8。一般情況下，參數(shù)重構(gòu)問題是不適定的，也就是說，即使測(cè)量數(shù)據(jù)的噪聲水平很小，也可能導(dǎo)致重構(gòu)結(jié)果嚴(yán)重偏離真實(shí)參數(shù)9。正則化方法主要就是為了克服參數(shù)重構(gòu)的不適定性，通過選取合適的正則化方法能夠抑制觀測(cè)數(shù)據(jù)中誤差對(duì)于參數(shù)重構(gòu)的不良影響，獲得較為準(zhǔn)確的重構(gòu)值。使用最為普遍的正則化方法是吉洪諾夫正規(guī)化方法，它的目標(biāo)函數(shù)是由擬合項(xiàng)和懲罰項(xiàng)構(gòu)成。吉洪諾夫正則化方法已被用于許多參數(shù)辨識(shí)問題中，很多學(xué)者對(duì)其數(shù)值計(jì)算方法進(jìn)行了研究，如Landweber方法10，高斯牛頓法11，Newton-Kaczmarz方法12和多尺度平滑方法13，這些方法能夠有效地重構(gòu)光滑參數(shù)。隨著經(jīng)

3、濟(jì)和金融理論的發(fā)展，波動(dòng)率和經(jīng)濟(jì)政策參數(shù)的重構(gòu)已經(jīng)廣泛應(yīng)用于許多實(shí)際問題中。在實(shí)際應(yīng)用中，很多需要重構(gòu)的參數(shù)都是稀疏的，即參數(shù)的非零元素的數(shù)目非常有限，遠(yuǎn)遠(yuǎn)小于零元素的數(shù)目。即使重構(gòu)參數(shù)稀疏程度不夠，也可以利用小波和曲波變換使參數(shù)稀疏化。本文先對(duì)該方法在經(jīng)濟(jì)和金融領(lǐng)域中的應(yīng)用進(jìn)行研究，因?yàn)檫@兩個(gè)領(lǐng)域中的參數(shù)通常可以分成已知部分和未知部分。已知部分通常與已有的經(jīng)濟(jì)和金融政策相對(duì)應(yīng)，而未知部分通過參數(shù)重構(gòu)，再將結(jié)果進(jìn)行經(jīng)濟(jì)學(xué)和金融學(xué)解釋，能夠?yàn)檎咧贫ú块T提供行之有效的對(duì)策建議。通過這兩類問題的研究表明該方法能夠應(yīng)用到參數(shù)分解成已知部分和未知部分的問題中。然后再將該方法應(yīng)用到較為復(fù)雜的山體表面重構(gòu)

4、問題中。傳統(tǒng)的吉洪諾夫正則化方法對(duì)于重構(gòu)稀疏參數(shù)效果很不好，而稀疏正則化方法卻能較好地重構(gòu)稀疏參數(shù)，但是稀疏正則化是不可微的，因此需要采用一些技巧來克服這一困難，典型的方法是Bregman迭代14, 15-17。本文構(gòu)造同倫攝動(dòng)稀疏正則化方法，達(dá)到提高算法精度和提高計(jì)算效率的目的。數(shù)值實(shí)驗(yàn)表明同倫攝動(dòng)稀疏正則化方法是收斂和穩(wěn)定的。2. 稀疏正則化方法在實(shí)際應(yīng)用中，對(duì)于光滑參數(shù)重構(gòu)，吉洪諾夫正則化方法具有良好的收斂性和穩(wěn)定性。但是，在稀疏參數(shù)重構(gòu)時(shí)，吉洪諾夫正則化方法的重構(gòu)效果很差，很難滿足工程實(shí)踐的要求。因此，需要采用稀疏正則化方法進(jìn)行參數(shù)重構(gòu)。2.1 吉洪諾夫正則化方法參數(shù)重構(gòu)問題可以歸結(jié)為

5、非線性算子方程(1)的形式：(1)其中， 分別代表非線性算子、需要重構(gòu)的參數(shù)和觀測(cè)數(shù)據(jù)。在實(shí)際問題中，往往還需要考慮觀測(cè)數(shù)據(jù)和理想數(shù)據(jù)之間的噪音水平，即：(2)其中，分別代表真實(shí)的觀測(cè)數(shù)據(jù)和噪音水平。 參數(shù)重構(gòu)的難點(diǎn)在于不適定性，很小的噪音水平 也會(huì)使得重構(gòu)結(jié)果嚴(yán)重背離真實(shí)的物理參數(shù)，從而造成結(jié)果的無意義。解決這一難點(diǎn)最重要的方法就是吉洪諾夫正則化方法，與之對(duì)應(yīng)的吉洪諾夫正則化目標(biāo)函數(shù)定義為：(3)其中， 是在范數(shù)意義下的數(shù)據(jù)擬合項(xiàng), 是起到穩(wěn)定作用的懲罰項(xiàng), 參數(shù)為正則化參數(shù), 該參數(shù)主要是起到平衡數(shù)據(jù)擬合項(xiàng)和懲罰項(xiàng)的作用。吉洪諾夫正則化方法主要是求解下面的優(yōu)化問題.(4)最小值滿足如下的歐

6、拉方程：(5)其中， 是F-導(dǎo)數(shù)。在求解方程(5)的時(shí)候，采用的最普遍的方法是Landweber方法，該方法可以寫成下面的表達(dá)式：,(6)其中，表示迭代次數(shù)。方程(6)是著名的Landweber迭代格式，該數(shù)值格式的顯著優(yōu)點(diǎn)是穩(wěn)定性特別好，但是，收斂速度很慢，不適合應(yīng)用到大型實(shí)際問題中。另一個(gè)重要的方法就是高斯-牛頓方法，該方法收斂速度較快，但是，不如Landweber方法穩(wěn)定。本文主要是在Landweber方法基礎(chǔ)上研究同倫攝動(dòng)稀疏正則化方法，整個(gè)過程可以推廣到高斯-牛頓方法上。2.2 稀疏正則化方法為了能夠有效地重構(gòu)稀疏變量，將標(biāo)準(zhǔn)的吉洪諾夫正則化方法進(jìn)行改進(jìn)，使得吉洪諾夫正則化目標(biāo)函數(shù)轉(zhuǎn)

7、換為如下形式:,(7)其中，表示向量的非零元素的個(gè)數(shù)。近些年，由于稀疏正則化方法能夠有效地重構(gòu)稀疏變量，使得反問題領(lǐng)域?qū)W者越來越多地重視該方法在實(shí)際問題中的應(yīng)用。該方法的難點(diǎn)在于，泛函(7)的懲罰項(xiàng)是不可微的，而且該問題還是NP不可解問題。為了解決NP不可解難點(diǎn)，采用下面的形式進(jìn)行替代：,(8)其中，范數(shù)表示. 在泛函(8)中, 利用范數(shù)代替了原有的范數(shù)，這樣的改進(jìn)使得計(jì)算效率得到了顯著的提高。盡管進(jìn)行了上述的改進(jìn)，但是(8)的懲罰項(xiàng)依然是不可微的，因此，還要在此基礎(chǔ)上進(jìn)一步改進(jìn)，將(8)中的懲罰項(xiàng)改為下面的帶有阻尼系數(shù)的范數(shù)：,(9)其中，利用代替范數(shù), 輔助參數(shù)是一個(gè)正實(shí)數(shù)。與吉洪諾夫正則

8、化方法相仿，將Landweber方法應(yīng)用到(9)中，得到下面的迭代格式：.(10)3. 同倫攝動(dòng)稀疏正則化方法 為了提高數(shù)值迭代格式(10)的收斂速度，將同倫攝動(dòng)反演算法與稀疏正則化方法有機(jī)地融合到一起，達(dá)到提高重構(gòu)效果和收斂速度的目的。構(gòu)造如下的同倫映射 ,(11)其中，代表嵌入?yún)?shù), 迭代算法的初始值。因此，有下面的關(guān)系式：(12)將寫成關(guān)于嵌入?yún)?shù)的冪級(jí)數(shù)形式：(13)目標(biāo)函數(shù)(9)的解可以表示為：.(14)將方程(11)中的在點(diǎn)寫成泰勒級(jí)數(shù)的形式：(15)可以得到下面的關(guān)系式：.(16)和.(17)根據(jù)(17)，對(duì)于含有噪音的觀測(cè)數(shù)據(jù)，采用前兩項(xiàng)構(gòu)造稀疏變量:,(18)其中，表示迭代次

9、數(shù)。 采用第一項(xiàng)構(gòu)造稀疏變量，則有下面的迭代格式：,(19)迭代格式(19)正是前面介紹的Landweber方法，而迭代格式(18)可以看做是在其上面進(jìn)行的修正，將更多的信息應(yīng)用到迭代算法中，可以提高計(jì)算效率。4 數(shù)值算例4.1 重構(gòu)隱含波動(dòng)率在金融領(lǐng)域中，基于布萊克一斯科爾斯期權(quán)定價(jià)模型(Black-Scholes (B-S) option pricing model)重構(gòu)隱含波動(dòng)率一直是重要的研究方向。利用有限差分方法求解該模型的正演問題，在重構(gòu)隱含波動(dòng)率的時(shí)候采用同倫攝動(dòng)稀疏正則化方法，達(dá)到提高重構(gòu)精度和計(jì)算效率的目的。布萊克一斯科爾斯期權(quán)定價(jià)模型表明金融衍生品的價(jià)格變化滿足標(biāo)準(zhǔn)幾何布朗

10、運(yùn)動(dòng)。該模型的邊界條件隨著金融衍生品的變化而變化。當(dāng)邊界條件給定的時(shí)候，金融衍生品的價(jià)格可以通過求解布萊克一斯科爾斯期權(quán)定價(jià)模型獲得。當(dāng)金融衍生品選定為期權(quán)時(shí)，本文定義正演問題為求解期權(quán)的價(jià)格，反演則是重構(gòu)隱含波動(dòng)率。根據(jù)期權(quán)的行權(quán)方式的不同，可以分為歐式期權(quán)和美式期權(quán)，本節(jié)以歐式期權(quán)為研究對(duì)象。 在時(shí)間上的歐式看漲期權(quán), 設(shè) 是歐式看漲期權(quán)在區(qū)域上的價(jià)格, 滿足下面的方程：,(20) 其中，: 股票價(jià)格, : 執(zhí)行價(jià)格, : 利率, : 紅利, : 到期日, : 時(shí)間, : 隱含波動(dòng)率。.反問題可以表述為: 從觀測(cè)數(shù)據(jù), 重構(gòu)隱含波動(dòng)率.定義非線性向量值函數(shù), 即：. 為了測(cè)試本文提出的方法

11、和標(biāo)準(zhǔn)的吉洪諾夫正則化方法之間的差別，設(shè)置如下參數(shù)：到期日, 股票價(jià)格，利率, 執(zhí)行價(jià)格. 為了測(cè)試算法的抑制噪音的能力，在觀測(cè)數(shù)據(jù)上添加1%的高斯噪音。真實(shí)的隱含波動(dòng)率如下：(21)將真實(shí)的波動(dòng)率分解為兩部分：,(22)其中，分別表示常值隱含波動(dòng)率和稀疏隱含波動(dòng)率。設(shè)是已知的, 需要重構(gòu)部分，該部分為：(23) 目標(biāo)函數(shù)(9) 可以寫成：,(24)采用同倫攝動(dòng)稀疏正則化方法重構(gòu)結(jié)果為： (25)采用標(biāo)準(zhǔn)吉洪諾夫正則化方法重構(gòu)結(jié)果為： (26)從上面的重構(gòu)結(jié)果可以看出，同倫攝動(dòng)稀疏正則化方法能夠提高標(biāo)準(zhǔn)吉洪諾夫正則化方法對(duì)于重構(gòu)稀疏隱含波動(dòng)率的精度。4.2 托達(dá)羅模型的政策參數(shù)重構(gòu)托達(dá)羅模型是

12、著名的描述人口流動(dòng)的經(jīng)濟(jì)學(xué)模型，在勞動(dòng)經(jīng)濟(jì)學(xué)領(lǐng)域具有重要的地位。托達(dá)羅模型將人口流動(dòng)的數(shù)量和不同地區(qū)的收入差距聯(lián)系到一起：(27)其中，分別表示人口流動(dòng)的數(shù)量和不同地區(qū)之間的收入差距。函數(shù)是單調(diào)遞增的，即. 隨著經(jīng)濟(jì)的發(fā)展，標(biāo)準(zhǔn)的托達(dá)羅模型(27)已經(jīng)不能準(zhǔn)確地描述人口流動(dòng)和收入差距之間的關(guān)系，還需要考慮政策參數(shù)的作用，將方程(27)改為下面的形式： (28)其中，表示政府部門的政策參數(shù)。 設(shè)有個(gè)農(nóng)村地區(qū)和個(gè)城市地區(qū)，第個(gè)農(nóng)村區(qū)域和第個(gè)城市區(qū)域之間的收入差距表示為(), 第個(gè)和 第個(gè)城市區(qū)域之間收入差距表示為(). 人口轉(zhuǎn)移到第個(gè)城市區(qū)域的數(shù)量用表示。 對(duì)于第城市區(qū)域, 政策參數(shù)可以分解為兩部

13、分, 其中分別表示政策參數(shù)對(duì)于農(nóng)村和城市的作用效果。為了簡(jiǎn)化討論，假設(shè)從城市區(qū)域轉(zhuǎn)移到農(nóng)村區(qū)域的人口數(shù)量為零。在實(shí)際應(yīng)用中，函數(shù)有多種表達(dá)方式，由于本文主要是測(cè)試同倫攝動(dòng)稀疏正則化方法的有效性，因此，函數(shù)選取為線性函數(shù)的形式。托達(dá)羅模型可以寫成下面的形式：(29)其中，、和 已知，政策參數(shù)向量是未知的。 令 ,(30)和,.(31) 方程 (30) 可以寫成下面的形式：.(32)值得注意的是方程(32)是欠定的，因此，需要正則化方法求解該方程。將政策參數(shù)進(jìn)行分解 , 其中分別表示以前的政策參數(shù)(已知)和稀疏政策參數(shù)(未知)。目標(biāo)函數(shù)改寫為：.(33)在數(shù)值算例中， 設(shè), 的每個(gè)元素都是. 為了

14、測(cè)試算法的抑制噪音的能力，在觀測(cè)數(shù)據(jù)上添加1%的高斯噪音。重構(gòu)5個(gè)稀疏政策參數(shù)，每一個(gè) 有且只有一個(gè)非零元素，大小為0.5. 同倫攝動(dòng)稀疏正則化方法能夠精確地將非零元素的位置判斷出來，但是標(biāo)準(zhǔn)的吉洪諾夫正則化方法就不具備這樣的能力。表1給出了兩種不同正則化方法重構(gòu)的均方差。表1：兩種不同方法的均方差Table 1: Mean Square Errors非零元素的位置同倫攝動(dòng)稀疏正則化方法吉洪諾夫正則化方法25.25%10.07%45.19%12.96%64.69%15.78%84.92%14.06%104.71%9.93%4.3重構(gòu)初始地形地貌 地形地貌演變模型將引起地球表面變化的內(nèi)因和外因有

15、機(jī)的結(jié)合到一起，影響地形地貌演變的內(nèi)因是地球內(nèi)部板塊速度場(chǎng), 而外因則非常多. 主要有以下三個(gè)因素：(I) 地球表面擴(kuò)散過程(diffusion of hillslope topography): 這里所說的擴(kuò)散過程是一個(gè)包含很多復(fù)雜因素的統(tǒng)稱，例如：風(fēng)化作用(weathering), 坡面沖刷(slope wash), 地面水流(overland flow), 土壤滑動(dòng)(soil creep)和基巖滑坡引起的質(zhì)量損失(mass wasting by bedrock-involved landsliding); (II) 河流沖刷對(duì)于基床的切割(bedrock incision)；(III)

16、冰川的融化影響(melting effect).建立地形地貌演變數(shù)學(xué)模型是屬于地形地貌學(xué)的范疇，描述上面因素的數(shù)學(xué)公式一般都是高度非線性的，將這些具體模型的模型進(jìn)行統(tǒng)一考慮，可以得到如下模型： (34)模型正演定義為： 給 定 地 球 表 面 的 速 度 場(chǎng)和 初 始 地 形 地 貌 分 布 函數(shù)求解地形地貌隨時(shí)間變化關(guān)系，寫成數(shù)學(xué)表達(dá)式為： (36)模型反演的主要目的是利用現(xiàn)在的地形地貌，重構(gòu)初始地形地貌，可以表示為： (37)值得注意的是我們這里假設(shè)山體的隨時(shí)間運(yùn)動(dòng)的速度場(chǎng)是已知的，將研究的重點(diǎn)放到了地形地貌運(yùn)動(dòng)模型上，并沒有放到地球內(nèi)部熱場(chǎng)分布模型上，該部分的研究已經(jīng)在上一年度完場(chǎng)。方程

17、 (37) 可以寫成下面的線性形式：.(38)值得注意的是方程(38)是不適定的，因此，需要正則化方法求解該方程。將初始地形地貌進(jìn)行分解 , 其中分別表示通過其他探測(cè)手段已知的地形地貌成分(尺度較大)和細(xì)節(jié)地形地貌成分(尺度較小)。目標(biāo)函數(shù)改寫為：.(39)在數(shù)值算例中，在小波基下表示為稀疏系數(shù)，重構(gòu)這些稀疏系數(shù)。 為了測(cè)試算法的抑制噪音的能力，在觀測(cè)數(shù)據(jù)上添加1%的高斯噪音。同倫攝動(dòng)稀疏正則化方法能夠精確地將非零元素的位置判斷出來，但是標(biāo)準(zhǔn)的吉洪諾夫正則化方法就不具備這樣的能力。表2給出了兩種不同正則化方法重構(gòu)的均方差。從表2中還可以看出新方法更加穩(wěn)定。表2：兩種不同方法的均方差Table2

18、: Mean Square Errors實(shí)驗(yàn)次數(shù)同倫攝動(dòng)稀疏正則化方法吉洪諾夫正則化方法12.21%8.07%22.16%8.16%32.23%7.98%42.39%9.01%52.31%8.58%5 結(jié)論本文構(gòu)造了同倫攝動(dòng)稀疏正則化方法，提高了標(biāo)準(zhǔn)吉洪諾夫正則化方法重構(gòu)稀疏變量的精度，同時(shí)也提高了計(jì)算效率，并將該方法應(yīng)用到重構(gòu)初始地形地貌問題中。對(duì)于實(shí)際情況而言，絕大部分參數(shù)變量都可以分解為已知部分和稀疏部分，因此只需要重構(gòu)對(duì)應(yīng)的稀疏部分即可，說明該方法具有較好的應(yīng)用前景。同時(shí)，數(shù)值算例表明本文構(gòu)造的算法具有較高的精度和收斂速度。參考文獻(xiàn)1 Egger H, Engl H W. Tikhon

19、ov regularization applied to the inverse problem of option pricing: convergence analysis and rateJ. Inverse Problems, 2005, 21: 10271045.2 Trong D D, Thanh D, Lan N N, Uyen P H. Calibration of the purely T-dependent Black-Scholes implied volatilityJ. Applicable Analysis, 2014, 93: 859-874.3 Wu C D.

20、A new denoising method based on multistage median filter and nonsubsampled contourlet transformJ. Journal of Computational Information Systems, 2012, 8(12): 4881-4888.4 Zhao Q, Ye B, Cao J. Image denoising based on improved non-local means and nonsubsampled contourlet transform wiener filteringJ. Jo

21、urnal of Computational Information Systems, 2010, 6 (2): 601-609.5 Baig A M, M. Campillo M, Brenguier F. Denoising seismic noise cross correlationsJ. Journal of Geophysical Research, 2009, 114: doi: 10.1029/2008JB006085. 6 Hennenfent G, Herrmann F J. Simply denoise: wavefield reconstruction via jitt

22、ered undersamplingJ. Geophysics, 2008, 73 (3): V19-V28.7 Istepanian R H, Petrocian A A. Optimal zonal wavelet-based ECG data compression for a mobile telecardiology systemJ. IEEE Transactions on Information Technology in Biomedicine, 2000, 4 (3): 200-211.8 Sameni R, Shamsollahi M B, Jutten C, Cliffo

23、rd G.D. A nonlinear Bayesian filtering framework for ECG denoisingJ. IEEE Transactions on Biomedical Engineering, 2007, 54 (12): 2172-2185.9 Engl H W, Hanke M, Neubauer A. Regularization of Inverse ProblemsM. Kluwer: Dordrecht, 1996.10 Hanke M, Neubauer A, Scherzer O. A convergence analysis of the Landweber iteration for nonlinear ill-posed problemsJ. Numer. Math., 1995, 72: 21-37.11 Blaschke B, Neubauer A, Scherzer O. On convergence rates for the iteratively regularized Gauss-Newton methodJ. IMA J. Numer. Anal., 1997, 17: 421-436.12 Burger M, Kaltenbacher B. R