基于GRA-GPR的滾石運動特征參數(shù)敏感性分析及預(yù)測模型

趙天豪,馬春輝,楊 杰,程 琳

(1. 西安理工大學(xué) 水利水電學(xué)院,陜西 西安 710048; 2. 西安理工大學(xué) 省部共建西北旱區(qū)生態(tài)水利國家重點實驗室,陜西 西安 710048)

0 引言

滾石是指個別塊石因某種原因從地質(zhì)體表面失穩(wěn)后經(jīng)過下落、回彈、跳躍、滾動或滑動等運動方式中的一種或幾種的組合沿著坡面向下快速運動,最后在較平緩的地帶或障礙物附近靜止的動力演化過程[1]。滾石災(zāi)害雖發(fā)生規(guī)模較小,但由于其具有泛生性、突發(fā)性和隨機性的特點,其所造成的生命、財產(chǎn)損失與大規(guī)模滑坡災(zāi)害基本相當(dāng)。隨著我國國民經(jīng)濟的高質(zhì)量發(fā)展和西部大開發(fā)戰(zhàn)略的實施,水利、公路和建筑等基礎(chǔ)建設(shè)范圍日益廣泛,所面臨的滾石災(zāi)害問題愈發(fā)凸顯[2-3]。因此,明確合理、有效的支擋結(jié)構(gòu)設(shè)計方案和位置,是防止?jié)L石災(zāi)害發(fā)生、減少災(zāi)害損失的有效手段之一[4],而滾石的運動特征參數(shù)是制定支擋結(jié)構(gòu)方案的關(guān)鍵指標,直接關(guān)系到滾石防護治理工程的安全性、經(jīng)濟性和有效性。由于邊坡滾石的影響因素眾多,滾石的運動特征參數(shù)具有明顯的隨機性和模糊性[5],目前缺少易用性高和精確度高的邊坡滾石運動特征參數(shù)預(yù)測方法[6]。

國內(nèi)外專家學(xué)者針對滾石的運動模型運動參數(shù)開展了大量的研究:在滾石的運動參數(shù)方面,黃潤秋等[6]通過現(xiàn)場滾石試驗,分析了滾石形狀對滾動摩擦系數(shù)的影響;AZZON等[7]應(yīng)用數(shù)學(xué)模型進行落石的研究;DUSSAUGE等[8]研究了3個實際落石的堆積,并采用冪律分布擬合落石體積分布數(shù)據(jù);秦志英等[9]通過單球碰撞系統(tǒng)的數(shù)值仿真,從精度、效率和微觀接觸過程等方面對各種碰撞模型進行了比較;何思明等[10]研究了滾石法向和切向碰撞恢復(fù)系數(shù)的計算模式與計算公式,并闡明了影響碰撞恢復(fù)系數(shù)的主要因素。在滾石運動預(yù)測方面,吳順川等[11]提出了落石參數(shù)的正交試驗確定方法,及基于概率分析的落石預(yù)測分析思路;吳琴等[12]通過分析崩塌落石的運動形式和受力情況,對災(zāi)害的危險性進行預(yù)測;王學(xué)良等[13]對山區(qū)輸變電工程中滾石災(zāi)害的識別與預(yù)測方法進行了研究;郭亮等[14]通過室內(nèi)模型試驗考察基-覆界面傾角、傍依坡角和堆積高度對震后崩塌堆石體穩(wěn)定性的影響,揭示失穩(wěn)前兆特征及變形破壞模式;何宇航等[15]對九寨溝景區(qū)懸溝危巖體的威脅范圍進行了預(yù)測;熊紹真等[16]采用3DEC離散元數(shù)值模擬方法,對高陡斜坡在地下開采作用下崩塌所產(chǎn)生的機理、失穩(wěn)模式和運動軌跡進行了全過程模擬;魏新容等[17]基于顆粒接觸理論,并考慮了斜坡滾石碰撞過程中的隨機因素,建立了斜坡滾石速度預(yù)測模型;HUANG等[18]利用K近鄰算法預(yù)測落石跳動范圍,根據(jù)與坡腳的距離將落石跳動范圍分類為不同的子區(qū)間;HUANG等[19]將邏輯模型樹與裝飾、隨機子空間和旋轉(zhuǎn)森林集成技術(shù)相結(jié)合,比較混合方法的預(yù)測性能以對不同類型的滑坡進行敏感性評估;ZHOU等[20]引入SVM支持向量機和LM-BP神經(jīng)網(wǎng)絡(luò)建立斜坡特征落石預(yù)測系統(tǒng),并通過工程實例對預(yù)測系統(tǒng)進行了驗證;劉福臻等[21]基于GIS環(huán)境下三維落石模擬軟件Rockfall Analyst,對物源區(qū)塊石進行了運動軌跡、速度和彈跳高度的數(shù)值模擬。上述研究在滾石問題中均取得了重要進展,但對于滾石運動特征參數(shù)的研究目前多針對單一因子,尚未有不同因素的敏感性研究,且缺少科學(xué)、快捷的滾石運動特征參數(shù)預(yù)測系統(tǒng)。

綜上,本文將灰色關(guān)聯(lián)法與滾石運動數(shù)值模擬相結(jié)合,以滾石的水平運動距離、彈跳高度和沖擊運動能量為目標指標,分析評價不同滾石運動特征影響因素的顯著性。在此基礎(chǔ)上,采用多種預(yù)測方法對滾石運動特征參數(shù)進行預(yù)測,并依據(jù)預(yù)測精度選取精度最優(yōu)的預(yù)測模型,以快速估算滾石運動參數(shù),為支擋結(jié)構(gòu)的設(shè)計施工提供科學(xué)合理的依據(jù)。

1 GRA-GPR模型原理簡介

1.1 GRA原理

灰色關(guān)聯(lián)法(grey relation analysis, GRA)用于判斷多因素之間的關(guān)聯(lián)性,可以確定各因素與參考變量之間的主次關(guān)系,進而掌握事件的主要特征,相較于其他多目標優(yōu)化方法,灰色關(guān)聯(lián)法具有搜索性強和效率高的特點[22-23]。本文采用鄧氏關(guān)聯(lián)度[23]評價滾石運動特征參數(shù)與其影響因素的顯著性關(guān)系。鄧氏灰色關(guān)聯(lián)法原理如下:

設(shè)有參考序列X0(t)和比較序列Xi(t);由于不同變量的量綱不同,故需進行無量綱處理:

(1)

求差序列Δi(t)及差序列的最大差Δmax和最小差Δmin公式為

(2)

(3)

(4)

可得灰色關(guān)聯(lián)系數(shù):

(5)

式中ρ為分辨系數(shù),取值區(qū)間為(0,1),具體取值要具體分析情況。當(dāng)ρ≤0.5463時,分辨力最好,通常取值0.5。

滾石災(zāi)害受到地形、地質(zhì)、環(huán)境和人為因素等各種因素的影響,有較強的不確定性,屬于灰色系統(tǒng)[24],因此,滾石運動特征參數(shù)及其影響因素的研究也屬于灰色關(guān)聯(lián)分析的范疇,灰色關(guān)聯(lián)法還具有原始信息量少、計算簡單及預(yù)測精度較高等優(yōu)點。綜上,灰色關(guān)聯(lián)法適合滾石運動特征參數(shù)影響因素的敏感性分析。

1.2 GPR預(yù)測模型

高斯過程回歸(Gaussian process regression, GPR)是基于貝葉斯網(wǎng)絡(luò)的新型機器學(xué)習(xí)算法[25],它適用于處理小樣本、隨機性強及多維復(fù)雜因素預(yù)測問題[26-27]。GPR模型本質(zhì)上是非參數(shù)回歸模型[28],是將函數(shù)的分布從有限維度空間推廣到無限維。對于數(shù)據(jù)集合:

D={(xi,yi)|i=1,…,n}

(6)

輸入矩陣xi∈Rd,輸出矩陣為yi∈R。其中f(x(1)),f(x(2)),…,f(x(n))可構(gòu)成輸入變量的聯(lián)合高斯分布集合,即

f(x)～GP(m(x),k(x,x′))

(7)

考慮到訓(xùn)練樣本目標值中可能存在噪聲,所以GPR的一般模型為

y=f(x)+ε

(8)

式中ε為獨立的高斯白噪聲,方差為σ2。

(9)

式中:I為單位矩陣;C(X,X)和K(X,X)分別為N×N的協(xié)方差矩陣與核矩陣。

訓(xùn)練數(shù)據(jù)的輸入向量X和輸出向量f,以及測試數(shù)據(jù)的輸入向量X*和輸出向量f*間的聯(lián)合高斯分布為

(10)

式中,K(X*,X)=KT(X,X*)=[k(x1,x*),k(x2,x*),…,k(xn,x*)]T為n×1維協(xié)方差函數(shù),測試數(shù)據(jù)的自協(xié)方差函數(shù)值為K(X*,X*)。由此可以得到高斯過程回歸方程,即

f*|X,y,X*～N(mx,cov(f*))

(11)

(12)

(13)

式中f*為預(yù)測均值向量。

2 基于GRA-GPR的滾石運動特征參數(shù)預(yù)測模型構(gòu)建

影響滾石運動的因素眾多,針對單一工程邊坡,目前缺少簡單易操作的滾石運動特征參數(shù)預(yù)測系統(tǒng)。因此,深入研究滾石影響因素是進行防護措施設(shè)計、施工和構(gòu)建滾石運動特征參數(shù)預(yù)測模型的關(guān)鍵。綜上,本文基于GRA計算分析滾石運動影響因素的敏感性,并根據(jù)敏感性分析結(jié)果,基于GPR等預(yù)測方法構(gòu)建了不同滾石運動特征參數(shù)的預(yù)測模型,其主要計算步驟如下:

1)灰色關(guān)聯(lián)分析:①確定參考序列和比較序列,即滾石運動特征參數(shù)和各個影響因子;②針對數(shù)據(jù)進行無量綱化處理;③求解參考序列和比較序列之間的灰色關(guān)聯(lián)系數(shù)值;④求解關(guān)聯(lián)度值;⑤對關(guān)聯(lián)度值進行排序,得出結(jié)論。并將灰色關(guān)聯(lián)的分析結(jié)果與極差方差的結(jié)果結(jié)合分析,得出不同影響因子對3個滾石運動特征參數(shù)的敏感性排序。將灰色關(guān)聯(lián)分析的結(jié)果與方差極差分析的結(jié)果相結(jié)合,對各個影響因子的敏感性進行綜合排序。

2)預(yù)測模型的構(gòu)建:基于不同影響因子的敏感性排序結(jié)果和現(xiàn)場實際的勘察情況構(gòu)建預(yù)測模型數(shù)據(jù)集,對敏感性較高的因子取較多的因素水平,對敏感性較小的因子取較少的因素水平。完成因素水平的選取后,將部分數(shù)據(jù)集按照一定比例歸為驗證集,利用二維滾石數(shù)值模擬軟件,計算工程邊坡的滾石運動特征參數(shù)。

3)數(shù)據(jù)預(yù)測:利用不同預(yù)測方法對3個滾石運動特征參數(shù)進行預(yù)測,并將預(yù)測結(jié)果與驗證集中滾石模擬軟件的計算結(jié)果、誤差進行對比分析,以此來驗證模型的精確性,為實際工程的工作人員提供科學(xué)精準方便易操作的滾石運動特征參數(shù)預(yù)測模型。模型構(gòu)建的流程如圖1所示。

圖1 模型構(gòu)建流程圖Fig.1 Model building flowchart

3 基于GRA的滾石運動影響因素敏感性分析

滾石運動特征參數(shù)的影響因素眾多,且均具有隨機性和模糊性的特點,對滾石防護措施的設(shè)計、施工及管理都造成了困難。因此,探究眾多影響因素對滾石運動特征參數(shù)的影響程度,對實際工程中防護方案的制定有重要的意義。GRA是分析因子間相互關(guān)系的有效方法,且對樣本量、樣本概率分布形式要求低,因此,能夠?qū)RA與滾石數(shù)值模擬結(jié)合,用于分析不同滾石運動影響因素的敏感性。

影響滾石運動特征參數(shù)因素眾多。本文將水平運動距離L、彈跳高度H和總動能J作為滾石運動特征參數(shù),其主要受初始速度V、邊坡摩擦角φ、滾石質(zhì)量M、法向恢復(fù)系數(shù)Rn和切向恢復(fù)系數(shù)Rt等5個滾石運動影響因素的影響。其中Rn和Rt的取值主要受邊坡和滾石巖體性質(zhì)的影響,因此,在以下滾石運動分析中,不將Rn和Rt作為獨立因素考慮。

本研究通過滾石運動數(shù)值模擬進行計算,將滾石簡化為質(zhì)點以分析崩落巖體的運動。數(shù)值模擬是基于滾石基本參數(shù)進行大量模擬和概率統(tǒng)計,能夠模擬出滾石運動過程中的運動軌跡、運動狀態(tài)及運動特征參數(shù)。

3.1 敏感性分析試驗方案設(shè)計

滾石運動影響因素的具體取值范圍如表1所示。對于初始速度V,通常認為滾石在啟動時V較小,本文

表1 滾石運動影響因素的取值范圍Table 1 Value range of factors affecting rolling stone movement

取值范圍為0.1～0.5 m/s;Rn和Rt參數(shù)的取值參考文獻[29];在此次研究中,假設(shè)邊坡為單一的無植被覆蓋的堆石邊坡。

作為研究多因素、多水平的高效和快速設(shè)計方法,采用正交試驗設(shè)計方法選取的代表點具有“均勻分散,齊整可比”的特點[30],既能減少試驗次數(shù),又能收到較好的代表效果。本研究采用正交試驗設(shè)計方法組成25組試驗方案,結(jié)果如表2所示,進行數(shù)值模擬計算。

3.2 滾石運動影響因素敏感性分析

通過對表2中的25組試驗方案進行計算,獲得相對應(yīng)的水平運動距離、彈跳高度和總動能等滾石運動特征參數(shù),使用灰色關(guān)聯(lián)方法對計算結(jié)果進行敏感性分析,結(jié)果如圖2所示。灰色關(guān)聯(lián)法的關(guān)聯(lián)度值介于0～1之間,該值越大代表其與滾石運動特征參數(shù)之間的相關(guān)性越強。

表2 數(shù)值模擬計算結(jié)果Table 2 Numerical simulation calculation results

由圖2可知,Rn、Rt、φ對于水平運動距離的影響較大,對V、M的影響結(jié)果較小;Rn、Rt、φ對于彈跳高度的影響較大,對V、M的影響結(jié)果較小;M對于總動能的影響較大,對其余因素的影響較小。此結(jié)果符合邊坡滾石的一般規(guī)律特性。

圖2 滾石運動影響因素關(guān)聯(lián)度計算結(jié)果Fig.2 Calculation results of correlation degree of influencing factors of rolling stone

在上述分析的基礎(chǔ)上,使用極差和方差方法進一步分析數(shù)值模擬結(jié)果。極差和方差的計算結(jié)果如表3和表4所示。極差分析與灰色關(guān)聯(lián)法的結(jié)果分析類似,根據(jù)極差分析的計算結(jié)果,可以對影響因素的敏感性進行排序;方差分析的貢獻值p是衡量3個滾石運動特征參數(shù)與影響因子差異大小的指標,*代表p值小于0.01,表示2組存在顯著差異。

表3 數(shù)值模擬計算結(jié)果極差分析Table 3 Range analysis of numerical simulation calculation results

綜合圖2、表4的計算結(jié)果進行分析,得出以下結(jié)論:

1)對于水平運動距離L,灰色關(guān)聯(lián)法和極差分析的滾石運動影響因素敏感性結(jié)果一致,結(jié)合方差分析結(jié)果可知:Rn、Rt、φ、V均對L影響較大,而M對L無顯著性影響。

2)對于彈跳高度H,灰色關(guān)聯(lián)法和極差分析的滾石運動影響因素敏感性結(jié)果一致,結(jié)合方差分析結(jié)果可知:Rn、Rt、φ均對H影響較大,而V和M對L無顯著性影響。

3)對于滾石的總動能J,灰色關(guān)聯(lián)度和極差分析對各因素顯著性的排序結(jié)果分別為M>V>Rn>Rt>φ和M>Rn>Rt>φ>V;方差分析認為僅滾石質(zhì)量對滾石的總動能有顯著性;因此結(jié)合方差分析結(jié)果可知:M對于J有較強的顯著性,其余因素對J均無顯著性。

4)由圖2、表3和表4可知,V對于3個滾石運動特征參數(shù)的影響均較小;φ對于L和H的影響較大,對于J的影響較小;M對于J的影響較大,對于L和H的影響較小;Rn和Rt對于水平運動距離L和H的影響較大,對于J的影響較小。

表4 數(shù)值模擬計算結(jié)果方差分析Table 4 Analysis of variance of numerical simulation results

綜上所述,基于GRA分析的結(jié)果可知:Rn、Rt、φ對水平運動距離L和彈跳高度H的影響較大,對總動能J無顯著性影響;M僅對J有較大影響,對水L和H無顯著性影響;V僅對L有較小影響,對H和J無顯著性影響。

4 基于GPR的滾石運動特征參數(shù)預(yù)測模型

滾石運動特征參數(shù)多需要通過數(shù)值模擬計算獲得,當(dāng)影響滾石運動影響因素發(fā)生較大變化時,需要重新進行數(shù)值模擬計算,制約著防護方案的制定。因此,本研究將采用支持向量機、神經(jīng)網(wǎng)絡(luò)和高斯過程回歸等預(yù)測方法,以某抽水蓄能電站的特大型堆石邊坡為例,建立滾石運動特征參數(shù)的預(yù)測模型,為防護結(jié)構(gòu)的設(shè)計、施工及時提供科學(xué)數(shù)據(jù)。

4.1 預(yù)測模型數(shù)據(jù)集建立

本研究以某一等大(1)型抽水蓄能電站的堆石邊坡為研究對象,滾石的邊坡和業(yè)主營地的位置關(guān)系如圖3所示。建成后的特大型堆石邊坡最大高差為199.80 m,一旦發(fā)生滾石災(zāi)害,極有可能對上下庫連接路、業(yè)主營地造成堵塞,存在極大的安全隱患。因此,建立快速、準確的滾石運動特征參數(shù)預(yù)測模型,對工程的安全管理與可靠運行具有重要的意義。

圖3 滾石邊坡與工程棄渣場的平面布置圖Fig.3 Floor plan of rolling stone slope and engineering waste slag yard

模型構(gòu)建方法如下:根據(jù)現(xiàn)場勘察的情況可知,邊坡上半段主要為矮植被的巖屑堆積層,邊坡下半段主要為無植被的塊石堆積層,邊坡的實際狀況及三維模型如圖4 所示。因此,在進行數(shù)值模擬建模時,兩段邊坡的恢復(fù)系數(shù)及摩擦角取值范圍也不同,具體的取值范圍參考文獻[31]。各參數(shù)具體取值范圍如表5所示。由基于GRA的滾石運動影響因素敏感性分析結(jié)果可知:恢復(fù)系數(shù)、質(zhì)量和摩擦角對滾石運動特征參數(shù)影響程度較大,對初始速度的影響程度較小,基于敏感性分析結(jié)果和實際調(diào)查情況,擬定該邊坡滾石運動特征參數(shù)的水平取值,共組成216組試驗組合,選取數(shù)據(jù)集后40組作為驗證集。

圖4 邊坡分段示意圖Fig.4 Schematic diagram of slope section

對表5的數(shù)據(jù)排列組合并進行數(shù)值模擬計算,得出預(yù)測模型所需數(shù)據(jù)集,如表6所示。

表5 預(yù)測模型各參數(shù)取值范圍Table 5 Value range of each parameter of the prediction model

根據(jù)表6的數(shù)據(jù),建立滾石運動特征參數(shù)的預(yù)測模型,具體步驟如下:①確定取數(shù)據(jù)集176組試驗數(shù)據(jù),采用SVM、GPR和神經(jīng)網(wǎng)絡(luò)方法進行回歸訓(xùn)練,用后40組樣本進行預(yù)測驗證;②通過貝葉斯優(yōu)化器進行參數(shù)尋優(yōu),確定最小MSE;③基于生成的預(yù)測模型,使用交叉驗證法進行驗證;④對驗證集的預(yù)測值和實際值進行誤差分析,確定模型的可行性。

表6 預(yù)測模型試驗數(shù)據(jù)集Table 6 Experiment dataset of prediction model

4.2 預(yù)測結(jié)果分析

各模型對不同特征參數(shù)的預(yù)測相對誤差、對角誤差圖如圖5所示,預(yù)測值和實際值的均方誤差(MSE)、R2如表7所示。

圖5 不同特征參數(shù)的預(yù)測相對誤差及對角誤差Fig.5 Prediction relative error and diagonal error of different feature parameters

對于模型精確度的評價通過兩方面進行:①MSE可以評價數(shù)據(jù)的變化程度,MSE的值越小,說明預(yù)測模型描述試驗數(shù)據(jù)具有更好的精確度;②R2也稱擬合優(yōu)度和決定系數(shù):反應(yīng)回歸模型擬合數(shù)據(jù)的優(yōu)良程度,R2越接近1,代表模型的擬合效果越好。由圖5和表7可知:

表7 不同預(yù)測模型對不同參數(shù)預(yù)測的MSE及R2Table 7 MSE and R2 of different prediction models for different parameters

1)對于滾石的水平距離,基于GPR的滾石運動特征參數(shù)預(yù)測模型的相對誤差及MSE最小,R2值最大,模型預(yù)測的精度最高;

2)對于滾石的彈跳高度,基于神經(jīng)網(wǎng)絡(luò)的相對誤差及MSE最小,R2值最大,模型預(yù)測的精度最高,基于GPR的滾石運動特征參數(shù)預(yù)測模型精度略差;

3)對于滾石的總動能,基于GPR的滾石運動特征參數(shù)預(yù)測模型的相對誤差及MSE最小,R2值最大,模型預(yù)測的精度最高;

4)從模型的訓(xùn)練時長看,GPR預(yù)測模型的訓(xùn)練時間要明顯優(yōu)于SVM預(yù)測及神經(jīng)網(wǎng)絡(luò)預(yù)測。

綜上所述,基于GPR的滾石運動特征參數(shù)預(yù)測模型雖存在一定的誤差,但未超過10%[32],模型精度能夠滿足工程要求,并且相較于SVM預(yù)測及神經(jīng)網(wǎng)絡(luò)預(yù)測有較為明顯的優(yōu)勢,在模型預(yù)測時間上也明顯優(yōu)于SVM預(yù)測及神經(jīng)網(wǎng)絡(luò)預(yù)測。基于GPR的滾石運動特征參數(shù)預(yù)測模型,省去了傳統(tǒng)滾石數(shù)值模擬繁瑣的計算過程,預(yù)測模型只需修改不同影響因子的取值即可快速獲得相應(yīng)的滾石運動特征參數(shù),有效提高了設(shè)計施工人員的工作效率。

5 結(jié)論

本文將灰色關(guān)聯(lián)法與滾石數(shù)值模擬相結(jié)合,以滾石的運動特征參數(shù)為目標參數(shù),對滾石運動影響因素的敏感性進行分析,得出結(jié)論如下:

1)Rn、Rt、φ對水平運動距離L和彈跳高度H的影響較大,對總動能J無顯著性影響。

2)M僅對J有較大影響,對L和H無顯著性影響。

3)V僅對L有較小影響,對H和J無顯著性影響。

根據(jù)影響因素的敏感性分析結(jié)果,基于GPR、SVM和神經(jīng)網(wǎng)絡(luò)方法對滾石運動特征參數(shù)建立預(yù)測模型,對不同模型的精度和訓(xùn)練耗時等進行對比,結(jié)果表明基于GPR的方法滾石運動特征參數(shù)預(yù)測模型能夠?qū)μ卣鲄?shù)進行精準快速估算,提高了工作人員設(shè)計施工管理的工作效率,對實際工程中支擋結(jié)構(gòu)的設(shè)計施工有重要的參考價值和意義。