李 橋, 巨能攀, 黃 健, 王昌明

(成都理工大學(xué)地質(zhì)災(zāi)害防治與地質(zhì)環(huán)境保護(hù)國家重點(diǎn)實(shí)驗(yàn)室，成都 610059)

泥石流是一種在山區(qū)非常普遍的地質(zhì)災(zāi)害。中國山地面積分布廣泛且人口聚集，決定了泥石流災(zāi)害的多發(fā)性和致災(zāi)性[1]。尤其是在高強(qiáng)度高烈度地震后，極震區(qū)產(chǎn)生了大量的滑坡和崩塌，為泥石流的形成提供了充足的物源條件。故而極震區(qū)內(nèi)泥石流活動(dòng)呈現(xiàn)出頻率增高、規(guī)模增大，并具有區(qū)域群發(fā)性、雨季突發(fā)性和周期復(fù)發(fā)性等特點(diǎn)，成為了眾多專家和學(xué)者們重點(diǎn)關(guān)注的焦點(diǎn)[2]。

泥石流動(dòng)儲(chǔ)量作為設(shè)計(jì)泥石流防治工程中不可或缺的參數(shù)，傳統(tǒng)的物源動(dòng)儲(chǔ)量計(jì)算方法有：人工經(jīng)驗(yàn)法、比例統(tǒng)計(jì)法和遙感解譯法等[3-5]。但與普通泥石流相比，極震區(qū)泥石流的物源類型在規(guī)模、數(shù)量和分布方式上都有極大差別，更為重要的是，極震區(qū)泥石流物源在后期強(qiáng)降雨的作用下啟動(dòng)方式發(fā)生了明顯改變，采用傳統(tǒng)方法獲取泥石流動(dòng)儲(chǔ)量可信度不高。為此，相關(guān)學(xué)者在該領(lǐng)域提出了多種方法來解決這一問題，并取得了諸多成果。Dong等[6]通過收集1999年“9·21”集集地震后中國臺(tái)灣地區(qū)臺(tái)中縣區(qū)域內(nèi)泥石流溝在兩次臺(tái)風(fēng)影響下的物源量變化數(shù)據(jù)，提出了基于判別因素的多元回歸模型。喬建平等[7]調(diào)查統(tǒng)計(jì)汶川地震極震區(qū)內(nèi)44條泥石流溝的物源信息，總結(jié)出了汶川地震災(zāi)區(qū)泥石流物源的主要類型和啟動(dòng)地質(zhì)模式，基于數(shù)學(xué)統(tǒng)計(jì)的方法發(fā)現(xiàn)總物源量與動(dòng)儲(chǔ)量呈線性相關(guān)。顧文韜等[8]以四川安縣高川鄉(xiāng)區(qū)域內(nèi)多條泥石流溝為研究對象，提出了“地震高程指數(shù)放大經(jīng)驗(yàn)?zāi)Ｐ汀?，通過多元統(tǒng)計(jì)擬合法提出了極震區(qū)震后泥石流的動(dòng)儲(chǔ)量計(jì)算公式。方群生等[9]基于震后泥石流調(diào)查踏勘資料，將泥石流流域內(nèi)物源分解成崩滑體、溝道物源和坡面物源，分別進(jìn)行單因子回歸分析再疊加，建立了新的泥石流動(dòng)儲(chǔ)量計(jì)算模型。以上計(jì)算模型均取得了較好的效果，但具有一定的局限性，不利于推廣。

極震區(qū)泥石流物源動(dòng)儲(chǔ)量影響因素眾多，屬于復(fù)雜的非線性問題，籠統(tǒng)的對其進(jìn)行數(shù)學(xué)分析，無法達(dá)到理想的計(jì)算精度。而神經(jīng)網(wǎng)絡(luò)具有很強(qiáng)的非線性信息處理能力、自適應(yīng)學(xué)習(xí)能力和容錯(cuò)性，可以很好地解決這一問題。極限學(xué)習(xí)機(jī)(extreme learning machine, ELM)作為一種單隱層前饋神經(jīng)網(wǎng)絡(luò)算法，具有結(jié)構(gòu)簡單、適應(yīng)性強(qiáng)和訓(xùn)練學(xué)習(xí)速度快的特點(diǎn)[10]，在各工程領(lǐng)域取得了不錯(cuò)的效果。廉城[11]運(yùn)用經(jīng)驗(yàn)?zāi)B(tài)分解(ensemble empirical mode decomposition, EEMD)將滑坡位移曲線分解成多個(gè)子序列，分別運(yùn)用ELM進(jìn)行預(yù)測分析，預(yù)測結(jié)果精度較高。Xu等[12]基于生存分析模型和ELM從定量的角度對文家溝泥石流治理工程效果進(jìn)行了分析評價(jià)。李驊錦等[13]通過巖移數(shù)據(jù)決策和ELM，提出了一中礦山開采最大下沉值的新方法。同時(shí)相關(guān)研究也表明使用智能優(yōu)化算法對ELM進(jìn)行參數(shù)選取可進(jìn)一步提高預(yù)測精度并提升網(wǎng)絡(luò)穩(wěn)定性[14-15]。因此，提出一種基于ELM的極震區(qū)泥石流物源動(dòng)儲(chǔ)量計(jì)算方法，并且使用粒子群算法(particle swarm optimization, PSO)對ELM進(jìn)行優(yōu)化，并與AdaBoost算法進(jìn)行耦合。

采集汶川極震區(qū)區(qū)內(nèi)60條泥石流溝的物源信息作為樣本數(shù)據(jù)，從泥石流物源形成與啟動(dòng)方式入手，提出了流域面積、相對高差、主溝長度、較發(fā)震斷裂帶距離、溝床平均縱比降和物源總儲(chǔ)量作為泥石流物源動(dòng)儲(chǔ)量的影響因子，運(yùn)用Person相關(guān)系數(shù)(Pearson correlation coefficient，PCC)、灰色關(guān)聯(lián)度(grey relational grade，GRG)和最大互信息系數(shù)(maximal information coefficient，MIC)對影響因子進(jìn)行了敏感性分析；為了避免信息冗余，基于主成分分析(principal component analysis,PCA)對樣本數(shù)據(jù)進(jìn)行處理，再采用AdaBoost算法和粒子群優(yōu)化的極限學(xué)習(xí)機(jī)(PSO-ELM)相結(jié)合的PSO-ELM_AdaBoost耦合模型進(jìn)行訓(xùn)練和預(yù)測，并將結(jié)果與BP(back propagation)、支持向量機(jī)(SVM)、ELM、PSO-ELM模型和傳統(tǒng)計(jì)算 模型計(jì)算值進(jìn)行比較；最后從每個(gè)子研究區(qū)中抽取一條泥石流溝和其他極震區(qū)的三條泥石流溝應(yīng)用PSO-ELM_AdaBoost模型進(jìn)行泥石流物源動(dòng)儲(chǔ)量預(yù)測，驗(yàn)證了本文模型的準(zhǔn)確性和適宜性。

1 理論與方法

1.1 最大互信息系數(shù)

現(xiàn)實(shí)世界數(shù)據(jù)之間的不一定總是呈現(xiàn)線性關(guān)系，采用余弦相似度(cosine similarity，CS)和Person相關(guān)系數(shù)等線性相關(guān)性系數(shù)，可能會(huì)造成關(guān)聯(lián)度誤判。為此Reshef等[16]基于互信息理論提出一種具有普適性和公平性的新型變量關(guān)聯(lián)評價(jià)指標(biāo)：最大互信息系數(shù)，有效地解決了變量間的非線性關(guān)聯(lián)分析。

對于給定的兩個(gè)變量A=(a1,a2,…,an)和B=(b1,b2,…,bn)，n為樣本個(gè)數(shù)，構(gòu)建一個(gè)二元數(shù)據(jù)集D，利用一個(gè)x×y的網(wǎng)格G將D進(jìn)行網(wǎng)格化，那么關(guān)于D的MIC可定義為

(1)

式(1)中：B(n)為網(wǎng)格劃分x×y的最大值，一般取B(n)=n0.6；M(D)為D的特征矩陣，其計(jì)算公式為

(2)

式(2)中：MI*是D雙變量間的最大互信息，即

MI*(D,x,y)=max[MI(D|G)]

(3)

式(3)中：D|G為每個(gè)單元網(wǎng)格的概率分布。

換言之，MIC是一種歸一化的最大互信息，其取值區(qū)間為[0，1]。MIC越大，表明兩個(gè)變量的相關(guān)性越強(qiáng)；反之，則表明兩個(gè)變量間的相關(guān)想越弱[17]。通過計(jì)算各因子與極震區(qū)泥石流動(dòng)儲(chǔ)量的MIC，綜合PCC和GRG，就是為了全面地了解不同因子的影響能力，進(jìn)而進(jìn)行敏感性分析。

1.2 主成分分析

主成分分析可以有效處理輸入因子間存在一定相關(guān)性，對問題的反應(yīng)存在一定的信息重疊問題[18]。經(jīng)PCA處理后生成的新變量可以包含原變量大部分的信息，且相互之間不存在相關(guān)性。對于有n個(gè)樣本，每個(gè)樣本含有p個(gè)變量的原始數(shù)據(jù)矩陣Xn×p，具體分析過程如下。

(1)為消除量綱影響，對樣本各數(shù)據(jù)進(jìn)行歸一化處理，并作為輸出樣本數(shù)據(jù)Yn×p。

(2)計(jì)算相關(guān)系數(shù)矩陣R，Rij(i,j=1,2,…,p)為原始變量間的相關(guān)性系數(shù)：

(4)

(3)根據(jù)矩陣R求出其特征值λi和特征向量ui，并按照從大到小排列。

(4)計(jì)算主成分貢獻(xiàn)率em和累計(jì)方差貢獻(xiàn)率Em，從而確定主成分的個(gè)數(shù)。

(5)

(6)

式中：m取值標(biāo)準(zhǔn)是使Em達(dá)到設(shè)定閾值，一般要求Em≥90%。

(5)輸出主成分樣本值Z：

(7)

用主成分樣本值Z代替原來數(shù)據(jù)樣本X，消除了原始數(shù)據(jù)間的相關(guān)性，從而達(dá)到了簡化結(jié)構(gòu)的效果。

1.3 極限學(xué)習(xí)機(jī)

極限學(xué)習(xí)機(jī)是由Huang等[19-20]提出的一種基于單隱層前向反饋型神經(jīng)網(wǎng)絡(luò)(SLGNs)的監(jiān)督型學(xué)習(xí)算法。相較于傳統(tǒng)的神經(jīng)網(wǎng)絡(luò)，該算法隨機(jī)產(chǎn)生輸入層與隱含層的連接權(quán)值及隱含層神經(jīng)元的閾值，且在訓(xùn)練過程中只需要設(shè)置隱含層神經(jīng)元的個(gè)數(shù)，便可以獲取唯一的最優(yōu)解。

設(shè)有n個(gè)任意的樣本(xi,ti)，其中xi=(xi1,xi2,…,xin)T∈Rn，ti=(ti1,ti2,…,tin)T∈Rn，Rn為n元矩陣。對于一個(gè)有L個(gè)隱層節(jié)點(diǎn)的單隱層神經(jīng)網(wǎng)絡(luò)可以表示為

(8)

式(8)中：wi=(wi1,wi2,…,win)T為輸入權(quán)重;βi為輸出權(quán)重;bi為第i個(gè)隱層單元的偏置;Oj為網(wǎng)絡(luò)輸出值;g(x)為激活函數(shù);激活函數(shù)均采用sigmoid方程，其形式如式(9)所示：

(9)

(10)

用矩陣表述為

Hβ=T

(11)

式(11)中：H為隱含層輸出矩陣；β為輸出權(quán)重矩陣；T為期望輸出矩陣。網(wǎng)絡(luò)訓(xùn)練中，由于鎖定了隨機(jī)選擇的wi和bi，H為固定矩陣。此時(shí)β可通過求解式(11)最小二乘解進(jìn)行求解，即：

(12)

此線性方程的最小二乘解為

(13)

式(13)中，H?是矩陣H的Moore-Penrose廣義逆[20]。

1.4 粒子群算法優(yōu)化極限學(xué)習(xí)機(jī)

粒子群算法是一種全局優(yōu)化算法[21]。其主要思想為將每個(gè)優(yōu)化問題的潛在解設(shè)為一個(gè)粒子，在初始化階段每個(gè)粒子都被賦予初始位置和速度，并且為了衡量每個(gè)粒子的優(yōu)越性，定義一個(gè)適應(yīng)度函數(shù)，并設(shè)定迭代次數(shù)。在每次迭代中，所有粒子向全局最優(yōu)解pbest(整個(gè)種群目前搜索到的最優(yōu)解)與個(gè)體最優(yōu)解qbest(個(gè)體自身所能達(dá)到的個(gè)體最優(yōu)解)進(jìn)行逼近，并比較適應(yīng)度，以更新自己的速度和位置，最終獲得全空間搜索最優(yōu)解。

因此在ELM預(yù)測模型中，基于PSO優(yōu)化網(wǎng)絡(luò)參數(shù)w、β和b，可以避免參數(shù)的盲目試算，提高了預(yù)測模型的準(zhǔn)確性。

1.5 PSO-ELM_AdaBoost耦合模型

AdaBoost算法是基于弱學(xué)習(xí)定理的一種Boosting應(yīng)用算法[22]，可以提高任意給定弱預(yù)測器的預(yù)測精度。模型采用PSO-ELM作為弱預(yù)測器，多次調(diào)用PSO-ELM，并根據(jù)每次訓(xùn)練樣本預(yù)測的優(yōu)劣，更新對應(yīng)的權(quán)重，再將改變權(quán)重后的樣本重新對弱預(yù)測器進(jìn)行訓(xùn)練，最后采用AdaBoost算法對這些弱預(yù)測器訓(xùn)練結(jié)果進(jìn)行集成，輸出最終結(jié)果。建立基于PCA的PSO-ELM_AdaBoost模型的步驟如下。

(1)數(shù)據(jù)選擇與網(wǎng)絡(luò)初始化。首先將原始數(shù)據(jù)樣本集X進(jìn)行歸一化處理得到數(shù)據(jù)集Y，再利用PCA將Y降維，得到了消除冗余信息的新數(shù)據(jù)樣本集Z。從Z中隨機(jī)選擇m組訓(xùn)練數(shù)據(jù)，初始化測試數(shù)據(jù)的分布權(quán)值Dt(i)=m-1，并確定預(yù)測誤差閾值。

(2)弱分類器預(yù)測。利用訓(xùn)練數(shù)據(jù)訓(xùn)練PSO-ELM并且預(yù)測訓(xùn)練數(shù)據(jù)輸出，得到預(yù)測序列g(shù)(t)的預(yù)測誤差和et：

(14)

式(14)中：Dt(i)為第t次迭代權(quán)值；g(t)為預(yù)測結(jié)果且g(t)≠y，其中y為預(yù)測期望誤差。

(3)計(jì)算預(yù)測序列權(quán)重。權(quán)重at的計(jì)算公式為

(15)

(4)更新樣本權(quán)重。根據(jù)權(quán)重at調(diào)整下一輪訓(xùn)練樣本權(quán)重：

(16)

式(16)中：Bt為歸一化因子；g(t)為預(yù)測結(jié)果。

(5)強(qiáng)預(yù)測器函數(shù)。將訓(xùn)練得到的多個(gè)弱分類函數(shù)集合成強(qiáng)預(yù)測函數(shù)h(x)。

(17)

式(17)中：ft(x)為弱分類函數(shù)。

根據(jù)上述步驟，基于PCA的PSO-ELM_AdaBoost預(yù)測模型算法流程如圖1所示。

圖1 基于主成分分析的粒子群優(yōu)化極限學(xué)習(xí)機(jī)和AdaBoost算法耦合模型流程圖Fig.1 Modeling flow chart of the PSO-ELM_AdaBoost based on PCA

2 研究區(qū)概況及影響因子選取

2.1 研究區(qū)概況

2008年“5·12”汶川地震發(fā)生后，極震區(qū)內(nèi)泥石流呈現(xiàn)出頻率高、規(guī)模大和周期復(fù)發(fā)等特征，泥石流溝間物源類型、啟動(dòng)方式也和發(fā)育特征各有不同，對該區(qū)域進(jìn)行震后泥石流物源動(dòng)儲(chǔ)量特征研究具有代表性。主要選取汶川地震極震區(qū)震后泥石流物源量較多、啟動(dòng)條件較低、危險(xiǎn)性較大的泥石流作為研究對象(部分樣本現(xiàn)場航拍圖如圖2所示)。泥石流樣本數(shù)據(jù)來源于汶川縣映秀鎮(zhèn)(13條)、G213公路沿線(5條)、都江堰市龍池鎮(zhèn)(10條)、綿竹市清平鄉(xiāng)(9條)、安州區(qū)高川鄉(xiāng)(20條)和北川縣縣城附近(3條)，共計(jì)60條樣本，樣本分布如圖3所示。

圖2 部分樣本現(xiàn)場航拍圖Fig.2 Aerial photos of some sample sites

圖3 研究區(qū)地形地貌、地震烈度及泥石流溝分布Fig.3 Geomorphologic,seismicity and distribution of debris flow in research area

2.2 影響因子的選取與敏感性和響應(yīng)分析

影響因子的選取既要考慮因子的是否具有代表性，也要確保各因子之間相互獨(dú)立且選取因子易量化。極震區(qū)震后泥石流動(dòng)儲(chǔ)量主要來源有：①由強(qiáng)震效應(yīng)誘發(fā)可參加泥石流活動(dòng)的溝道崩滑堆積體；②先前溝內(nèi)潛在物源及震后破碎不穩(wěn)定山體經(jīng)強(qiáng)降雨條件下導(dǎo)致的洪水沖刷以及侵蝕等作用累積產(chǎn)生的補(bǔ)給性和次生性物源。因此，綜合考慮泥石流溝的流域規(guī)模、地質(zhì)環(huán)境背景和地震效應(yīng)等條件下初步選取流域面積(A)、相對高差(h)、主溝長度(l)、距發(fā)震斷裂帶距離(d)、溝床平均縱比降(J)和物源靜儲(chǔ)量(V)作為泥石流物源動(dòng)儲(chǔ)量(V0)的影響因子。為了進(jìn)一步研究各影響因子對極震區(qū)震后泥石流物源動(dòng)儲(chǔ)量的敏感程度，分別采用Person相關(guān)系數(shù)(PCC)、灰色關(guān)聯(lián)度(GRG)和最大信息系數(shù)(MIC)，以收集到的汶川地震極震區(qū)區(qū)內(nèi)60條泥石流溝的物源信息為樣本(表1)，對各因子與極震區(qū)泥石流物源動(dòng)儲(chǔ)量進(jìn)行相關(guān)性分析。為消除量綱影響，對樣本各數(shù)據(jù)進(jìn)行歸一化處理，灰色關(guān)聯(lián)度[23]計(jì)算中分辨系數(shù)取ζ=0.5，最大信息系數(shù)基礎(chǔ)參數(shù)取?=0.6，c=15，得到的計(jì)算結(jié)果，如表2所示。

表1 汶川極震區(qū)泥石流溝物源信息樣本Table 1 Material source parameters of the debris flow gully sample in the Wenchuan meizoseismal area

續(xù)表1

No.A/m3h/ml/md/mJ/‰V/104 m3V0/104 m3No.A/m3h/ml/md/mJ/‰V/104 m3V0/104 m3298.439634.023.80239.55964.86292.10305.248762.087.60421.15234.6451.42316.518864.9810.50177.91287.5267.93326.271 8705.6016.89333.93358.2674.12335.351 2883.600.13357.78358.1498.403410.71 8425.820.42316.491 151.41218.72352.181 8202.720.72669.12322.3080.403654.22 90014.20.47204.232 180.57505.34370.069841.580.65622.78222.0077.30387.509355.204.39179.81647.48150.5392.181 2202.684.93455.22122.3092.51405.211 6783.402.51493.53321.3289.16411.211 5961.142.071 400.0013.803.624210.391 4535.512.61263.70727.04108.91430.466501.301.68500.0031.638.97448.631 6054.451.49360.67858.94250.06451.98964.52.440.92395.2981.4028.89461.541 0022.351.53426.3877.6016.3470.689521.951.71488.21191.2946.67488.321 6684.763.86350.42136.0229.25490.204341.210.91358.6821.466.29500.214401.261.82349.2119.172.54510.294601.350.85340.7426.137.94520.646601.984.80333.3359.9015.60531.551 1204.010.45279.30270.1694.50549.801 1624.519.58257.651 754.64162.645536.771 20311.927.04100.921 200.26311.935616.492 3826.20.68384.19742.68194.525571.121 0003.251.36307.6935.938.21580.971 6001.8752.10853.3348.8415.3590.291 5802.11.97752.3817.224.386021.702 9528.90.57331.69366.67136.07

注：1～20為安州區(qū)高川鄉(xiāng)內(nèi)泥石流溝；21～29為綿竹市清平鄉(xiāng)內(nèi)泥石流溝；30～42為汶川縣映秀鎮(zhèn)內(nèi)泥石流溝；43～52為都江堰市龍池鎮(zhèn)內(nèi)泥石流溝；53～55為北川縣境內(nèi)泥石流溝；56～60為汶川縣G213公路沿線泥石流溝。

表2 影響因子與極震區(qū)泥石流物源動(dòng)儲(chǔ)量敏感性分析Table 2 Sensibility analysis between impact factors and dynamic reserve of the debris flows

注：**表示在0.01水平(雙側(cè))上顯著相關(guān)。

由表2可知，基于PCC的相關(guān)性分析中，只能判定泥石流物源總儲(chǔ)量和流域面積兩個(gè)因子與泥石流物源動(dòng)儲(chǔ)量有明顯的線性關(guān)系。而在灰色關(guān)聯(lián)度計(jì)算中，各因子的GRG均大于0.6，可認(rèn)為各影響因子與泥石流物源動(dòng)儲(chǔ)量密切相關(guān)[23]。MIC相較于GRG區(qū)分度十分明顯，更能有效地反映極震區(qū)泥石流物源動(dòng)儲(chǔ)量對各影響因子的敏感性。最終判定，泥石流物源總儲(chǔ)量因子最敏感，而距發(fā)震斷裂帶距離因子最末，選取的各因子均會(huì)在不同程度上影響極震區(qū)泥石流物源動(dòng)儲(chǔ)量，驗(yàn)證了影響因子選取的合理性。需要指出的是，通過對5個(gè)影響因子之間進(jìn)行相關(guān)性系數(shù)計(jì)算，結(jié)果表明：流域面積與主溝長度、流域面積與相對高差和主溝長度與相對高差Person相關(guān)系數(shù)值分別為0.624、0.814和0.715，相關(guān)性較高，其他各影響因子之間Person相關(guān)系數(shù)值均小于0.419，相關(guān)性較弱。在此考慮使用PCA，對輸入因子進(jìn)行降維，避免信息冗余。

3 泥石流物源動(dòng)儲(chǔ)量計(jì)算

以汶川地震極震區(qū)區(qū)內(nèi)60條泥石流溝的物源信息為樣本，從每個(gè)子研究區(qū)中抽取一條泥石流溝作為測試樣本，另外的所有樣本作為訓(xùn)練樣本和驗(yàn)證樣本，其中驗(yàn)證樣本采用“留一法”交叉驗(yàn)證方式得到。對樣本數(shù)據(jù)進(jìn)行PCA降維處理后分別基于五種神經(jīng)網(wǎng)絡(luò)算法對研究區(qū)內(nèi)泥石流動(dòng)儲(chǔ)量分別建立計(jì)算模型，樣本數(shù)據(jù)如表1所示。

3.1 PCA對原始數(shù)據(jù)進(jìn)行預(yù)處理

利用PCA求得的原始樣本數(shù)據(jù)得出的特征值(從大到小)和累計(jì)貢獻(xiàn)率如表3所示。

表3 各主成分因子的特征值及累積貢獻(xiàn)率Table 3 Eigenvalues of each principal factor and its cumulative contribution rate

計(jì)算得到前4個(gè)主成分的累積貢獻(xiàn)率達(dá)到94.197%，因此可以選取第1主成分(PCA-1)、第2主成分(PCA-2)、第3主成分(PCA-3)和第4主成分(PCA- 4)作為神經(jīng)網(wǎng)絡(luò)模型的輸入，這4個(gè)主成分因子特征根對應(yīng)的特征向量如表4所示。將上述特征向量與原始樣本數(shù)據(jù)對應(yīng)相乘，即可得到PCA處理后的樣本數(shù)據(jù)。

表4 前4個(gè)主成分因子特征值對應(yīng)的特征向量Table 4 Eigenvectors of 4 big eigenvalues

3.2 各神經(jīng)網(wǎng)絡(luò)模型參數(shù)的設(shè)置

3.2.1 BP神經(jīng)網(wǎng)絡(luò)模型

網(wǎng)絡(luò)設(shè)置為4- 4-1三層結(jié)構(gòu)，即輸入層節(jié)點(diǎn)為4，隱含層節(jié)點(diǎn)數(shù)按照經(jīng)驗(yàn)公式[24]確定取值范圍后經(jīng)試算，設(shè)置為4，輸出層節(jié)點(diǎn)為1。在訓(xùn)練過程中，學(xué)習(xí)率Lr=0.05，訓(xùn)練精度Ggoal=0.01，最大訓(xùn)練次數(shù)為5 000。

3.2.2 SVM模型

網(wǎng)絡(luò)核函數(shù)設(shè)置為RBF(radial basis function)，基于五折交叉驗(yàn)證法(5-cross validation)的參數(shù)選取結(jié)果如圖4所示。最終選取懲罰函參數(shù)c和核函數(shù)參數(shù)g分別為5.656 9、0.062 5。

CvMES表示尋優(yōu)過程中的均方誤差圖4 基于交叉驗(yàn)證法的SVM參數(shù)選取結(jié)果Fig.4 The results of SVM parameters selection based on cross validation

3.2.3 ELM模型

網(wǎng)絡(luò)設(shè)置sigmoid為激活函數(shù)，設(shè)置隱含層節(jié)點(diǎn)個(gè)數(shù)為一個(gè)循環(huán)數(shù)列，在[1，100]中尋找最優(yōu)隱含層節(jié)點(diǎn)數(shù)。選擇不同的隱含層節(jié)點(diǎn)數(shù)對ELM模型的計(jì)算準(zhǔn)確性有較大影響，最終模型計(jì)算得到不同隱含層節(jié)點(diǎn)數(shù)下均方誤差(MSE)歸一化后的變化曲線。由圖5(a)可知，ELM模型隱含層節(jié)點(diǎn)數(shù)到達(dá)57后MSE(歸一化后)收斂到0.03以下，且隨著節(jié)點(diǎn)數(shù)的增多，節(jié)點(diǎn)數(shù)達(dá)到76后，MSE(歸一化后)呈振蕩趨勢。

3.2.4 PSO-ELM模型

同樣建立以sigmoid方程為激活函數(shù)的PSO-ELM計(jì)算模型，并在[1，100]中搜索最優(yōu)隱含層節(jié)點(diǎn)數(shù)。模型MSE(歸一化后)隨著隱含層節(jié)點(diǎn)數(shù)變化曲線如圖5(b)所示。在PSO算法中：學(xué)習(xí)因子c1=1.5、c2=1.7，慣性權(quán)重w=1，種群規(guī)模Zsizepop=20，最大進(jìn)化代數(shù)為Gmaxgan=100。由圖5(b)可知，PSO-ELM模型隱含層節(jié)點(diǎn)數(shù)達(dá)到40后MSE(歸一化后)收斂到0.03以下，之后趨于平穩(wěn)。說明PSO算法減少了“無用的”隱含層節(jié)點(diǎn)，獲得了更為緊湊的網(wǎng)絡(luò)體系結(jié)構(gòu)，提升了模型的穩(wěn)定性。

圖5 MSE與隱含層節(jié)點(diǎn)數(shù)關(guān)系圖Fig.5 Relationship diagram of the number of hidden neurons with MSE

3.2.5 PSO-ELM_Adaboost耦合模型

PSO-ELM參數(shù)設(shè)置與模型4相同，基分類器數(shù)經(jīng)試算設(shè)置為10。

3.3 預(yù)測結(jié)果與誤差分析

各模型的訓(xùn)練樣本擬合結(jié)果如圖6所示，測試樣本預(yù)測結(jié)果如圖7所示，各樣本間誤差比較如表5 所示。

圖6 訓(xùn)練樣本預(yù)測結(jié)果對比Fig.6 Comparison of test sample training results

為比較各模型之間的可靠性與準(zhǔn)確性，模型誤差采用均方根誤差(RMSE)和平均絕對百分誤差(MAPE)，計(jì)算公式為

(18)

(19)

式中：di為實(shí)際值；Di為預(yù)測值；n為樣本個(gè)數(shù)。

圖7 測試樣本預(yù)測結(jié)果對比Fig.7 Comparison of test sample prediction results

表5 各神經(jīng)網(wǎng)絡(luò)模型誤差對比Table 5 Comparison of error of each neural network model

注：模型計(jì)算用時(shí)基于Intel i7-7700K、16GB DDR4內(nèi)存和windows 7 64位系統(tǒng)的MATLAB R2016a 平臺(tái)。

由表5可知，在單一神經(jīng)網(wǎng)絡(luò)模型中，ELM計(jì)算效果最好，優(yōu)于BP和SVM。并且在模型用時(shí)方面，ELM計(jì)算用時(shí)為0.042 s，顯著低于BP和SVM(用時(shí)分別為3.024、1.431 s)，說明了ELM方法具有耗時(shí)短，效率高的優(yōu)點(diǎn)。接著將兩種耦合模型與單一神經(jīng)網(wǎng)絡(luò)模型進(jìn)行對比，發(fā)現(xiàn)在計(jì)算精度上PSO-ELM和PSO-ELM_AdaBoost具有顯著優(yōu)勢，訓(xùn)練樣本和測試樣本的誤差評價(jià)參數(shù)RMSE、MAPE均小于前兩者。在兩種耦合模型中，由于AdaBoost算法“能夠提高任意給定弱預(yù)測器的預(yù)測精度”的特點(diǎn)， 訓(xùn)練樣本和測試樣本RMSE分別降低了4.17和4.43，MAPE分別降低了1.20%和1.86%，模型誤差有較大程度的減小，所以PSO-ELM_AdaBoost的計(jì)算精度優(yōu)于PSO-ELM。但由于網(wǎng)絡(luò)結(jié)構(gòu)復(fù)雜程度的提升，模型用時(shí)有不可避免的增加。最后將PSO-ELM_AdaBoost模型的精度與其他計(jì)算模型[4-9]相比，精度也相對更高，可以有效地對汶川極震區(qū)內(nèi)泥石流物源動(dòng)儲(chǔ)量進(jìn)行計(jì)算。

對模型產(chǎn)生的誤差進(jìn)行分析，可能的原因有：①由于儀器、方法或人為原因造成野外采集的樣本數(shù)據(jù)不準(zhǔn)確，樣本數(shù)據(jù)的準(zhǔn)確性直接影響模型的精度；②流域地質(zhì)環(huán)境背景有所差異。例如，測試樣本中No.26和No.53在震時(shí)分別誘發(fā)了兩個(gè)巨型滑坡，從而影響了模型預(yù)測精度，總體上各模型計(jì)算值均小于實(shí)際動(dòng)儲(chǔ)量，剔除后測試樣本的各模型的RMSE為23.02、21.59、16.79、13.63和9.91，MAPE為16.32%、14.31%、13.94%、9.73%和6.69%；③坡體結(jié)構(gòu)和巖性組合方式的區(qū)別：地層巖性對形成崩滑體具有重要的作用，其主要影響基巖和堆積體的物理力學(xué)性質(zhì);同時(shí)，不同地層巖性中崩滑體發(fā)育的規(guī)模、類型特征也不盡相同，而崩滑體多寡正是動(dòng)儲(chǔ)量形成的重要內(nèi)部因素[25]；④其他影響因素：物源組成顆粒的級配、密實(shí)度和泥石流溝坡向等。

3.4 模型在其他極震區(qū)的適宜性檢驗(yàn)

為檢驗(yàn)?zāi)Ｐ驮谄渌麡O震區(qū)的適宜性，隨機(jī)選擇玉樹地震極震區(qū)結(jié)古鎮(zhèn)布慶隆溝、蘆山地震極震區(qū)冷木溝和中崗溝作為通用適宜性檢驗(yàn)樣本，將數(shù)據(jù)進(jìn)行PCA降維處理后，根據(jù)所得的PSO-ELM_AdaBoost耦合模型，進(jìn)行動(dòng)儲(chǔ)量計(jì)算及誤差分析。通用適宜性檢驗(yàn)樣本具體參數(shù)及預(yù)測值如表6所示。

表6 通用適宜性樣本預(yù)測結(jié)果及誤差Table 6 Prediction results and errors of suitability samples outside the study area

由表6可知，預(yù)測精度依然較好，能夠滿足實(shí)際要求。但與汶川地震極震區(qū)內(nèi)預(yù)測誤差相比，誤差相對較大，且計(jì)算動(dòng)儲(chǔ)量均大于實(shí)際值。分析考慮是：玉樹地震和蘆山地震強(qiáng)度不及汶川地震且區(qū)域內(nèi)平均降雨量較汶川地震極震區(qū)少[26-27]。因此，此處產(chǎn)生的誤差不僅與3.3節(jié)提出的誤差相關(guān)可能還與地震強(qiáng)度的高低和區(qū)域降雨條件的不同等因素有關(guān)。

4 結(jié)論

(1)以泥石流溝的流域面積、相對高差、主溝長度、溝床平均縱比降、較發(fā)震斷裂帶距離和物源總儲(chǔ)量作為泥石流物源動(dòng)儲(chǔ)量的影響因素，對采集的汶川極震區(qū)區(qū)內(nèi)的60條泥石流溝的物源信息數(shù)據(jù)，經(jīng)PCA降維后，進(jìn)行訓(xùn)練與預(yù)測，建立的PSO-ELM_AdaBoost耦合模型具有精度高、可控性強(qiáng)和穩(wěn)定性好的特點(diǎn)。計(jì)算精度顯著優(yōu)于BP、SVM、ELM和PSO-ELM，與其他計(jì)算模型相比這些誤差也相對較小，并且對其他極震區(qū)泥石流物源動(dòng)儲(chǔ)量的預(yù)測精度依然較高，滿足實(shí)際要求。因此該模型可以在極震區(qū)泥石流動(dòng)儲(chǔ)量計(jì)算中發(fā)揮一定的作用，為設(shè)計(jì)泥石流防治工程提供有價(jià)值的參考。

(2)通過對模型誤差的分析，在同一極震區(qū)內(nèi)考慮是地質(zhì)環(huán)境背景、坡體結(jié)構(gòu)、巖性組合、物源組成等因素的影響。 而在不同極震區(qū)內(nèi)誤差的產(chǎn)生可能還與地震自身強(qiáng)度和區(qū)域降雨條件等因素有關(guān)。因此，下一步工作的重點(diǎn)就是如何對以上因素進(jìn)行量化再分析，以繼續(xù)提升本模型的適宜性和精度。

(3)隨著神經(jīng)網(wǎng)絡(luò)等機(jī)器學(xué)習(xí)技術(shù)的不斷發(fā)展，算法預(yù)測模型精度越來越高，但準(zhǔn)確的泥石流實(shí)測數(shù)據(jù)是建立好的算法預(yù)測模型的基礎(chǔ)，因此建立泥石流多發(fā)區(qū)域的泥石流溝數(shù)據(jù)庫具有重要意義。