劉超 ，喬圣揚(yáng)

河北地質(zhì)大學(xué)a.城市地質(zhì)與工程學(xué)院b.河北省高校生態(tài)環(huán)境地質(zhì)應(yīng)用技術(shù)研發(fā)中心，河北石家莊050031

0 引言

泥石流是山區(qū)常見的一種地質(zhì)災(zāi)害，嚴(yán)重威脅到了人類的生命與財產(chǎn)安全[1-3]，因此對泥石流開展危險性預(yù)測尤為重要。泥石流危險性的預(yù)測方法很多，如專家經(jīng)驗(yàn)法[4]、經(jīng)驗(yàn)公式計算法[5]、軟件模擬法[6]、GIS法[7]等。

近年來，隨著機(jī)器學(xué)習(xí)算法熱度的逐漸升高，許多學(xué)者將機(jī)器學(xué)習(xí)的方法引入到泥石流危險性等級預(yù)測中，常見的機(jī)器學(xué)習(xí)方法有隨機(jī)森林[8]、神經(jīng)網(wǎng)絡(luò)[9-11]、支持向量機(jī)[12-13]等。劉永垚[14]等（2018）建立基于交叉驗(yàn)證的隨機(jī)森林模型對泥石流易發(fā)性進(jìn)行評價。張永宏[15]等（2018）使用平均影響值算法篩選主要影響指標(biāo)，結(jié)合神經(jīng)網(wǎng)絡(luò)算法預(yù)測泥石流易發(fā)程度。Qian X[16]等（2016）使用粒子群算法的支持向量機(jī)模型對泥石流危險性有較高的辨識精度。Shiuan W[17]等（2009）建立決策樹與支持向量機(jī)的混合模型對樣區(qū)發(fā)生泥石流的危險性類別進(jìn)行預(yù)測，模型預(yù)測準(zhǔn)確率為73%。寧志杰與周愛紅[18]（2020）提出空間變異性會影響支持向量機(jī)模型的適用性，需要合理地選擇評價指標(biāo)。

本文基于白龍江流域泥石流，通過主成分分析（Principle Component Analysis，PCA）提取3個互為獨(dú)立的主成分，進(jìn)一步使用5折交叉驗(yàn)證的方式抽取訓(xùn)練樣本與預(yù)測樣本。聯(lián)合使用多分類支持向量機(jī)（Multi-classification Support Vector Machine，MSVM）對泥石流危險性類別進(jìn)行分類，建立主成分分析與多分類支持向量機(jī)預(yù)測模型（PCAMSVM），可為預(yù)測白龍江流域泥石流危險性提供科學(xué)依據(jù)。

1 基本理論

1.1 主成分分析

泥石流危險性通常被眾多指標(biāo)所影響，而指標(biāo)與指標(biāo)中的信息會存在部分程度的相關(guān)性，使各指標(biāo)間出現(xiàn)信息重疊問題，因此利用主成分分析把存在相關(guān)性的指標(biāo)線性組合后，形成互為獨(dú)立的新綜合指標(biāo)，從而達(dá)到降維目的。具體步驟如下：

（1）標(biāo)準(zhǔn)化原始數(shù)據(jù)

采集到的數(shù)據(jù)構(gòu)成矩陣X，即：

式中，n為樣本個數(shù)，X1，X2，…，XP為涉及的P個指標(biāo)。其中：

（2）計算相關(guān)系數(shù)矩陣

計算X*的相關(guān)系數(shù)矩陣：

（3）計算特征值與特征向量

計算相關(guān)系數(shù)矩陣R的p個特征值分別為λ1，λ2，…，λp，以及每個特征值λk（k=1,2,…,p）對應(yīng)的特征向量為：

（4）選取主成分

主成分是原始指標(biāo)的線性組合。通過主成分分析得到p個主成分，前m個主成分的累計方差貢獻(xiàn)率為：

在解決實(shí)際問題時，累計方差貢獻(xiàn)率通常選取大于85%時所確定的前m個主成分。

（5）計算主成分

主成分是通過選擇前m個特征值所對應(yīng)的特征向量來計算的。計算的主成分為：

式中，把Yk稱為第k主成分。主成分的線性表達(dá)式寫為：

1.2 多分類支持向量機(jī)

支持向量機(jī)是一種分類方法，常用的解決多分類問題（分類數(shù)大于2）的總體思路是把多分類轉(zhuǎn)化為若干個二分類。這樣可將具有N個分類的數(shù)據(jù)模型轉(zhuǎn)化為個二分類模型，即每兩個類別之間用二分類方法進(jìn)行一次分類[19]。

在MSVM中將樣本分為訓(xùn)練樣本與測試樣本。設(shè)訓(xùn)練樣本為：

其中，xt為p維向量，p為指標(biāo)數(shù)，yt為分類標(biāo)簽，l為樣本個數(shù)。二分類的具體步驟如下：

（1）構(gòu)造最優(yōu)超平面

設(shè)二分類模型的兩個分類為i和j，構(gòu)造最優(yōu)超平面如圖1所示：

圖1 最優(yōu)超平面示意圖Fig.1 Hyperplane diagram of optimal classification

圖1中，黑色小球表示分類i的樣本（分類標(biāo)簽用1來表示），白色小球表示分類j的樣本（分類標(biāo)簽為用-1來表示）。設(shè)H為最優(yōu)超平面，H1、H2分別平行于H，是經(jīng)過兩類樣本并且離H最近的面，H1與H2的平分面為H。最大分類間隔是指H1與H2之間的距離。因此H的方程不妨設(shè)為：

式中，ω和b為待求解的參數(shù)。圖1中，H1與H2的方程可以相應(yīng)的設(shè)為：

有時兩類樣本中個別數(shù)據(jù)點(diǎn)不能精確地被最優(yōu)超平面正確劃分，則引入松弛項tξ≥0，放寬約束條件，即：

式中，C為懲罰因子，是對錯誤分類樣本增加懲罰因素的參數(shù)。

為了求解上述凸二次規(guī)劃問題，構(gòu)建拉格朗日函數(shù)，并對此函數(shù)進(jìn)行求偏導(dǎo)，再轉(zhuǎn)化為對偶問題進(jìn)行求解，可求得最優(yōu)解w*、b*，確定分類函數(shù)[20]。對于非線性問題引入核函數(shù)，最終得到的分類函數(shù)為：

式中，sgn（·）為符號函數(shù)，x為預(yù)測樣本，yt為分類標(biāo)簽，at*為對偶變量最優(yōu)解，K（xt ,x）為核函數(shù)。

（3）判別分類

經(jīng)過對相鄰玻璃折彎角度多次進(jìn)行計算機(jī)三維實(shí)體建模及實(shí)體渲染。總結(jié)出玻璃折彎角度應(yīng)控制在不小于170°，使玻璃之間光線折射效果可以得到有效的緩解，且觀感效果極佳（見圖2）。

通過分類函數(shù)判斷類別，f（x）=1時，判斷類別為分類i，f（x）=-1時，判斷類別為分類j。對于多分類問題，應(yīng)用MSVM進(jìn)行分類，對于每一個樣本數(shù)據(jù)進(jìn)行CN2次二分類，采用投票法判別最終分類。即在i、j分類時，分類函數(shù)判斷為i類，就在i類上加一票，判斷為j類，在j類加一票，然后進(jìn)行其他二分類，直到所有分類函數(shù)都完成判別后，再判別樣本所屬分類，樣本所屬分類是根據(jù)最后統(tǒng)計類別票數(shù)最多而確定的。

2 預(yù)測模型的建立

以白龍江流域泥石流為例[20]，聯(lián)合使用主成分分析和多分類支持向量機(jī)的方法，利用Matlab軟件編譯程序，進(jìn)行泥石流危險性預(yù)測。

白龍江流域位于青藏高原與四川盆地的過渡區(qū)，該地區(qū)高差大，縱坡陡峭，山坡坡度較大，降雨具有強(qiáng)度大、歷時短、突發(fā)性等特點(diǎn)[21]。此流域極易發(fā)生泥石流，給當(dāng)?shù)卦斐闪藝?yán)重的危害，特別是2010年舟曲發(fā)生了特大泥石流災(zāi)害，縣城幾乎被泥石流摧毀殆盡，死亡和失蹤人數(shù)超過兩千余人[3]。

2.1 評價指標(biāo)確定

以文獻(xiàn)[20]中白龍江流域泥石流30組數(shù)據(jù)為例，選取7個指標(biāo)作為評價因子,如表1所示。泥石流危險性等級分為四類，低度危險用標(biāo)簽“1”表示、中度危險用標(biāo)簽“2”表示、高度危險用標(biāo)簽“3”表示、極高危險用標(biāo)簽“4”表示[22]。

表1 白龍江流域泥石流原始數(shù)據(jù)[20]Table 1 Primary data of debris flow in Bailong River basin

2.2 主成分分析提取

對原始數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化，進(jìn)一步根據(jù)式（6）可得到相關(guān)系數(shù)矩陣。

式中，相關(guān)系數(shù)越接近1表明指標(biāo)與指標(biāo)之間越相關(guān)，表明各指標(biāo)間存在大量重疊信息。因此通過PCA提取主要信息。

然后計算相關(guān)系數(shù)矩陣的特征值及相應(yīng)的貢獻(xiàn)率，如表2所示。從表2中可以看出，前3個主成分累計貢獻(xiàn)率達(dá)到85.512%，因此提取前3個主成分。

表2 特征值、貢獻(xiàn)率及累計貢獻(xiàn)率Table 2 Eigenvalue, contribution rate and cumulative contribution rate

其中，X1*、X2*、X3*、X4*、X5*、X6*、X7*分別表示流域面積、松散固體物質(zhì)總量、爆發(fā)規(guī)模、工程治理攔砂量、溝床平均比降、流域切割密度及不穩(wěn)定溝床比這7個指標(biāo)標(biāo)準(zhǔn)化之后的數(shù)值。主成分的計算結(jié)果如表3所示。

表3 主成分的計算結(jié)果Table 3 Principal component calculation results

2.3 MSVM模型的建立

MSVM模型選取5折交叉驗(yàn)證方法抽取訓(xùn)練樣本與測試樣本，再進(jìn)行MSVM預(yù)測。將30個樣本按序號等分為5個子集，每個子集6個樣本。在每次實(shí)驗(yàn)中，按照順序選擇其中1個子集作為預(yù)測樣本，另外4個子集作為訓(xùn)練樣本，進(jìn)行模型的訓(xùn)練與預(yù)測，得到預(yù)測的準(zhǔn)確率。同樣的實(shí)驗(yàn)再進(jìn)行4次，模型準(zhǔn)確率取5次預(yù)測的準(zhǔn)確率平均值。

其具體過程如圖2所示。

圖2 5折交叉驗(yàn)證示意圖Fig.2 5-fold cross validation schematic diagram

MSVM模型中選取高斯核函數(shù)，參數(shù)尋優(yōu)算法選取粒子群算法。

2.4 模型預(yù)測結(jié)果

按照5折交叉驗(yàn)證方法進(jìn)行PCA-MSVM預(yù)測，預(yù)測結(jié)果如圖3和表4所示。由圖3可知，在第2次實(shí)驗(yàn)、第3次實(shí)驗(yàn)以及第5次實(shí)驗(yàn)各有1次誤判。由表4可知，PCA-MSVM模型準(zhǔn)確率達(dá)到了90%，能夠滿足實(shí)際工程的需求。

表4 基于交叉驗(yàn)證的PCA-MSVM模型預(yù)測結(jié)果統(tǒng)計Table 4 Prediction result statistics of PCA-MSVM model based on cross validation

圖3 基于5折交叉驗(yàn)證的PCA-MSVM模型預(yù)測結(jié)果Fig.3 Prediction results of PCA-MSVM model based on 5-fold cross validation

3 討論

3.1 PCA對MSVM模型準(zhǔn)確率的影響

基于白龍江流域泥石流30組數(shù)據(jù)，對比單獨(dú)MSVM和PCA-MSVM模型預(yù)測的準(zhǔn)確率，預(yù)測結(jié)果如圖4所示。

圖4 PCA-MSVM模型與MSVM模型預(yù)測結(jié)果對比Fig.4 Comparison of prediction results between PCAMSVM model and MSVM model

由圖4可知，PCA-MSVM模型準(zhǔn)確率可達(dá)90%，MSVM模型準(zhǔn)確率僅為76.67%，PCA-MSVM模型明顯優(yōu)于MSVM模型。在建立MSVM模型之前，對原始數(shù)據(jù)進(jìn)行主成分分析，可以去除原始變量之間的相關(guān)性，達(dá)到一定的降維效果，建立互不相關(guān)且保留原有信息的綜合指標(biāo)，可有效提高預(yù)測性能。

3.2 交叉驗(yàn)證折數(shù)K對模型準(zhǔn)確率的影響

在基于PCA-MSVM模型下分析模型預(yù)測準(zhǔn)確率隨K值的變化，如圖5所示。

圖5 模型準(zhǔn)確率隨K值的變化Fig.5 Model accuracy varies with K value

由圖5可知，模型準(zhǔn)確率隨著K值增大而增加，當(dāng)K≥5時，模型準(zhǔn)確率可以達(dá)到最大值90%。K值較小時，訓(xùn)練樣本數(shù)量較少，訓(xùn)練出的模型泛化能力較差，導(dǎo)致模型準(zhǔn)確率較低。隨著K值的增大，模型泛化能力增強(qiáng)，模型準(zhǔn)確率也隨之提高。K值大于5時，模型準(zhǔn)確率隨K值增大變化并不明顯。K值增大會加大程序運(yùn)行時間與運(yùn)算復(fù)雜程度，因此上文預(yù)測模型中選取了K為5，即5折交叉驗(yàn)證。

4 結(jié)論

本文選取白龍江流域泥石流為例，通過主成分分析對原始數(shù)據(jù)進(jìn)行處理，提取主成分，基于交叉驗(yàn)證提取訓(xùn)練樣本與預(yù)測樣本，聯(lián)合多分類支持向量機(jī)算法對泥石流危險性進(jìn)行預(yù)測，可以得到以下結(jié)論：

（1）PCA-MSVM模型能夠很好的預(yù)測泥石流危險性，以白龍江流域泥石流30組數(shù)據(jù)為例，使用5折交叉驗(yàn)證的PCA-MSVM模型進(jìn)行預(yù)測，模型準(zhǔn)確率可達(dá)90%。

（2）通過MSVM模型與PCA-MSVM模型的對比分析，PCA-MSVM模型準(zhǔn)確率較高。引入PCA可以起到降維作用，減少相關(guān)性，提高M(jìn)SVM模型準(zhǔn)確率。

（3）模型準(zhǔn)確率隨交叉驗(yàn)證折數(shù)K增大而增加，當(dāng)K達(dá)到某一值后，模型準(zhǔn)確率趨于平緩。本文中，K<5時，模型準(zhǔn)確率較低，并隨著K增大而增加，K≥5時，K值達(dá)到最大為90%，趨于穩(wěn)定。在程序計算過程中，隨著K的增大運(yùn)算時間增加。