林逸暉,李廣濤,楊天雨,喬登攀,王 俊,張 希,趙懷軍
(1.昆明理工大學(xué) 國(guó)土資源工程學(xué)院,云南 昆明 650093;2.涼山礦業(yè)股份有限公司,四川 涼山 615141)
近年來(lái),隨著現(xiàn)代工業(yè)對(duì)礦產(chǎn)資源的需求量日益增加,我國(guó)露天礦山的開(kāi)采規(guī)模不斷擴(kuò)大,開(kāi)采深度越來(lái)越大,高陡邊坡數(shù)量越來(lái)越多,邊坡災(zāi)害發(fā)生頻次逐年上升。邊坡災(zāi)害會(huì)對(duì)人民生命財(cái)產(chǎn)安全造成嚴(yán)重威脅,并使生態(tài)環(huán)境遭受巨大破壞[1-2]。因此,為減少邊坡災(zāi)害的發(fā)生,開(kāi)展邊坡穩(wěn)定性預(yù)測(cè)研究具有重要的現(xiàn)實(shí)意義。
露天礦山邊坡穩(wěn)定性分析方法主要有定性分析法、定量分析法和數(shù)值分析法等[3]。隨著計(jì)算機(jī)技術(shù)的發(fā)展,機(jī)器學(xué)習(xí)算法被引入邊坡穩(wěn)定性預(yù)測(cè)研究中。目前用于邊坡穩(wěn)定性預(yù)測(cè)的機(jī)器學(xué)習(xí)算法有粗糙集理論、模糊邏輯理論、支持向量機(jī)、蟻群算法、貝葉斯算法、信息量法、多元邏輯回歸、決策樹(shù)法等[4]。瞿婧晶等[5]采用粗糙集屬性約簡(jiǎn)和相互關(guān)系矩陣法進(jìn)行綜合分析,建立了邊坡穩(wěn)定性評(píng)價(jià)體系,并將其應(yīng)用于鎮(zhèn)江地區(qū)下蜀土邊坡的穩(wěn)定性預(yù)測(cè)中。張?jiān)蒲鉡6]采用網(wǎng)格搜索算法對(duì)支持向量機(jī)模型進(jìn)行參數(shù)優(yōu)化,并根據(jù)邊坡特征參數(shù),采用優(yōu)化后的支持向量機(jī)模型對(duì)邊坡穩(wěn)定性進(jìn)行了預(yù)測(cè)。趙允坤等[7]利用改進(jìn)粒子群算法(IPSO)對(duì)極限學(xué)習(xí)機(jī)的輸入權(quán)值和隱層偏置進(jìn)行了優(yōu)化,建立了基于IPSO-ELM的邊坡穩(wěn)定性預(yù)測(cè)模型,實(shí)現(xiàn)了對(duì)邊坡穩(wěn)定性的有效預(yù)測(cè)。
但機(jī)器學(xué)習(xí)算法也存在一定缺陷,如:粗糙集理論需結(jié)合其他算法使用,否則難以有效反映不確定性問(wèn)題;支持向量機(jī)模型對(duì)于缺失數(shù)據(jù)以及參數(shù)和核函數(shù)的選擇較敏感,且核函數(shù)及其參數(shù)是根據(jù)經(jīng)驗(yàn)選取的,帶有一定的主觀性;極限學(xué)習(xí)機(jī)存在易過(guò)擬合以及可控性差等缺點(diǎn)。
隨機(jī)森林(Random Forest,RF)是由美國(guó)科學(xué)家Leo Breiman[8]于2001年提出的,這是一種將Bagging集成學(xué)習(xí)理論與隨機(jī)子空間方法相結(jié)合的機(jī)器學(xué)習(xí)算法。隨機(jī)森林以決策樹(shù)為基本分類器,相對(duì)于決策樹(shù)法,隨機(jī)森林克服了其易出現(xiàn)過(guò)擬合的缺點(diǎn),降低了異常值可能帶來(lái)的影響,但計(jì)算量較大[9-10]。此外,隨機(jī)森林對(duì)高維數(shù)據(jù)分類問(wèn)題有著良好的可擴(kuò)展性和并行性,并且作為一種由數(shù)據(jù)驅(qū)動(dòng)的非參數(shù)分類方法,使用隨機(jī)森林算法時(shí)無(wú)需調(diào)參以及分類先驗(yàn)知識(shí)[11]。針對(duì)隨機(jī)森林算法計(jì)算性能開(kāi)銷較大的缺點(diǎn),本文利用主成分分析法(Principal Component Analysis,PCA)進(jìn)行數(shù)據(jù)處理,過(guò)濾冗雜數(shù)據(jù),降低隨機(jī)森林算法的計(jì)算量,構(gòu)建基于PCA-RF的邊坡穩(wěn)定性預(yù)測(cè)模型,并將其應(yīng)用于露天邊坡的穩(wěn)定性預(yù)測(cè)。
主成分分析是利用降維思想,在信息損失較少的前提下,將原來(lái)具有一定相關(guān)性的指標(biāo)重新組合成一組新的信息互不重疊的綜合指標(biāo),來(lái)反映大部分原指標(biāo)中所攜帶的信息量的多元統(tǒng)計(jì)方法[12]。主成分分析計(jì)算步驟敘述如下。
(1)建立m×n矩陣,m為樣本數(shù)量,n為每個(gè)樣本的指標(biāo)個(gè)數(shù)。
(1)
(2)將原始數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化,生成標(biāo)準(zhǔn)化矩陣。
(2)
(3)進(jìn)行KMO(Kaiser-Meyer-Olkin)和巴特利特球形檢驗(yàn)。若KMO值大于0.5,表明變量間的相關(guān)程度差異較小,數(shù)據(jù)很適合作因子分析[13];若巴特利特球形檢驗(yàn)的結(jié)果小于0.05,球形假設(shè)被拒絕,表明原始變量之間存在相關(guān)性,數(shù)據(jù)適合作因子分析[14]。
(4)根據(jù)標(biāo)準(zhǔn)化矩陣建立協(xié)方差矩陣。
(3)
式中,i=1,2,…,n,j=1,2,…,n。
協(xié)方差矩陣由標(biāo)準(zhǔn)化矩陣的兩兩列變量的協(xié)方差組成,因此cov(zi,zj)為列變量zi與zj的相關(guān)系數(shù),即協(xié)方差矩陣為相關(guān)系數(shù)矩陣。
(5)計(jì)算協(xié)方差矩陣R的非負(fù)特征根λi(i=1,2,…,p)。
|λE-R|=0。
(4)
(6)由計(jì)算得到的p個(gè)非負(fù)特征根計(jì)算前q個(gè)主成分的累計(jì)貢獻(xiàn)率。
(5)
(7)根據(jù)主成分的累計(jì)貢獻(xiàn)率選取主成分個(gè)數(shù),并生成新變量指標(biāo)。前q個(gè)主成分的累計(jì)貢獻(xiàn)率達(dá)到85%時(shí),即可認(rèn)為這些主成分中包含了絕大部分原指標(biāo)所攜帶的信息量[15]。將生成的新變量指標(biāo)表示原變量指標(biāo):
(6)
式中,x1,x2,…,xn為原變量指標(biāo),y1,y2,…,yp為新變量指標(biāo),cij與yi互不相關(guān)且cij滿足ci12+ci22+…+cin2=1。由此完成了數(shù)據(jù)降維,在信息損失較少的前提下減少了變量數(shù)。
隨機(jī)森林{h(X,θk),k=1,2,…,K}是以K個(gè)決策樹(shù)模型作為基本分類器,通過(guò)集成學(xué)習(xí)形成的組合分類器[16]。隨機(jī)森林算法以待分類樣本作為輸入,基于決策樹(shù)模型,采用自上而下的遞歸方式,從樹(shù)的根節(jié)點(diǎn)開(kāi)始,在其內(nèi)部節(jié)點(diǎn)上進(jìn)行屬性值的測(cè)試比較,然后按照給定的屬性值確定對(duì)應(yīng)分支,在決策樹(shù)的葉子節(jié)點(diǎn)得到結(jié)果,最后根據(jù)每個(gè)決策樹(shù)模型的分類結(jié)果投票決定最終輸出的分類結(jié)果。
隨機(jī)森林中的{θk,k=1,2,…,K}為隨機(jī)變量序列,其體現(xiàn)在兩個(gè)方面:
(1)數(shù)據(jù)的隨機(jī)選取。在Bagging集成算法的基礎(chǔ)上,利用Bootstrap方法從原始樣本集中采取有放回抽樣,隨機(jī)創(chuàng)建K個(gè)與原樣本集數(shù)據(jù)量相同的子數(shù)據(jù)集{Tk,k=1,2,…,K},其中不同子數(shù)據(jù)集的元素可以重復(fù),同一個(gè)子數(shù)據(jù)集中的元素也可以重復(fù)。
(2)待選屬性的隨機(jī)選取。采用隨機(jī)子空間方法,在隨機(jī)森林中的子決策樹(shù)的各節(jié)點(diǎn)進(jìn)行分裂時(shí),從所有的待選屬性中等概率隨機(jī)選取一定的屬性[通常選取(log2M+1)個(gè)屬性,M為屬性總個(gè)數(shù)[17]],再?gòu)乃x屬性中選取最優(yōu)屬性進(jìn)行節(jié)點(diǎn)分裂。
隨機(jī)森林的分類過(guò)程為:
(1)利用Bootstrap方法隨機(jī)選取數(shù)據(jù),生成子數(shù)據(jù)集Tk,將Tk作為訓(xùn)練集構(gòu)建K棵決策樹(shù)。
(2)每棵決策樹(shù)生長(zhǎng)過(guò)程中在所有的待選屬性中等概率隨機(jī)選取一定的屬性,并從中選取最優(yōu)屬性對(duì)每棵決策樹(shù)的各個(gè)節(jié)點(diǎn)進(jìn)行分裂,并且讓每棵決策樹(shù)均不進(jìn)行剪枝地完整生長(zhǎng)。
(3)由完整生長(zhǎng)的K棵決策樹(shù)構(gòu)成隨機(jī)森林。
(4)輸入測(cè)試集樣本進(jìn)行預(yù)測(cè),由每棵決策樹(shù)的輸出結(jié)果投票決定隨機(jī)森林輸出的分類結(jié)果(即每棵決策樹(shù)輸出結(jié)果的眾數(shù))。
影響邊坡穩(wěn)定性的因素主要有邊坡角度、邊坡高度、巖土體物理力學(xué)性質(zhì)、巖石強(qiáng)度、巖石結(jié)構(gòu)、孔隙水壓、水文條件及地震等[18]。但巖石強(qiáng)度、巖石結(jié)構(gòu)等因素對(duì)于邊坡穩(wěn)定性的影響較小,地震等自然災(zāi)害出現(xiàn)的概率較低,樣本較少。因此,為減少隨機(jī)森林算法的計(jì)算量,提高計(jì)算速度,選取巖石容重γ、黏聚力c、內(nèi)摩擦角φ、邊坡角度α、邊坡高度H和孔隙水壓μ等6個(gè)主要影響因素進(jìn)行邊坡穩(wěn)定性預(yù)測(cè)[19]。本文基于此6項(xiàng)指標(biāo)構(gòu)建隨機(jī)森林邊坡穩(wěn)定性預(yù)測(cè)模型,將使用主成分分析降維后的指標(biāo)作為實(shí)際輸入、邊坡穩(wěn)定狀態(tài)“穩(wěn)定”和“失穩(wěn)”作為輸出進(jìn)行邊坡穩(wěn)定性預(yù)測(cè)。
本文選取文獻(xiàn)[20-21]中的54組不同且具有明確穩(wěn)定狀態(tài)的邊坡實(shí)例作為研究對(duì)象,構(gòu)建PCA-RF邊坡穩(wěn)定性預(yù)測(cè)模型。邊坡實(shí)例數(shù)據(jù)見(jiàn)表 1,其中:邊坡?tīng)顟B(tài)“1”表示邊坡穩(wěn)定,“2”表示邊坡失穩(wěn);穩(wěn)定邊坡為34組,失穩(wěn)邊坡為20組,第55-第58組數(shù)據(jù)為拉拉銅礦東部露天采場(chǎng)終了邊坡數(shù)據(jù),以此作為預(yù)測(cè)實(shí)例。
表1 邊坡實(shí)例數(shù)據(jù)
對(duì)表 1中的58組數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理,計(jì)算得到KMO值為0.665,大于0.5;巴特利特球形檢驗(yàn)值為0.000,小于0.05,表明樣本數(shù)據(jù)適合進(jìn)行主成分分析。對(duì)表 1中的數(shù)據(jù)進(jìn)行主成分分析,得到各影響因素間的相關(guān)矩陣熱圖(見(jiàn)圖 1)。將相關(guān)性系數(shù)取絕對(duì)值后,如在0~0.1,表示無(wú)相關(guān);如在0.1~0.3,表示弱相關(guān);如在0.3~0.5,表示中等相關(guān);如在0.5~1.0,表示強(qiáng)相關(guān)[22]。
圖 1 相關(guān)矩陣熱圖
由圖 1可知,某些影響因素間存在較強(qiáng)的相關(guān)性,如巖石容重與邊坡角度、邊坡高度間的相關(guān)性系數(shù)分別為0.535、0.679,均大于0.5,表明具有強(qiáng)相關(guān)性。主成分分析能夠?qū)⑦@些指標(biāo)重新組合成一組新的信息互不重疊的綜合指標(biāo),從而避免相關(guān)性過(guò)強(qiáng)可能帶來(lái)的多重共線性問(wèn)題。
為確保選取的指標(biāo)中攜帶原樣本數(shù)據(jù)的絕大部分信息量,需計(jì)算各因素對(duì)邊坡穩(wěn)定性影響的貢獻(xiàn)率及累計(jì)貢獻(xiàn)率,結(jié)果見(jiàn)表2。
表2 各影響因素貢獻(xiàn)率
表3 主成分系數(shù)矩陣
由表 2可知,前4項(xiàng)主成分的累計(jì)貢獻(xiàn)率達(dá)到了91.283%,大于85%,表明前4項(xiàng)主成分已攜帶原樣本數(shù)據(jù)的絕大部分信息量,因此選取前4項(xiàng)主成分生成新變量指標(biāo)。各主成分系數(shù)矩陣見(jiàn)表 3。
將各主成分系數(shù)分別除以相應(yīng)的主成分特征根的平方根,得到因子載荷矩陣(見(jiàn)表4)。
表4 因子載荷矩陣
由因子載荷矩陣可得出新變量指標(biāo)的表達(dá)式:
(7)
經(jīng)主成分分析后,在降低變量維數(shù)的同時(shí),還保證了原樣本數(shù)據(jù)攜帶的信息量盡可能少丟失。經(jīng)主成分分析降維后的邊坡數(shù)據(jù)見(jiàn)表5。
表5 經(jīng)主成分分析降維后的邊坡數(shù)據(jù)
表6 PCA-RF模型預(yù)測(cè)結(jié)果
從降維后的前54組邊坡數(shù)據(jù)中隨機(jī)選取45組數(shù)據(jù)作為輸入訓(xùn)練樣本建立PCA-RF邊坡穩(wěn)定性預(yù)測(cè)模型,將剩下的9組數(shù)據(jù)作為測(cè)試樣本,用來(lái)檢驗(yàn)該預(yù)測(cè)模型的準(zhǔn)確率。
在使用隨機(jī)森林進(jìn)行預(yù)測(cè)之前,需先設(shè)置2個(gè)重要參數(shù):隨機(jī)森林中包含的決策樹(shù)個(gè)數(shù)ntree和每棵樹(shù)節(jié)點(diǎn)預(yù)選的變量指標(biāo)個(gè)數(shù)mtry。若預(yù)設(shè)的ntree值太小會(huì)導(dǎo)致隨機(jī)森林模型訓(xùn)練不夠充分,太大則會(huì)導(dǎo)致模型運(yùn)行速度緩慢;若預(yù)設(shè)的mtry值太小可能會(huì)導(dǎo)致分類器過(guò)擬合,使預(yù)測(cè)精度降低,太大則會(huì)增加模型的運(yùn)算量。綜合考慮,本文將ntree值設(shè)為500,mtry值設(shè)為4。PCA-RF模型預(yù)測(cè)結(jié)果見(jiàn)表 6。
由表 6可知,PCA-RF模型的預(yù)測(cè)準(zhǔn)確率達(dá)到了100%。但因隨機(jī)森林生成訓(xùn)練集時(shí)采用Bootstrap有放回自助采樣法,輸入訓(xùn)練樣本中有(1-1/N)N的樣本不會(huì)被抽取(N為輸入訓(xùn)練樣本中的樣本個(gè)數(shù)),當(dāng)N足夠大時(shí),(1-1/N)N將收斂于1/e,約為0.368,即有約36.8%的樣本數(shù)據(jù)不會(huì)被抽取[23]。因此,該預(yù)測(cè)模型有可能存在誤差,為此需要對(duì)該模型進(jìn)行性能評(píng)估。
模型性能的評(píng)估方法有留一法和交叉驗(yàn)證法等,本文選用6折交叉驗(yàn)證法對(duì)PCA-RF模型進(jìn)行性能評(píng)估。6折交叉驗(yàn)證法驗(yàn)證步驟為:
(1)將54組樣本數(shù)據(jù)隨機(jī)分為數(shù)據(jù)量相同的6組數(shù)據(jù)子集D1-D6。
(2)取其中的第i組數(shù)據(jù)子集Di作為測(cè)試樣本,將其余5組數(shù)據(jù)子集作為輸入訓(xùn)練樣本,構(gòu)建6個(gè)隨機(jī)森林預(yù)測(cè)模型。
(3)匯總預(yù)測(cè)結(jié)果,開(kāi)展PCA-RF模型的性能評(píng)估。
交叉驗(yàn)證過(guò)程中預(yù)測(cè)模型的參數(shù)設(shè)置與前文相同,驗(yàn)證結(jié)果見(jiàn)表7。由表7可知,模型的預(yù)測(cè)精度較高,達(dá)到了94.44%,表明該預(yù)測(cè)模型可以用于邊坡穩(wěn)定性預(yù)測(cè)。
表 7 交叉驗(yàn)證結(jié)果
拉拉銅礦礦體主要部分埋藏不深,礦山原設(shè)計(jì)采用露天開(kāi)采,露天開(kāi)采區(qū)域分為東部露天采場(chǎng)和西部露天采場(chǎng),東部露天采場(chǎng)尚未開(kāi)采至終了境界,其終了邊坡穩(wěn)定狀態(tài)未知,因此本文選取東部露天采場(chǎng)終了邊坡作為預(yù)測(cè)對(duì)象,現(xiàn)場(chǎng)照片見(jiàn)圖2。該礦東部露天采場(chǎng)終了邊坡東幫設(shè)計(jì)高度為252 m,坡角為45.57°;西幫設(shè)計(jì)高度為216 m,坡角為41.75°;南幫設(shè)計(jì)高度為296 m,坡角為45.89°;北幫設(shè)計(jì)高度為312 m,坡角為30.05°。露天采場(chǎng)區(qū)域巖性主要為片巖,巖石力學(xué)參數(shù)見(jiàn)表 1中第55-第58組數(shù)據(jù)。東部露天采場(chǎng)邊坡現(xiàn)場(chǎng)照片見(jiàn)圖3。
圖2 拉拉銅礦東部露天采場(chǎng)
圖3 拉拉銅礦東部露天采場(chǎng)邊坡
將經(jīng)主成分分析降維后的拉拉銅礦東部露天采場(chǎng)終了邊坡數(shù)據(jù)作為預(yù)測(cè)數(shù)據(jù),輸入已訓(xùn)練好的隨機(jī)森林預(yù)測(cè)模型中,輸出的預(yù)測(cè)結(jié)果見(jiàn)表8。由表8可知,拉拉銅礦東部露天采場(chǎng)終了邊坡的東幫、西幫、南幫、北幫均處于穩(wěn)定狀態(tài)。
表8 拉拉銅礦東部露天采場(chǎng)終了邊坡穩(wěn)定性預(yù)測(cè)結(jié)果
a.本文采用主成分分析對(duì)原始邊坡數(shù)據(jù)進(jìn)行降維,減少了隨機(jī)森林算法的計(jì)算量,建立了PCA-RF邊坡穩(wěn)定性預(yù)測(cè)模型。
b.邊坡實(shí)例預(yù)測(cè)結(jié)果表明,PCA-RF模型預(yù)測(cè)準(zhǔn)確率達(dá)到100%,并且該模型經(jīng)6折交叉驗(yàn)證法評(píng)估的預(yù)測(cè)精度高達(dá)94.44%,表明該模型可以用于邊坡穩(wěn)定性預(yù)測(cè)。
c.將PCA-RF模型應(yīng)用于拉拉銅礦東部露天采場(chǎng)邊坡穩(wěn)定性預(yù)測(cè),結(jié)果表明,該礦東部露天采場(chǎng)東幫、西幫、南幫、北幫終了邊坡均處于穩(wěn)定狀態(tài),該預(yù)測(cè)結(jié)果可為礦山生產(chǎn)和管理提供參考。