基于二值Logistic回歸模型的砂土液化指標敏感性分析

潘 翔

（周口師范學(xué)院基建處，河南 周口 466001）

砂土液化是地震災(zāi)害的一種主要形式，常伴有噴砂、冒水、地面沉陷、誘發(fā)高速滑坡等現(xiàn)象，是地基失效并造成結(jié)構(gòu)物損傷的重要原因之一［1］。20 世紀60 年代以來，幾次大的地震，如日本新瀉地震（1964）、美國阿拉斯加地震（1964）、中國唐山大地震（1976）、日本神戶大地震（1995）以及中國汶川特大地震（2008）等均因砂土液化導(dǎo)致堤防或建筑物大規(guī)模的破壞，造成了巨大的經(jīng)濟損失，因此砂土液化問題已引起了工程界的廣泛關(guān)注［2］。目前，對砂土液化已經(jīng)發(fā)展出多種預(yù)測方法，如標準貫入法（SPT）、靜力觸探法（CPT）、剪切波速法（vs），它們是最常用的砂土液化評價方法，也是許多砂土液化預(yù)測方法建立的基礎(chǔ)［3］。其中，標準貫入法較為成熟，已成為國內(nèi)外工程界普遍接受的方法，我國《建筑抗震設(shè)計規(guī)范》和美國Seed都建議采用以SPT-N 值為指標的經(jīng)驗方法進行砂土液化預(yù)測［4］。

由于砂土液化受眾多因素影響，同時各因素與砂土液化之間存在著高度的非線性關(guān)系，因此建立多因素綜合模型對砂土液化進行準確預(yù)測對抗震減災(zāi)具有重要意義。據(jù)此，國內(nèi)外的學(xué)者提出了多種綜合預(yù)測模型和方法，如神經(jīng)網(wǎng)絡(luò)模型［4-5］、模糊綜合評價法［6］、距離判別法［7］等，但他們在建立預(yù)測模型的同時，缺乏對變量的重要性和敏感性等方面的研究。由于不同的變量具有不同的敏感性，它們對于砂土液化的影響是不同的，因此可以利用現(xiàn)有數(shù)據(jù)提供的信息去分析變量的敏感性，在此基礎(chǔ)上選擇變量建立更為合理的預(yù)測模型。砂土液化的預(yù)測結(jié)果一般分為兩種情況：液化（1）和未液化（0），而砂土液化指標一般都為連續(xù)變量，在這種情況下，二值Logistic回歸模型可以用來作為砂土液化預(yù)測的基本模型。本文根據(jù)國內(nèi)外173組樣本數(shù)據(jù)［8］，建立了砂土液化預(yù)測的二值Logistic回歸模型，在此基礎(chǔ)上利用優(yōu)勢比的概念對砂土液化指標的敏感性進行了分析。

1 二值Logistic回歸模型

1.1 數(shù)學(xué)模型

現(xiàn)實中的很多現(xiàn)象可以劃分為兩種可能，或者歸結(jié)為兩種狀態(tài)，這兩種狀態(tài)分別用0和1表示，如砂土液化（1）和未液化（0）兩種狀態(tài)。如果采用多個因素對0和1表示的砂土液化現(xiàn)象進行因果關(guān)系解釋，則可以利用二值logistic回歸的數(shù)學(xué)模型來表示，其表達式為［9］

式中：x1，…，xm為自變量；β0，…，βm 為待定參數(shù)，根據(jù)已有觀測值進行估計。

對于給定的自變量x1，…，xm的值，其概率表達式為

由（2）式可以預(yù)測事件“0”和“1”的概率。一般情況下，若p＞0.5，則確定為事件“1”；若p＜0.5，則確定為事件“0”。

1.2 優(yōu)勢比

在Logistic回歸模型中，設(shè)p 表示砂土液化的概率，則砂土液化的概率p 與砂土未液化的概率1-p之比為優(yōu)勢（Odds），有

固定其他自變量，比較xi與xi+1的ln（Odds）變化。

對于xi，有

對于xi+1，有

由（4）、（5）兩式可得

經(jīng)反對數(shù)變換得到

式中：ORi為第i個指標對應(yīng)的優(yōu)勢比。

由此，系數(shù)βi 取值大小及符號說明變量xi對砂土液化發(fā)生的影響大小和方向，具體地說：若βi＞0，則ORi=exp（βi）＞1，即變量xi的水平每增加一個單位，砂土液化發(fā)生的相對危險度便提高exp（βi）倍，且xi處于高水平時，砂土液化發(fā)生的概率更大；若βi＜0，則ORi=exp（βi）＜1，即變量xi的水平每增加一個單位，砂土液化發(fā)生的相對危險度便降低exp（βi）倍，且xi處于低水平時，砂土液化發(fā)生的概率更大；若βi=0，則ORi=exp（βi）=1，即變量xi的各水平下砂土液化發(fā)生的相對危險度一致，表明該變量對砂土液化發(fā)生沒有影響。若x1，…，xm的量綱不同，則可將數(shù)據(jù)標準化，得到標準化系數(shù)，標準化系數(shù)的絕對值越大，表示其對事件發(fā)生的影響也越大。

1.3 模型評價

評價二值Logistic回歸方程的擬合效果一般采用以下幾個指標［9］：

（1）模型擬合顯著性檢驗的p 值，如果p＜0.05，說明回歸方程的擬合效果是顯著的；

（2）擬合優(yōu)度檢驗中對應(yīng)于Pearson、偏差、Hosmer-Lemeshow 3 種方法 的p 值，如 果p ＞0.05，說明回歸方程的擬合效果是好的；

（3）自變量系數(shù)顯著性檢驗的p 值，如果p＜0.05，說明自變量系數(shù)是顯著的；

（4）事件預(yù)測準確率，根據(jù)回歸方程可以預(yù)測事件的概率，一般情況下，若p＞0.5，則確定為事件“1”，若p＜0.5，則確定為事件“0”，在此基礎(chǔ)上與實際情況對比即可計算事件預(yù)測準確率作為擬合效果評估的依據(jù)。

2 砂土液化指標敏感性分析

2.1 指標和數(shù)據(jù)

本文選擇Stark等［10］研究的180組砂土液化與未液化現(xiàn)場實例中的173組數(shù)據(jù)作為分析樣本，其中7組由于缺失平均粒徑數(shù)據(jù)未入選。173組數(shù)據(jù)中，107組發(fā)生砂土液化，66組未發(fā)生砂土液化。在這些實例中，沙涌、建筑物破壞或下沉、地面?zhèn)认虼笞冃蔚缺灰曌魃巴烈夯l(fā)生的標志，缺乏上述標志則認為砂土液化未發(fā)生。限于篇幅，本文并未列出具體數(shù)據(jù)，可參見文獻［10］。

原始數(shù)據(jù)中統(tǒng)計了震級M、深度ds、地下水位dw、總應(yīng)力σ0、有效應(yīng)力σ′0、平均粒徑D50、修正的CPT 錐尖阻力qc1、地面峰值加速度amax和地震剪應(yīng)力比SSR。其中，qc1和SSR 的計算公式如下：

上式中：Cq為有效應(yīng)力修正系數(shù)，根據(jù)Kayen等［11］提供的公式進行計算；qc為CPT 錐尖阻力值；amax為地面峰值加速度；g 為重力加速度（9.81 m／s2）；σ0為總應(yīng)力（MPa）；σ′0為有效應(yīng)力（MPa）；rd為深度折減因子，可由下式計算［11］：

其中，z為土層埋深（m）。

在上述這些指標中，部分指標之間具有高度相關(guān)性，如σ0和σ′0，amax和SSR。對于二值Logistic回歸模型而言，如果指標之間高度相關(guān)，則可能導(dǎo)致模型不穩(wěn)定，解釋上存在矛盾等［12］。例如σ0和σ′0、amax和SSR 同時存在于二值Logistic回歸模型中，采用上述優(yōu)勢比分析，所反映的結(jié)果是σ0和σ′0或者amax和SSR 與砂土液化潛勢之間的趨勢關(guān)系恰好相反，這是不符合實際的。因此，需對原始數(shù)據(jù)進行相關(guān)性分析，并選擇相對獨立且具有代表性的指標做進一步的敏感性分析。本文在對原始數(shù)據(jù)進行統(tǒng)計分析時，對于原始數(shù)據(jù)中的區(qū)間值，取所在區(qū)間的中值代替。原始數(shù)據(jù)相關(guān)性分析結(jié)果見表1。

表1 原始數(shù)據(jù)相關(guān)性分析Table 1 Correlation analysis of data

由表1 可以看出：M、D50、qc13 個指標 相對獨立，而ds、σ0與σ′0以及amax與SSR 之間高度相關(guān)，相關(guān)性系數(shù)大于0.9，而dw與σ′0之間的相關(guān)性較好，相關(guān)性系數(shù)為0.773。根據(jù)各指標的代表性和獨立性，本文從上述9個指標中選擇M、σ′0、D50、qc1和SRR 作為基本分析指標。

2.2 模型建立

以173組數(shù)據(jù)建立二值Logistic回歸模型，具體操作可通過SPSS統(tǒng)計分析軟件實現(xiàn)，得到二值Logistic回歸的數(shù)學(xué)模型如下：

由此概率表達式為

2.3 模型評價

對建立的二值Logistic回歸方程的擬合效果進行評價，其評價結(jié)果如下：

（1）模型擬合顯著性檢驗的p=0.000＜0.001，說明回歸方程的擬合效果非常顯著；

（2）擬合優(yōu)度檢驗中對應(yīng)于Pearson、偏差、Hosmer-Lemeshow 3種方法的p 值分別為0.840、1.000和0.068，三者都滿足p＞0.05，說明回歸方程的擬合效果是好的；

（3）自變量系數(shù)顯著性檢驗的p 值中，M 對應(yīng)的p=0.017＜0.05，SRR 對應(yīng)的p=0.005＜0.05，說明兩者的系數(shù)是顯著的，其余3個指標對應(yīng)的p=0.000＜0.001，說明三者的系數(shù)非常顯著；

（4）由式（12）預(yù)測事件“1”砂土液化發(fā)生的準確率為95.3%，預(yù)測事件“0”砂土液化未發(fā)生的準確率為83.3%，平均預(yù)測準確率為90.8%，預(yù)測效果良好。

綜上所述，本文所建立的二值Logistic回歸模型能夠較好地解釋砂土液化現(xiàn)象之間的因果關(guān)系，因此可用二值Logistic回歸模型進行指標敏感性分析，也可用于砂土液化預(yù)測。

2.4 指標敏感性分析

本文利用優(yōu)勢比分析各指標的敏感性以及其與砂土液化潛勢的趨勢性關(guān)系。由于各指標的意義和量綱不同，為了使優(yōu)勢比具有可比性，需對指標數(shù)據(jù)進行標準化處理后再建立二值Logistic回歸模型，而指標數(shù)據(jù)的標準化并不影響模型的評價和預(yù)測結(jié)果。

首先對原始數(shù)據(jù)進行均值化處理，這樣不僅可以消除數(shù)據(jù)量綱和數(shù)量級的影響，同時也可以保持數(shù)據(jù)內(nèi)的變異信息［13］。但由于均值化后的數(shù)據(jù)過小，因此將所有數(shù)據(jù)增大10倍處理，均值化公式如下：

式中：i為指標編號（i=1，2，…，m）；j為樣本編號（j=1，2，…，n）。

然后通過SPSS統(tǒng)計分析軟件，得到標準化處理后的二值Logistic回歸的數(shù)學(xué)模型如下：

模型的系數(shù)和對應(yīng)的優(yōu)勢比匯總于表2中。

表2 模型的系數(shù)和對應(yīng)的優(yōu)勢比Table 2 Coefficients of each model and corresponding odds ratios

由表2可以看出：未標準化的模型系數(shù)和對應(yīng)的優(yōu)勢比差異相當(dāng)大，根本不具有可比性，無法說明指標的敏感性；而標準化后的模型系數(shù)和對應(yīng)的優(yōu)勢比處于同一個數(shù)量級，具有可比性，可根據(jù)各指標的優(yōu)勢比來分析其對砂土液化的影響大小和方向。σ′0和qc1的系數(shù)小于零，對應(yīng)的優(yōu)勢比小于1，說明在其他指標不變的情況下，這兩個指標任意一個增大，將會導(dǎo)致砂土液化的相對危險性減小，以qc1敏感性較大；M、D50、SRR 的系數(shù)大于零，對應(yīng)的優(yōu)勢比大于1，說明在其他指標不變的情況下，這三個指標任意一個增大，將會導(dǎo)致砂土液化的相對危險性增大，以M 敏感性最大。對于σ′0、qc1、M 和SRR 與砂土液化的趨勢性關(guān)系符合一般認識，但是一般認為平均粒徑D50越大，則抗液化能力越強［14］，而優(yōu)勢比給出的是D50越大，砂土液化的相對危險性越大，因此這兩者的統(tǒng)計關(guān)系還有待進一步驗證。

將σ′0和qc1優(yōu)勢比求倒數(shù)（表2中括號內(nèi)數(shù)值）即可與其他指標進行對比，由此可以得到各指標敏感性排序為：M＞qc1＞σ′0＞D50≈SRR，也就是說砂土液化對M 指標最為敏感，對qc1指標次之，說明qc1指標是除M 以外最為有效的預(yù)測指標，其次σ′0、D50和SRR 三者之間的敏感性相差不大。敏感性排序結(jié)果對應(yīng)于指標的優(yōu)選順序。由于5個指標的優(yōu)勢比或者優(yōu)勢比倒數(shù)都與1具有一定的差距，說明它們都是預(yù)測砂土液化不可或缺的指標。

3 結(jié)論

本文將二值Logistic回歸模型應(yīng)用到砂土液化預(yù)測中，以國內(nèi)外173組樣本數(shù)據(jù)作為研究對象，對砂土液化指標的敏感性進行了分析，所得結(jié)論如下。

（1）通過原始數(shù)據(jù)相關(guān)性分析選擇M、σ′0、D50、qc1和SRR5個指標建立二值Logistic回歸模型，模型評價結(jié)果都滿足要求，預(yù)測準確率達到了90.8%，說明二值Logistic回歸模型可用于描述指標與砂土液化現(xiàn)象之間的因果關(guān)系，其概率表達式可用于砂土液化預(yù)測。

（2）砂土液化指標敏感性分析結(jié)果表明：各指標的敏感性排序為M＞qc1＞σ′0＞D50≈SRR；σ′0和qc1指標與砂土液化的相對危險性呈負相關(guān)趨勢，以qc1指標敏感性較大，M、D50、SRR 指標與砂土液化的相對危險性呈正相關(guān)趨勢，以M 指標敏感性最大。

［1］陳建偉，蘇經(jīng)宇，蘇幼坡.場地地震液化的主要影響因素分析［J］.安全與環(huán)境工程，2013，20（4）：18-22.

［2］劉勇健.基于聚類-二叉樹支持向量機的砂土液化預(yù)測模型［J］.巖土力學(xué)，2008，29（10）：2764-2768.

［3］Lai S Y，Hsu S C，Hsieh M J.Discriminant model for evaluating soil liquefaction potential using cone penetration test data［J］.Journal of Geotechnical and Geoenvironmental Engineering，2004，130（12）：1271-1282.

［4］陳國興，李方明.基于徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)模型的砂土液化概率判別方法［J］.巖土工程學(xué)報，2006，28（3）：301-305.

［5］林志紅，項偉.基于貝葉斯正則化BP神經(jīng)網(wǎng)絡(luò)的砂土地震液化研究［J］.安全與環(huán)境工程，2011，18（2）：23-27.

［6］薛新華，張我華，劉紅軍.砂土地震液化的模糊綜合評判法［J］.重慶建筑大學(xué)學(xué)報，2006，28（1）：55-58.

［7］金志仁.基于距離判別分析方法的砂土液化預(yù)測模型及應(yīng)用［J］.巖土工程學(xué)報，2008，30（5）：776-780.

［8］Xue，X H，Yang，X G.Application of the adaptive neuro-fuzzy inference system for prediction of soil liquefaction［J］.Natural Hazards，2013，67（2）：901-917.

［9］吳令云，吳家祺，吳誠鷗，等.MINITAB軟件入門：最易學(xué)實用的統(tǒng)計分析教程［M］.北京：高等教育出版社，2012.

［10］Stark T D，Olson S M.Liquefaction resistance using CPT and field case histories［J］.Journal of Geotechnical Engineering，1995，121（12）：856-869.

［11］Kayan R E，Mitchell J K，Seed R B，et al.Evaluation of SPT-，CPT-，and Shear Wave-basedMethodsfor Liquefaction Potential Assessments Using Lorna Prieta Data［R］.Buffalo，USA：National Center for Earthquake Engineering Research，1992：177-192.

［12］襲炯良，郟劍寧，張揚.主成分改進的logistic回歸模型方法在流行病學(xué)分析中的應(yīng)用［J］.中國熱帶醫(yī)學(xué)，2005，5（2）：207-209.

［13］葉宗裕.關(guān)于多指標綜合評價中指標正向化和無量綱化方法的選擇［J］.浙江統(tǒng)計，2003（4）：24-25.

［14］羅戰(zhàn)友，龔曉南.基于經(jīng)驗的砂土液化灰色關(guān)聯(lián)系統(tǒng)分析與評價［J］.工業(yè)建筑，2002，32（11）：36-39.