董 雪,張德平
(南京航空航天大學計算機科學與技術學院,江蘇 南京 210016)
作戰(zhàn)效能是指武器裝備系統(tǒng)在作戰(zhàn)中發(fā)揮作用的有效程度[1],是衡量武器價值的重要指標。在有限的預算和時間內提升武器的作戰(zhàn)效能,是武器裝備發(fā)展的重點。作戰(zhàn)效能是多指標綜合作用的結果,由于預算和時間的限制,提高武器作戰(zhàn)效能最有效的途徑是找到對武器作戰(zhàn)效能影響較大的影響因素,即對武器作戰(zhàn)效能指標進行敏感性分析與篩選,進而找到與其相關聯(lián)的作戰(zhàn)設備,進行功能的完善與提高,從而提高作戰(zhàn)武器的整體作戰(zhàn)效能。
敏感性分析分為局部敏感性分析和全局敏感性分析[2]。局部敏感性分析通過改變模型中某個參數的取值,同時固定其他參數的取值,來探索模型響應的變化。但是局部敏感性分析不能計算變量間的交互作用對輸出結果的影響?;诖?,克服這些問題的全局敏感性分析得到廣泛應用。全局敏感性分析[3-4]方法包括回歸分析法、傅里葉振幅敏感性檢驗法、響應曲面法、互信息指數法、Sobol指數法等,其中Sobol指數法綜合性能更好。
全局敏感性分析方法多用于氣候、動力學方向的研究,Nossent等人[5]為評估復雜環(huán)境模型中參數的重要度,采用Sobol指數全局敏感性分析,產生復雜環(huán)境模型參數使用的關鍵信息。Zhang等人[6]使用基于方差的Sobol靈敏度分析,以中國河流流域為例,分析了干旱、正常、潮濕年份SWAT中參數的敏感性,從分析獲得的結果中,加深了對土壤水分評估工具(SWAT)中敏感性參數的理解及不同指標和氣候條件下的潛在水文過程。Mathieu等人[7]則將Sobol指數法應用在簡化功能性結構植物模型中的參數,用于增強冬季油菜生長的研究。
傳統(tǒng)Sobol指數敏感性分析采用蒙特卡洛方法計算敏感系數,需要大量的數據樣本。如果只單純依靠原復雜模型生成這些樣本,計算成本過高。為此,可以采用代理模型替代原模型生成樣本,節(jié)省計算時間,降低計算成本。目前,很多學者對代理模型進行了相關研究。Zhan等人[8]將多變量選擇適應回歸樣條(MARS)作為分布式時變增益模型(DTVGM)的代理模型,進行參數敏感性分析。孔凡哲等人[9]將支持向量回歸(SVR)引入水文模型構建代理模型,快速有效地定量評估參數敏感性,識別水文模型的關鍵參數。Luo等人[10]分別用徑向基人工神經網絡(RBFANN)與Kriging模型作為代理函數,在三氯乙烯(TCE)污染含水層上進行Sobol敏感性分析,以評估水井的修復時間、表面活性劑濃度和注射速率等設計變量對修復效率的敏感性。
當前基于代理模型的Sobol指數全局敏感性分析方法仍存在以下幾個問題:1)代理模型存在擬合時間緩慢、擬合精度不高的問題。2)實驗樣本質量問題,首先是生成代理模型,需要大量的訓練樣本進行訓練才能達到預訂擬合精度;其次是Sobol指數法敏感系數的計算,傳統(tǒng)蒙特卡洛抽樣方法,收斂速度慢,計算量大。為此,本文針對第一個問題,引入擬合精度高、學習速度快的極限學習機[11]作為代理模型,構建一種基于極限學習機的全局靈敏度分析模型,并將該模型應用到武器裝備作戰(zhàn)效能的靈敏度分析。針對第二個問題,為減少代理模型所需樣本數量,選擇采樣質量高的拉丁超立方體采樣方法,提高擬合精度;同時為提高Sobol指數法計算效率,引入低差異、分布均勻的擬蒙特卡洛抽樣方法代替蒙特卡洛抽樣方法,減少抽樣次數,提高收斂速度,降低計算成本。
武器作戰(zhàn)效能評估方法[12]繁多,主要包括試驗統(tǒng)計法、解析法、作戰(zhàn)模擬法和多指標綜合評估法等。其中,解析法中的ADC模型[13]是美國工業(yè)界武器系統(tǒng)效能咨詢委員會為美國空軍而建立的,其評估模型為:
E(t)=A×D×C
(1)
其中,A=[a1,a2,a3,…,an]為系統(tǒng)的可用性向量,表示系統(tǒng)開始執(zhí)行任務瞬間處于不同狀態(tài)的概率。ai為開始執(zhí)行任務時處于i狀態(tài)的概率;D=(dij)n×n為可信性矩陣,dij表示開始瞬間系統(tǒng)處于i狀態(tài)而在使用過程中轉移到j狀態(tài)的概率;C=(cjk)n×m為系統(tǒng)的能力矩陣,cjk表示在最后的可能狀態(tài)j中達到的第k項效能指標值。利用ADC模型,能夠合理地評估潛艇的作戰(zhàn)效能,但是模型參數較多,計算復雜。
代理模型可以對一組輸入輸出數據之間的關系用具體的數學表達式或數學模型表示。為精確擬合輸入輸出變量之間的關系,提高計算效率,本文引入極限學習機(ELM)作為代理模型,代替復雜的效能評估模型。ELM是一種單隱層的神經網絡算法,ELM可以隨機初始化輸入權重和偏置并得到相應的輸出權重。
假設有N個不同樣本的訓練集{(xi,oi)|xi∈Rn,ti∈Rm,i=1,…,N},其中輸入樣本xi=(xi1,xi2,…,xin)T,真實值oi。假設隱藏層含有K個節(jié)點,極限學習機[14]可以表示為:
(2)
式(2)中,xj為第j個樣本輸入;βi是第i個隱層節(jié)點與輸出節(jié)點的輸出權重;g(x)是隱層激活函數;ωi=(ωi1,ωi2,…,ωin)是第i隱層節(jié)點與輸入節(jié)點之間的權重;bi是第i個隱層單元的偏置;tj為模型輸出值。該網絡模型可以用矩陣形式表示:
Hβ=T
(3)
其中,H是隱藏層節(jié)點輸出矩陣:
H(ω1,…,ωL,b1,…,bL,x1,…,xL)=
該方法的學習目標是尋找最優(yōu)權值,使得輸出的誤差最小。數學模型可表示為:
(4)
式(4)中,e表示模型輸出值與真實值之間的誤差,ej是第j個樣本的模型輸出值與真實值之間的誤差。
單隱層的神經網絡可以轉化為求解一個線性模型的問題。由此確定輸出權重:
(5)
實驗設計方法有很多,常見的經典采樣法[15]有中心復合實驗、全因子設計實驗、正交設計實驗等。這些采樣方法,在多變量情況下,需要大量采樣點才能獲得可靠實驗結果,計算成本高。拉丁超立方體采樣法[16]采用等概率隨機正交分布的原則,樣本點的個數可以靈活設定,適用于多變量問題。所以本文在構建代理模型時,采用拉丁超立方體采樣方法生成訓練樣本的輸入變量集。拉丁超立方體采樣點生成策略如下:
Sobol指數法敏感系數的計算,傳統(tǒng)方法是通過蒙特卡洛抽樣[17]實現,需要大量的實驗樣本,收斂速度慢,計算時間長。為有效減少試驗次數,提高計算準確性,采用低差異、分布均勻的擬蒙特卡洛方法。常見的擬蒙特卡洛采樣方法有Halton序列、Hammersley序列、Sobol序列等。其中Sobol序列分布均勻且不受樣本數量限制,可替代蒙特卡洛抽樣進行敏感性分析[18-20]。Sobol序列采樣點生成策略如下:
Sobol序列是基于一組直接數di構造的隨機序列,設qi是小于2i的正奇數,則di=qi/2i。
di以及qi的生成需要借助系數只為0或1的簡單多項式,可表示為:
f(x)=xp+a1xp-1+…+ap-1x+ap
(6)
式(6)中,p為多項式的度數,a1,a2,…,ap為多項式系數。對于i>p,由遞推公式求得di:
di=a1di-1⊕a2di-2⊕…⊕apdi-p⊕?di-p/2p」
(7)
式(7)中,⊕表示二進制按位異或。對于qi,遞推公式為:
qi=2a1qi-1⊕22a2qi-2⊕…⊕2papqi-p⊕qi-p
(8)
綜合以上推理,可以利用公式(9)生成Sobol序列x1,x2,x3,…
xn=b1d1⊕b2d2⊕…
(9)
式(9)中,bi是n的二進制形式。
Sobol指數法是由俄羅斯學者Sobol[21]提出,并以他的名字命名的一種基于方差分解的全局敏感性分析方法。其核心思想是方差分解,把模型以單參數及參數之間組合的方式表示,通過計算單個輸入參數或輸入參數集的方差對總輸出方差的影響來分析參數的重要性以及參數之間的交互效應。
假設數學模型為Y=f(X),平方可積,分解為單個模型參數及參數之間相互作用的子項函數之和:
(10)
其中,X=(x1,x2,…,xn),xi屬于n維單位立方體Hn,式(10)中一共含有2n個子項。如式(10)滿足:
(11)
(12)
公式(10)兩邊對除xi以外的其他各項求積分得到:
(13)
公式(10)兩邊對除xi、xj以外的其他各項求積分可得:
(14)
以此類推,可以得到公式(10)中等式右邊的各個分解函數。
基于以上條件,Sobol指數敏感性分析方法定義了偏方差和總方差,并通過偏方差占總方差的比率來表示模型參數及其交互作用對目標響應的影響程度,其中模型f(X)的總方差為:
(15)
各子項的偏方差Di1,i2,…,is為:
(16)
變量的敏感性指數Si1,i2,…,is為:
(17)
式(17)中,Si∈(0,1)表示參數xi的一階敏感性指數,描述了參數xi對輸出的貢獻度。參數的一階敏感性指數越大,它對輸出值的影響越大。為描述參數的整體影響,即參數的一階敏感性影響及其與其他所有參數的交互影響對輸出值的貢獻度,引入了參數總敏感性影響指數[22]STi,總敏感性指數包含了變量之間的交互效應。若一個輸入變量的全效應指數很小,表明該變量的變化不僅對輸出影響小,而且與其他變量之間的交互效應也很小。因此,可以考慮對全效應指數小的變量取固定值,減少可變變量個數,從而簡化模型。
根據文獻[22],上述一階敏感性指數可用公式(18)計算:
(18)
總敏感性指數:
(19)
為解決基于ADC模型的作戰(zhàn)效能評估存在的計算成本高、計算時間長的問題,本文提出基于極限學習機的武器裝備作戰(zhàn)效能全局敏感性分析模型,目的是找到影響作戰(zhàn)效能的關鍵指標,進行功能完善與提高。模型計算流程如下:
算法1基于極限學習機的Sobol指數法
Step1明確作戰(zhàn)任務,構建指標體系。
Step2采樣設計,生成訓練樣本輸入集,輸入原效能模型,計算并輸出訓練樣本結果集。
Step3訓練樣本集歸一化處理。
Step4輸入訓練樣本集,構建代理模型。
Step5采樣設計,生成計算樣本輸入集并歸一化,運行代理模型,輸出計算樣本結果集。
Step6計算樣本集反歸一化處理,利用Sobol指數敏感性分析,篩選敏感指標。
Step7根據敏感指標找到對應作戰(zhàn)武器裝備,對其進行功能優(yōu)化。
算法1各步驟說明如下:
Step1針對潛艇具體作戰(zhàn)任務,構建潛艇作戰(zhàn)效能評估模型及其指標體系,效能評估模型為ADC模型。
Step2確定各指標的取值范圍,采用拉丁超立方體采樣方法生成N組訓練樣本輸入變量,通過公式(1)計算出相應的效能值,得到N組訓練樣本集。
Step3進行樣本預處理,對效能指標進行無量綱化處理,將樣本值歸一化到[0,1]。歸一化公式為:
(20)
Step4生成代理模型,利用生成的訓練樣本集代入公式(4)訓練極限學習機,擬合精度達到預先設定值,訓練結束。代理模型擬合優(yōu)度用均方誤差(MES)和平均絕對誤差(MAE)來衡量,值越小,預測結果越好。其中yi表示真實值,yi′為預測值,n為樣本個數。
均方誤差計算公式:
(21)
平均絕對誤差計算公式:
(22)
Step5進行敏感性分析,利用擬蒙特卡洛方法中的Sobol序列,根據公式(6)~公式(9)生成計算樣本輸入集,由公式(20)歸一化,通過代理模型計算效能值,生成計算樣本集,進行反歸一化處理,得到敏感性分析數據。
Step6根據公式(18)~公式(19)計算每個指標的一階敏感性系數及總體敏感性系數。
Step7篩選出影響武器裝備作戰(zhàn)效能的敏感性指標,找到該指標對應的性能指標,進而得到影響武器效能的關鍵裝備。
本章包含2個實驗,實例1為潛艇攻擊海面艦艇的實例,將極限學習機作為代理模型進行敏感性分析,并與BP代理模型、SVR代理模型對比,證明方法的可行性和有效性。實例2為方法的應用,將本文提出的方法應用到新的案例中,驗證方法的可用性。
實例1將潛艇作戰(zhàn)系統(tǒng)劃分為潛艇平臺系統(tǒng)和武器系統(tǒng)這2部分,每個系統(tǒng)的初始狀態(tài)分為故障和正常2種。本例構建潛艇作戰(zhàn)能力評估指標如圖1所示。
圖1 潛艇作戰(zhàn)效能評估指標
潛艇效能評估ADC模型中,能力矩陣C是各指標綜合計算的結果,本例以潛艇各系統(tǒng)正常狀態(tài)工作下的作戰(zhàn)效能為例進行敏感性分析。由于效能指標多,本實驗選取火力打擊能力為敏感性分析的目標,其中,武器類型取值為離散型數值,分別為1、2、3,對應3種不同的武器類型。剩余5個子指標,如表1所示。根據每個指標的屬性,確定取值范圍。
表1 火力打擊能力子指標
指標取值范圍最小值最大值發(fā)射速度/(km·h-1)50120射程距離/km60110制導精度/°140武器數量/個310毀傷半徑/m1050
實驗采用基于拉丁超立方體采樣的ELM代理模型(L-ELM)、BP代理模型(L-BP)、SVR代理模型(L-SVR)進行對比試驗。隱層神經元個數為20,激活函數設置為sigmoid函數。訓練樣本選取樣本集的90%,測試樣本為樣本集的10%。樣本集設定為300個,分別運行3種代理模型。
表2 實例1代理模型實驗結果與真實值對比
衡量指標模型名稱L-ELML-BPL-SVRMSE0.0030.0100.020MAE0.0420.0660.135
由表2可知,3種代理模型分別對測試樣本運算,結果與真實結果相比,其中L-ELM模型的MSE值、MAE值最小,說明L-ELM代理模型的擬合度最優(yōu)。
圖2為3種代理模型測試樣本值與真實樣本值擬合效果圖。ADC曲線為真實效能值,從圖中可以看出,L-ELM代理模型擬合效果優(yōu)于L-BP模型、L-SVR模型,擬合效果最好。
圖2 代理模型與真實模型效能值對比
圖3為3種代理模型計算不同樣本數量所用的時間,樣本數量取[10,100,1000,10000,100000]進行測試,當樣本數量超過10000時,L-BP模型與L-SVR模型運行時間變慢,L-ELM模型速度保持基本不變。
潛艇效能敏感性分析,設置3組對比試驗與真實敏感系數進行對比,抽樣方法為低差異的Sobol序列,設置樣本數為10000。對比試驗為:1)基于L-BP模型的敏感性分析;2)基于L-ELM模型的敏感性分析;3)基于L-SVR模型的敏感性分析。實驗結果如表3所示。
圖3 測試時間對比圖
表3 火力打擊能力指標敏感系數
模型名稱S1S2S3S4S5S6真實值0.0910.0450.4830.00013.95e-050.234L-ELM0.0850.0480.5160.00040.00020.245L-BP0.0690.0430.4830.02430.00430.188L-SVR0.0660.0380.6290.0009-2.28e-050.231
由表3可知,基于L-ELM模型計算的敏感性系數整體與真實的敏感性系數非常接近且穩(wěn)定,而基于L-BP模型與L-SVR模型計算的敏感性系數不穩(wěn)定,有些與真實值非常接近,如L-BP模型中敏感系數S3,L-SVR模型中S6與真實值非常接近,但有些偏離過大,如L-BP模型中S4、S5,L-SVR模型中S1、S3。
表4為火力打擊能力指標的總敏感系數,可以看出,基于L-ELM代理模型計算得到的總敏感系數與真實值更接近,整體效果較好,而基于L-BP模型、L-SVR模型計算得到的總敏感系數與真實值差別較大,尤其是St4和St5與真實值差別明顯。
表4 火力打擊能力指標總敏感系數
模型名稱St1St2St3St4St5St6真實值0.13910.07300.60043.09e-059.70e-050.3346L-ELM0.11430.08190.59320.00790.00650.2999L-BP0.14940.10940.57710.11790.08750.2770L-SVR0.07790.04440.65600.00380.00250.2500
根據表3敏感系數Si從大到小排序S3>S6>S1>S2>S4>S5,得到影響火力打擊能力的指標從強到弱依次為制導精度、毀傷半徑、發(fā)射速度、射程距離、武器類型、武器數量。根據表4總敏感系數Sti從大到小排序St3>St6>St1>St2>St5>St4,得到火力打擊能力指標間相互影響程度從強到弱依次為制導精度、毀傷半徑、發(fā)射速度、射程距離、武器數量、武器類型。分析發(fā)現,制導精度、毀傷半徑、發(fā)射速度是火力打擊能力的關鍵因素,可以通過提高其相應的設備來提高火力打擊能力。而射程距離、武器數量、武器類型自身的變化對效能值的影響較小,且對其他效能指標影響較小,所以在效能優(yōu)化中,其值可以設定為固定值,減少計算成本及計算時間。
潛艇裝備防空導彈對于空中威脅目標是十分必要的,在潛艇對空作戰(zhàn)中,根據防空導彈裝備體系的構成,作戰(zhàn)能力劃分為戰(zhàn)場感知能力、指揮控制能力、火力打擊能力,防空導彈效能值計算所需指標及取值范圍如表5所示。
表5 防空導彈武器作戰(zhàn)指標
指標名稱最小值最大值雷達平均發(fā)射功率/kw410脈沖寬度/μs15脈沖重復頻率/Hz5001000雷達電磁波長/cm3.757.5雷達天線增益/dBi1020雷達目標截面積/m215距離/km20004000接收機帶寬/MHz30500接收機等效噪聲系數12決策時間/min410戰(zhàn)斗準備時間/min68目標停留時間/min25制導精度/m50200殺傷半徑/m100200
訓練樣本集設為1000,通過拉丁超立方體取樣,將3種代理模型與真實模型仿真結果進行對比,結果如表6所示。由表6可知,L-EML模型計算結果與真實結果更加相近,擬合效果優(yōu)于L-BP,L-SVR模型。
表6 實例2代理模型實驗結果與真實值對比
評價指標模型名稱L-ELML-BPL-SVRMSE0.00680.01070.0133MAE0.06350.08400.0981
防空導彈作戰(zhàn)效能指標一階敏感性分析結果、總敏感性分析結果分別如圖4、圖5所示。
圖4 防空導彈作戰(zhàn)效能指標一階敏感性分析結果
圖5 防空導彈作戰(zhàn)效能指標總敏感性分析結果
從圖4可以看出,影響防空導彈作戰(zhàn)效能的指標敏感系數最大為雷達目標截面積S7,最小為決策時間S10、戰(zhàn)斗準備時間S11。根據圖4和圖5分析對比,可以將總敏感系數小的指標設置為常數,表示它們的變化與其他指標間的相互影響較小,如決策時間St10、戰(zhàn)斗準備時間St11、目標停留時間St3等,而雷達目標截面積S6、敵我距離S7、雷達電磁波長指標一階敏感性指數與總敏感指數都很高,對提高防空導彈的效能影響較大,可以針對性地改進雷達中與其相關的裝置,提高作戰(zhàn)效能。
為提高武器裝備的作戰(zhàn)效能,通過全局敏感性分析找到影響作戰(zhàn)效能的關鍵指標,引入ELM學習機可有效解決效能評估模型及仿真計算模型中計算成本高、時間長的問題。將本文方法與BP神經網、SVR支持向量回歸模型進行對比,在計算精度相同的前提下,本文方法可節(jié)省大量計算時間。通過實例分析,基于ELM的全局敏感性分析方法分析得到的敏感系數合理,同時通過與真實結果的比對,驗證了本方法的可行性,并且在樣本數量極大的情況下,基于ELM的全局敏感性分析方法的時間消耗與其他算法相比,有明顯的降低,使得在有限時間和預算內提升潛艇作戰(zhàn)效能評估成為可能。
參考文獻:
[1] 付東,方程,王震雷. 作戰(zhàn)能力與作戰(zhàn)效能評估方法研究[J]. 軍事運籌與系統(tǒng)工程, 2006,20(4):35-39.
[2] Iooss B, Lemaitre P. A review on global sensitivity analysis methods[M]// Operations Research/ Computer Science Interfaces Series(Vol.59). Springer, 2014:101-122.
[3] 羅鵬程,傅攀峰. 武器裝備敏感性分析方法綜述[J]. 計算機工程與設計, 2008,29(21):5546-5549.
[4] 張曉航. 防空導彈武器裝備體系作戰(zhàn)效能全局敏感性分析方法研究[D]. 長沙:國防科學技術大學, 2010.
[5] Nossent J, Elsen P, Bauwens W. Sobol’ sensitivity analysis of a complex environmental model[J]. Environmental Modelling and Software, 2011,26(12):1515-1525.
[6] Zhang Chi, Chu Jinggang, Fu Guangtao. Sobol’s sensitivity analysis for a distributed hydrological model of Yichun River Basin, China[J]. Journal of Hydrology, 2013,480:58-68.
[7] Mathieu A, Vidal T, Jullien A, et al. Sensitivity analysis to help individual plant model parameterization for winter oilseed rape[C]// IEEE International Conference on Functional-Structural Plant Growth Modeling, Simulation, Visualization and Applications. 2017:133-139.
[8] Zhan Che-sheng, Song Xiao-meng, Xia Jun, et al. An efficient integrated approach for global sensitivity analysis of hydrological model parameters[J]. Environmental Modelling & Software, 2013,41(41):39-52.
[9] 孔凡哲,宋曉猛,占車生,等. 水文模型參數敏感性快速定量評估的RSMSobol方法[J]. 地理學報, 2011,66(9):1270-1280.
[10] Luo Jiannan, Lu Wenxi. Sobol’ sensitivity analysis of NAPL-contaminated aquifer remediation process based on multiple surrogates[J]. Computers & Geosciences, 2014,67:110-116.
[11] Huang Guang-bin, Zhu Qin-yu, Siew C K. Extreme learning machine: A new learning scheme of feedforward neural networks[C]// Proceedings of IEEE International Joint Conference on Neural Networks. 2004:985-990.
[12] 高尚,婁壽春. 武器系統(tǒng)效能評定方法綜述[J]. 系統(tǒng)工程理論與實踐, 1998,18(7):109-114.
[13] 吳曉鋒,錢東. 用于系統(tǒng)效能分析的WSEIAC模型及其擴展[J]. 系統(tǒng)工程理論與實踐, 2000,20(8):1-6.
[14] Ding Shifei, Zhao Han, Zhang Yanan, et al. Extreme learning machine: Algorithm, theory and applications[J]. Artificial Intelligence Review, 2015,44(1):103-115.
[15] 魏昕. 基于元模型的全局優(yōu)化算法研究[D]. 武漢:華中科技大學, 2012.
[16] Mckay M D, Beckman R J, Conover W J. A comparison of three methods for selecting values of input variables in the analysis of output from a computer code[J]. Technometrics, 2000,42(1):55-61.
[17] Buslenko N P, Golenko D I, Sobol I M, et al.The Monte Carlo method[J]. Journal of the American Statistical Association, 1949,44(247):335-341.
[18] Sobol I M. On the distribution of points in a cube and the approximate evaluation of integrals[J]. USSR Computational Mathematics & Mathematical Physics, 1976,7(4):86-112.
[19] 周心蓮. 幾個常用隨機數及其性質的比較[J]. 鄖陽師范高等??茖W校學報, 2010,30(6):13-17.
[20] Wang Chen, Duan Qingyun, Gong Wei, et al. An evaluation of adaptive surrogate modeling based optimization with two benchmark problems[J]. Environmental Modelling & Software, 2014,60(76):167-179.
[21] Sobol I M. Sensitivity estimates for nonlinear mathematical models[J]. Matem Mod, 1993,2(1):112-118.
[22] Paruggia M. Sensitivity analysis in practice: A guide to assessing scientific models[J]. Journal of the Royal Statistical Society, 2005,168(2):466.