鄧 欽 宣, 張 韶 鵬, 邵 偉, 楊 洋
(1.四川足木足河流域水電開發(fā)有限公司,四川 成都 610041;2.四川大學(xué)水利水電學(xué)院,四川 成都 610200)
在水利水電工程中,大壩壩基及邊坡穩(wěn)定是保證水利水電工程安全至關(guān)重要的影響因素之一[1]。針對(duì)水利水電工程結(jié)構(gòu),《水利水電工程結(jié)構(gòu)可靠度設(shè)計(jì)統(tǒng)一標(biāo)準(zhǔn)》[2]要求采用基于可靠度理論的分項(xiàng)系數(shù)設(shè)計(jì)方法進(jìn)行結(jié)構(gòu)設(shè)計(jì),而可靠度設(shè)計(jì)要求提供巖石物理力學(xué)參數(shù)的分布信息。巖石在形成過(guò)程中經(jīng)歷了長(zhǎng)時(shí)間跨度的復(fù)雜地質(zhì)過(guò)程,巖體的材料和工程參數(shù)具有不同的分布特征且通常具有相關(guān)性,其聯(lián)合分布具有復(fù)雜特征[3-4]。若各參數(shù)的聯(lián)合概率分布估計(jì)不當(dāng),或者簡(jiǎn)單地忽略它們的相關(guān)性,會(huì)導(dǎo)致錯(cuò)估結(jié)構(gòu)失效概率,導(dǎo)致有偏差的設(shè)計(jì)結(jié)果,從而直接影響工程安全[5-6]。
以準(zhǔn)確估計(jì)巖石物理力學(xué)概率分布為目標(biāo),國(guó)內(nèi)外諸多學(xué)者進(jìn)行了大量研究,經(jīng)歷了從確定性分析到不確定性分析,從簡(jiǎn)單概率統(tǒng)計(jì)到概率分布擬合,從單參數(shù)邊緣分布到多參數(shù)聯(lián)合分布的研究歷程,發(fā)展到了考慮不確定性的各參數(shù)不同邊緣分布類型及不同相關(guān)結(jié)構(gòu)的概率分布構(gòu)建階段[7-12]。然而,因巖石土體材料的場(chǎng)地或區(qū)域特異性、材料天然變異性和形成過(guò)程的復(fù)雜性,其聯(lián)合分布常有悖于經(jīng)典概率分布,且需要事先進(jìn)行分布類型和相關(guān)結(jié)構(gòu)的假設(shè),而現(xiàn)有方法無(wú)法恰當(dāng)解決該難題[13]。為此,本文將基于貝葉斯機(jī)器學(xué)習(xí)框架的高斯混合模型(Gaussian Mixture Model, GMM)應(yīng)用于大渡河上游巖石物理力學(xué)參數(shù)概率分布構(gòu)建任務(wù)中,在考慮統(tǒng)計(jì)不確定條件下精準(zhǔn)刻畫聯(lián)合分布特征及各參數(shù)間的相關(guān)性,構(gòu)建區(qū)域性多元巖石物理力學(xué)參數(shù)概率分布,為該區(qū)域的后續(xù)工程建設(shè)和設(shè)計(jì)提供參考。
高斯混合模型是混合分量為高斯分布的一種特殊混合模型,于1894年由生物統(tǒng)計(jì)學(xué)家Karl Person首次提出并應(yīng)用于生物統(tǒng)計(jì)學(xué)中的偏態(tài)數(shù)據(jù)分析[14-15],并在之后的百年間獲得了持續(xù)而深入的發(fā)展,已成為機(jī)器學(xué)習(xí)領(lǐng)域不可或缺的模型之一。高斯混合模型不僅具有極強(qiáng)的靈活性,而且繼承了高斯分布的數(shù)學(xué)便利性,已被成功應(yīng)用于多學(xué)科、多領(lǐng)域的研究,包括航空航天、醫(yī)學(xué)、信號(hào)處理與分析、經(jīng)濟(jì)學(xué)以及社會(huì)科學(xué)等諸多方面[16]。
p(X|ω)=π1N(X|μ1,Σ1)+…+πkN(X|μk,Σk)+πKN(X|μK,ΣK)
(1)
一般形式為:
(2)
式中K為高斯分量的個(gè)數(shù);N(X|μk,Σk)為1到K個(gè)高斯分量中第k個(gè)分量的概率密度函數(shù);πk為第k個(gè)高斯分量的權(quán)重參數(shù);向量π=[π1,…,πK]也被稱為權(quán)重分布。為使GMM模型的概率密度函數(shù)在其定義域上積分為1,權(quán)重πk必須滿足以下約束條件:
(3)
式中μk為第k個(gè)高斯分量的均值向量;Σk為第k個(gè)高斯分量的協(xié)方差矩陣,由各維度變量的方差和各變量之間的相關(guān)系數(shù)矩陣表示如下:
(4)
而在研究模型復(fù)雜度時(shí),協(xié)方差矩陣Σk作為d(d+1)/2個(gè)參數(shù)進(jìn)行計(jì)算。因此,多維高斯混合模型共有K(d+d(d+1)/2+1)-1個(gè)模型參數(shù),模型參數(shù)求解涉及高維求解問(wèn)題。為直觀展示其參數(shù)規(guī)模,圖1展示了不同變量維度條件下模型參數(shù)個(gè)數(shù)Np與高斯分量個(gè)數(shù)K的對(duì)應(yīng)關(guān)系。
圖1 不同變量維度條件下模型參數(shù)個(gè)數(shù)Np與高斯分量個(gè)數(shù)K的對(duì)應(yīng)關(guān)系
由圖1可知,高斯混合模型的模型參數(shù)個(gè)數(shù)Np隨著高斯分量個(gè)數(shù)K的增加而線性增加,10維高斯混合分布在高斯分量個(gè)數(shù)為3時(shí),即有高達(dá)198個(gè)參數(shù)。高斯模型復(fù)雜度在為概率分布擬合帶來(lái)極大便利的同時(shí),也給模型學(xué)習(xí)帶來(lái)了挑戰(zhàn),因此,本文采用貝葉斯機(jī)器學(xué)習(xí)方法解決模型學(xué)習(xí)難題。
除去上文提到的模型復(fù)雜度難題外,高斯混合模型的學(xué)習(xí)還因其為隱變量模型的特性存在標(biāo)簽切換問(wèn)題,而在考慮不確定性分析的前提下,被廣泛應(yīng)用于求解GMM模型的EM算法失效[17]。因此,本節(jié)應(yīng)用筆者之前提出的貝葉斯學(xué)習(xí)框架解決模型學(xué)習(xí)問(wèn)題。
(5)
(6)
(7)
式中p(xi|Si,ωK,MK)為觀測(cè)數(shù)據(jù)xi屬于第Si個(gè)高斯分量的似然函數(shù);p(Si|ωK,Mg)為在給定模型參數(shù)ωK和模型MK的情況下,觀測(cè)數(shù)據(jù)xi屬于第Si個(gè)高斯分量的概率。一旦隱變量SK被確定,則觀測(cè)數(shù)據(jù)xi所歸屬的高斯分量即被確定。為定量描述觀測(cè)數(shù)據(jù)xi的歸屬情況,引入變量nk,k= 1, 2, …,K,其代表了屬于第k個(gè)高斯分量的觀測(cè)數(shù)據(jù)樣本的個(gè)數(shù)。
(8)
式中xi,k代表了屬于第k個(gè)高斯分量的第i個(gè)觀測(cè)數(shù)據(jù);p(xi,k|μk,∑k為第i個(gè)觀測(cè)數(shù)據(jù)歸屬于第k個(gè)高斯分量概率密度函數(shù)(PDF)值,可表示為:
(9)
GMM參數(shù)ωK(例如:π,μk, Σk,k= 1, 2, …,K)的共軛先驗(yàn)分布p(ωk|Mk)的設(shè)置如下[16-19]:
π~Dirichlet(α1,…,αk)
(10)
μk~N(bk,∑k/Bk)
(11)
(12)
式中 權(quán)重參數(shù)π= [π1,π2, …,πK] 服從狄利克雷分布(Dirichlet distribution),其概率密度函數(shù)PDF為:
(13)
式中a1,a2, …,aK為迪利克雷分布的分布參數(shù)。
均值參數(shù)μk的共軛先驗(yàn)分布為以bk為均值,∑k/Bk為方差的多維高斯分布,其中Bk為事先給定的比例因子。
(14)
式中 Гd(·)為d維的Gamma函數(shù),其表達(dá)式為:
(15)
在多維高斯分布的求解算法中,Wishart分布常被用來(lái)構(gòu)造精度矩陣的共軛先驗(yàn)分布。概率密度函數(shù)表達(dá)式中Ck為Wishart分布的尺度參數(shù),其為一個(gè)d×d的對(duì)稱非奇異方陣。ck為Wishart分布的自由度參數(shù),其與分布的自由度vk存在線性關(guān)系(即ck=vk/2)。
至此,表征GMM參數(shù)和隱變量的聯(lián)合后驗(yàn)分布公式(5)的分子部分已經(jīng)完全給出。但是,由于歸一化常數(shù)仍然未知,仍舊無(wú)法得到參數(shù)和隱變量的后驗(yàn)分布樣本。為了在考慮標(biāo)簽切換的前提下解決此問(wèn)題,本文采用隨機(jī)吉布斯抽樣方法進(jìn)行模型參數(shù)求解,此處不贅述。
在應(yīng)用數(shù)學(xué)模型進(jìn)行數(shù)據(jù)分析建模時(shí),模型不確定性是無(wú)法避免的,而對(duì)于高斯混合模型而言,其模型不確定性的最大來(lái)源就是高斯分量個(gè)數(shù)K的不確定性。確定高斯混合模型中的高斯分量個(gè)數(shù)會(huì)不可避免地涉及過(guò)擬合混合模型的學(xué)習(xí)問(wèn)題,而由于標(biāo)簽切換和過(guò)擬合問(wèn)題的存在,此時(shí)的高斯混合模型經(jīng)常是不可識(shí)別的[18]。
隨著計(jì)算科學(xué)的發(fā)展,科研工作者們提出了許多不同的模型選擇方法來(lái)嘗試考慮模型選擇的不確定性。其中基于模型證據(jù)的貝葉斯模型比選方法同信息準(zhǔn)則等方法比較具有優(yōu)越性,因此,本文應(yīng)用貝葉斯模型比選方法,權(quán)衡模型復(fù)雜度和擬合優(yōu)度,選擇最優(yōu)備選模型。
(16)
(17)
然而,當(dāng)高斯分量個(gè)數(shù)大于1時(shí),高斯混合模型的模型證據(jù)并沒(méi)有解析解。在算法實(shí)踐中,模型證據(jù)需要通過(guò)合適的后處理手段(postprocessing manner)求得。在后處理算法中,基于隨機(jī)模擬方法的模型證據(jù)求解方法(Simulation-based approximations)的應(yīng)用最為成功。本研究主要應(yīng)用隨機(jī)模擬方法中的橋采樣(Bridge Sampling, BS)方法求解GMM的模型證據(jù)。
大渡河位于四川西部,是岷江水系最大的支流,年徑流量470億m3,干流全長(zhǎng)1 062 km,天然落差4 175.0 m,為我國(guó)重要的水電能源基地之一。如今,大渡河干流形成了以下爾呷為龍頭的28級(jí)開發(fā)方案,其中上游共規(guī)劃有3級(jí)水電站,自上而下分別為下爾呷、巴拉和達(dá)維水電站,其中巴拉水電站在建,下爾呷、達(dá)維正處于項(xiàng)目前期階段。
大渡河上游巖石主要以變質(zhì)砂巖和板巖為主,巴拉水電站區(qū)域分布有花崗巖侵入?yún)^(qū),本文收集了大渡河上游3級(jí)水電站的巖石物理力學(xué)參數(shù)試驗(yàn)數(shù)據(jù),以巖體干密度ρd、飽和吸水率w、飽和抗壓強(qiáng)度Rw為主要研究對(duì)象,考慮統(tǒng)計(jì)不確定性和參數(shù)相關(guān)性,通過(guò)所提方法構(gòu)建大渡河上游巖石物理力學(xué)參數(shù)概率分布,驗(yàn)證所提方法有效性,為大渡河上游區(qū)域水電工程可靠度設(shè)計(jì)提供依據(jù)。大渡河上游巖石物理力學(xué)參數(shù)基本統(tǒng)計(jì)信息見表1。
表1 大渡河上游巖石物理力學(xué)參數(shù)基本統(tǒng)計(jì)信息
如表1所示,所選的三個(gè)參數(shù)的統(tǒng)計(jì)特征之間存在較大差異,其中干密度ρd的分布較集中,變異系數(shù)僅為0.01,而飽和吸水率w及飽和抗壓強(qiáng)度Rw的變異性較大,其中飽和吸水率的變異系數(shù)甚至高達(dá)0.91,這無(wú)疑為參數(shù)概率分布特征的準(zhǔn)確表征帶來(lái)了困難。
大渡河上游巖石物理力學(xué)參數(shù)試驗(yàn)數(shù)據(jù)的二維散點(diǎn)矩陣(圖2),展示了各參數(shù)的頻率統(tǒng)計(jì)直方圖以及各參數(shù)之間的二維數(shù)據(jù)分布散點(diǎn)。
圖2 大渡河上游巖石物理力學(xué)參數(shù)試驗(yàn)數(shù)據(jù)的二維散點(diǎn)圖矩陣(n=222)
如圖2所示,參數(shù)并不嚴(yán)格服從正態(tài)分布且具有多模態(tài)特征,且各參數(shù)之間明顯具有非線性相關(guān)關(guān)系,很難通過(guò)現(xiàn)有的方法對(duì)這種復(fù)雜相關(guān)的特征進(jìn)行描述。為了恰當(dāng)表述巖石力學(xué)參數(shù)的多模態(tài)特征,同時(shí)處理參數(shù)之間的非線性相關(guān)關(guān)系,接下來(lái)應(yīng)用所提方法對(duì)大渡河上游巖石物理力學(xué)參數(shù)試驗(yàn)數(shù)據(jù)進(jìn)行分析。
首先假設(shè)共有5個(gè)備選模型M1,M2,M3,M4,M5,即高斯分量個(gè)數(shù)K的取值范圍為從1到5的正整數(shù)(Kmax= 5),模型下標(biāo)的數(shù)字即代表了其具有的高斯分量個(gè)數(shù)。備選模型的權(quán)重參數(shù),各高斯分量的均值和精度矩陣的共軛先驗(yàn)分布,分別取為狄利克雷分布、正態(tài)分布和Wishart分布。然后應(yīng)用RGS-GS方法學(xué)習(xí)高斯混合模型,設(shè)置RGS方法后驗(yàn)樣本數(shù)、burn-in樣本數(shù)和BS方法中的重要性抽樣樣本數(shù)均為10 000。
貝葉斯模型選擇結(jié)果見圖3。圖中展示了5個(gè)備選模型的模型概率和模型證據(jù)的對(duì)數(shù)值,其中模型證據(jù)對(duì)數(shù)值用帶有實(shí)心方形標(biāo)記的實(shí)線表示,模型概率用直方圖表示。由圖3可知,隨著模型分量個(gè)數(shù)的增加,模型證據(jù)的對(duì)數(shù)值在K= 4時(shí)達(dá)到最大值-489.08,之后下降至-489.68。因此,最優(yōu)模型高斯分量個(gè)數(shù)為4的高斯混合模型GMM4,其對(duì)應(yīng)的模型概率為0.59。
圖3 貝葉斯模型選擇結(jié)果
在進(jìn)行最優(yōu)模型結(jié)構(gòu)進(jìn)行選擇后,應(yīng)用所提貝葉斯學(xué)習(xí)框架對(duì)模型參數(shù)進(jìn)行學(xué)習(xí)。圖4(a)到(c)分別展示了各巖石力學(xué)參數(shù)的邊緣概率分布擬合結(jié)果。
圖4 各巖石力學(xué)參數(shù)的邊緣概率分布擬合結(jié)果
圖4中,黑色實(shí)線代表了PDF的最可能值,由黑色虛線圍成的區(qū)域代表了PDF的95%置信區(qū)間,反映了模型統(tǒng)計(jì)不確定性的大小,直方圖為歸一化的頻率分布直方圖,展示了數(shù)據(jù)本身的分布特征及模型統(tǒng)計(jì)不確定性的大小。由圖4可知,歸一化的頻率分布直方圖分布于95%置信區(qū)間內(nèi),且與MPV值十分接近,這表明高斯混合模型GMM4的邊緣概率分布不僅能很好地?cái)M合觀測(cè)數(shù)據(jù)的實(shí)際邊緣分布特征,還能正確地表征其不確定性。
通過(guò)二維變量聯(lián)合概率密度對(duì)數(shù)值等值線圖,展示多維高斯混合模型對(duì)多維巖土體參數(shù)的聯(lián)合分布的表征能力和對(duì)參數(shù)間相關(guān)性及相關(guān)結(jié)構(gòu)的刻畫能力,各巖石力學(xué)參數(shù)的聯(lián)合概率分布擬合結(jié)果見圖5。
圖5 各巖石力學(xué)參數(shù)的邊緣概率分布擬合結(jié)果
圖5展示了各變量間二維聯(lián)合概率密度函數(shù)對(duì)數(shù)值的等值線圖,圖中用黑色空心圓形標(biāo)記代表實(shí)測(cè)巖石物理力學(xué)參數(shù)數(shù)據(jù);黑色實(shí)心正方形標(biāo)記代表了識(shí)別出的各高斯分量的均值,也可認(rèn)為是各高斯分量的中心;黑色虛線代表了各參數(shù)的聯(lián)合概率密度對(duì)數(shù)值等值線,為增強(qiáng)結(jié)果的可視性,填充顏色從淺到深的演變對(duì)應(yīng)概率密度函數(shù)對(duì)數(shù)值的從大到小變化。由各子圖可知,觀測(cè)樣本點(diǎn)集中分布在顏色較淺的高概率密度區(qū)域,隨著填充顏色從淺到深,觀測(cè)樣本點(diǎn)的分布密度也逐漸下降。此現(xiàn)象說(shuō)明學(xué)習(xí)出的高斯混合模型能夠有效刻畫參數(shù)空間中的概率密度變化情況。各子圖中的概率密度函數(shù)等值線并不是規(guī)則的橢圓形,而是呈現(xiàn)出隨著數(shù)據(jù)點(diǎn)密度變化的不規(guī)則圖形。這說(shuō)明GMM模型表征的相關(guān)結(jié)構(gòu)并不是高斯型或者其他傳統(tǒng)類型的相關(guān)結(jié)構(gòu),而是由數(shù)據(jù)的特征決定的“數(shù)據(jù)驅(qū)動(dòng)”相關(guān)結(jié)構(gòu),可以有效表述巖石物理力學(xué)參數(shù)分布特征。
GMM4的參數(shù)的最可能值如表2所示。
表2 GMM4參數(shù)最可能值
本文將基于貝葉斯機(jī)器學(xué)習(xí)框架的高斯混合模型應(yīng)用于大渡河上游巖石物理力學(xué)參數(shù)概率分布構(gòu)建任務(wù)中,在考慮統(tǒng)計(jì)不確定條件下精準(zhǔn)刻畫聯(lián)合分布特征及各參數(shù)間相關(guān)性,構(gòu)建了區(qū)域性多元巖石物理力學(xué)參數(shù)概率分布,得出以下結(jié)論:
(1)大渡河上游巖體干密度ρd、飽和吸水率w及飽和抗壓強(qiáng)度Rw具有較大變異性,且分布具有多峰、多模態(tài)特征,具有復(fù)雜相關(guān)結(jié)構(gòu)。
(2)所提方法和模型打破了現(xiàn)有方法必須事先假設(shè)各參數(shù)概率分布類型及相關(guān)結(jié)構(gòu)類型的假設(shè),能夠在考慮統(tǒng)計(jì)不確定條件下精準(zhǔn)刻畫聯(lián)合分布特征及各參數(shù)間相關(guān)性,有效表述了大渡河上游巖石物理力學(xué)參數(shù)分布特征。
(3)本文給出了所選參數(shù)聯(lián)合概率分布的GMM模型,明確了模型參數(shù),可直接用于后續(xù)工程設(shè)計(jì)工作,為后續(xù)工程提供了參考。