對(duì)數(shù)正態(tài)分布貝葉斯更新方法比較研究

楊振中

摘 要設(shè)備可靠性數(shù)據(jù)是PSA定量化分析的基礎(chǔ)，對(duì)于可靠性數(shù)據(jù)，一般分為通用數(shù)據(jù)和電廠的特定數(shù)據(jù)。由于兩種數(shù)據(jù)在應(yīng)用中均存在一定的缺陷，現(xiàn)在通用的做法是將通用數(shù)據(jù)作為先驗(yàn)數(shù)據(jù)，然后結(jié)合電廠的特定數(shù)據(jù)貝葉斯更新得到后驗(yàn)分布在PSA模型中使用。在操作實(shí)踐中，發(fā)現(xiàn)對(duì)數(shù)正態(tài)分布的通用數(shù)據(jù)，在進(jìn)行后期的貝葉斯擬合時(shí)存在一定的困難。在工程實(shí)踐上，一般將其轉(zhuǎn)化為Gamma或Beta分布再進(jìn)行處理，并形成了幾種成熟的解決方案。本文對(duì)對(duì)數(shù)正態(tài)分布轉(zhuǎn)化為Gamma分布的幾種方法做了比較研究，并簡(jiǎn)述其差異和特性。

關(guān)鍵詞貝葉斯方法;PSA;Gamma分布;Beta分布;對(duì)數(shù)正態(tài)分布

0 前言

隨著PSA技術(shù)在核電廠廣泛的推廣，PSA的應(yīng)用也逐漸地深入，這就對(duì)作為PSA技術(shù)基礎(chǔ)的可靠性數(shù)據(jù)工作提出了更高的要求?，F(xiàn)階段，國(guó)內(nèi)各個(gè)電廠一級(jí)PSA模型基本已開(kāi)發(fā)完畢，隨著電廠運(yùn)行時(shí)間的積累，勢(shì)必伴隨著電廠的變更改造和設(shè)備可靠性數(shù)據(jù)的更新。電廠的變更改造是對(duì)模型結(jié)構(gòu)的修改，而設(shè)備可靠性數(shù)據(jù)的更新即是對(duì)模型中數(shù)據(jù)的修改，通過(guò)這一過(guò)程，保證模型反映電廠實(shí)際，這也就是所謂的living PSA。PSA技術(shù)發(fā)展至今，針對(duì)數(shù)據(jù)的處理方法基本上有三種思路。第一種是采用國(guó)際上已發(fā)布的同型電廠通用數(shù)據(jù)，現(xiàn)在主流的就是美系的NUREG-6928、法系的EPS-900數(shù)據(jù)源，此外我國(guó)國(guó)家核安全局在國(guó)內(nèi)機(jī)組運(yùn)行數(shù)據(jù)統(tǒng)計(jì)的基礎(chǔ)上，于2016年發(fā)布了我們自己的通用數(shù)據(jù);第二種是使用核電廠自己的運(yùn)行數(shù)據(jù)做經(jīng)典估計(jì)。對(duì)于第一種方法，同型電廠與自身電廠不可避免的存在差異，所以直接使用通用數(shù)據(jù)并不能完全反映特定電廠的所有特征。對(duì)于第二種方法，由于在最初階段電廠的運(yùn)行時(shí)間并不長(zhǎng)，并且核電廠設(shè)備的可靠性一般較高，所以收集到的可靠性數(shù)據(jù)較少，往往是沒(méi)有數(shù)據(jù)，使用經(jīng)典估計(jì)存在著困難。最好的方法就是將第一種與第二種方法結(jié)合起來(lái)，這也就是在可靠性數(shù)據(jù)處理中最常用的貝葉斯方法。在這種方法中，將通用數(shù)據(jù)作為先驗(yàn)數(shù)據(jù)，將電廠收集到數(shù)據(jù)作為特定數(shù)據(jù)，然后應(yīng)用貝葉斯公式進(jìn)行擬合，使得到的后驗(yàn)數(shù)據(jù)同時(shí)包括了兩方面的特征。

使用貝葉斯方法時(shí)需將數(shù)據(jù)的先驗(yàn)分布通過(guò)貝葉斯定理轉(zhuǎn)化為后驗(yàn)分布。在貝葉斯處理中，存在一種共軛分布的現(xiàn)象，就是對(duì)于既定的樣本數(shù)據(jù)的分布，存在一種先驗(yàn)分布函數(shù)，使得后驗(yàn)分布函數(shù)和先驗(yàn)分布函數(shù)同屬于一個(gè)分布族。對(duì)于運(yùn)行失效來(lái)說(shuō)，這個(gè)分布是Gamma分布，對(duì)于需求失效而言，為Beta分布。對(duì)于Gamma和Beta分布，其通過(guò)貝葉斯轉(zhuǎn)化之后仍為Gamma和Beta分布，可以說(shuō)處理起來(lái)是非常方便的。但很多的數(shù)據(jù)源，如法國(guó)EDF、EPRI的數(shù)據(jù)，假設(shè)失效率的先驗(yàn)分布是對(duì)數(shù)正態(tài)分布，對(duì)于這樣的數(shù)據(jù)，在貝葉斯轉(zhuǎn)化時(shí)存在一定困難。

本文在總結(jié)貝葉斯處理方法和原理的基礎(chǔ)上，對(duì)對(duì)數(shù)正態(tài)分布的先驗(yàn)數(shù)據(jù)貝葉斯擬合的幾種主要方法做比較研究，并分析其優(yōu)劣。

1 可靠性數(shù)據(jù)的三種典型分布

由定義可以看出，部件運(yùn)行失效的概率服從泊松分布，λ是泊松分布的參數(shù)。而泊松分布的共軛分布是Gamma分布，所以對(duì)于運(yùn)行失效率，常以Gamma分布來(lái)表示。

由定義可以看出，需求失效概率p服從二項(xiàng)分布，而二項(xiàng)分布的共軛分布是Beta分布，所以需求失效常以Beta分布來(lái)表示。

此外，在可靠性和維修性方面，隨機(jī)變量的對(duì)數(shù)可能符合正態(tài)分布，對(duì)此情況在可靠性數(shù)據(jù)的統(tǒng)計(jì)分析中也常應(yīng)用對(duì)數(shù)正態(tài)分布。

1.1 Gamma分布

Gamma分布因其概率密度函數(shù)很像Gamma函數(shù)定義中積分號(hào)內(nèi)的被積部分而得名。對(duì)于X服從形狀參數(shù)為α，尺度參數(shù)為β的Gamma分布，我們記為：

因其在數(shù)學(xué)上處理方便，且有兩個(gè)參數(shù)α，β的調(diào)節(jié)，適應(yīng)范圍很廣，所以工程上常常將運(yùn)行失效率的分布假設(shè)為Gamma分布，如美國(guó)所做的柴油發(fā)電機(jī)可靠性分析、IAEA-TE CDOC -478中收集的瑞典可靠性數(shù)據(jù)[1]。在意義上，Gamma分布可看作是泊松分布參數(shù)（單位時(shí)間內(nèi)隨機(jī)事件發(fā)生率）的概率分布，且泊松分布的共軛分布即為Gamma分布。

1.2 Beta分布

1.3 對(duì)數(shù)正態(tài)分布

對(duì)數(shù)正態(tài)分布有如下特點(diǎn)。

對(duì)數(shù)正態(tài)分布是自變量取對(duì)數(shù)時(shí)，其故障密度函數(shù)符合正態(tài)分布的一種偏態(tài)性概率分布。它的故障率基本屬于遞增型的，但遞增的速度是變化的，先快后慢然后趨于平緩。對(duì)數(shù)變換可將較大的數(shù)縮小為較小的數(shù)，且越大的數(shù)縮小的越多，這一特性可以使較為分散的數(shù)據(jù)通過(guò)對(duì)數(shù)變換相對(duì)的集中起來(lái)，所以常把跨幾個(gè)數(shù)量級(jí)的數(shù)據(jù)用對(duì)數(shù)分布去擬合。在機(jī)械零件及材料的疲勞壽命中，對(duì)數(shù)正態(tài)分布應(yīng)用的較多。在核電廠設(shè)備可靠性參數(shù)、共因參數(shù)方面，對(duì)數(shù)正態(tài)分布同樣有著大量的應(yīng)用。

2 貝葉斯方法

Bayes方法在20個(gè)世紀(jì)50年代提出，用于解決復(fù)雜系統(tǒng)的可靠性評(píng)估問(wèn)題，在處理小樣本數(shù)據(jù)時(shí)具有明顯優(yōu)勢(shì)，在武器和航天等領(lǐng)域取得很好的評(píng)估效果和經(jīng)濟(jì)效益。經(jīng)過(guò)半個(gè)多世紀(jì)的發(fā)展，Bayes小樣本理論已經(jīng)形成了比較成熟的研究體系。

Bayes方法同經(jīng)典的統(tǒng)計(jì)理論的區(qū)別是所求參數(shù)的本身被看成隨機(jī)變量，如單元的壽命，求得的后驗(yàn)分布為參數(shù)的概率分布區(qū)間。Bayes方法利用了先驗(yàn)分布，不需要很大的樣本也可以得到較好的概率估計(jì)值。相比經(jīng)典方法而言，Bayes方法的優(yōu)勢(shì)在于能夠充分利用各類先驗(yàn)信息，包括可靠性試驗(yàn)數(shù)據(jù)、歷史數(shù)據(jù)、專家信息和仿真試驗(yàn)信息等，降低了經(jīng)典方法對(duì)現(xiàn)場(chǎng)試驗(yàn)樣本容量的依賴程度，使得在相同評(píng)估精度要求條件下，現(xiàn)場(chǎng)樣本容量可以相對(duì)減少。需要強(qiáng)調(diào)的是，Bayes方法不是少用信息，而是充分運(yùn)用產(chǎn)品試驗(yàn)過(guò)程中的各類信息（可靠性試驗(yàn)數(shù)據(jù)、歷史數(shù)據(jù)、專家信息和仿真信息等），因而在實(shí)際工程可以得到較好的應(yīng)用[2]。

連續(xù)分布的貝葉斯公式是可靠性數(shù)據(jù)分析中經(jīng)常用到的，因?yàn)閷?duì)于部件的運(yùn)行失效率λ和需求失效率γ（這里統(tǒng)一用θ表示）而言，在其可能的區(qū)間里，其分布顯然是連續(xù)的。設(shè){x1，…，xn}代表一個(gè)樣本，它服從一個(gè)具有似然函數(shù)L、連續(xù)參數(shù)θ和先驗(yàn)密度g（θ）的一個(gè)分布。那么，θ的后驗(yàn)分布密度為：

對(duì)于對(duì)數(shù)正態(tài)分布作為先驗(yàn)分布的情況，如果直接進(jìn)行貝葉斯估計(jì)，需要進(jìn)行煩瑣的數(shù)值積分，且得到的后驗(yàn)分布不一定是正態(tài)分布，這樣會(huì)對(duì)今后進(jìn)行PSA定量化的不確定性分析帶來(lái)困難[2]。在工程上，對(duì)于對(duì)數(shù)正態(tài)分布的情況，一般的解決方法是將其轉(zhuǎn)化為Gamma分布和Beta分布，由于在設(shè)備可靠性較高的情況下，Beta分布可近似看作Gamma分布處理，所以本文著重討論對(duì)數(shù)正態(tài)分布轉(zhuǎn)化為Gamma分布的情形。

3 對(duì)數(shù)正態(tài)分布為先驗(yàn)分布的貝葉斯處理

工程上，如果先驗(yàn)分布是對(duì)數(shù)正態(tài)分布，一般將其轉(zhuǎn)化為Gamma分布，然后進(jìn)行貝葉斯處理。這一轉(zhuǎn)化過(guò)程，采用的是一種數(shù)學(xué)上的近似，在實(shí)踐領(lǐng)域也有了多種成熟的方法。

3.1 轉(zhuǎn)化方法

1）EF值轉(zhuǎn)化法：在這種方法中，令轉(zhuǎn)化前后的mean值和EF值相等。

由于先驗(yàn)分布的mean值和EF值已知，所以問(wèn)題轉(zhuǎn)化為已知Gamma分布的mean值和EF值，求取α和β。

Gamma分布均值mean，誤差因子EF，由于mean已知，只需通過(guò)EF求取出α，即可確定β。這里我們可以選用查表的方法，也可以選用數(shù)值計(jì)算對(duì)α的數(shù)值進(jìn)行搜索。

清華大學(xué)就EF值轉(zhuǎn)化法進(jìn)行過(guò)成功應(yīng)用，并且已經(jīng)在大亞灣核電廠的數(shù)據(jù)處理中進(jìn)行了實(shí)踐[3]。

2）95%分位點(diǎn)轉(zhuǎn)化法：在轉(zhuǎn)化過(guò)程中，使這兩種分布在兩個(gè)關(guān)鍵點(diǎn)（均值和95%分位點(diǎn)）上重合，找到滿足上述條件的Gamma分布參數(shù)α，β，即完成了將先驗(yàn)分布的對(duì)數(shù)正態(tài)分布轉(zhuǎn)換為Gamma分布的工作。

在具體轉(zhuǎn)化時(shí)，同樣已知均值mean和EF，利用對(duì)數(shù)正態(tài)分布的函數(shù)關(guān)系，求出其x0.95，利用參考文獻(xiàn)[4]中與α的關(guān)系曲線，可以求得α值。再根據(jù)已知的mean關(guān)系，即可得出Gamma分布。同方法1），也可以利用數(shù)值計(jì)算的方法完成對(duì)α的搜索。

在國(guó)內(nèi)，機(jī)械科學(xué)研究院與秦山核電壓水堆電廠即秦山一二期與方家山電廠的數(shù)據(jù)處理采用該種方法。

3）方差轉(zhuǎn)化法：這種方法是令轉(zhuǎn)化前后的均值和方差相等。由已知LN分布的均值和EF值，可以方便地求出其方差。再利用Gamma分布均值mean，可以順利地求出α和β。通過(guò)三種方法的比較，可以看到此種方法最為簡(jiǎn)便。在實(shí)際應(yīng)用中，秦山三期重水堆采用的是該種方法。

3.2 三種轉(zhuǎn)化方法的比較

1）三種方法的準(zhǔn)確性分析

假設(shè)一均值mean=0.5，誤差因子EF=3的對(duì)數(shù)正態(tài)分布，當(dāng)累計(jì)失效次數(shù)為2，累計(jì)運(yùn)行小時(shí)數(shù)為2500時(shí)，得到如下后驗(yàn)分布圖像。

由圖1可以看出，當(dāng)失效數(shù)據(jù)較小時(shí)，三種轉(zhuǎn)化方法差異不大，與LN分布經(jīng)數(shù)值積分后的貝葉斯轉(zhuǎn)化結(jié)果也無(wú)明顯偏差。

假設(shè)對(duì)數(shù)正態(tài)分布的均值和誤差因子不變，當(dāng)累計(jì)失效次數(shù)為6，累計(jì)運(yùn)行小時(shí)數(shù)為4500時(shí)，得到如下圖像。

由圖2可知，當(dāng)失效數(shù)據(jù)較大時(shí)，EF值轉(zhuǎn)化法與LN經(jīng)數(shù)值積分后的結(jié)果最為接近，95%分位點(diǎn)法次之，總體上這兩種方法差異不大，而方差轉(zhuǎn)化法的誤差偏大。

2）LN分布EF值與Gamma分布α值相關(guān)性比較

設(shè)先驗(yàn)數(shù)據(jù)的分布類型為對(duì)數(shù)正態(tài)分布，其先驗(yàn)均值為0.5，利用EF值轉(zhuǎn)化法，可得EF關(guān)于α的關(guān)系曲線如圖3所示。

由圖3可以觀察到，當(dāng)EF值較小時(shí)，α值是關(guān)于EF值的減函數(shù)，隨著EF值的增大，α值的變化有一個(gè)相對(duì)平緩的區(qū)間，此時(shí)的α相對(duì)于EF的變化敏感性不強(qiáng)。隨著EF值繼續(xù)增大，α值有一個(gè)迅速增大的過(guò)程，從形態(tài)上來(lái)看，呈現(xiàn)出指數(shù)增長(zhǎng)的特征。

4 結(jié)論

在針對(duì)核電廠小樣本失效概率分析的情況下，貝葉斯方法有著其明顯的優(yōu)越性，所以在核電廠PSA分析中有著廣泛的應(yīng)用。但是當(dāng)先驗(yàn)分布為對(duì)數(shù)正態(tài)分布的情況，在求取后驗(yàn)分布時(shí)又面臨著數(shù)學(xué)上的困難。為此工程上形成了三種行之有效的解決方案。經(jīng)過(guò)對(duì)三種方法進(jìn)行對(duì)比，發(fā)現(xiàn)EF值轉(zhuǎn)化法較為復(fù)雜，但所得結(jié)果準(zhǔn)確性較高;95%分位點(diǎn)法次之;方差轉(zhuǎn)化法最為簡(jiǎn)單，但形成的偏差也最大。此外，研究發(fā)現(xiàn)，LN分布的EF值與轉(zhuǎn)化后Gamma分布的α值并不呈單調(diào)關(guān)系，前部隨著EF值增大會(huì)有一個(gè)遞減的過(guò)程，且當(dāng)EF值增大超過(guò)一定數(shù)值后，α將會(huì)呈現(xiàn)一個(gè)指數(shù)增大的過(guò)程。

參考文獻(xiàn)

[1]馬靜嫻，閆國(guó)奎，劉志軍.用貝葉斯方法處理核電站PSA分析中的設(shè)備可靠性數(shù)據(jù)[J].2004年全國(guó)機(jī)械可靠性學(xué)術(shù)交流會(huì)論文集.

[2]劉方亮.核電廠小樣本數(shù)據(jù)Bayes處理方法應(yīng)用研究[D].清華大學(xué)碩士學(xué)位論文，2010年.

[3]茆定遠(yuǎn)，薛大知.核電站PSA分析中可靠性數(shù)據(jù)處理的貝葉斯方法，核動(dòng)力工程Vol.21.No.5 Oct.2000.

[4]Some Properties of Distributions Useful in the Study of Rare Events，IEEE Transaction on Reliability Vol.R-31，No.1，APRIL 1982.