基于Copula模型的多維地震動參數相關性分析

劉明旭, 姚國文, 彭剛輝, 姜云木, 喻宣瑞, 宋安祥, 靳紅華

(1.重慶交通大學土木工程學院, 重慶 400074; 2.重慶市黔江區城市管理局, 重慶 409000;3.大連理工大學海岸和近海工程國家重點實驗室, 大連 116024)

地震動是典型的多維運動,工程上將其劃分為含兩個水平向分量和一個豎向分量。地震動的多維特性易對地震作用下的結構動力響應和破壞機理產生顯著影響,因此,細化研究地震動多維特性是合理開展結構抗震設計與設防的基礎[1]。

諸多學者針對多維地震動參數相關性領域開展了研究。黃玉平等[2]研究發現,多維地震動水平向分量的相關性較高,且平均自功率譜的形狀大致相同,互動率譜的形狀也與自功率譜相近。劉章軍等[3]基于實測記錄給出了多維地震動各方向分量間的調制函數和功率譜參數比值。Bozorgnia等[4]提出了一種可以預測地震動豎向分量與水平向分量加速度反應譜比值的模型。薛素鐸等[5]根據地震動波動理論總結了水平向分量與豎直向分量的場地卓越圓頻率以及譜強度因子之間的關系。劉子心等[6]引入隨機方法,實現了多維地震動的降維模擬。雖然現階段關于多維地震動各向分量間的相關性研究結果較多,但少有考慮多維地震動數據的概率信息和特征。

Copula理論最早應用于金融界,近年來被逐漸應用于土木工程領域。姜云木等[7]通過Copula理論分析了主余震的地震動演變功率譜相關參數,推導了主余震參數間的計算公式,并構建了主余震型地震動的降維模型。陳建兵等[8]借助Vine Copula理論對99條混凝土受壓本構試驗曲線展開研究,分析了混凝土抗壓強度、峰值應變、彈性模量的邊緣分布情況和參數間的相關性,為混凝土結構非線性研究提供了重要基礎。唐小松[9]將Copula函數應用于水利水電巖土工程,提出巖土體參數聯合概率分布模型構建方法,建立了有限數據條件下相關非正態巖土體參數二維和多維分布模型。朱瑞廣等[10]利用Copula函數對662條主余震序列地震動展開研究,建立了主余震地震動強度參數間的聯合分布,得到Copula條件均值可以用來預測余震地震動強度參數取值的結論。陳將宏等[11]利用Copula函數研究了抗剪強度參數相關變量的相關性,并給出了聯合概率模型。李幸鈺等[12]提出一種新型的優化藤Copula,并用于橋梁結構系統易損性方法研究,避免了大量聯合概率密度函數的求解以及信息判定準則的試算。孟祥成等[13]提出了高維 Gumbel Copula 參數的擬合估計公式,并將其嵌入推廣的遞推分解算法,形成大型生命線網絡地震動力可靠度的計算方法。

基于此,為了在考慮隨機特性和概率特性影響下,研究u-v(u、v為地震動兩個水平向分量和一個豎向分量中的任意兩個分量,u、v=x,y,z)向分量間地震動參數的相關性,首先從太平洋地震動工程研究中心篩選出1 500組實測地震動,并從強度、持時和頻譜3個方面篩選出12組地震動參數;然后,通過計算12組地震動參數的Pearson線性相關系數、Kendall秩相關系數和Spearman秩相關系數,分析u-v(u、v=x,y,z)向分量間地震動參數的相關性。其次,以對數正態分布(LOGN)、耿貝爾分布(G)、廣義極值分布(GEV)、威布爾分布(WEI)為備選概率模型,并借助柯爾莫哥洛夫-斯米爾諾夫(Kolmogorov-Smirnov,K-S)檢驗和貝葉斯信息準則(the Bayesian information criteria,BIC準則),建立12組地震動參數在x、y、z向分量上的最優概率模型;最后,通過Copula理論和貝葉斯公式,確定u-v向分量間12組地震動參數的最優Copula函數,建立它們之間的聯合概率函數,并在給定u向分量地震動參數條件下,得到v向分量地震動參數預測值。以期為地震動多維特性研究分析、結構抗震設防提供參考價值。

1 基于Copula模型的多維地震動參數相關性分析方法

1.1 歸一化處理

在開展相關性分析前,為增強研究結果可靠性,先將1 500組實測地震動的12組地震動參數進行了歸一化處理,如式(1)所示。

(1)

1.2 計算相關系數

選取Pearson線性相關系數、Kendall秩相關系數和Spearman秩相關系數用于判斷u-v(u、v=x,y,z)向分量地震動參數間的強弱關系。所用的λ′u和λ′v分別為實測地震動u向分量和實測地震動v向分量的歸一化地震動參數。

Pearson線性相關系數r又稱積差相關[14],可反映兩組變量λ′u和λ′v間的相關程度,如式(2)所示。

(2)

式(2)中：cov(λ′u,λ′v)為變量λ′u和λ′v之間的協方差;σλ′u、σλ′v分別為變量λ′u和λ′v的方差。

Kendall秩相關系數τ又稱為等級相關系數[15],如式(3)所示。兩組變量λ′u和λ′v的集合為D=[λ′uλ′v],當λ′u,i>λ′u,j且λ′v,i>λ′v,j或λ′u,i<λ′u,j且λ′v,i<λ′v,j時,稱為同序對,其數量用A表示;當λ′u,i<λ′u,j且λ′v,i>λ′v,j或λ′u,i>λ′u,j且λ′v,i<λ′v,j時,稱為異序對,其數量用B表示,其中,i,j=1,2,…,n。

(3)

式(3)中：n為樣本數量,即n=1 500。

Spearman秩相關系數ρ是一個非參數性質(與分布無關)的秩統計參數[16],主要用于解決稱名數據和順序數據問題,如式(4)所示。

(4)

1.3 選取備選概率模型

得到u-v(u、v=x,y,z)向分量間地震動參數的相關系數后,為進一步分析多維地震動參數的概率信息和特征,選取對數正態分布(LOGN)、耿貝爾分布(G)、廣義極值分布(GEV)、威布爾分布(WEI)作為備選概率模型。

對數正態分布在短期上與正態分布相似,對數正態分布具有更大向上波動的可能和更小向下波動的可能,其地震動參數概率密度函數為

(5)

式(5)中：x為地震動參數;μ1、λ1分別為對數正態分布下地震動參數的位置參數、尺度參數。

Gumbel分布可表示為極值I型分布[17],其地震動參數概率密度函數為

(6)

式(6)中：μ2、λ2分別為Gumbel分布下地震動參數的位置參數與尺度參數。

GEV分布是具有三參數的極值分布通式,對各種數據均有較強的適應能力[18],其地震動參數概率密度函數為

(7)

式(7)中：μ3、λ3、k3分別為廣義極值分布下地震動參數的位置參數、尺度參數和形狀參數。

Weibull分布最常用于可靠性數據建模[19],其地震動參數概率密度函數為

(8)

式(8)中：k4、λ4分別為Weibull分布下地震動參數的形狀參數和尺度參數。

1.4 確定最優概率模型

為從上述4種備選概率模型中篩選出最優概率模型,將利用K-S檢驗對其進行驗證。若有兩種及兩種以上備選概率模型通過K-S檢驗時,則采用BIC準則[20]選出最優概率模型,如式(9)所示。

BIC=ln(n)b-2lnL

(9)

式(9)中：BIC為備選概率模型的后驗概率值,一般后驗概率值越小,代表模型擬合度越好;b為分布模型中參數的個數;L為模型最大似然函數;n為樣本容量。

1.5 建立聯合分布函數

Copula理論最先由Sklar[21]提出,可用來描述任意兩組隨機變量間的相關性。假設有兩組任意尺寸相等的一維隨機變量X和Y,則X和Y的聯合累積分布函數為

F(x,y)=C[FX(x),FY(y)]=C(u1,u2)

(10)

式(10)中：C(·)為Copula分布函數,其中u1=FX(x)、u2=FY(y);FX(x)、FY(y)分別為隨機變量X和Y的邊緣累積分布函數。

因此,X和Y的聯合概率密度函數為

f(x,y)=c(u1,u2)fX(x)fY(y)

(11)

式(11)中：fX(x)、fY(y)分別為隨機變量X和Y的邊緣概率密度函數;c(·)為Copula密度函數,可表示為

(12)

1.6 選取Copula函數備選模型

采用Gumbel Copula、Clayton Copula、Frank Copula和Plackett Copula對u-v(u、v=x,y,z)向分量間的地震動參數進行分析。這4種阿基米德Copula函數均具有不同的尾部相關性。其中Gumbel Copula具有上尾相關性,適合描述兩個變量同時上漲的情況;Clayton Copula具有下尾相關性,適合描述兩個變量同時下降的情況;Frank Copula可以同時擬合上、下尾相關且可以同時描述變量正相關性與負相關性;Plackett Copula對變量間相關性的正負無要求。4種Copula函數備選模型如表1[22]所示。

表1 Copula函數備選模型[22]Table 1 Alternative model of Copula function[22]

表1中,對以上4種阿基米德Copula函數模型采用BIC準則檢驗,便可選出u-v(u、v=x,y,z)向分量間同一組參數的最優Copula模型。

1.7 得到Copula條件均值和條件隨機數

在前文基礎上,根據貝葉斯條件概率公式,若給定隨機變量X的條件下,則可得到隨機變量Y的條件概率密度函數為

(13)

條件均值可進一步表示為

(14)

同時,對式(13)進行積分,并對其進行概率反變換,便可得到給定X條件下Y的條件隨機數。變換過程如式(15)所示。

(15)

式(15)中：θ為區間[0,1]的隨機變量。

因此,Y的條件隨機數如式(16)所示。

(16)

2 選取實測多維地震動和參數

2.1 實測多維地震動

以斷層距離大于10 km和實測震級大于5為選取原則,從太平洋地震動工程研究中心篩選1 500組實測地震動,每一組地震動均有3個方向分量地震動。并根據《中國地震動參數區劃圖》(GB 18306—2015)[23]建議的5種場地類型(I0～IV),對實測地震動記錄進行了分類。

《中國地震動參數區劃圖》(GB 18306—2015)[24]與Eurocode8：DesignofStructuresforEarthquakeResistance[25]、MinimumDesignLoadsforBuildingsandOtherStructures(ASCE/SEI 7-10)[26]的場地劃分標準差異較大,為了將外國實測地震動與中國場地類型相匹配,結合文獻[27-28]給出5種場地類型(I0～IV)與地下30 m平均剪切波速Vs30之間的對應關系和每類場地對應的實測地震動數量,如表2所示。在分析實測地震動前,對實測地震動進行四階Butterworth濾波處理,并在1%～99%能量范圍內截取,以增強研究可信度。

表2 場地分類與Vs30的對應關系Table 2 Corresponding relationship between site classification and Vs30

2.2 多維地震動參數

從強度、持時、頻譜3個方面選取12組地震動參數,如表3[29]所示。

表3 多維地震動參數[29]Table 3 Multi-dimensional ground motion parameters[29]

3 結果與分析

3.1 不同向分量間地震動參數的相關系數

以II類場地為例,通過式(1)～式(4),得到u-v(u、v=x,y,z)向分量間12組地震動參數的3種相關系數變化趨勢,如圖1所示。

圖1 不同向分量間地震動參數的相關系數變化趨勢

由圖1可知,除反應譜峰值對應周期地震動參數的相關系數約為0.5,其余u-v(u、v=x,y,z)向分量間地震動參數相關性均較好,基本維持在0.8～1.0。但x-z和y-z向分量間的阿里亞斯強度地震動參數相關系數較低,平均保持在0.5。

工程上常見的II類、III類場地不同向分量間地震動參數的相關系數計算結果,如表4所示。

表4 不同向分量間地震動參數的相關系數Table 4 Correlation coefficient of ground motion parameters between different directional components

由表4可知,III類場地與II類場地的地震動參數相關性結論基本保持一致。即u-v(u、v=x,y,z)向分量間地震動參數的相關性均較好;水平向分量間反應譜峰值對應周期地震動參數的相關性較差;豎直向分量間阿里亞斯強度和反應譜峰值對應周期地震動參數的相關性較差。

3.2 多維地震動參數概率模型

以II類場地峰值加速度(PGA)和5%～95%有效持時為例,通過式(5)～式(8),得到這兩個地震動參數的概率密度函數(probability density function,PDF),如圖2、圖3所示。

圖2 峰值加速度概率密度模型

圖3 5%～95%有效持時概率密度模型

由圖2、圖3可知,從4種備選概率模型中無法直接篩選出最優概率模型。因此需采用K-S檢驗對其進行驗證,若有多種備選概率模型通過K-S檢驗時,則采用BIC準則確定出最優概率模型。因此得到x、y、z向地震動參數最優概率模型的分布情況,如表5所示。鑒于數據結果較多,只給出了工程上常見的II類、III類場地的多維地震動參數最優概率模型分布情況。

表5 多維地震動參數最優概率模型的分布情況Table 5 Distribution of the probabilistic model of optimal multi-dimensional ground motion parameters

由表5可知,多維地震動參數最優概率模型以對數正態分布、廣義極值分布和韋布爾分布為主;x向和y向大多數地震動參數的最優概率模型均為同一種概率模型,且分布參數也大致相同,但是z向參數的最優分布與水平向有較大差異。

上述結果表明,多維地震動兩個水平向分量具有高度相似性,而豎向分量與水平向分量相比,豎向分量的能量更小、持時更長、頻率成分更豐富。

3.3 不同向分量間地震動參數的最優Copula模型

同樣,以II類和III類場地為例,u-v(u,v=x,y,z)向分量間12組地震動參數的最優Copula模型及其參數信息,如表6所示。

表6 不同向分量間地震動參數的最優Copula模型Table 6 Optimum Copula model of ground motion parameters between different directional components

由表6可知,不同向分量間地震動參數的最優Copula模型大多為Plackett Copula和Gumbel Copula為主。根據表6結果,得到II類場地x-z向分量間的PGA和5%～95%有效持時這兩個地震動參數的條件概率密度函數,如圖4所示。

Td為x向地震動分量條件下5%～95%能量的有效持時

由圖4可知,在x向分量地震動參數不同取值的條件下,z向分量對應的地震動參數條件概率密度函數差異較大,這充分說明有必要考慮條件概率。

3.4 不同向分量間地震動參數的預測

通過式(14)得到II類場地峰值加速度(PGA)和5%～95%有效持時這兩個地震動參數的Copula條件均值,如圖5、圖6所示。

圖5 峰值加速度的Copula條件均值

圖6 5%～95%有效持時的Copula條件均值

由圖5、圖6可知,在給定x向或y向地震動參數的情況下,Copula條件均值與y向或z向地震動參數觀測值較為一致。即Copula條件均值可實現多維地震動參數的預測。

應用Monte Carlo隨機數法,通過式(16)得到II類場地PGA和5%～95%有效持時這兩個地震動參數的Copula條件隨機數,如圖7、圖8所示。

圖7 峰值加速度的Copula條件隨機數

圖8 5%～95%有效持時的Copula條件隨機數

由圖7、圖8知,Copula條件隨機數也可以較好的預測不同向分量間地震動參數。但與Copula條件均值相比,Copula條件隨機數可以更好體現地震動參數的隨機性。

4 結論

從太平洋工程抗震研究中心選取了1 500組實測地震動,并基于Copula模型對不同向分量間地震動參數的相關性展開研究,得到以下主要結論。

(1)整體來說,12組地震動參數在u-v(u、v=x,y,z)向分量間的相關性均較好,只有反應譜峰值對應周期地震動參數的相關性較差,其相關系數約為0.5;x-z向和y-z向分量間的阿里亞斯強度和反應譜峰值對應周期地震動參數相關性也較差,其相關系數約為0.5。

(2)K-S檢驗和BIC準則可用于擬合多維地震動參數的最優概率模型,其擬合結果精度高。

(3)在u分量地震動參數條件下,v分量對應地震動參數的條件概率密度函數與單獨v分量對應的地震動參數條件概率密度函數差異顯著,該現象表明考慮條件概率是有必要的。

(4)Copula條件均值與Copula條件隨機數,可以在給定u向分量地震動參數條件下,實現v向分量地震動參數的預測。