基于WLs的二維地球化學景觀多重分形建模

高 歆

(中國地質(zhì)大學(北京)地球科學與資源學院， 北京 100083)

利用多重分形理論建模構(gòu)成了過去數(shù)十年來客觀自然現(xiàn)象比較流行的一種尺度分析手段，廣泛的被應(yīng)用在包括地學、湍流和因特網(wǎng)等自然科學領(lǐng)域中。為了有效的指導(dǎo)地球化學勘探和礦產(chǎn)資源評價，Mandelbrot和Cheng[1-2]提出了多種多重分形模型來描述地球化學數(shù)據(jù)的空間尺度結(jié)構(gòu)，這些新的模型不僅考慮地殼中礦石元素的分布幾何形狀，而且將傳統(tǒng)的地質(zhì)統(tǒng)計學納入到了計算分維譜的框架之下，其中尤其是數(shù)據(jù)域的奇異性分析，即計算空間中每點的奇異值，一種高通濾波方式，可以有效的從由各種地質(zhì)過程生成的混合模式和特征中分解背景值和異常值。

對多重分形譜的估計是多重分形分析的重要任務(wù)。針對多重分形譜建模，有兩種最常用的建模方法，其中一種是基于分維譜f(α)，另外一種是基于余維數(shù)c(γ)。兩種模型之間的內(nèi)在聯(lián)系現(xiàn)已得到證明[3]。多重分形理論已廣泛用于各種學科領(lǐng)域，并因此產(chǎn)生了各種各樣的計算多重分形譜函數(shù)的方法，如矩方法[4]、直方圖法、小波方法[5],以及二次維矩法[6]等，其中應(yīng)用最廣泛的是將空間分成許多格子的矩統(tǒng)計方法。另外，多重分形模型將服從正態(tài)或者對數(shù)正態(tài)分布的微量元素的背景值與服從分形分布(冪律分布)的高低異常值統(tǒng)一到多重分形譜中，多重分形分析體現(xiàn)了統(tǒng)計手段上的進步,目前地學數(shù)據(jù)處理中，最常用到的是基于盒子的多重分析方法。Cheng等[7]通過頻率-濃度、數(shù)量-濃度、濃度-濃度、濃度-周長、S-A等分形模型及矩分析等多重分形方法，將地球化學異常從背景中分離出來；郭科等[8]運用多重分形模型，通過研究Pb、Zn、Sb、Ag、Cu共五種元素的廣義分形維數(shù)Dq，認為這些元素在空間上的分布服從多重分形分布特征，同時指出多重分形方法相比于傳統(tǒng)的化探方法來說在判斷元素的共生組合方面具有一定的優(yōu)越性。Zheng等[9]運用周長-面積等分形學研究方法研究了中國3個金礦床中黃鐵礦、雄黃等礦物表面As、S、Fe、Si等元素在成礦過程中的遷移和轉(zhuǎn)換規(guī)律。

自從Parish和Frisch[10]提出用湍流域中的速度增量來刻畫尺度函數(shù)：Sp(l)=<(δfl)p>～lξp，δfl(x)=f(x+l)-f(x)是尺度域中的增量。不幸的是，Muzy[11]發(fā)現(xiàn)了該種方法只能刻畫信號中比較強的奇異性，對于一些弱的奇異性和譜值左半部分顯得無能為力。Mallat[12]提出了利用小波極大模值捕捉信號的奇異性以來，Arneodo等人提出了WTMM代替直接的速度增量來計算尺度函數(shù)，從而避免了上述缺點。然而，WTMM需要連續(xù)小波變換，以及利用不同的半徑在不同的尺度上進行搜索鄰域中的模極大值，復(fù)雜的計算過程，高昂的計算量都使的這種方法變的不很實用。為了避免上述所述WTMM涉及到的繁瑣的工作量，法國數(shù)據(jù)學Jaffard[13]等人提出了WLs的計算方法，不僅能夠描述譜值的整個區(qū)域，還避免了連續(xù)小波變換以及尋值等繁瑣步驟。多次研究證明，計算分形譜的方法嚴重依賴方法和數(shù)據(jù)本身，WLs的引入為這種情況提供了一種補充。Dojnow[14]通過基于小波系數(shù)的多重分形模型，考察了含有噪聲信號的多條腦電圖曲線，指出分析信號的相位攜帶了最多的尺度不變性信息；Khalil等[15]利用基于WTMM的多重分形模型考察了間歇染色體結(jié)構(gòu)，將染色體版圖從背景中分離出來，提出了一種新的邊界分割方法，闡述了染色體結(jié)構(gòu)、位置和不穩(wěn)定性之間的潛在關(guān)系。

1 基于盒子的多重分形分析

(1)

這里的非整數(shù)αi稱為coarse holder指數(shù)[16]。把在分形上具有相同α值的小盒子數(shù)目記為Nα(ε)，它與ε的大小有關(guān)，并且可以寫成：

(2)

將上式與N(ε)～ε-D的簡單分形公式對比，可見f(α)的物理意義是表示具有相同α值的子集的分形維數(shù)，或者說它描述了ε→0+直方圖N(ε)的變化，f(α)定義為具有相同α值的子集的Hausdorff維數(shù)，稱為多重分形譜或奇異性譜。在多重分形譜的眾多計算方法當中，矩方法是最常用的方法之一。為了解f(α)的分布特性，定義分割函數(shù)：

(3)

可以獲得質(zhì)量指數(shù)函數(shù)τ(q)，τ(q)是關(guān)于q的凸函數(shù)。當q=1時達到最大曲率，該點的曲率可以用于判斷是否具有多重分形性。通過勒讓德變換，可以得到α(q)和f(α)如下的兩個函數(shù)：

(4)

三個函數(shù)τ(q)、α(q)和f(α)是相互關(guān)聯(lián)的。

2 基于WLs的多重分形分析

(5)

繼續(xù)設(shè)H0(k)和G0(k)為兩個具有有限沖擊響應(yīng)的濾波器，引入二進區(qū)間以及并集：

λj,k1,k2={[k12j,(k1+1)2j],[k22j,

(6)

(7)

理論上的WLs定義為(圖1)

(8)

圖1 WL的定義示意圖

WLs正是通過其局部奇異性刻畫能力構(gòu)成了了尺度分析模型，如果X(t)在t0點存在奇異性指數(shù)α(t0)≥0，消失矩Nψ≥α(t0)，WLs與尺度之間的冪律關(guān)系為：

(9)

通過計算位于尺度2j上的WLs的q階矩的平均值，可以得到結(jié)構(gòu)函數(shù)和質(zhì)量指數(shù)函數(shù)。

(10)

(11)

通過勒讓德變換，奇異值函數(shù)和分維譜函數(shù)定義為：

(12)

同樣，這三個函數(shù)τL(q)、αL(q)和fL(α)是相互聯(lián)系的。

3 方法驗證

3.1 De Wijs模擬數(shù)據(jù)

為了對比這兩種方法的分析效果，本文采用De Wijs模擬數(shù)據(jù)[18]驗證包括N=64和N=1028兩種不同尺度數(shù)據(jù)，生成算子為1.96、0.84、0.84、0.36。圖2是尺度指數(shù)和多重分形譜的估計值，對于兩種不同數(shù)據(jù)盒子估計法均獲得了比較的擬合，但是在譜值的首末端有點偏差，主要是由于邊緣效應(yīng)造成的。結(jié)果表明，對于小尺度的模擬數(shù)據(jù)，圖2A中與理論曲線相比，質(zhì)量指數(shù)曲線在左邊和右邊分別出現(xiàn)了上移和下移的情況，相應(yīng)的譜值左邊和右邊偏離標準的譜域，因為基于WLs的計算中對原小波系數(shù)進行挑選導(dǎo)致的發(fā)散，而基于盒子法的估計比較接近理論值，也顯示了盒子法對于小尺度的測度數(shù)據(jù)的穩(wěn)定性和準確性。同樣，值得注意的是，對于n=1的消失矩WLs估計，原理上等同于盒子法(盒子法：測度上的求和或者積分，一階消失矩：加權(quán)求和或者積分)，譜值的左邊現(xiàn)實了良好的逼近，即模擬數(shù)據(jù)中的小值模擬準確，而對于大值出現(xiàn)了偏差。對于大尺度多重分形分析，盒子估計法和WLs的基本都和理論曲線實現(xiàn)了基本重疊，表明對于大尺度數(shù)據(jù)，WLs和盒子估計都可以用來多重分形分析。總之，對于小尺度的逼近，由于分辨率有限，邊緣效應(yīng)，以及高階消失矩的過分平滑都會導(dǎo)致估計發(fā)散，而對于大尺度的數(shù)據(jù)，兩者都可以作為分析手段，對于消失矩，可以去2、3、4階的均值來近似估計。

3.2 地球化學數(shù)據(jù)驗證

本文地球化學數(shù)據(jù)驗證來自云南個舊錫礦水系沉積物元素Sn、Cu和Pb(圖3)，研究區(qū)域為29400 km2，，采樣間隔為2km，分辨率為2km×2km，N=76×50，分析小波為Daubechies，消失矩為Nψ=1,2,3,4。如圖4所示，對于Sn元素譜值估計，基于WLs不同的消失矩的計算，差別明顯，顯示了WLs估計得不穩(wěn)定性，而Cu和Pb的估計，在左邊還表示了部分的重疊，正則性方面強于Sn元素。另外，利用經(jīng)典級聯(lián)過程驗證多重分形數(shù)據(jù)的正則性，也可獲得同樣的結(jié)果，本文限于篇幅，另文闡述。通過地球化學元素數(shù)據(jù)的驗證，盒子估計法展現(xiàn)出來了優(yōu)良的穩(wěn)定性和準確性，對于小尺度數(shù)據(jù)WLs估計需慎重。從另外的角度來說，利用基于函數(shù)奇異性的這套框架某種程度上不符合以測度刻畫的地學數(shù)據(jù)，或許WLs優(yōu)勢在于震蕩的一維信號。

圖2 De Wijs模擬數(shù)據(jù)的質(zhì)量指數(shù)和分形譜(盒子法: “+”；WL:消失矩N=1,2,3,4，符號：“o”;“×”; “*”; “矩形”)

圖3 云南個舊水系沉積物Sn、Cu和Pb元素IDW插值圖(A-C: 分別為Sn、Cu和Pb元素濃度值)

圖4 地球化學元素數(shù)據(jù)的質(zhì)量指數(shù)和分形譜(盒子法: “+”；WL:消失矩N=1,2,3,4，符號：“o”; “×”; “*”; “矩形”)

4 結(jié)論

對于分辨率有限、大小且不規(guī)則的現(xiàn)實地學數(shù)據(jù)，基于WLs的尺度分析，在準確性和穩(wěn)定性方面都弱于盒子法，但是有時候尺度分析并不涉及到其多重分形譜的計算，僅僅利用系數(shù)與尺度之間的冪律關(guān)系，需要找到一個拐點，因此，小波系數(shù)不失為找不同尺度閾值一個好方法。對于小波消失矩的選擇，消失矩過小，捕捉奇異性的強度不夠，過大會導(dǎo)致分析的數(shù)據(jù)濾波過于平滑，從而導(dǎo)致譜值的兩端發(fā)散，因此，一般取消失矩的比較大的加權(quán)平均較合適。另外，基于WLs的估計涉及到小波變換的多分辨率分析，在圖像分割過程中，利用小波刻畫的奇異性進行分割有著內(nèi)在的優(yōu)越性，延伸繼續(xù)刻畫多重分形譜，正如前文引用，是一個有效的方法。

這篇文章主要討論了以積分為主的測度和以累積量為主的尺度分析，以經(jīng)典礦物分割模型de Wijs模型模擬數(shù)據(jù)為例，結(jié)果表明，對于盒子法的尺度分析，不管是小尺度還是大尺度，都表現(xiàn)出了更優(yōu)的結(jié)果，基于WLs的尺度分析對于大尺度的數(shù)據(jù)也表現(xiàn)出了良好的結(jié)果。此外，基于小波系數(shù)的WLs不僅可以為分析數(shù)據(jù)提供多分辨率分析，還可以檢測信號的振蕩器異性等行為?？傊环N是基于概率測度下的奇異性，一種基于連續(xù)函數(shù)的連續(xù)性，分別來源于不同的數(shù)據(jù)應(yīng)用，對于大分辨率的數(shù)據(jù)，兩者的刻畫能力相同，而對于小尺度的地學數(shù)據(jù)，盒子法的效果更優(yōu)。

[1] Mandelbrot BB. Intermittent turbulence in self similar cascades: Divergence of high moments and dimension of the carrier [J]. Fluid. Mech., 1974, 62: 331.

[2] Cheng Q. Multifractality and spatial statistics [J]. Computers & Geosciences, 1999, 25(9): 949-961.

[3] Cheng Q, Agterberg FP. Comparison between two types of multifractal modeling [J]. Mathmatical Geoglogy, 1996, 28(8):1001-1016.

[4] Halsey TC, Jensen MH, Kadanoff LP. Fractal measures and their singularities: The characterization of strange sets [J]. Phys. Rev., 1986, 33(2): 1141-1151

[5] Arneodo A, Bacry E, Muzy JF. The thermodynamics of fractals revisited with wavelets [J]. Physica A, 1995, 213: 232-275.

[6] Schertzer D, Lovejoy S, eds. Non-linear variability in geophysics [M]. The Netherlands: Kluwcr Acad. Publ., Dordrecht, 1991.

[7] Parisi G., Frisch U. On the singularity structure of fully developed turbulence, in Turbulence and Predictability in Geophysical Fluid Dynamics [M], edited by M. Ghil, R. Benzi, and G. Parisi, New York: Elsevier,1985.

[8] Muzy JF, Bacry E, Arneodo A.. Wavelets and multifractal formalism for singular signals-application to turbulent data [J]. Physical Review Letters, 1991, 67(25):3515-3518.

[9] Mallat S, Hwang WL. Singularity detection and processing with wavelets [J]. IEEE Trans. Inf. Theory, 1992, 38, 617-643.

[10] Jaffard S. Multifractal formalism for functions [J]. SIAM J. Math. Anal., 1997, 28(4): 944-998.

[11] Cheng QM, Agterberg FP, Ballantyne SB. the separation of geochemical anomalies from background by fractal methods [J]. Journal of geochemical exploration, 1994, 51(2):109-130.

[12] 郭科，施澤興，唐菊興等.用多重分形研究元素的共生組合.電子科技大學學報,2004,33(2)：221-224.

[13] Zheng Z, Mao H, Cheng QM. Fractal geometry of element distribution on mineral surfaces [J]. mathematical of geology, 2001, 33(2): 217-228.

[14] Dojnow P. Multifractal analysis of amplitude and phase components of EEG-based analytical signals [J]. Comptes Rendus De L Academic Bulgare Des Sciences, 2007, 60(10):1071-1076.

[15] Khalil A, Grant GL, Caddle Lb, et al. Chromosome territories have a highly nonspherical morphology and nonrandom positioning [J]. Chromosome Research, 2007, 15(7):899-916.

[16] Evertsz CJG., Mandelbrot BB. Multifractal measures, In: Peitgen, H.O., Jurgens, H., Saupe, D. (Eds.), Chaos and Fractals. Springer [M], New York: Springer, 1992.

[17] Wendt H, Roux SG., Jaffard S, et al. Wavelet leaders and bootstrap for multifractal analysis of images [J]. Signal Processing, 2008, 89(6): 1100-1114.

[18] De Wijs, HJ, Statistics of ore distribution [J], Geologic en Mijnbouw. 1951, 13 (8), 365-375.