田 曉，鄭洪艷，張 超，王世進(jìn)

(中國(guó)地震局第一監(jiān)測(cè)中心，天津 300180)

多面函數(shù)擬合區(qū)域垂直形變場(chǎng)中核函數(shù)的優(yōu)化選擇

田 曉，鄭洪艷，張 超，王世進(jìn)

(中國(guó)地震局第一監(jiān)測(cè)中心，天津 300180)

應(yīng)用4種核函數(shù)對(duì)3個(gè)不同地區(qū)的垂直形變場(chǎng)進(jìn)行多面函數(shù)擬合比較，并通過(guò)外符合精度以及得到的實(shí)際等值線圖來(lái)綜合選擇較優(yōu)核函數(shù)。通過(guò)對(duì)比發(fā)現(xiàn)，錐面函數(shù)是一種很好的核函數(shù)，適應(yīng)性較強(qiáng)，能較好地?cái)M合出區(qū)域垂直形變場(chǎng)的升降變化，具有一定的實(shí)際應(yīng)用價(jià)值。

多面函數(shù)；垂直形變場(chǎng)；核函數(shù)；優(yōu)化選擇

0 引言

目前，利用重復(fù)高精度水準(zhǔn)測(cè)量仍然是研究地殼垂直運(yùn)動(dòng)的主要方法。但由于水準(zhǔn)重合點(diǎn)是離散的，且水準(zhǔn)網(wǎng)內(nèi)部沒(méi)有數(shù)據(jù)，導(dǎo)致動(dòng)態(tài)平差得到的垂直運(yùn)動(dòng)速率在時(shí)間域和空間域都是不連續(xù)的，實(shí)際上地殼垂直運(yùn)動(dòng)是一種連續(xù)的形態(tài)變化。為了實(shí)現(xiàn)在時(shí)空域內(nèi)客觀準(zhǔn)確地描述形變分布特征，往往需要基于離散的測(cè)點(diǎn)數(shù)據(jù)構(gòu)建合適的數(shù)學(xué)模型擬合整個(gè)地區(qū)的地殼垂直形變速率面。

在地殼垂直形變場(chǎng)的擬合中，多面函數(shù)法是一種傳統(tǒng)的、得到廣泛應(yīng)用的方法。多面函數(shù)法中核函數(shù)的選擇是一個(gè)至關(guān)重要的方面，會(huì)嚴(yán)重影響整個(gè)地區(qū)的擬合效果。在我國(guó)，20世紀(jì)80年代末期開(kāi)始，多位學(xué)者研究了多面函數(shù)法在擬合地殼垂直形變場(chǎng)中的應(yīng)用，這其中也對(duì)核函數(shù)的選擇做了相關(guān)研究，并得出了不同的觀點(diǎn)。楊國(guó)華等[1]對(duì)多面函數(shù)的特性進(jìn)行了初步數(shù)值分析，指出就地殼垂直形變分析來(lái)講，選用冪次不高的“缽”形函數(shù)作為核函數(shù)為佳；趙承坤等[2]研究了速率面擬合法中核函數(shù)中心點(diǎn)的選擇，提出用逐步回歸理論來(lái)篩選核函數(shù)中心點(diǎn)，收到了較好的效果，并指出雙曲型核函數(shù)一般優(yōu)于三次曲面核函數(shù)，其中又以正雙曲型核函數(shù)為佳；陶本藻等[3-4]對(duì)多面函數(shù)擬合法做了相關(guān)試驗(yàn)研究，建議在地殼垂直形變分析中，使用正雙曲面型核函數(shù)，并使平滑因子為0，應(yīng)用此核函數(shù)研究了鮮水河地區(qū)地殼垂直運(yùn)動(dòng)，獲得了該地區(qū)地殼垂直形變速率面的預(yù)測(cè)結(jié)果，較好地反映了地殼垂直形變及應(yīng)變積累的特征；黃立人等[5]在多面函數(shù)用于擬合地殼垂直運(yùn)動(dòng)的研究中，指出雙曲型函數(shù)優(yōu)于三次函數(shù)，且在大多數(shù)情況下，倒雙曲型更好一些；董鴻聞等[6-7]在中國(guó)大陸現(xiàn)今地殼垂直運(yùn)動(dòng)的研究中，認(rèn)為以全國(guó)范圍采用倒雙曲型核函數(shù)為宜，并且核函數(shù)中心點(diǎn)取所有重合點(diǎn)的1/3較好。另外，郝明在研究青藏高原東緣地區(qū)的垂直運(yùn)動(dòng)時(shí)采用了正雙曲面核函數(shù)，平滑因子取了0.01，核函數(shù)中心點(diǎn)取了1/3水準(zhǔn)重合點(diǎn)[8]；秦姍蘭和劉同文在研究西秦嶺地區(qū)的垂直形變時(shí)分別采用了倒雙曲面核函數(shù)和正雙曲面核函數(shù)[9-10]。

總之，學(xué)者們對(duì)多面函數(shù)擬合地殼垂直形變場(chǎng)中核函數(shù)的選擇主要集中在選擇正雙曲面還是倒雙曲面核函數(shù)，并且大多學(xué)者是采用回代點(diǎn)(擬合點(diǎn))或者插值點(diǎn)(檢核點(diǎn))的擬合精度這一指標(biāo)來(lái)判別選擇哪種核函數(shù)。實(shí)際應(yīng)用中這樣選擇核函數(shù)存在不妥之處，因?yàn)橐粋€(gè)好的核函數(shù)應(yīng)該對(duì)于已測(cè)點(diǎn)上的運(yùn)動(dòng)速率有良好的擬合效果，對(duì)于未測(cè)點(diǎn)的預(yù)測(cè)值較為平穩(wěn)，不至于在預(yù)測(cè)點(diǎn)偏離數(shù)據(jù)點(diǎn)時(shí)，內(nèi)插值有急劇的變化。但是對(duì)于實(shí)測(cè)水準(zhǔn)網(wǎng)不管采用擬合點(diǎn)(內(nèi)符合)精度還是檢核點(diǎn)(外符合)精度都只是檢核水準(zhǔn)網(wǎng)上是否有良好的擬合效果，對(duì)于水準(zhǔn)網(wǎng)內(nèi)部的擬合效果得不到驗(yàn)證。因此，本文采用檢核點(diǎn)(外符合)精度以及利用該核函數(shù)得到的實(shí)際等值線圖(驗(yàn)證水準(zhǔn)網(wǎng)內(nèi)部是否有急劇的變化)來(lái)選擇較優(yōu)的核函數(shù)。

1 多面函數(shù)擬合法

多面函數(shù)擬合法是美國(guó)Hardy教授于1977年提出的，并建議用于地殼垂直形變分析。該方法的基本思想是：任何數(shù)學(xué)表面和任何不規(guī)則的圓滑表面，總可以用一系列有規(guī)則的數(shù)學(xué)表面的和以任意精度逼近[10]。根據(jù)這一思想，假設(shè)地球表面上有一點(diǎn)(x,y)，它的垂直速率值ξ(x,y)可以用下式來(lái)表示：

(1)

式中，aj為待定系數(shù)，Q(x,y,xj,yj)是多面函數(shù)的核函數(shù)。核函數(shù)有多種類型，常用的是距離型核函數(shù)和錐面函數(shù)，其中距離型核函數(shù)表達(dá)式為：

Q(x,y,xj,yj)=[(x-xj)2+(y-yj)2+δ]μ

(2)

式中：δ為平滑因子，通?？扇∫恍≌龜?shù)或0；μ一般取1/2、-1/2、3/2。當(dāng)μ=1/2時(shí)，核函數(shù)稱為正雙曲面函數(shù)；當(dāng)μ=-1/2時(shí)，核函數(shù)稱為倒雙曲面函數(shù)；當(dāng)μ=3/2時(shí)，核函數(shù)稱為三次曲面函數(shù)。

錐面函數(shù)表達(dá)式為：

(3)

式中：c為待定參數(shù)。

設(shè)有m個(gè)已知水準(zhǔn)點(diǎn)(xi,yi)，選取其中n(n≤m)個(gè)點(diǎn)作為核函數(shù)的中心點(diǎn)(xj,yj)，令Qij=Q(xi,yi,xj,yj)，則各數(shù)據(jù)點(diǎn)應(yīng)滿足：

(4)

由此可列出誤差方程：

V=QX-ξ

(5)

根據(jù)最小二乘原理可求得待定系數(shù)X，即

X=(QTQ)-1QTξ

(6)

待定系數(shù)求出后，根據(jù)式(4)可計(jì)算測(cè)區(qū)各未知點(diǎn)的垂直形變速率值。

2 算例與分析

2.1 算例介紹

本文選取3個(gè)不同地區(qū)的各2期一等水準(zhǔn)網(wǎng)數(shù)據(jù)，分別采用自由網(wǎng)動(dòng)態(tài)平差計(jì)算得到各個(gè)地區(qū)的垂直形變速率值，并剔除平差后速率中個(gè)別量值過(guò)大的突變點(diǎn)，然后均勻選取1/5點(diǎn)作為外部檢核點(diǎn)，其余點(diǎn)用于擬合多面函數(shù)模型參數(shù)。應(yīng)用多面函數(shù)中4種核函數(shù)對(duì)3個(gè)地區(qū)進(jìn)行擬合比較，其中距離型核函數(shù)的平滑因子是從0到5(約100D，D為水準(zhǔn)點(diǎn)平均間距)間隔0.01共501次計(jì)算中根據(jù)外符合精度選取的最優(yōu)值，錐面函數(shù)參數(shù)是從0到100共101次計(jì)算中根據(jù)外符合精度選取的最優(yōu)值，核函數(shù)中心點(diǎn)是均勻選取1/3點(diǎn)。平差情況如下：

算例1：晉冀蒙地區(qū)2002年和2006年2期數(shù)據(jù)，以大同基巖點(diǎn)為參考點(diǎn)，平差后的單位權(quán)中誤差為1.15 mm/a，共有627個(gè)具有平差后速率的數(shù)據(jù)點(diǎn)(圖1)。

圖1 晉冀蒙地區(qū)垂直形變速率矢量圖

算例2：鄂爾多斯地區(qū)1990年和2013年2期數(shù)據(jù)，以棋東49基為參考點(diǎn)，平差后的單位權(quán)中誤差為0.84 mm/a，共有343個(gè)具有平差后速率的數(shù)據(jù)點(diǎn)(圖2)。

算例3：蘇魯皖地區(qū)1985年和1995年2期數(shù)據(jù)，以青島原點(diǎn)為參考點(diǎn)，平差后的單位權(quán)中誤差為1.10 mm/a，共有796個(gè)具有平差后速率的數(shù)據(jù)點(diǎn)(圖3)。

圖2 鄂爾多斯地區(qū)垂直形變速率矢量圖

圖3 蘇魯皖地區(qū)垂直形變速率矢量圖

為了選出較優(yōu)的核函數(shù)，本文采用外符合精度檢核水準(zhǔn)網(wǎng)上的擬合效果，以及繪制等值線圖檢核水準(zhǔn)網(wǎng)內(nèi)部的擬合效果，進(jìn)行綜合選取。外符合精度σ的表達(dá)式為：

(7)

式中：Δi表示檢核點(diǎn)速率擬合值與實(shí)際觀測(cè)值之差，n為檢核點(diǎn)個(gè)數(shù)。

2.2 不同核函數(shù)擬合比較

分別應(yīng)用多面函數(shù)法中的4種核函數(shù)對(duì)3個(gè)不同地區(qū)的垂直形變場(chǎng)進(jìn)行擬合，得到的最優(yōu)外符合精度以及對(duì)應(yīng)的參數(shù)(表1～3)。

表1 晉冀蒙地區(qū)擬合比較表

表2 鄂爾多斯地區(qū)擬合比較表

表3 蘇魯皖地區(qū)擬合比較表

對(duì)3個(gè)地區(qū)應(yīng)用各自得到的4種核函數(shù)進(jìn)行格網(wǎng)化插值，插值間隔設(shè)為1′×1′，并繪制3個(gè)地區(qū)的垂直形變速率等值線圖(圖4～6)，其中4種核函數(shù)分別對(duì)應(yīng)每幅圖中的a、b、c、d四幅圖。

a 正雙曲面；b 倒雙曲面； c 三次曲面；d 錐面圖4 晉冀蒙地區(qū)垂直形變速率等值線圖

a 正雙曲面；b 倒雙曲面；c 三次曲面；d 錐面圖5 鄂爾多斯地區(qū)垂直形變速率等值線圖

a 正雙曲面；b 倒雙曲面；c 三次曲面；d 錐面圖6 蘇魯皖地區(qū)垂直形變速率等值線圖

2.3 算例結(jié)果分析

1)從表1和圖4可以看出，當(dāng)錐面函數(shù)作為核函數(shù)時(shí)擬合精度優(yōu)于其他3種核函數(shù)，并且用錐面函數(shù)繪制的等值線圖與圖1反映的升降情況一致，不僅水準(zhǔn)網(wǎng)上具有一致的速率值，而且水準(zhǔn)網(wǎng)內(nèi)部也變化平穩(wěn)，較好地反映了整個(gè)地區(qū)的升降情況。前3種核函數(shù)雖然對(duì)水準(zhǔn)網(wǎng)上具有較好的擬合效果，但是水準(zhǔn)網(wǎng)內(nèi)部偏離數(shù)據(jù)點(diǎn)的位置產(chǎn)生急劇的升降變化，不符合地殼垂直形變場(chǎng)擬合的原則，失真度較高，不應(yīng)采用。綜合來(lái)看，晉冀蒙地區(qū)應(yīng)該選擇錐面函數(shù)作為核函數(shù)。

(2)從表2和圖5可以看出，正雙曲面函數(shù)和錐面函數(shù)的擬合精度一樣，都是最高的，而且參數(shù)也相同，實(shí)際上此情況下，2種核函數(shù)完全一致，具有相同的表達(dá)式，說(shuō)明當(dāng)平滑因子為0時(shí)的正雙曲面函數(shù)是錐面函數(shù)的一個(gè)特例。從等值線圖看，用這2個(gè)核函數(shù)繪制的等值線與圖2反映的升降情況一致，在水準(zhǔn)網(wǎng)上和內(nèi)部都表現(xiàn)出較好的擬合效果。但是倒雙曲面函數(shù)同樣能反映圖2的升降情況，并且水準(zhǔn)網(wǎng)內(nèi)部也變化比較平穩(wěn)，同樣適用該地區(qū)，不過(guò)從擬合精度來(lái)看，比正雙曲面和錐面函數(shù)略低。三次曲面為核函數(shù)時(shí)插值效果不好，失真度較高，不應(yīng)采用。綜合來(lái)看，鄂爾多斯地區(qū)應(yīng)該選擇正雙曲面和錐面函數(shù)為核函數(shù)。

(3)從表3和圖6可以看出，三次曲面為核函數(shù)時(shí)的擬合精度最高，說(shuō)明此例中以三次曲面為核函數(shù)時(shí)對(duì)路線上有較好的擬合效果，由于擬合精度不能對(duì)環(huán)內(nèi)部給出真實(shí)的反映，在此例中對(duì)比等值線圖可以發(fā)現(xiàn)水準(zhǔn)網(wǎng)內(nèi)部的插值效果很差(繪圖時(shí)對(duì)于小于-10 mm/a的區(qū)域一律用-10 mm/a顏色代替，對(duì)于大于10 mm/a的區(qū)域一律用10 mm/a顏色代替)，偏離數(shù)據(jù)點(diǎn)時(shí)產(chǎn)生了急劇的升降變化，失真度很高，不是好的核函數(shù)。正雙曲面和倒雙曲面存在同樣的問(wèn)題。雖然錐面函數(shù)為核函數(shù)時(shí)擬合精度最低，但是從等值線圖可以發(fā)現(xiàn)其能反映圖3的升降情況，并且水準(zhǔn)網(wǎng)內(nèi)部變化比較平穩(wěn)，能較好地反映該地區(qū)垂直形變場(chǎng)的升降變化。綜合來(lái)看，蘇魯皖地區(qū)應(yīng)該選擇錐面函數(shù)作為核函數(shù)。

3 結(jié)語(yǔ)

1)多面函數(shù)擬合區(qū)域垂直形變場(chǎng)時(shí)，核函數(shù)選擇的好壞能產(chǎn)生至關(guān)重要的作用。并且不能單靠擬合精度來(lái)判別選擇哪種核函數(shù)最好，應(yīng)該結(jié)合等值線圖中水準(zhǔn)網(wǎng)內(nèi)部的擬合情況來(lái)選擇最優(yōu)核函數(shù)。

2)正雙曲面函數(shù)和倒雙曲面函數(shù)作為核函數(shù)有時(shí)也能有較好的擬合效果，但是這兩種核函數(shù)具有反映特別靈敏的平滑因子需要選取，如果選取不當(dāng)，擬合效果會(huì)很差，實(shí)際選取比較麻煩。但是錐面函數(shù)此時(shí)同樣具有較好的擬合效果，而且其參數(shù)的選取對(duì)擬合效果靈敏度較低。

3)從本文3個(gè)具體算例來(lái)看，錐面函數(shù)是一種很好的核函數(shù)，適應(yīng)性較強(qiáng)，能較好地?cái)M合出區(qū)域垂直形變場(chǎng)的升降變化。綜合來(lái)看，建議多面函數(shù)擬合垂直形變場(chǎng)時(shí)選用錐面函數(shù)作為核函數(shù)。

[1] 楊國(guó)華, 黃立人. 速率面擬合法中多面函數(shù)幾個(gè)特性的初步數(shù)值研究[J]. 地殼形變與地震, 1990,10(4):70-82.

[2] 趙承坤, 黃立人. 速率面擬合法中核函數(shù)中心點(diǎn)的選擇[J]. 地殼形變與地震, 1991,11(2):48-54.

[3] 陶本藻, 王新洲, 于正林,等. 用于垂直形變模型的多面函數(shù)擬合法的試驗(yàn)研究[J]. 地殼形變與地震, 1992,12(1):1-13.

[4] 陶本藻, 杜方. 多面函數(shù)法研究鮮水河地區(qū)現(xiàn)今地殼垂直運(yùn)動(dòng)[J]. 武漢測(cè)繪科技大學(xué)學(xué)報(bào), 1991, 16(3):12-21.

[5] 黃立人, 陶本藻, 趙承坤. 多面函數(shù)擬合在地殼垂直運(yùn)動(dòng)研究中的應(yīng)用[J]. 測(cè)繪學(xué)報(bào), 1993，22(1):25-32.

[6] 董鴻聞, 顧旦生, 李國(guó)智,等. 中國(guó)大陸現(xiàn)今地殼垂直運(yùn)動(dòng)研究[J]. 測(cè)繪學(xué)報(bào), 2002, 31(2):100-103.

[7] 王文利, 陳士銀, 董鴻聞,等. 利用多面函數(shù)擬合中國(guó)陸地垂直運(yùn)動(dòng)速率圖[J]. 測(cè)繪通報(bào), 2002(8):6-8，11.

[8] 郝明. 基于精密水準(zhǔn)數(shù)據(jù)的青藏高原東緣現(xiàn)今地殼垂直運(yùn)動(dòng)與典型地震同震及震后垂直形變研究[D]. 北京:中國(guó)地震局地質(zhì)研究所, 2012.

[9] 秦姍蘭, 王慶良, 季靈運(yùn),等. 利用水準(zhǔn)資料研究西秦嶺地區(qū)的垂直形變[J]. 大地測(cè)量與地球動(dòng)力學(xué), 2012, 32(2):16-19.

[10] 劉同文, 楊志強(qiáng), 王慧敏. 西秦嶺地區(qū)垂直形變的分析研究[J]. 測(cè)繪科學(xué), 2014, 39(4):61-63.

OptimizationSelectionofKernelFunctioninFittingRegionVerticalDeformationFieldUsingMulti-surfaceFunction

TIAN Xiao, ZHENG Hong-yan, ZHANG Chao, WANG Shi-jin

(The First Monitoring and Application Center, CEA, Tianjin 300180, China)

Selecting kernel function in fitting region vertical deformation field using multi-surface function is an important aspect. Four kernel functions are compared in fitting vertical deformation field of three different regions in this paper, and the best kernel function is selected by the external accuracy and the contour map. Tapered surface function is a good kernel function with a strong adaptability by comparison, and it can well fit the rise and fall of region vertical deformation field. Therefore it has some practical values.

multi-surface function; vertical deformation field; kernel function; optimization selection

田曉，鄭洪艷，張超，等. 多面函數(shù)擬合區(qū)域垂直形變場(chǎng)中核函數(shù)的優(yōu)化選擇[J]．華北地震科學(xué)，2017，35(4):11-15.

2017-04-27

中國(guó)綜合地球物理場(chǎng)觀測(cè)——大華北地區(qū)(201508009)；中國(guó)地震局監(jiān)測(cè)、預(yù)測(cè)、科研三結(jié)合課題(CEA-JC/3JH-163711)

田 曉(1988—)，男，山東人，助理工程師，主要從事精密水準(zhǔn)和跨斷層數(shù)據(jù)處理．E-mail：tx_yice@163.com

P315.7

A

1003-1375(2017)04-0011-05

10.3969/j.issn.1003-1375.2017.04.002