徐廣波, 胡甚平
(1. 上海海事大學 商船學院,上海 201306; 2. 江陰海事局,江蘇 江陰 214431)
隨著經濟建設發展和科技水平的不斷進步,航運業得到快速發展,船舶逐漸向大型化、高速化、綜合化、智能化方向邁進.大量頻繁活躍的船舶使港口、航道等水域交通流量急增,航行環境日益復雜.[1]因此,有必要對港口水域的交通風險進行分析和評估,為港航安全管理部門制定有針對性的預防措施提供依據.目前,國內外已有很多學者采用定性、定量以及定性定量相結合的方法從不同角度研究水上交通風險問題,如綜合安全評估(FSA)法[2]、灰色系統理論法[3]、灰色馬爾可夫鏈方法[4]等.這些方法在一定程度上提供定性和定量的決策依據,然而,對動態、線性的數據進行實時評估,需要更多的數據推理及模擬實驗分析發展趨勢.
本文在港口交通風險定量化評估的基礎上給出交通事故率和事故后果的貝葉斯概率統計,接著建立港口交通系統風險的蒙特卡羅仿真模型,并通過馬爾可夫鏈蒙特卡羅(Markov Chain Monte Carlo,MCMC)方法對風險模型的參數進行推斷與優化,最后進行仿真實驗.
港口交通系統風險[5]是指在某一特定的客觀狀態下港口交通系統中人、船(貨)、環境受到傷害的可能性和這種傷害的嚴重程度,可表示為
R|S=f(P,C)|S
(1)
式中:P為某一時間事故發生的可能性;C為事件發生的后果;S代表某一特定的客觀狀態;R|S為在S狀態下分析對象的風險度;f為關于P和C的實數函數.可能性指事故發生的機會,用于描述概率或頻率的性質.概率是理論值,由事件的本質決定,只能取唯一值,它能精確反映事件發生可能性的大小.稱隨機事件在一定時間內統計取得的發生次數為頻率. 頻率是試驗值或使用時的統計值,具有隨機性,可能取多個數值,因此只能近似反映事件發生可能性的大小.可能性
(2)
式中:i∈[1,n],i∈Z;P的物理意義是指在選取的樣本中發生事故的頻率;ai為年事故數;ai為每年船舶活動量.
后果是描述有害事件或非正常事件發生所造成損害的程度.在實際事故后果分析中,由于考慮不同的風險,常采用“事故等效后果”衡量.“事故等效后果”是經歸一化以后事故的各種后果總和.
(3)
MCMC方法是最近發展起來的一種簡單且行之有效的貝葉斯計算方法.其基本思想是通過建立一個平穩分布為π(z)的馬爾可夫鏈得到π(z)的樣本,基于這些樣本作各種統計推斷.其主要目的是借助貝葉斯概率估計,通過頻率數據獲取風險事件的概率.
將MCMC算法概括為3步:(1)在z上選一個“合適”的馬爾可夫鏈,使其轉移核為p(·,·),這里“合適”的含義主要指π(z)應是其相應的平穩分布;(2)由z中某一點Z(0)出發,用(1)中的馬爾可夫鏈產生點序列Z1,…,Zn;(3)對某個m和n,任一函數f(z)的期望估計
(4)
采用文獻[6]中所建立的模型作為風險仿真模型,即把要解決的目標問題抽象成一個概率模型R=f(P,C)(式中R表示風險度;P表示某一事件發生的概率;C表示事件發生的后果).分別對相關資料進行統計,獲取P和C的樣本信息;通過結合樣本信息和模擬風險值確定概率模型;最后選定適當的模擬值個數N與次數M,對所獲得的M·N個樣本值進行統計分析,得到分布曲線和宏觀風險的特征.
主要步驟[7]如下:(1)收集、分析主客觀先驗信息,確定合適的先驗分布形式以及先驗參數.(2)結合試驗數據確定第i個模型,利用Gibbs抽樣對模型的后驗進行MCMC模擬.(3)判斷馬爾可夫鏈是否已收斂、MC誤差是否足夠小.如果馬爾可夫鏈已收斂、MC誤差足夠小則轉入下一步,否則需進一步調試模型.重新確定抽樣迭代次數及抽樣方法等,若效果依然不盡如人意,則返回(2),重新考慮修改先驗參數和模型.(4)修改(3)中的模型、選擇更高一級的i+1個模型并返回(3),比較模型、選擇相對更優秀的后驗模型,進行模型的貝葉斯推斷,并根據有關準則得出正確的結論.
基于MCMC方法的模型運行流程見圖1.
構建風險的貝葉斯多層對數正態模型:
R[i]~d lnorm(μ,σ)//風險服從對數正態分布
R[i]<-P[i]*C[i]//風險由概率和后果確定
P[i]~d lnorm(μ,σ) //可能性服從對數正態分布
C[i]~d lnorm(μ,σ) //后果服從對數正態分布
μ~d norm(0,1.0E-6) //超參數服從正態分布
σ~d gamma(0.01,0.001) //超參數服從伽馬分布
sigma<-1/sqrt(σ) //符號函數

圖1 基于MCMC方法的模型運行流程
利用WinBUGS軟件對各參數進行迭代計算,從而獲得其后驗分布和數學特征值.
仿真數據檢驗方法主要采用文獻[8]中的變異因數,即通過判斷估計值與待定值的收斂性檢驗仿真數據的可靠性,
V=σ/m
(5)
式中:σ為仿真樣本的風險標準差;m為仿真樣本的風險平均值.
用我國某港口的交通事故數據[2]進行應用測試.由于事故數據樣本量較少,屬于小樣本事件,為獲取港口交通系統風險事件發生的可能性往往只能用頻率進行估計,通過蒙特卡羅方法獲取概率.后果采用因數進行換算.文獻[8]用基于對數正態的蒙特卡羅模型進行風險推理.本文則通過MCMC算法先修正概率,進而對風險事件數據進行仿真.
WinBUGS是英國劍橋公共衛生研究所推出的用MCMC方法來分析復雜統計模型的貝葉斯推理軟件.其基本原理就是通過Gibbs sampling和Metropolis算法,從完全條件概率分布中抽樣,從而生成馬爾可夫鏈,通過迭代最終估計出優化后的概率參數[9].
在對風險事件樣本數據整理之后,得到先驗分布數字特征,然后利用WinBUGS軟件進行MCMC參數計算[10],獲得所建模型各參數的后驗分布數字特征,見表1.

表1 各參數的后驗分布數字特征
利用表1中的參數進行蒙特卡洛仿真,結果見圖2.圖中實心點為樣本數據,圓圈表示仿真獲得的隨機數據.
從圖2(b)中可以看出,樣本量得到有效增加,新增加樣本對原樣本有較好的覆蓋.從圖2(c)和2(d)中可見高風險出現的概率明顯較小.通過擴充樣本信息量,在風險合成后可得到整體風險.
當模擬運行1 000次后(每次抽取樣本100個),仿真結果相差極小.經統計可得出該港口水域頻率均值為1.279 3×10-4,后果均值為0.764 2,后果均值與實際情況也較為符合,風險平均值為0.975 1×10-4,風險標準差小,為2.518 2×10-4,說明仿真的風險值穩定性很好,仿真結果可信.

圖2經優化后的港口水域樣本及模擬樣本仿真風險分布
利用上面得到的模擬風險值(1 000個),可以畫出一個統計直方圖,見圖3.

(a)直方圖 (b)三維離散圖
圖3多次模擬樣本的仿真風險
按照文獻[8]的模型,對數正態下的概率參數與MCMC優化參數后的概率參數的風險統計特征值見表2.

表2 風險統計特征值
從圖3可得到以下結論:(1)該港口不同水域不同事故類型的平均風險為0.975 1×10-4,風險值標準差較小,穩定性較好.(2)風險值主要集中在0.4×10-4~1.5×10-4之間,峰值主要集中在0.95×10-4處,與風險均值接近重合.該水域的風險值是較為典型的對數正態分布的隨機變量.(3)利用統計數據,對風險進行分級.根據二八定律,可將低風險、合理可行風險、高風險等3類風險度的閾值進行重新劃分,進一步保證相對風險分級的可行性.(4)與參數優化前MC仿真的結果相比較,參數優化對仿真結果起到一定的優化作用,主要表現在所獲得的參數本身具有更強的穩定性,仿真所得的頻率均值、風險均值及后果均值與實際結果更為接近,同時風險值區間上下限值有所擴大,較高風險可能性增大,更能提醒相關職能部門和船舶駕駛人員防范交通風險.
從第3.1節風險信息(N=100,M=1 000)中隨機抽取100次仿真數據進行變異因數計算,得到圖4(其中穩定數據線為樣本的變異因數,波動數據線為仿真樣本變異因數).

圖4 基于WinBUGS軟件的仿真風險的變異因數
圖4中,雖然仿真樣本的數據曲線振蕩幅度較為顯著,但是樣本的統計結果基本落在仿真樣本的變異因數的變化范圍之內,兩種樣本的關聯性較強.
應用MC方法建立的仿真模型能夠有效擴展小樣本問題分析數據,為確立風險分級標準提供依據.利用WinBUGS軟件通過MCMC方法能夠實現對風險模型中概率參數的優化,更加接近風險事件發生的可能性.實例表明,參數優化后的風險模型具有更強的穩定性,與實際結果具有較好的一致性.如何利用MCMC方法實現風險的動態轉移是后續的研究工作.
參考文獻:
[1] 徐廣波, 軒少永, 尤慶華. 基于熵權的模糊集對模型在港口水域通航風險評價中的應用[J].上海海事大學學報, 2012, 33(1): 7-11.
[2] 方泉根, 胡甚平. FSA 在船舶引航風險評估中的應用[J]. 哈爾濱工程大學學報, 2006, 27(3): 329-334.
[3] 周麗麗, 胡甚平. 船舶引航風險成因灰色綜合評價模型[J]. 上海海事大學學報, 2008, 29(3): 21-25.
[4] 趙佳妮, 吳兆麟. 基于灰色馬爾可夫模型的水上交通事故預測[J]. 大連海事大學學報, 2005, 31(4): 15-18.
[5] 胡甚平. 不確定性信息下的海上交通風險評估方法及應用研究[D]. 上海: 上海海事大學, 2010.
[6] 高帥. 基于蒙特卡羅方法的沿海水上交通風險仿真[D]. 上海: 上海海事大學, 2011.
[7] 林靜. 基于MCMC的貝葉斯生存分析理論及其在可靠性評估中的應用[D]. 南京: 南京理工大學, 2008.
[8] 胡甚平. 海上交通系統風險蒙特卡洛仿真[J].上海海事大學學報, 2011, 32(4): 7-11.
[9] 余芳. 基于 MCMC 方法的 WinBUGS 軟件應用[J]. 經營管理者, 2010(15): 15.
[10] LUNN D J, THOMAS A, BEST N,etal. WinBUGS-a Bayesian modelling framework: concepts, structure, and extensibility[J]. Stat & Computing, 2000, 10(4): 325-337.