基于s變換的presley非平穩(wěn)隨機過程演變功率譜密度

基于s變換的presley非平穩(wěn)隨機過程演變功率譜密度

1時間-尺度分析工具之二:st果變換隨機噪聲的研究通常涉及對隨機過程的描述。對于非平穩(wěn)隨機過程而言,其非平穩(wěn)性不僅依賴于時間的幅值,而且其頻率成分也是隨著時間變化的。觀察典型的地震動過程,其幅值表現為上升-平穩(wěn)-衰減的趨勢,此現象稱之為時域非平穩(wěn);而由于地震波在傳播過程中的行波效應或高頻成分不斷被吸收導致地震動時程越零率不斷減少,此現象稱之為頻域非平穩(wěn)。獲取非平穩(wěn)隨機過程的時間-頻率特性,成為描述和模擬隨機過程的重要環(huán)節(jié)。最早的平穩(wěn)隨機過程頻域描述為Fourier變換。因其基函數為諧和函數,決定了Fourier系數在物理意義上能解釋為頻率。由于借助諧和函數可為諸多物理和數學方程提供解析解,決定了Fourier變換在自然科學領域數百年時間內的重要地位。但Fourier系數卻不具有時間定位特征,或更具體地講,其時間特性湮滅在相位信息之中。1946年,Gabor提出了首個時間-頻率分析工具,即Gabor展開(Gaborexpansion),Gabor展開的系數可以通過短時Fourier變換STFT(ShortTimeFourierTransform)得到。其基本思想是,為待分析信號加窗,然后,對每個加窗的子信號進行Fourier變換。這樣,所得到的加窗分段Fourier變換便具有了時間信息。其缺點在于,當窗函數給定后,對于每個時間點和頻帶,其時間和頻率分辨率是無法改變的。Wigner和Ville將經典的平穩(wěn)隨機過程譜理論延伸于非平穩(wěn)隨機過程方面,提出了Wigner-Ville分布WVD(Wigner-VilleDistribution)。在這一時間-頻率分析工具中,時變譜是通過時變相關函數的Fourier變換得到的。一般來講,由于WVD的變換核仍為時間無限長的諧和函數,因此,從本質上無法合理反映局部的時間-頻率信息;而且,在某些情況下會有負值出現,因此很難在物理上給出合理的解釋。小波變換WT(WaveletTransform)的提出,克服了上述時間-頻率分析工具的缺點。小波變換是一種時間-尺度分析工具,依賴于尺度因子,每個小波相對于“母小波”進行拉伸或壓縮;依賴于平移因子,每個小波相對于“母小波”進行平移。因此,尺度因子提供了局部的頻率信息,而平移因子提供了局部的時間信息。但是,在特定尺度下,并非所有小波基函數都能與頻率顯示聯系,比如正交低階Daubechies小波,其在頻域內表現為非緊支且多峰的情況。因此,雖然小波分析能很自然地進行時間-尺度分析,但在利用其進行時間-頻率分析時,必須注意尺度與頻率之間的聯系。最近提出的一種時間-頻率分析工具,即Stockwell變換ST(StockwellTransform)可以視為WT與STFT的結合。其優(yōu)點在于:(1)提供了類似于WT的多尺度的恒Q(ConstantQuality)分析;(2)同時提供了STFT才具有的絕對相位信息。近年來,S變換已廣泛應用于地球物理,地震和醫(yī)學信號處理,電子和機械工程以及頻譜分析。關于S變換的研究,Gibson等分析了S變換與小波變換之間的區(qū)別及聯系;Pinnegar等對S變換提出了一些改進;Wong等提出了廣義S變換并對S變換的相位信息進行了進一步考察。在非平穩(wěn)隨機過程的演變功率譜估計方面,Iyama等和Basu等建立了不同尺度小波變換均方值與隨機過程時變譜之間的關系。Spanos等得到了諧和小波變換與非平穩(wěn)隨機過程演變功率譜密度之間的聯系;與此同時,Spanos等結合Priestley非平穩(wěn)隨機過程,利用一般非正交小波估計了非平穩(wěn)隨機過程的演變功率密度,此文的根本在于,相對于一般小波基函數,Priestley非平穩(wěn)隨機過程的調制函數為相對慢變函數。2009年,Huang等將這種方法推廣到多變量非平穩(wěn)隨機過程的情況。受文獻啟發(fā),本文提出了一種利用S變換估計Priestley非平穩(wěn)隨機過程演變功率譜密度的方法。方法的根本在于,相對于S變換的“變換核”,Priestley非平穩(wěn)隨機過程的調制函數為相對慢變函數。因此,非平穩(wěn)隨機過程的S變換可視為相位修正后的另一非平穩(wěn)隨機過程。進而推導出了以非平穩(wěn)隨機過程演變功率譜密度表達的在特定頻率點處的S變換瞬時均方值。考慮所有頻率點的表達式,且通過將前者寫為有限個頻率點的級數展開,則以S變換瞬時均方值顯式地表達非平穩(wěn)隨機過程演變功率譜密度的問題,轉化為了求解一組代數方程。此代數方程組的解為級數展開中的依賴于時間的未知項,非平穩(wěn)隨機過程的演變功率譜密度隨之而得。2速度s變換與stft的關系為了與S變換的原始文獻保持一致,本文中所述推導式中的頻率除非明顯標出,均為自然頻率。如果確定性過程為x(t)∈L2(R),則其S變換可以表達為式中f為頻率?？梢?S變換的核函數是經過高斯函數調制后的諧和函數,即與小波變換的相同之處在于,其高斯窗函數的寬度是隨著頻率而變化的,因此導致了其具有多尺度分辨率;與STFT的相同之處在于,其核函數是隨著平移因子τ而變化的,因而使其具有絕對相位信息。由于S變換可以看作為卷積形式,利用Fourier變換的卷積定理,可以得到S變換的頻域表達式為式中X(α)為確定性過程x(t)的Fourier變換?；谑?3),其離散S變換可表達為式中T為確定性過程x(t)的時間采樣間隔,N為其采樣總點數,Sj[T,n/NT]為變換的采樣點。研究表明,如采用快速Fourier變換(FFT)技術,能很大程度上加快S變換的計算速度?？梢宰C明,如果Sxφ(τ,f)為確定性過程x(t)的S變換,則其平移過程x(t-r)的S變換可以表達為S變換還具有其他重要性質,限于篇幅限制,這里只給出對于估計隨機過程功率譜重要的性質,其余性質可參見文獻。3子問題的估計考慮一個形如的非平穩(wěn)隨機過程,其中A(ω,t)為依賴于時間的慢變函數,Z(ω)為復值正交增量的隨機過程,即式中E為期望算子,SX-X-(ω)為相應的平穩(wěn)隨機過程的雙邊演變功率譜密度函數。因此,當調制項A(ω,t)為時間慢變函數時,非平穩(wěn)隨機過程X(t)的演變功率譜密度函數可以近似表達為下文中常把非平穩(wěn)隨機過程的功率譜近似理論表達(式(9))稱為目標功率譜,并與估計功率譜對比以驗證所提方法的合理性?？紤]某特定頻率fj處的S變換,并將式(6)代入式(1)可得根據Priestley非平穩(wěn)隨機過程模型(式(6))的定義,A(ω,t)相對于S變換核為時間慢變函數。具體而言,圖1所示為典型地震工程中使用的某調制函數與S變換核對比?？梢钥闯?當S變換頻率點選定后,在相應的S變換核有效積分區(qū)間內,調制函數可視為常值,因此可近似將此慢變函數移至積分號外,同時其時間變量取為高斯調制函數的對稱軸t=τ,即如將式(11)內部積分單獨考慮,顯見,式(12)即為h(t)=exp[i2πω(t-τ)]的S變換。利用S變換的平移特性(式(5))可知,如果u(t)=exp(i2πωt)的S變換已知,則h(t)的S變換亦可求得。利用S變換的頻域表達,u(t)的S變換可寫為進而,h(t)的S變換可以表示為將式(14)代入式(11)可得式中保留了dZ(ω)的正交特性,亦為正交隨機過程,即若令W(τ,fj)=S(τ,fj)exp(i2πfjπ),則W(τ,fj)亦為Priestley非平穩(wěn)隨機過程。因此,原非平穩(wěn)隨機過程在某特定頻率處的S變換的相位修正可認為是連續(xù)變量的非平穩(wěn)隨機過程,且此隨機過程的功率譜密度函數為因此,這一隨機過程的瞬時均方值為4基于時間條件的變分迭代式(18)建立了在某頻率點fj處原非平穩(wěn)隨機過程的S變換的均方值與該隨機過程演變功率譜密度的關系?？紤]若干頻率點上的S變換的相位修正項,式(18)可以寫為式中N為S變換的頻率點取樣總點數。如果將原非平穩(wěn)隨機過程寫為的級數展開。其中cj(τ)為某頻率點的時變系數。結合式(20)與式(19),可得式中E[|S(τ,f)|2]為N維列向量,c(τ)為N維列向量,Q為N×N維矩陣,即式中每個元素Qi,j為由此可知,通過解矩陣方程(21),可得時變系數c(τ)。因此,可由式(20)給出的非平穩(wěn)隨機過程的演變功率譜估計。在矩陣方程(21)中,形狀函數矩陣Q是不依賴時間τ的。因此,對于所有時間取樣點,此矩陣只需求取一次逆,且由于基于FFT技術的S變換的引入,故本文建議方法具有較高的計算效率。值得特別指出的是,上述求解過程并不需要隨機過程的時-頻功率譜函數的具體形式,因此,特別適合于利用實際觀測記錄估計非平穩(wěn)隨機過程的時-頻譜。當然,若預設某類時-頻譜形式,利用上述譜估計方法,可以通過求解反問題方式識別譜參數。有關研究工作,將另文發(fā)表。5機過程估計算法本節(jié)以兩種具有演變功率譜的非平穩(wěn)隨機過程,即均勻調制和非均勻調制的非平穩(wěn)隨機過程為例,利用本文所建議的算法,估計其演變功率譜密度。5.1平穩(wěn)過程/功率譜密度的定義考慮如下均勻調制的非平穩(wěn)隨機過程,即式中g(t)為只依賴于時間t的調制函數,為平穩(wěn)隨機過程,其雙邊功率譜密度為Kanai-Tajimi譜式中S0=1.0cm2s-3rad-1為基巖白噪聲輸入強度,ζs=0.4和ωs=20rads-1分別為場地土阻尼與卓越頻率系數。調制函數g(t)定義為如圖1所示,其中a=0.25,b=0.5,k為使gmax=1的正規(guī)化系數。據式(9)可知,均勻調制非平穩(wěn)隨機過程的目標演變功率譜密度可以表達為可以看出,不論在樣本層次上,還是在功率譜意義上,均勻調制非平穩(wěn)隨機過程的樣本/功率譜都能由相應平穩(wěn)過程的樣本/功率譜乘以僅依賴于時間的調制函數(/的平方)而得。因此,文獻[20,23,24]也將調制項A(ω,t)=g(t)的均勻調制隨機過程稱為可分離(Separable)非平穩(wěn)隨機過程或“準平穩(wěn)”(Quasi-stationary)隨機過程。在本文建議方法中,確定性樣本函數的生成為基于S變換的譜估計的基礎。在此,平穩(wěn)過程的樣本函數用譜表達方法生成,即式中Δω=ωu/n為頻率取樣間隔,ωu為上限截止頻率,φn為[0,2π]之間均勻分布的隨機變量。在譜表達方法中,樣本函數數目為500個,時間取樣點數目為1024,取樣間隔為Δt=0.02s,截止頻率取為Nyquist頻率ωu=π/Δtrad/s,諧和項數目為512。據第3節(jié)所述估計功率譜密度的方法(式(21)和式(20)),圖2所示為該非平穩(wěn)隨機過程在6s處的瞬時目標功率譜密度與瞬時估計功率譜密度對比,可以看出,二者吻合良好。為了更清楚地表現在其他時刻處的瞬時功率譜密度,圖3和圖4分別給出了目標演變功率譜密度和估計演變功率譜密度曲面圖。5.2非平穩(wěn)隨機過程的演化過程算例1為均勻調制的非平穩(wěn)隨機過程的功率譜密度函數,其頻譜特征并不隨時間變化,只是幅值隨著調制函數變化。因此,這種均勻調制的隨機過程并不符合實際地震動的時變頻譜特征。正如引言中所述,由于行波效應或傳播路徑中高頻成分不斷被吸收等因素,隨機地震動的瞬時頻譜特征是依賴于時間的。因此,利用本文建議方法估計具有演變譜隨機過程的功率譜密度函數,更具有實際意義。為此,考慮具有目標演變頻譜密度的非平穩(wěn)隨機過程為其演變功率譜密度如圖5所示?？梢钥闯?在不同的時間段,其功率譜密度不僅在幅值上,且其卓越頻率所在區(qū)間也隨時間變化,式(29)較為充分地反應了地震動的幅值及頻譜特點,即非平穩(wěn)不僅存在于幅值方面,而且存在于頻譜方面。由式(29)亦知其并不能表示為如式(27)所示的時間-頻率分離形式,因此,這種同時反應時間-頻率非平穩(wěn)特征的隨機過程,在諸多文獻中稱之為不可分離(Non-separable)非平穩(wěn)隨機過程。非平穩(wěn)樣本函數由譜表達方法合成,即式中SXX(t,kΔω)為時變的目標演變功率譜密度函數,其他諸參數與式(28)相同。圖6為利用本文建議方法所估計的演變功率譜密度的三維視圖;圖7為在4s和12s時刻的目標瞬時功率譜密度與估計瞬時功率譜密度對比?？梢钥闯?估計功率譜密度與目標功率譜密度吻合良好。5.3非平穩(wěn)隨機過程隨時間的變化顯然,調制函數相對于S變換核的慢變特性對于演變功率譜估計值有較大影響。但本文所建議方法,對于大多數工程實際中幅值非快變的非平穩(wěn)隨機過程均適用。另一方面,合適地選擇S變換頻率點使其變換核較之調制函數為快變亦至關重要,如圖1所示。如何定量描述由S變換核和調制函數相對變化快慢而導致的估計值偏差,可由后續(xù)研究考察。其次,由式(20)可知,非平穩(wěn)隨機過程瞬時功率譜密度展開為若干高斯函數平方的加權和。參與計算的頻率點fj決定了高斯函數的對稱軸位置及寬度。數值計算表明,合理選擇頻率點fj決定了估計瞬時功率譜的光滑程度。圖2所示估計功率譜中出現的若干小幅抖動實為式(20)中不同高斯函數峰值;圖7中4s瞬時功率譜峰值附近處,目標與估計功率譜的差別亦由此導致。最后值得注意的是,理論上,式(6)之非平穩(wěn)隨機過程的精確時變譜密度并非為式(9)之演變譜密度。考察Priestley非平穩(wěn)隨機過程可知,只有當式(6)中的調制項A(ω,t)為時間慢變函數時,其時變功率譜密度才能近似以式(9)表示。因此,Priestley定義的非平穩(wěn)隨機過程為一類幅值慢變的振蕩(Oscillatory)非平穩(wěn)隨機過程。只有當調制項A(ω,t)依時間慢變時,才不會對非平穩(wěn)隨機過程功率譜近似值(式(9))帶來時間-頻率域內的偏移(Shift)。換言之,當調制項A(ω,t)依時間非慢變時(其自身帶有頻率成分),須考慮由其導致的時間-頻率域內的偏移。對演變功率譜分析感興趣的讀者可參見文獻。6非平穩(wěn)隨機過程的s變換本文首先總結和回顧了各種時間-頻率分析工具,包括Fourier變換、短時Fourier變換、WignerVille分布以及小波變換。針對最近出現的S變換,做了詳細描述。結合非平穩(wěn)隨機過程的Priestley模型,提出