基于S型增長(zhǎng)曲線(xiàn)組合預(yù)測(cè)模型的滑坡變形預(yù)測(cè)研究

彭 鑫，賀小黑，賀鑫焱，彭必建

(1. 東華理工大學(xué)水資源與環(huán)境工程學(xué)院，330013，南昌；2. 中國(guó)人民武裝警察部隊(duì)研究院 工程設(shè)計(jì)研究所，100012，北京；3.北京國(guó)信華源科技有限公司，100055, 北京；4.云南地質(zhì)工程勘察設(shè)計(jì)研究院，650041，昆明)

0 引言

斜坡巖土體在構(gòu)造、地應(yīng)力等內(nèi)在因素控制作用及降雨、人類(lèi)工程活動(dòng)等外在誘發(fā)因素共同影響下，沿著某一軟弱結(jié)構(gòu)面產(chǎn)生成塊的滑移現(xiàn)象，稱(chēng)為滑坡。一旦下滑的滑體達(dá)到一定規(guī)模時(shí)，就可能造成其周邊影響范圍內(nèi)的重大經(jīng)濟(jì)損失以及危害人員生命安全[1]。為了盡可能避免滑坡的危害，根據(jù)監(jiān)測(cè)信息進(jìn)行滑坡的預(yù)測(cè)預(yù)報(bào)成為滑坡防治工程的重要一環(huán)。

滑坡的預(yù)測(cè)預(yù)報(bào)研究包括時(shí)間尺度上的失穩(wěn)時(shí)間預(yù)報(bào)和空間尺度上的變形趨勢(shì)預(yù)測(cè)等內(nèi)容[2]。學(xué)者們基于不同的理論對(duì)滑坡變形監(jiān)測(cè)數(shù)據(jù)分析采用了多種方法，其中，在滑坡變形預(yù)測(cè)方面，學(xué)者們?cè)诨诨依碚?、生長(zhǎng)曲線(xiàn)理論的統(tǒng)計(jì)型模型以及依據(jù)突變理論、協(xié)同理論等非線(xiàn)性理論的非線(xiàn)性模型投入了非常大的精力，取得了豐碩的研究成果，比較典型的變形預(yù)測(cè)模型有基于生物生長(zhǎng)曲線(xiàn)的S型生長(zhǎng)曲線(xiàn)擬合模型(Verhulst、Pearl)[3-4]；基于灰理論的GM(1,1)模型[5-6]；基于機(jī)器學(xué)理論的神經(jīng)網(wǎng)絡(luò)模型[7-8]、支持向量機(jī)等模型[9-10]等。

然而，有些模型方法過(guò)程復(fù)雜，可能需要較深的數(shù)學(xué)理論知識(shí)儲(chǔ)備(如協(xié)同理論[2,11]、混沌理論[12])，當(dāng)模型的掌握不夠深時(shí)容易造成模型預(yù)測(cè)精度不可控。Verhulst和Pearl模型在滑坡的變形預(yù)測(cè)和時(shí)間預(yù)報(bào)方面都有應(yīng)用，其數(shù)學(xué)理論相對(duì)簡(jiǎn)單，有高等數(shù)學(xué)基礎(chǔ)的人就能熟練運(yùn)用。因此，針對(duì)Verhulst和Pearl單一預(yù)測(cè)模型預(yù)測(cè)精度的不足，本文基于加權(quán)組合預(yù)測(cè)的思想[13]，采用最優(yōu)加權(quán)組合Verhulst和Pearl子模型對(duì)其進(jìn)行優(yōu)勢(shì)互補(bǔ)，提高滑坡預(yù)測(cè)模型的變形預(yù)測(cè)精度，并通過(guò)已有滑坡實(shí)例對(duì)其進(jìn)行驗(yàn)證。

1 滑坡的變形演化特征

通過(guò)總結(jié)大量的監(jiān)測(cè)位移-時(shí)間曲線(xiàn)的形態(tài)特征，許強(qiáng)[14]等人將滑坡分為：穩(wěn)定型滑坡、漸進(jìn)型滑坡、突發(fā)型滑坡(圖1)。對(duì)于突變型滑坡，有學(xué)者研究發(fā)現(xiàn)其一般具有較深的滑動(dòng)面，且一般會(huì)發(fā)生液化現(xiàn)象，如修德皓[15]、許強(qiáng)[14]等人通過(guò)對(duì)甘肅黑方臺(tái)滑坡研究發(fā)現(xiàn)，由于黃土固有的濕陷性，在降雨、地下水、灌水等因素下發(fā)生液化流動(dòng)現(xiàn)象，是典型的突變型滑坡。

大量的研究表明，滑坡大多都有一定的變形

圖1 不同滑坡演進(jìn)類(lèi)型的位移隨時(shí)間的變化

演化過(guò)程[16-17]，現(xiàn)有的預(yù)測(cè)預(yù)報(bào)模型也大多針對(duì)漸進(jìn)型這一類(lèi)滑坡。漸進(jìn)型滑坡的變形演化特征可分為緩慢變形、等速變形、加速變形及失穩(wěn)破壞4個(gè)階段[4]。晏同珍[18]、孫景恒[19]等學(xué)者認(rèn)為滑坡的孕育、生長(zhǎng)、成熟及消亡過(guò)程具有生物生長(zhǎng)曲線(xiàn)類(lèi)似的機(jī)制：滑坡的緩慢變形及等速變形階段，滑坡主要以蠕滑為主，相當(dāng)于生物生長(zhǎng)的孕育階段，加速變形階段，滑帶土的內(nèi)摩擦角、抗剪強(qiáng)度不斷弱化，曲線(xiàn)斜率逐漸增大，相當(dāng)于生長(zhǎng)模型曲線(xiàn)的生長(zhǎng)階段，滑坡進(jìn)一步發(fā)展，變形急劇加速，曲線(xiàn)斜率呈陡崖式發(fā)展，然后趨于穩(wěn)定，相當(dāng)于生長(zhǎng)模型曲線(xiàn)的成熟階段[4,19]。用生物增長(zhǎng)曲線(xiàn)取擬合滑坡歷史監(jiān)測(cè)數(shù)據(jù)對(duì)滑坡進(jìn)行預(yù)測(cè)預(yù)報(bào)是可行的。

2 基于滑坡演化過(guò)程的預(yù)測(cè)模型

2.1 灰色Verhulst模型基本原理

Verhulst模型是1987年德國(guó)物理學(xué)家發(fā)現(xiàn)的一種生物生長(zhǎng)模型，由于滑坡的位移特征與生物的生長(zhǎng)規(guī)律類(lèi)似，很多學(xué)者用Verhulst 模型對(duì)滑坡進(jìn)行預(yù)報(bào)研究[3]。Verhuslt模型的白化微分方程形式為：

(1)

式中：a1，b1為系數(shù)。

對(duì)于一組等時(shí)距非負(fù)增量位移監(jiān)測(cè)序列x(1)=[x(0)(1),x(0)(2),…,x(0)(n)]，經(jīng)一次累加后得到原始累計(jì)位移監(jiān)測(cè)序列x(1)=[x(1)(1),x(1)(2),…,x(1)(n)]。

對(duì)累加序列x(1)作緊鄰均值生成序列z(1)，

式中

(2)

構(gòu)造數(shù)據(jù)矩陣，用最小二乘法求解系數(shù)a1，b1，計(jì)算表達(dá)式為

(3)

(4)

2.2 以累積位移為參量的Pearl預(yù)測(cè)預(yù)報(bào)模型

Pearl預(yù)測(cè)預(yù)報(bào)模型最早由孫景恒提出并應(yīng)用在新灘滑坡和意大利Vaiont滑坡監(jiān)測(cè)數(shù)據(jù)上，取得了好的成果，少量學(xué)者利用Pearl模型開(kāi)展了滑坡變形預(yù)測(cè)預(yù)報(bào)的研究，均取得不錯(cuò)的進(jìn)展。Pearl生長(zhǎng)曲線(xiàn)的數(shù)學(xué)模型表達(dá)式為：

(5)

式中：k為常數(shù)；f(t)為自變量t的多項(xiàng)式。

f(t)=a0+a1t+…

(6)

一般多項(xiàng)式階數(shù)取1，則其數(shù)學(xué)表達(dá)式轉(zhuǎn)換為一般logistic函數(shù)的表達(dá)式：

(7)

式中：a、b、c為待擬合系數(shù)；y為t時(shí)刻的位移擬合值。

因此，只需確定好a、b、c3個(gè)參數(shù)即可對(duì)滑坡進(jìn)行擬合回歸預(yù)測(cè)。采取非線(xiàn)性擬合求參的方法[4]，通過(guò)matlab中的Levenberg-Marquardt算法進(jìn)行非線(xiàn)性擬合來(lái)確定a、b、c3個(gè)參數(shù)。

2.3 基于Verhulst增長(zhǎng)曲線(xiàn)和Pearl增長(zhǎng)曲線(xiàn)的組合預(yù)測(cè)模型

由于這些曲線(xiàn)各有不同，對(duì)某一滑坡的預(yù)測(cè)結(jié)果可能偏高或者偏低。因此，本文引入組合預(yù)測(cè)的思想，基于最小二乘原理將最優(yōu)權(quán)重與前述2種單一模型組合在一起，通過(guò)對(duì)單一模型的取長(zhǎng)補(bǔ)短以提高預(yù)測(cè)的準(zhǔn)確性。

基于最優(yōu)加權(quán)算法的組合模型計(jì)算流程如下。

令一組原始監(jiān)測(cè)數(shù)據(jù)表示為xt,可用m個(gè)模型來(lái)擬合，擬合結(jié)果表示為xit，其中i=1,2,…,m；t=1,2,…,N。單個(gè)模型的權(quán)重可表示為wi，其滿(mǎn)足下述要求：

(8)

組合預(yù)測(cè)模型的擬合值可表示為：

(9)

令eit為預(yù)測(cè)模型i在時(shí)刻t的擬合殘差，則組合預(yù)測(cè)模型的擬合殘差可表示為：

(10)

基于殘差平方和最小的原則，通過(guò)最小二乘方法，最優(yōu)權(quán)重wi可在下列約束條件下獲取。

(11)

3 實(shí)例分析

3.1 模型輸入數(shù)據(jù)的選取

依據(jù)文獻(xiàn)[21]提供的臥龍寺新滑坡位移監(jiān)測(cè)資料進(jìn)行預(yù)測(cè)預(yù)報(bào)，監(jiān)測(cè)數(shù)據(jù)見(jiàn)表1。新滑坡發(fā)生時(shí)間為1971年5月5日，因此，取最后2 d監(jiān)測(cè)數(shù)據(jù)為驗(yàn)證數(shù)據(jù)用以預(yù)測(cè)，選取3月15日至5月3日經(jīng)過(guò)平滑處理提取趨勢(shì)項(xiàng)位移，其具有類(lèi)似齋藤室內(nèi)實(shí)驗(yàn)?zāi)M滑坡累計(jì)位移的“三段式”特征。通過(guò)正態(tài)檢驗(yàn)及相關(guān)性分析方法[5]綜合判斷得出1971年4月22號(hào)為滑坡體從等速變形階段進(jìn)入加速變形階段的臨界點(diǎn)，故選取1971年4月22日至1971年5月3日的監(jiān)測(cè)數(shù)據(jù)進(jìn)行預(yù)測(cè)預(yù)報(bào)。

3.2 單模型的預(yù)測(cè)分析

利用Matlab將前述Verhulst模型擬合選取的監(jiān)測(cè)數(shù)據(jù)，通過(guò)式(2)～(3)得出模型參數(shù)a1=-0.204 6、b1=-0.009 5；將其代入式(4)，發(fā)現(xiàn)

表1 臥龍寺新滑坡監(jiān)測(cè)數(shù)據(jù)

當(dāng)i=2時(shí)，其相對(duì)誤差最小。利用Matlab將前述Pearl模型擬合選取的監(jiān)測(cè)數(shù)據(jù)，依據(jù)擬合優(yōu)度最優(yōu)原則得出模型參數(shù)a=0.004 974、b=0.398 9、c=-0.349 9，計(jì)算結(jié)果見(jiàn)表2，擬合曲線(xiàn)與監(jiān)測(cè)數(shù)據(jù)曲線(xiàn)見(jiàn)圖2。

表2 單預(yù)測(cè)模型擬合結(jié)果

圖2 單模型擬合曲線(xiàn)與原始監(jiān)測(cè)曲線(xiàn)的關(guān)系

圖2虛線(xiàn)左側(cè)是根據(jù)選取的建模數(shù)據(jù)得到的擬合值，虛線(xiàn)右側(cè)是根據(jù)建模數(shù)據(jù)計(jì)算的預(yù)測(cè)值。從圖2可以看出，根據(jù)選取的監(jiān)測(cè)數(shù)據(jù)(39—50 d)進(jìn)行擬合得到曲線(xiàn)與實(shí)測(cè)曲線(xiàn)相差不大，Pearl模型擬合的曲線(xiàn)較為精確，但兩個(gè)模型在外推最后2 d得到的預(yù)測(cè)結(jié)果相差較大(一個(gè)較實(shí)測(cè)值偏大而另一個(gè)則偏小)，可以看出，往外預(yù)測(cè)期數(shù)越多，誤差越大。由第52天預(yù)測(cè)結(jié)果可知，Verhulst模型的預(yù)測(cè)結(jié)果可能會(huì)起到一個(gè)提前觸發(fā)臨滑預(yù)警預(yù)報(bào)的效果，而Pearl模型則可能不會(huì)觸發(fā)臨界失穩(wěn)預(yù)警預(yù)報(bào)。

3.3 組合預(yù)測(cè)模型預(yù)測(cè)分析

在得到前述單預(yù)測(cè)模型的擬合數(shù)據(jù)的基礎(chǔ)上，根據(jù)式(8)～(12)構(gòu)造基于Verhulst子模型和Pearl子模型的組合預(yù)測(cè)模型，由于本文只組合2個(gè)子模型，可令Verhust子模型的權(quán)重為w，則Pearl子模型的權(quán)重為1-w，計(jì)算得到的單預(yù)測(cè)模型權(quán)重見(jiàn)表3，單模型、組合模型的預(yù)測(cè)值及原始累計(jì)位移值之間的關(guān)系見(jiàn)圖3。

表3 S型增長(zhǎng)曲線(xiàn)子模型的權(quán)重

圖3 組合模型與單模型預(yù)測(cè)結(jié)果對(duì)比

從圖3可以看出，組合模型通過(guò)最優(yōu)權(quán)重組合各子模型的優(yōu)勢(shì)，使得擬合值最大限度“逼近”實(shí)測(cè)值，尤其從第50期數(shù)據(jù)之后(模型外推預(yù)測(cè)值)更為明顯，其預(yù)測(cè)值相較單模型預(yù)測(cè)值更加符合實(shí)際。

為了檢驗(yàn)單模型及組合模型的預(yù)測(cè)精度，采取擬合優(yōu)度R2及均方根誤差RMSE 2個(gè)指標(biāo)進(jìn)行分析。分析結(jié)果見(jiàn)表4。

表4中的數(shù)據(jù)從擬合區(qū)精度、預(yù)測(cè)區(qū)精度以及綜合精度3個(gè)方面對(duì)比分析，發(fā)現(xiàn)在選取建模數(shù)據(jù)基礎(chǔ)上Pearl模型的擬合效果最好(擬合區(qū))，但是在預(yù)測(cè)區(qū)域根據(jù)建模數(shù)據(jù)建立的擬合回歸方程計(jì)算的預(yù)測(cè)值效果最差，而Verhulst模型在擬合區(qū)的擬合效果比Pearl模型差，但在預(yù)測(cè)區(qū)域卻比Pearl模型的要好，雖然它的均方根誤差也好大(表4，預(yù)測(cè)區(qū)90.456 1)，但是從圖3中可以看出Verhulst模型能起到提前預(yù)警預(yù)報(bào)的作用(假設(shè)第51天或第52天的數(shù)據(jù)值為滑坡失穩(wěn)時(shí)設(shè)定的預(yù)警閾值)。這表明單一預(yù)測(cè)模型根據(jù)滑坡監(jiān)測(cè)歷史數(shù)據(jù)擬合得到的回歸方程可能不能很好地反映未來(lái)的滑坡變形趨勢(shì)，有些模型可能會(huì)造成預(yù)警設(shè)備漏報(bào)或不報(bào)而滑坡已經(jīng)發(fā)生的現(xiàn)象(如本文中得到的Pearl預(yù)測(cè)模型)。而通過(guò)最優(yōu)加權(quán)組合單一模型的方法可以避免這一不足，擬合優(yōu)度R2從單一模型的0.934、0.923提高到0.985，均方根誤差RMSE從單一模型的4.69、4.72降低到1.29(表4中綜合精度數(shù)據(jù))。表明該組合預(yù)測(cè)模型優(yōu)于上述的單一預(yù)測(cè)模型，基于最優(yōu)加權(quán)的組合預(yù)測(cè)模型改進(jìn)效果非常明顯。

表4 單模型及組合模型的精度性能對(duì)比

4 結(jié)論

為了更好地服務(wù)于滑坡監(jiān)測(cè)預(yù)警工作，基于提高預(yù)測(cè)精度的目的，本文采取基于最優(yōu)加權(quán)的組合預(yù)測(cè)方法，通過(guò)最優(yōu)權(quán)重值組合S型增長(zhǎng)曲線(xiàn)模型方法中用于滑坡預(yù)警預(yù)報(bào)的Verhulst模型、Pearl模型2種常用模型。將其應(yīng)用于臥龍寺新滑坡監(jiān)測(cè)資料，可以得出以下結(jié)論。

1)單一預(yù)測(cè)模型根據(jù)建模數(shù)據(jù)擬合出來(lái)的回歸方程，其擬合效果與外推的預(yù)測(cè)效果并不一致，擬合效果好的Pearl預(yù)測(cè)模型其外推的預(yù)測(cè)值的效果相較Verhulst模型的要差，且其變形趨勢(shì)呈近似穩(wěn)定的線(xiàn)性緩慢上升趨勢(shì)，而Verhulst模型的外推變形趨勢(shì)呈線(xiàn)和實(shí)測(cè)曲線(xiàn)近似的急劇上升趨勢(shì)。

2)基于最優(yōu)加權(quán)算法的組合預(yù)測(cè)模型充分利用了Pearl模型、Verhulst模型的優(yōu)勢(shì)，彌補(bǔ)了Verhulst模型在數(shù)據(jù)擬合、Pearl模型在外推預(yù)測(cè)方面的不足，且在外推預(yù)測(cè)區(qū)的變形趨勢(shì)跟實(shí)測(cè)曲線(xiàn)接近一致。組合預(yù)測(cè)模型能起到提高預(yù)測(cè)精度的作用。