基于BOX-COX變換的橋梁結構地震易損性分析

張鵬輝， 郭軍軍， 周連緒， 袁萬城

(土木工程防災國家重點實驗室(同濟大學)， 上海 200092)

理論易損性曲線作為確定橋梁結構潛在地震風險的重要工具，建立了地震動強度指標(Intensity Measure, IM)與結構在某一損傷極限狀態(Limit State, LS)下失效概率之間的聯系[1]，而云圖法由于其計算效率上的優越性已經成為獲取結構理論易損性曲線的主要途徑之一[2-3]。然而以往的研究發現云圖法在建立概率地震需求模型(Probabilistic Seismic Demand Model, PSDM)時，為簡化計算所作出的三個基本假設常常與現實相違背，存在以下三個問題，即：(1) 特定地震動強度水平下的工程需求參數(Engineering Demand Parameter, EDP)不服從對數正態分布；(2) ln(IM)與ln(EDP)間不滿足線性相關；(3) ln(EDP)的方差隨IM取值變化[4-6]。Karamlou等[5]考慮結構不確定性，通過大量數值模擬研究了上述三個問題對理論易損性曲線計算結果的影響，研究表明與問題(1)相比后兩者對易損性曲線計算結果的影響更嚴重。為改善問題(2)，Padgett等[7]考慮高效性、實用性、充分性以及災害的可計算性對IM的選取進行優化，Pan等[8]采用二次函數對ln(IM)與ln(EDP)進行回歸，Mackie等[9]采用雙線性函數對ln(IM)與ln(EDP)進行回歸，均取得了較好的效果。此外，為了避免云圖法三個基本假設對易損性曲線計算結果的影響，研究人員通過增量動力分析(Incremental Dynamic Analysis, IDA)方法獲取易損性曲線，但該方法存在對地震動時程的放縮且需要的時程分析次數劇增等問題，故應用受到限制[10]。如何改進理論易損性分析方法，使其能夠不受上述三個基本假設的限制，同時又不增加時程分析次數是當前理論易損性分析亟待解決的問題。

本文旨在不顯著增加計算量的前提下提高易損性分析結果的可靠性，通過引入BOX-COX變換使概率地震需求模型滿足對數正態分布、線性和同方差性假設，結合蒙特卡洛抽樣得到結構易損性曲線，并以一座三跨連續梁橋為例，驗證該方法的有效性。

1 基于BOX-COX變換的地震易損性分析方法

1.1 云圖法的基本假設和BOX-COX變換

Cornell等[11]假設EDP的均值(SD)與IM之間滿足對數線性關系，并建立普通最小二乘(Ordinary Least Square, OLS)回歸模型

(1)

式中，α1、α2均為回歸系數。因此式(1)需服從OLS回歸中的對e作出的正態性和同方差性假設，由此得到概率地震需求模型

(2)

式中：D為構件的需求;Φ(·)為標準正態分布函數；βD為構件需求的對數標準差。

(3)

若進一步假定構件在某一損傷極限狀態(LSj)的抗震能力C服從均值為SC對數標準差為βC的對數正態分布，則云圖法得到的構件在LSj下的失效概率為

(4)

由此可見，云圖法中所作的三個基本假設，本質是為OLS回歸模型式(1)服務的。當ln(SD)與ln(IM)之間關系非線性，或回歸的殘差e不再服從正態性和同方差假設時，回歸系數α1、α2的估計不再是最小方差線性無偏估計，式(1)的應用效果將不理想，進而影響構件易損性計算的可靠性[12]。

BOX-COX變換由Box等[13]提出，在回歸分析領域應用廣泛，通過對因變量y=ln(EDP)進行適當的變換，可在不丟失任何信息的前提下，使回歸模型滿足線性、正態性和同方差假設

(5)

其中a為使y+a的各個分量都大于0的系數，當min(y)>0時a=0；當min(y)≤0時a=roundup(|min(y)|)，roundup(·)為向上取整函數[14]。引入BOX-COX變換后，式(1)更新為BOX-COX回歸模型

(6)

(7)

記x=ln(IM)，那么y=(y1,y2,…,yn)的聯合概率密度函數可表示為

f(y)=

(8)

對式(8)兩邊取對數可得極大似然函數

(9)

考慮到式(9)求偏導較為困難，無法得到λ的解析解，因此本文以0.01作為搜索步長，找到使L(λ)達到最大值的λ作為其數值最優解。

1.2 基于BOX-COX變換和蒙特卡洛抽樣的地震易損性分析

(10)

則構件在LSj下的失效概率為

(11)

值得指出的是，由于蒙特卡洛抽樣所增加的計算時間是有限的，相較于時程分析來說可忽略不計。

2 算例分析

2.1 有限元模型

本文以一座預應力混凝土連續梁橋為例，如圖1所示。跨徑布置為60 m+100 m+60 m，主梁采用C50混凝土的預應力箱梁；支座為球鋼支座，支座布置如圖1所示。滑動支座摩擦因數取0.02，屈服位移取2 mm[17]，屈服力為最大靜摩擦力；橋墩采用C40混凝土，鋼筋種類為HRB335，縱向鋼筋配筋率0.8%，體積配箍率1%；以梅花形布置的直徑2 m鉆孔灌注樁作為基礎；兩側橋臺均為樁接蓋梁橋臺。使用OpenSees程序建立全橋非線性有限元模型(如圖2)。主梁采用彈性梁柱單元模擬；橋墩采用彈塑性纖維單元模擬，混凝土本構采用Kent-Park-Scott模型[18]，不考慮混凝土的抗拉性能，鋼筋采用考慮強化的雙線性本構。支座、樁基、碰撞均采用零長單元模擬，采用M法計算土彈簧六個方向的剛度，碰撞單元采用Muthukumar提出的本構模型[19]，伸縮縫采用碰撞單元模擬，碰撞單元的初始縫隙根據橋臺處伸縮縫縫寬取為15 cm。由于橋臺對構件的地震響應有顯著影響[20]，本文參考Nielson給出的橋臺模擬方法[21]，被動土壓力方向的剛度由樁基和橋臺共同提供，主動土壓力方向剛度僅由樁基提供，以零長單元進行模擬。

圖1 橋梁結構總體布置圖(m)

圖2 有限元模型示意圖

2.2 地震動時程

地震動時程的選取是易損性分析的關鍵之一，Baker等[22]提出的非條件譜選波方法能夠根據地震動預測模型(Ground Motion Prediction Model, GMPM)得到的場地均值譜和方差選取符合目標反應譜概率分布的一組波。采用該選波方法，首先按焦馳宇[23]給出的中美剪切波速換算公式得到Ⅲ類場地30 m深度范圍內的土層平均剪切波速V30；然后根據Boore等[24]提出的地震動預測模型生成適用于中國區域震級M=7，斷層距Rjb=15 km，V30=205 m/s的反應譜概率分布；最后從NGA-West2強震數據庫中按M=6～8，Rjb=0～100 km，V30=150～260 m/s選出100條符合目標反應譜分布的原場單向地震動分量。圖3畫出了它們的加速度反應譜，圖4給出了這些地震動時程的震級和震源距分布直方圖。將這100條地震動時程按縱橋向輸入，得到結構各構件的響應。

圖3 地震波的加速度反應譜(阻尼比5%)

圖4 震級和震源距分布直方圖

2.3 概率地震需求模型建立與比較

Riddell[25]指出對于短周期結構選用加速度相關型的地震動強度指標(IM)是合適的，Moschonas等[26-27]采用PGA(Peak Ground Acceleration)作為IM對一系列中小跨徑梁橋和大跨度鋼構橋進行易損性分析，并取得了較好的效果。本文模型橋一階周期0.95 s，結構周期較短；同時PGA的地震危險性曲線目前已經能夠得到，即具有良好的災害可計算性，而結構相關型的地震動強度指標Sa(T1)的地震危險性曲線目前還未給出，故選用PGA(單位：g)作為地震動強度指標。連續梁橋的易損部位為橋臺、支座和墩柱，故選用被動土壓力方向橋臺位移δp、主動土壓力方向橋臺位移δa、滑動支座位移δb、固定墩墩頂漂移率Dr作為工程需求參數[28]。將時程分析結果，分別按式(1)和式(6)進行OLS回歸和BOX-COX回歸得到概率地震需求模型。

為評價殘差e的正態性，作出回歸殘差的頻率分布直方圖并得到其核密度估計與按正態分布擬合的概率密度函數，核密度估計曲線(Kernel Density Estimation Curve, KDE)與正態分布概率密度曲線(Norm)的偏差情況反映了殘差e的分布與正態分布的偏離程度。自變量與因變量的線性相關性強弱通過回歸決定系數(R2)來衡量，它是兩變量線性相關系數的平方。異方差性的診斷采用等級相關系數法，計算等級相關系數(rs)首先需要對回歸模型(1)和(6)中的自變量xi和殘差絕對值|ei|分別進行排序，則

(12)

其中di對應于xi和|ei|等級的差數。對rs進行顯著性檢驗得顯著性水平p，一般認為p<0.05時認為回歸模型存在異方差性。表1列出了OLS回歸和BOX-COX回歸的結果。

表1 OLS回歸和BOX-COX回歸的結果

由表1結果可知，當以δp、δa和δb作為EDP建立概率地震需求模型時，均存在顯著的異方差性，通過BOX-COX變換后，異方差問題得到消除。同時BOX-COX回歸的決定系數較OLS回歸平均提升了3.26%，δp的決定系數增加達到8.2%，說明自變量和因變量的線性相關程度有所提高。為直觀顯示殘差的非正態和異方差問題，圖5作出了δp回歸殘差的分布圖和頻率分布直方圖以及KDE曲線和Norm曲線。由圖像可看出OLS回歸的殘差在PGA較小時方差較小，隨著PGA的增大逐漸增加，即OLS回歸將方差作為一個定值過程中在PGA較小時對方差進行了高估，在PGA較大時低估了方差，且殘差的分布與正態分布具有一定程度的偏離，進行BOX-COX變換后，殘差的正態性得到改善，對δa和δb分析所得結果相同。需要說明的是，當應用云圖法以Dr作為EDP時，檢驗結果發現其概率地震需求模型滿足自身的三個基本假設。

(a) 基于OLS回歸的PSDM

(b) OLS回歸殘差圖

(c) 基于BOX-COX回歸的PSDM

(d) BOX-COX回歸殘差圖

2.4 構件易損性曲線比較

Hwang等[29]建議將橋梁的損傷描述為四個極限狀態，即輕微損傷(LS1)、中等損傷(LS2)、嚴重損傷(LS3)、完全損傷(LS4)。參考前人的研究結果，對應不同損傷極限狀態，本文采用的能力均值和對數標準差見表2。

表2 不同損傷極限狀態的能力均值和標準差

根據式(11)可求得各構件在不同損傷極限狀態下的易損性曲線，并與云圖法得到的構件易損性曲線進行比較，如圖6所示。由圖6(a)～6(c)可知，因為OLS回歸在PGA較小時高估了方差，導致在該PGA范圍內云圖法得到的失效概率偏大，最大偏差2.21%；而當PGA較大時OLS回歸低估了方差，故云圖法得到的失效概率偏小，最大偏差27.2%。且極限狀態對應的損傷程度越高，采用云圖法得到該極限狀態下的易損性曲線偏離實際越明顯。當PGA=1g時，δp達到LS3狀態、δa達到LS2狀態的和δb達到LS2狀態的概率分別被低估18.1%、524.7%、30.9%，若仍然采用云圖法的三個基本假設，將導致易損性計算結果失真。由圖6(d)可知，當采用云圖法能夠滿足自身基本假設時，基于BOX-COX變換的易損性分析方法(采用BOX-COX回歸建立需求模型)所得易損性曲線與云圖法(采用OLS回歸建立需求模型)計算結果相同，即BOX-COX變換不會導致信息的丟失。

3 結 論

(1) 提出了結合BOX-COX變換和蒙特卡洛抽樣的地震易損性分析方法，消除由于云圖法的三個基本假設帶來的易損性計算誤差，并通過算例驗證了其應用于連續梁橋地震易損性分析的有效性。

(a) 被動土壓力方向橋臺位移

(b) 主動土壓力方向橋臺位移

(c) 滑動支座位移

(d) 固定墩墩頂漂移率

(2) 采用BOX-COX變換建立的概率地震需求模型，無須遵循云圖法所特有的三個基本假設，且改善了云圖法特有的三個基本假設帶來的易損性計算誤差，使模型的線性程度提高，正態性改善，異方差性消除；同時，保持了云圖法非線性時程分析次數少的優點。

(3) 在云圖法適用的情況下，采用本文提出的方法計算的易損性曲線與云圖法結果一致；當概率地震需求模型存在異方差性時，使用云圖法將導致結構易損性曲線計算結果的失真，且極限狀態對應的損傷程度越高，該極限狀態下的易損性曲線偏離實際越明顯。