胡甚平，黃常海，鄧華，黃道正

(1.上海海事大學(xué) 商船學(xué)院，上海 200135; 2.江蘇海事職業(yè)技術(shù)學(xué)院 航海技術(shù)系，江蘇 南京 211170; 3.上海海事大學(xué) 交通運(yùn)輸學(xué)院， 上海 200135)

船舶引航過程風(fēng)險(xiǎn)馬爾科夫鏈模型動(dòng)態(tài)仿真

胡甚平1，黃常海1，鄧華2，黃道正3

(1.上海海事大學(xué) 商船學(xué)院，上海 200135; 2.江蘇海事職業(yè)技術(shù)學(xué)院 航海技術(shù)系，江蘇 南京 211170; 3.上海海事大學(xué) 交通運(yùn)輸學(xué)院， 上海 200135)

為定量研究特定水域船舶引航過程風(fēng)險(xiǎn)演化規(guī)律，有效進(jìn)行船舶進(jìn)出港安全預(yù)控，本文通過對(duì)港區(qū)水域船舶引航作業(yè)工作流程分析，確立不同船舶引航作業(yè)任務(wù)的風(fēng)險(xiǎn)程度，提出船舶引航狀態(tài)處于港內(nèi)外錨地、航道和泊位前沿等水域之間的轉(zhuǎn)移方程。計(jì)算結(jié)果表明：船舶引航狀態(tài)的轉(zhuǎn)移符合馬氏穩(wěn)態(tài)特性，經(jīng)潮汐型波動(dòng)后走向穩(wěn)定；船舶航行狀態(tài)轉(zhuǎn)移的過程變量對(duì)港區(qū)水域總體船舶秩序影響明顯?；隈R爾科夫鏈蒙特卡洛方法的過程風(fēng)險(xiǎn)仿真模型適用于海上交通風(fēng)險(xiǎn)的模擬分析，可為船舶安全管理提供風(fēng)險(xiǎn)演化的基本規(guī)律。

風(fēng)險(xiǎn)評(píng)估； 風(fēng)險(xiǎn)仿真； 船舶引航； 馬爾科夫鏈蒙特卡洛； 船舶流量； 港口交通； 馬爾科夫鏈； 動(dòng)態(tài)系統(tǒng)仿真

Abstract：The dynamic system simulation of effective ship embarking or disembarking at harbors is necessary to quantitatively analyze and identify risk evolution in ship pilotage in specific waters for safety control. In this study, the working process of the pilotage operation of a ship in port waters is analyzed and the risk degrees of different ship pilotage operation tasks are determined. The calculation result shows that the transfer of pilotage state complies with Markov steady- state characteristics and remains stable after tide- type fluctuations. In addition, the transferring variable remarkably affects the whole ship flow in the harbor. The dynamic system simulation model with the Markov Chain Monte Carlo algorithm is suitable for the risk analysis of maritime traffic risk and can provide the basic law for the evolution of risk in the safety management of marine traffic.

Keywords： risk assessment; risk simulation; ship pilotage; Markov Chain Monte Carlo; ship flow; harbour traffic; Markov Chain; dynamic system simulation

隨著航運(yùn)業(yè)的發(fā)展，大量的船舶頻繁活躍于港口水域，使港口水域通航密度大大增加，交通環(huán)境更加復(fù)雜，導(dǎo)致船舶水上交通風(fēng)險(xiǎn)逐漸增加[1]。目前，國內(nèi)進(jìn)出港口船舶中有50%左右的船舶是在引航員的引領(lǐng)下完成的，部分對(duì)外港口甚至達(dá)到90%，船舶引航是港口生產(chǎn)和快速發(fā)展中不可替代的環(huán)節(jié)，保證了港區(qū)的安全和水域的清潔，提升了港口綜合服務(wù)能力[2]。船舶引航在保障地方經(jīng)濟(jì)發(fā)展、港口生產(chǎn)能力和港口水域交通安全方面的作用不可替代[3]。

目前，國內(nèi)外已有很多學(xué)者采用定性、定量以及定性定量相結(jié)合的方法從不同角度研究港區(qū)水域船舶航行風(fēng)險(xiǎn)問題。方泉根等采用綜合安全評(píng)估方法對(duì)船舶引航風(fēng)險(xiǎn)形成原因進(jìn)行分析[3-4]。方誠等引入貝葉斯理論對(duì)船舶引航風(fēng)險(xiǎn)分別提出預(yù)測(cè)模型和動(dòng)態(tài)網(wǎng)絡(luò)推理模型[5]。Cameron提出引航安全決策研究[6]。馬飛等提出船舶引航中的風(fēng)險(xiǎn)脆弱性量化問題[7]。Wu等引入模糊綜合評(píng)判法對(duì)引航系統(tǒng)進(jìn)行評(píng)價(jià)[8]。G?r?ün等采用綜合安全評(píng)估方法對(duì)海峽水域船舶航行風(fēng)險(xiǎn)進(jìn)行量化研究[9]。Pak等采用模糊層次分析法對(duì)港口水域船舶安全航行風(fēng)險(xiǎn)進(jìn)行量化研究[10]。Praetorius等提出如何適應(yīng)不斷變化的操作條件的船舶交通系統(tǒng)(vessel traffic system,VTS)操作人員適應(yīng)力模型[11]。這些文獻(xiàn)在船舶引航靜態(tài)風(fēng)險(xiǎn)量化方面取得進(jìn)展，然而對(duì)于時(shí)序下的船舶引航動(dòng)態(tài)風(fēng)險(xiǎn)研究有待深入。

量化隨機(jī)過程下的船舶港區(qū)航行風(fēng)險(xiǎn)是進(jìn)行大數(shù)據(jù)條件下風(fēng)險(xiǎn)預(yù)防和控制的基礎(chǔ)問題。胡甚平等提出基于人工智能云模型的蒙特卡洛方法研究風(fēng)險(xiǎn)耦合機(jī)理[12]。Faghih Roohi引入馬爾科夫鏈蒙特卡洛方法(Markov chain Monte Carlo，MCMC)對(duì)海上交通風(fēng)險(xiǎn)進(jìn)行動(dòng)態(tài)仿真，通過船舶交通事件、事故與嚴(yán)重事故之間的風(fēng)險(xiǎn)轉(zhuǎn)移的動(dòng)態(tài)系統(tǒng)仿真得出風(fēng)險(xiǎn)系統(tǒng)動(dòng)態(tài)演化規(guī)律，分析安全風(fēng)險(xiǎn)的分布以及現(xiàn)有交通條件下的發(fā)展趨勢(shì)[13]。動(dòng)態(tài)風(fēng)險(xiǎn)納入隨機(jī)過程研究是一種新嘗試。

本文通過港區(qū)水域船舶引航作業(yè)工作流程分析，確立船舶引航狀態(tài)處于港內(nèi)外錨地、航道和泊位前沿等水域之間的狀態(tài)轉(zhuǎn)移方程。采用馬爾科夫鏈蒙特卡洛方法，提出港區(qū)水域船舶引航風(fēng)險(xiǎn)動(dòng)態(tài)系統(tǒng)仿真模型。結(jié)合港區(qū)水域船舶引航狀態(tài)轉(zhuǎn)移矩陣，進(jìn)行馬爾科夫鏈蒙特卡洛仿真實(shí)驗(yàn)與分析。

1 問題的描述

1．1 船舶引航過程風(fēng)險(xiǎn)分析

引航是指在一定的水域內(nèi)，由專門性的從業(yè)人員登上船舶，并就船舶駛抵目的地的有關(guān)問題，向船長(zhǎng)提出有關(guān)航行問題的建議，或者在不解除船長(zhǎng)對(duì)于全船駕駛責(zé)任的情況下，把船舶安全地從錨地或泊位引進(jìn)、帶出港口或在港內(nèi)移泊。因此，引航任務(wù)是由作業(yè)屬性和初始船舶狀態(tài)組成，作業(yè)屬性包括進(jìn)出港或移泊，初始船舶狀態(tài)包括內(nèi)外錨地或泊位。引航過程則是由引航任務(wù)的設(shè)定、人員與船舶內(nèi)部因素相關(guān)聯(lián)或相互作用的活動(dòng)以及預(yù)期的船舶安全營運(yùn)的狀態(tài)等組成，即船舶通過進(jìn)出港航行(含調(diào)頭)、靠離泊(含拋起錨、系離浮等)和移泊等，最終實(shí)現(xiàn)錨地與泊位之間的空間狀態(tài)轉(zhuǎn)移。船舶引航作業(yè)是船舶在不同區(qū)域間進(jìn)行狀態(tài)轉(zhuǎn)移，從而完成進(jìn)出港航行、靠離泊和移泊作業(yè)任務(wù)(見表1)。

表1 船舶引航作業(yè)流程

在港區(qū)水域，船舶引航過程需通過港外錨地、港內(nèi)錨地、港區(qū)航道和泊位(含前沿水域等附近水域)等改變而完成。鑒于船舶引航作業(yè)是船舶在狀態(tài)轉(zhuǎn)移過程中完成進(jìn)出港航行、靠離泊和移泊任務(wù)，船舶引航過程風(fēng)險(xiǎn)包含船舶引航中航行、靠泊、離泊、移泊等行為的風(fēng)險(xiǎn)[4]。

因而，船舶進(jìn)出港口(或移泊)的風(fēng)險(xiǎn)最終表現(xiàn)為不同區(qū)域不同作業(yè)風(fēng)險(xiǎn)的組合。比如，進(jìn)港船舶引航風(fēng)險(xiǎn)為外錨地起錨作業(yè)、航道航行作業(yè)和泊位靠泊作業(yè)等風(fēng)險(xiǎn)的組合(見圖1)。

圖1 船舶港區(qū)水域狀態(tài)轉(zhuǎn)移圖Fig.1 State transition diagram on ship pilotage

1.2 船舶引航過程風(fēng)險(xiǎn)動(dòng)態(tài)模型的建立

船舶引航過程動(dòng)態(tài)風(fēng)險(xiǎn)可以定義為在設(shè)定的時(shí)間或過程中船舶所處隨機(jī)狀態(tài)下安全引航作業(yè)風(fēng)險(xiǎn)的組合，數(shù)學(xué)表達(dá)為

(1)

式中：R(t)為設(shè)定時(shí)間t下的風(fēng)險(xiǎn)值，r為船舶不同引航作業(yè)狀態(tài)xi(i=1,2,3,4)非正常事件發(fā)生可能性和后果程度組合的風(fēng)險(xiǎn)值，S(t)為時(shí)間t時(shí)刻船舶暴露著的引航狀態(tài)，S(0)為船舶t=0的引航初始狀態(tài)，T為船舶引航狀態(tài)轉(zhuǎn)移變量，滿足：

(2)

式中:μij為不同狀態(tài)之間變化xi→xj的轉(zhuǎn)移概率。

2 船舶引航過程風(fēng)險(xiǎn)的動(dòng)態(tài)仿真

2.1 馬爾科夫鏈蒙特卡洛方法

船舶引航安全系統(tǒng)仿真是對(duì)船舶引航過程的一種模擬，利用數(shù)學(xué)模型來模擬船舶引航安全系統(tǒng)發(fā)展變化的規(guī)律。MCMC作為統(tǒng)計(jì)計(jì)算方法，也是動(dòng)態(tài)系統(tǒng)仿真的方法之一，通過設(shè)定隨機(jī)過程，重復(fù)生成時(shí)間序列，估計(jì)分布參數(shù)和統(tǒng)計(jì)數(shù)值，進(jìn)而研究其分布特征[14]。首先通過計(jì)算機(jī)進(jìn)行靜態(tài)模擬，然后將離散時(shí)間的馬爾科夫鏈或連續(xù)時(shí)間的馬爾科夫過程引入到蒙特卡洛模擬中，實(shí)現(xiàn)抽樣分布隨模擬條件進(jìn)而改變的動(dòng)態(tài)系統(tǒng)模擬。由于涉及時(shí)間序列反復(fù)生成，彌補(bǔ)了傳統(tǒng)的蒙特卡洛積分只能進(jìn)行靜態(tài)模擬的缺陷。該方法在機(jī)器學(xué)習(xí)、人工智能等領(lǐng)域具有廣泛的應(yīng)用[15]。

MCMC用一個(gè)復(fù)雜分布進(jìn)行采樣，馬爾科夫鏈存在唯一的平穩(wěn)分布，使得樣本點(diǎn)服從該復(fù)雜分布。有限狀態(tài)是不可約非周期的，當(dāng)轉(zhuǎn)移步數(shù)趨向于無窮時(shí)，每個(gè)狀態(tài)服從指定的概率分布。因此，從概率分布中進(jìn)行采樣，可以構(gòu)造一個(gè)馬爾科夫鏈過程，使得它的平穩(wěn)分布是復(fù)雜分布[16]。

2.2 船舶引航過程風(fēng)險(xiǎn)動(dòng)態(tài)系統(tǒng)仿真

馬爾科夫鏈體現(xiàn)的是離散狀態(tài)在離散時(shí)間中或連續(xù)時(shí)間轉(zhuǎn)換關(guān)系，下一個(gè)狀態(tài)只決定與當(dāng)前的狀態(tài)。就船舶引航任務(wù)而言，每一個(gè)引航任務(wù)只描述當(dāng)前狀態(tài)與下一個(gè)狀態(tài)，下一個(gè)狀態(tài)與當(dāng)前狀態(tài)相關(guān)。

狀態(tài)圖的轉(zhuǎn)換關(guān)系可以用轉(zhuǎn)換矩陣T來表示。考慮一般的情況，滿足條件下經(jīng)過一定的馬爾科夫鏈迭代后系統(tǒng)分布會(huì)趨近一個(gè)穩(wěn)定分布，即最后的μ(t)服從目標(biāo)分布S(x)(t)采樣。在采用 MCMC方法時(shí)，馬爾科夫鏈轉(zhuǎn)移核的構(gòu)造至關(guān)重要，不同的轉(zhuǎn)移核構(gòu)造方法，將產(chǎn)生不同的 MCMC 方法和目標(biāo)狀態(tài)轉(zhuǎn)移的結(jié)果。目前常用的 MCMC 中抽樣方法有三種，即MH算法(metropolis- hasting)、切片抽樣(slice sampling)和吉布斯抽樣(gibbs)，這三種方法各有優(yōu)缺點(diǎn)[11]：MCMC最一般的抽樣方法是MH方法。切片抽樣是一種通過一個(gè)定義好的不變分布來構(gòu)造可逆馬爾科夫轉(zhuǎn)移核的方法。Gibbs抽樣是一種特殊的MH抽樣，其隨機(jī)數(shù)值總是被接收(即接收概率為1)，見圖2。

MH方法是：取xi、xj兩狀態(tài)之間最小轉(zhuǎn)移概率作為二者之間轉(zhuǎn)入和轉(zhuǎn)出的概率，那么，轉(zhuǎn)出去概率大的狀態(tài)增加了自我轉(zhuǎn)移的概率。如圖3所示，假設(shè)狀態(tài)xi轉(zhuǎn)出到狀態(tài)xj的概率大。該算法能構(gòu)造最快收斂，通過調(diào)整每個(gè)狀態(tài)的自我轉(zhuǎn)移概率，使它最接近穩(wěn)態(tài)下的自我轉(zhuǎn)移概率，多余的概率都轉(zhuǎn)移出去，保證細(xì)致平穩(wěn)條件。對(duì)于引航作業(yè)而言，使用MH算法可以快速收斂到目標(biāo)狀態(tài)，減少中間狀態(tài)的過渡時(shí)間[14]。因此，使用該算法可以獲得引航初始時(shí)和結(jié)束時(shí)的船舶狀態(tài)。

圖2 船舶引航的狀態(tài)轉(zhuǎn)移變量Fig.2 Ship Pilotage state transition variable

圖3 MH算法的狀態(tài)概率轉(zhuǎn)移圖Fig.3 State probability transition diagram on MH Algorithms

馬爾科夫鏈的收斂性質(zhì)主要由轉(zhuǎn)移矩陣T決定，基于馬爾科夫鏈做采樣的關(guān)鍵問題是構(gòu)造轉(zhuǎn)移矩陣，使得平穩(wěn)分布恰好是需要的分布S(x)。依照?qǐng)D3的算法，可以根據(jù)轉(zhuǎn)移矩陣T?和 ?T來得到S(x)(t)和S(x)(t-1)的比率，進(jìn)而按照一定的概率對(duì)兩個(gè)樣本進(jìn)行選擇。通過反復(fù)大量類似處理，最終得到樣本符合原始的S(x)分布。

從一個(gè)高概率狀態(tài)x(t)向一個(gè)低概率狀態(tài)x(t-1)轉(zhuǎn)移的概率與從這個(gè)低概率狀態(tài)向高概率狀態(tài)轉(zhuǎn)移的概率相同。本文假設(shè)S(x)是一維標(biāo)準(zhǔn)高斯分布，xi是根據(jù)蒙特卡洛方法得到的一個(gè)樣本，對(duì)此數(shù)據(jù)進(jìn)行灰信息處理，使轉(zhuǎn)移矩陣中行概率和為1。

本文選用 MH 算法對(duì)標(biāo)準(zhǔn)高斯分布進(jìn)行采樣，轉(zhuǎn)移函數(shù)(非對(duì)稱)是船舶進(jìn)出港作業(yè)參數(shù)的高斯矩陣[16]。算法過程如下：

1)選取一個(gè)隨機(jī)點(diǎn)S(0)，作為一個(gè)采樣點(diǎn)。

2)選取隨機(jī)點(diǎn)，確定轉(zhuǎn)移函數(shù)為T，運(yùn)用MH算法確定新的矩陣，進(jìn)行灰性白化。

3)獲取S(1)=S(0)T，然后對(duì)S(1)進(jìn)行灰信息白化處理，進(jìn)行S(1)中數(shù)據(jù)的歸一化計(jì)算，獲得新的S(1)。

4)求取R(1)=r·S(1)。

5)重復(fù)第1)～4)，獲取S(n)=S(0)Tn，直至穩(wěn)態(tài)R(n)=r·S(n)。

3 基于MCMC動(dòng)態(tài)模型的船舶引航過程風(fēng)險(xiǎn)仿真

選用我國東部沿海某港區(qū)水域的船舶進(jìn)出港數(shù)據(jù)來進(jìn)行船舶引航風(fēng)險(xiǎn)仿真。選取的數(shù)據(jù)是2010-2015年間72個(gè)月的171 432艘船舶的引航作業(yè)。船舶種類包括集裝箱船(60.07%)，化工品船舶(6.31%)，干散貨船(30.55%)，其他(3.07%)。

3.1 樣本的采集

表2為72個(gè)月中船舶引航任務(wù)按照進(jìn)港、出港、移泊三種引航任務(wù)和外錨地、泊位、內(nèi)錨地三個(gè)船舶初始狀態(tài)統(tǒng)計(jì)基礎(chǔ)數(shù)據(jù)。

表2 港口船舶引航任務(wù)2010-2015年全年統(tǒng)計(jì)量

按每月的日均引航統(tǒng)計(jì)，得到31個(gè)不同日均引航進(jìn)港、出港和移泊作業(yè)量分布。日均船舶量在44與117艘次之間，日均79.6艘次。然后，通過灰信息處理，獲得歸一化(一天的船舶量為100%)的不同日期的平均引航任務(wù)組成比率，如圖4所示。

圖4 某港區(qū)每日船舶引航不同作業(yè)屬性的日均比例Fig.4 Ship daily activities at a certain harbour

3.2 動(dòng)態(tài)仿真初始條件的設(shè)定

1)日均引航數(shù)據(jù)顯示，進(jìn)港、出港、移泊三種引航任務(wù)還是波動(dòng)性強(qiáng)，具有隨機(jī)性和不穩(wěn)定性。假設(shè)不同任務(wù)下的比率符合高斯分布，那么確定t=0時(shí)船舶引航所處的初始狀態(tài)為

(3)

2)確定船舶引航作業(yè)的風(fēng)險(xiǎn)值。依據(jù)文獻(xiàn)[15]，對(duì)港口水域船舶航行、靠泊、離泊和移泊作業(yè)時(shí)發(fā)生的交通事故和作業(yè)量進(jìn)行分別計(jì)算。經(jīng)風(fēng)險(xiǎn)測(cè)量，再通過蒙特卡洛方法的仿真[15]，提出航行、靠泊、離泊和移泊的相對(duì)風(fēng)險(xiǎn)比值為6.9∶5.3∶5.1∶1.8。假設(shè)樣本符合威布爾分布，可以計(jì)算出各種風(fēng)險(xiǎn)的分布參數(shù)。

按照表1的作業(yè)流程，可推理出內(nèi)錨地、外錨地、泊位及前沿、港區(qū)航道四種狀態(tài)下的風(fēng)險(xiǎn)值，然后通過歸一化(相對(duì)風(fēng)險(xiǎn)量為100%)運(yùn)算，得出風(fēng)險(xiǎn)狀態(tài)矩陣r：

(4)

3)確定轉(zhuǎn)移矩陣。在數(shù)據(jù)中統(tǒng)計(jì)同一水域的引航作業(yè)任務(wù)的比例，從而確定轉(zhuǎn)移矩陣的相應(yīng)變量值。這里先選定一個(gè)進(jìn)出港平衡(對(duì)稱)的轉(zhuǎn)移矩陣T，數(shù)據(jù)已經(jīng)完成了灰信息白化處理，獲得下一個(gè)狀態(tài)轉(zhuǎn)移概率和為1。

在船舶進(jìn)出港口中，狀態(tài)轉(zhuǎn)移是波動(dòng)的，具有不穩(wěn)定特性。若假定船舶在內(nèi)外錨地、泊位及前沿、港區(qū)航道間的轉(zhuǎn)移變量也符合高斯分布，那么對(duì)原始數(shù)據(jù)進(jìn)行統(tǒng)計(jì)后，可確定隨機(jī)轉(zhuǎn)移矩陣變量值如表3所示，其他狀態(tài)轉(zhuǎn)移變量為0?？紤]到轉(zhuǎn)移矩陣中需要下一個(gè)狀態(tài)轉(zhuǎn)移概率和為1，因而對(duì)表3所確定的變量進(jìn)行行歸一化計(jì)算：

(5)

3.3 船舶引航過程馬爾科夫鏈仿真

選用標(biāo)準(zhǔn)高斯分布進(jìn)行采樣，通過確定的船舶進(jìn)出港作業(yè)參數(shù)的高斯矩陣(對(duì)稱或者非對(duì)稱)作為狀態(tài)轉(zhuǎn)移函數(shù)。因此，選定三個(gè)條件進(jìn)行分析。

1)條件1：出多進(jìn)少。進(jìn)出港船舶數(shù)量不平衡，出港船舶多于進(jìn)港船舶；

2)條件2：進(jìn)出相當(dāng)。進(jìn)出港船舶數(shù)量平衡；

3)條件3：進(jìn)多出少。進(jìn)出港船舶數(shù)量不平衡，進(jìn)港船舶多于出港船舶。

針對(duì)三種條件進(jìn)行MC動(dòng)態(tài)系統(tǒng)仿真，可以得到船舶港區(qū)水域流量分布和船舶引航動(dòng)態(tài)風(fēng)險(xiǎn)數(shù)值。

表3 引航船舶狀態(tài)轉(zhuǎn)移變量

圖5是在S(0)隨機(jī)狀態(tài)下進(jìn)行不同狀態(tài)轉(zhuǎn)移條件下的船舶流量值，分別表示不同水域中在條件1、2、3下船舶流量分布。數(shù)據(jù)顯示出：

1)在船舶引航狀態(tài)的轉(zhuǎn)移之初，船舶狀態(tài)變化呈現(xiàn)很強(qiáng)的潮汐型波動(dòng)，反映船舶引航作業(yè)的基本規(guī)律。船舶引航狀態(tài)分布最終走向穩(wěn)定，符合馬氏穩(wěn)態(tài)特性，同時(shí)港口船舶引航狀態(tài)保持穩(wěn)定。一般地，船舶狀態(tài)在經(jīng)歷多次遷移，港區(qū)水域船舶的狀態(tài)趨于穩(wěn)定，航道流量占主體，約50%。即，對(duì)于港區(qū)運(yùn)行而言，一般會(huì)保持在總流量的約50%數(shù)量的船舶在航道上航行。

2)船舶初始狀態(tài)與穩(wěn)定狀態(tài)的關(guān)聯(lián)不大，初始狀態(tài)不影響最終的穩(wěn)態(tài)。也就是，船舶引航任務(wù)的不確定不影響港口船舶流的分布。

3)船舶狀態(tài)轉(zhuǎn)移變量對(duì)港區(qū)水域總體船舶秩序影響明顯。若船舶進(jìn)港多于出港，那么最終錨地和泊位的船舶量趨于平衡，一般會(huì)均保持在總流量的約25%；若船舶出港多于進(jìn)港，那么最終內(nèi)錨地與泊位和港外錨地的船舶量趨于平衡，一般會(huì)均保持在總流量的約25%；若船舶進(jìn)出港平衡，那么最終港區(qū)泊位與錨地和航道的船舶量趨于平衡，一般會(huì)均保持在總流量的約50%。

4)內(nèi)錨地船舶流量不穩(wěn)定，在10%以內(nèi)波動(dòng)。

圖5 三個(gè)隨機(jī)條件下不同狀態(tài)轉(zhuǎn)移后的船舶流量值Fig.5 Traffic flow under various ship activities with 3 random conditions

3.4 船舶引航過程風(fēng)險(xiǎn)的動(dòng)態(tài)系統(tǒng)仿真

3.4.1 船舶引航過程風(fēng)險(xiǎn)蒙特卡羅仿真

若設(shè)定在不同水域船舶引航風(fēng)險(xiǎn)占總體風(fēng)險(xiǎn)值情況為r，如式(4)所示。設(shè)定航道水域有一定流量(中間態(tài))，根據(jù)蒙特卡羅方法進(jìn)行仿真，結(jié)果如圖6所示。

1)由于在船舶引航狀態(tài)的轉(zhuǎn)移之初，船舶狀態(tài)變化呈現(xiàn)很強(qiáng)的波動(dòng)性，風(fēng)險(xiǎn)的波動(dòng)性也很強(qiáng)，但趨勢(shì)最終走向穩(wěn)定。一般地，船舶狀態(tài)在經(jīng)歷多次遷移，港區(qū)水域船舶的狀態(tài)趨于穩(wěn)定，航道引航風(fēng)險(xiǎn)(相對(duì))數(shù)值約0.23。對(duì)于港區(qū)運(yùn)行而言，不論進(jìn)出港船舶比例如何，在航道上航行的船舶引航風(fēng)險(xiǎn)一般占主體。

2)船舶狀態(tài)轉(zhuǎn)移變量對(duì)港區(qū)水域船舶引航影響還是存在的。與進(jìn)出港口船舶流量相當(dāng)相比，若船舶出港多于進(jìn)港，外錨地引航風(fēng)險(xiǎn)增加40%，泊位引航風(fēng)險(xiǎn)降低38%，內(nèi)錨地引航風(fēng)險(xiǎn)降低50%；若船舶進(jìn)港多于出港，外錨地引航風(fēng)險(xiǎn)降低40%，泊位引航風(fēng)險(xiǎn)增加25%，內(nèi)錨地引航風(fēng)險(xiǎn)增加 1倍。

3)泊位水域引航作業(yè)的風(fēng)險(xiǎn)相對(duì)較高，尤其是進(jìn)港船舶多于出港船舶時(shí)，風(fēng)險(xiǎn)(相對(duì))數(shù)值達(dá)到0.091 3。

4)船舶引航作業(yè)受到船舶流結(jié)構(gòu)影響。進(jìn)港船舶比例較低時(shí)，內(nèi)錨地和泊位前沿水域引航作業(yè)的風(fēng)險(xiǎn)相對(duì)低。隨著結(jié)構(gòu)比例增加，風(fēng)險(xiǎn)也上升；外錨地風(fēng)險(xiǎn)略有降低。

圖6 不同狀態(tài)轉(zhuǎn)移條件下不同水域風(fēng)險(xiǎn)值Fig.6 Risk outcome under various ship activities and random condition

3.4.2 船舶引航過程風(fēng)險(xiǎn)MCMC仿真

引入MH算法，利用切片抽樣快速收斂，對(duì)不同船舶進(jìn)出港的流量引發(fā)船舶引航風(fēng)險(xiǎn)進(jìn)行MCMC仿真，獲得船舶引航作業(yè)初始狀態(tài)和最終狀態(tài)之間轉(zhuǎn)移產(chǎn)生的風(fēng)險(xiǎn)，如圖7所示。總體上，不同船舶進(jìn)出港的流量引發(fā)船舶引航風(fēng)險(xiǎn)波動(dòng)明顯，但是趨勢(shì)相同。作業(yè)量增加，風(fēng)險(xiǎn)趨于穩(wěn)態(tài)。船舶進(jìn)出港的流量不同，風(fēng)險(xiǎn)數(shù)值差異不大，船舶港區(qū)引航風(fēng)險(xiǎn)(相對(duì))值接近0.32。

圖7 不同船舶狀態(tài)轉(zhuǎn)移條件下風(fēng)險(xiǎn)值Fig.7 Risk outcome under random transition condition

3.5 風(fēng)險(xiǎn)動(dòng)態(tài)系統(tǒng)仿真靈敏性分析

為測(cè)試MCMC的靈敏性，這里將初始狀態(tài)S(0)，轉(zhuǎn)移矩陣T和各狀態(tài)下的風(fēng)險(xiǎn)值r都進(jìn)行隨機(jī)抽樣，初始狀態(tài)設(shè)定錨地和泊位前沿均有作業(yè)任務(wù)，轉(zhuǎn)移矩陣按照表2來獲取，運(yùn)用MH算法進(jìn)行快速收斂，可得到圖8。

圖8 船舶隨機(jī)轉(zhuǎn)移時(shí)船舶狀態(tài)和風(fēng)險(xiǎn)值的靈敏性(引航任務(wù)分配時(shí))Fig.8 Sensitivity on ship status and risk outcome under random transition condition

數(shù)據(jù)顯示：總體上，船舶進(jìn)港的流量稍大于出港流量，船舶流在50%左右波動(dòng)。船舶進(jìn)出港的日均流量變化引發(fā)船舶引航風(fēng)險(xiǎn)有波動(dòng)，但波動(dòng)幅度趨于收斂。

以上仿真情況表明，引入MCMC進(jìn)行動(dòng)態(tài)系統(tǒng)仿真，可以得到船舶動(dòng)態(tài)風(fēng)險(xiǎn)趨于穩(wěn)態(tài)的規(guī)律，港口水域船舶引航的風(fēng)險(xiǎn)(相對(duì))大致在30%總體風(fēng)險(xiǎn)。鑒于目前港區(qū)船舶作業(yè)(含航行、靠離泊和移泊等作業(yè))總體風(fēng)險(xiǎn)大致在2‰，那么單艘船引航風(fēng)險(xiǎn)維持在0.2‰～0.6‰的水平[13]。

引航實(shí)施過程中，將有部分船舶在航道水域，那么在四個(gè)水域均有一定比例的船舶流量，進(jìn)行上述步驟的仿真，將得到港區(qū)水域整體的動(dòng)態(tài)風(fēng)險(xiǎn)。經(jīng)計(jì)算，得出以下結(jié)論：

1)盡管轉(zhuǎn)移變量和作業(yè)次數(shù)變化，引航初始時(shí)，單艘次船舶引航風(fēng)險(xiǎn)值接近0.41‰的穩(wěn)態(tài)；

2)引航實(shí)施過程中，單艘次船舶引航風(fēng)險(xiǎn)值接近0.67‰的穩(wěn)態(tài)，結(jié)果如圖9所示。圖9(a)表示引航初始狀態(tài)的風(fēng)險(xiǎn)(最終狀態(tài)，收斂之后)，圖9(b)表示引航實(shí)施中的風(fēng)險(xiǎn)(考慮中間狀態(tài)，收斂之前)，橫軸是風(fēng)險(xiǎn)值。

圖9 不同時(shí)期風(fēng)險(xiǎn)值的靈敏性Fig.9 Sensitivity on risk outcome during initial and piloting

風(fēng)險(xiǎn)數(shù)值擴(kuò)大是由于存在航道航行的風(fēng)險(xiǎn)和過程的實(shí)施風(fēng)險(xiǎn)。這一結(jié)論與船員及引航員的實(shí)踐認(rèn)識(shí)一致。

4 結(jié)論

1)提出了船舶引航過程的風(fēng)險(xiǎn)分析, 采用基于馬爾科夫鏈蒙特卡洛方法進(jìn)行動(dòng)態(tài)系統(tǒng)仿真得到了港區(qū)水域船舶引航過程的動(dòng)態(tài)風(fēng)險(xiǎn)宏觀特征，特別是風(fēng)險(xiǎn)的穩(wěn)定性表現(xiàn)形式。

2)相對(duì)于以往的風(fēng)險(xiǎn)評(píng)價(jià)，馬爾科夫鏈蒙特卡洛模型適用于水上交通過程風(fēng)險(xiǎn)仿真與決策，克服事故特征和誘因獲取的數(shù)據(jù)困難，在確定狀態(tài)轉(zhuǎn)移變量后，就能量化論證基于內(nèi)部微觀結(jié)構(gòu)的水上交通宏觀動(dòng)態(tài)風(fēng)險(xiǎn)的變化趨勢(shì)。

3)基于隨機(jī)過程的分析方法可以研究港口船舶引航演化的基本規(guī)律，揭示港區(qū)水域船舶流量量化的分布狀態(tài)，闡明港口船舶流量演化模式。

當(dāng)然，下一步研究可以是具體船舶在狀態(tài)轉(zhuǎn)移條件下過程風(fēng)險(xiǎn)動(dòng)態(tài)系統(tǒng)仿真，也可以是研究船舶連續(xù)時(shí)間和連續(xù)狀態(tài)下馬爾科夫過程的風(fēng)險(xiǎn)動(dòng)態(tài)系統(tǒng)仿真，從而進(jìn)行引航風(fēng)險(xiǎn)預(yù)警預(yù)控。

[1] COOK T, SHIPLEY P. Human factors studies of the working hours of UK ship′s pilots：A field study of fatigue[J]. Applied ergonomics, 1980, 11(2): 85-92.

[2] HU S, FANG Q, XIA H, et al. Formal safety assessment based on relative risks model in ship navigation[J]. Reliability engineering & system safety, 2007, 2(3): 369-377.

[3] FANG Q, HU S. Study on risk control of ship pilotage in Shanghai harbor[C]//Proceedings of the 4th IEEE International Conference on Management of Innovation and Technology, ICMIT., 2008, 1187-1192.

[4] 方泉根, 胡甚平. FSA在船舶引航風(fēng)險(xiǎn)評(píng)估中的應(yīng)用[J].哈爾濱工程大學(xué)學(xué)報(bào), 2006, 27(3): 229-334. FANG Quangen, HU Shenping. Application of formal satety assessment to the risk assessment of the ship- pilotage[J].Journal of Harbin Engineering University， 2006, 27(3): 329-334

[5] 方誠,胡甚平,方泉根. 港口船舶引航風(fēng)險(xiǎn)預(yù)測(cè)[J]. 中國航海, 2008, 31(4): 388-391. FANG Cheng, HU Shenping, FANG Quangen.Risk prediction of ship pilotage in harbor[J].Navigation of China, 2008, 31(4): 388-391.

[6] CAMERON J R. Improving the safety of marine pilotage[C]//2005 International Oil Spill Conference, IOSC， 2005： 3767-3771.

[7] 馬飛,趙月林. 基于資產(chǎn)脆弱性風(fēng)險(xiǎn)評(píng)估模型的船舶引航安全評(píng)估[J]. 大連海事大學(xué)學(xué)報(bào), 2008(S1): 17-19. MA Fei, ZHAO Yuelin. Ship pilotage safety assessment with the risk modal based on assets vulnerability[J]. Journal of Dalian Maritime University, 2008(S1): 17-19.

[8] WU Y, PENG G. Evaluation system of the portable pilot system based on fuzzy comprehensive evaluation[J]. Advanced materials research:manufacturing science and technology III., 2013, 622: 112-116.

[9] ?MER F，SELMIN Z, BURAK. Formal safety assessment for ship traffic in the Istanbul Straits[J]. Procedia—social and behavioral sciences, 2015, 207: 252-261.

[10] PAK J, YEO G, OH S, et al. Port safety evaluation from a captain′s perspective: The Korean experience[J], Safety science, 2015, 72: 172-181.

[11] PRAETORIUS G, HOLLNAGEL E, DAHLMAN J. Modelling vessel traffic service to understand resilience in everyday operations[J]. Reliability engineering & system safety， 2015, 141: 10-21.

[12] 胡甚平, 黎法明, 席永濤,等. 海上交通系統(tǒng)風(fēng)險(xiǎn)成因耦合機(jī)理仿真[J]. 應(yīng)用基礎(chǔ)與工程科學(xué)學(xué)報(bào), 2015(2): 409-419. HU Shenping, LI Faming, XI Yongtao,et al. Novel simulation on coupling mechanism of risk formation segments for marine traffic system [J]. Journal of basic and engineering, 2015(2): 409-419.

[13] FAGHIH R S, XIE M, NG K M. Accident risk assessment in marine transportation via markov modelling and markov Chain monte carlo simulation[J]. Ocean engineering, 2014, 91: 363-370.

[14] SCOLLNIK D P M. Actuarial modeling with MCMC and BUGS[J]. North American actuarial journal, 2001， (20)： 29-34.

[15] 胡甚平, 黃常海, 張浩. 基于云模型的海上交通系統(tǒng)風(fēng)險(xiǎn)蒙特卡洛仿真[J]. 中國安全科學(xué)學(xué)報(bào), 2012, 22(4): 20-26. HU Shenping, HUANG Changhai, ZHANG Hao. Cloud model- based simulation of system risk of marine traffic by monte carlo algorithm[J]. China safety science journal, 2012, 22(4): 20-26.

[16] 王春峰, 萬海輝, 李剛. 基于MCMC的金融市場(chǎng)風(fēng)險(xiǎn)VaR的估計(jì)[J]. 管理科學(xué)學(xué)報(bào), 2000(2): 54-61. WANG Chunfeng, WAN Haihui, LI Gang. Estimation of value- at- risk using MCMC[J]. Journal of management sciences in China, 2000(2): 54-61.

Markovchainmodelforthedynamicsimulationofprocessriskinshippilotageatharbor

HU Shenping1, HUANG Changhai1, DENG Hua2, HUANG Daozheng3

(1.Merchant Marine College, Shanghai Maritime University, Shanghai 200135, China; 2.Marine Technology Department, Jiangsu Maritime Vocational Technology School, Nanjing 211170, China; 3.Transportation College, Shanghai Maritime University, Shanghai 200135, China)

10.11990/jheu.201605025

http://www.cnki.net/kcms/detail/23.1390.u.20170427.1510.082.html

X951

A

1006- 7043(2017)09- 1391- 08

2016-05-09. < class="emphasis_bold">網(wǎng)絡(luò)出版日期

日期：2017-04-27.

中國博士后科學(xué)基金項(xiàng)目(2016M591651)；國家自然科學(xué)基金項(xiàng)目(51279099)；上海市自然科學(xué)基金項(xiàng)目(12ZR1412500).

胡甚平(1974-)，男，教授.

胡甚平，Email:Sphu@shmtu.edu.cn.

本文引用格式：胡甚平，黃常海，鄧華，等. 船舶引航過程風(fēng)險(xiǎn)馬爾科夫鏈模型動(dòng)態(tài)仿真[J]. 哈爾濱工程大學(xué)學(xué)報(bào)， 2017， 38(9): 1391-1398.

HU Shenping, HUANG Changhai, DENG Hua, et al. Markov chain model for the dynamic simulation of process risk in ship pilotage at harbor[J]. Journal of Harbin Engineering University， 2017， 38(9): 1391-1398.