基于主成分分析與多分類(lèi)支持向量機(jī)的單溝泥石流危險(xiǎn)性預(yù)測(cè)

劉超 ，喬圣揚(yáng)

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

0 引言

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

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

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

1 基本理論

1.1 主成分分析

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

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

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

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

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

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

（3）計(jì)算特征值與特征向量

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

（4）選取主成分

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

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

（5）計(jì)算主成分

主成分是通過(guò)選擇前m個(gè)特征值所對(duì)應(yīng)的特征向量來(lái)計(jì)算的。計(jì)算的主成分為：

式中，把Yk稱(chēng)為第k主成分。主成分的線(xiàn)性表達(dá)式寫(xiě)為：

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

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

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

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

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

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

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

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

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

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

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

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

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

（3）判別分類(lèi)

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

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

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

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

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

2.1 評(píng)價(jià)指標(biāo)確定

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

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

2.2 主成分分析提取

對(duì)原始數(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)間存在大量重疊信息。因此通過(guò)PCA提取主要信息。

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

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

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

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

2.3 MSVM模型的建立

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

其具體過(guò)程如圖2所示。

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

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

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

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

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

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

3 討論

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

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

圖4 PCA-MSVM模型與MSVM模型預(yù)測(cè)結(jié)果對(duì)比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模型之前，對(duì)原始數(shù)據(jù)進(jìn)行主成分分析，可以去除原始變量之間的相關(guān)性，達(dá)到一定的降維效果，建立互不相關(guān)且保留原有信息的綜合指標(biāo)，可有效提高預(yù)測(cè)性能。

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

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

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

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

4 結(jié)論

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

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

（2）通過(guò)MSVM模型與PCA-MSVM模型的對(duì)比分析，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時(shí)，模型準(zhǔn)確率較低，并隨著K增大而增加，K≥5時(shí)，K值達(dá)到最大為90%，趨于穩(wěn)定。在程序計(jì)算過(guò)程中，隨著K的增大運(yùn)算時(shí)間增加。