基于CU變換的非齊次泊松過程的統(tǒng)計驗證模型

范朝霞，趙 明，楊劍鋒

(1.貴州大學(xué) 理學(xué)院，貴州 貴陽 550025;2.耶夫勒大學(xué) 技術(shù)與可持續(xù)發(fā)展學(xué)院，瑞典 80176;3.貴州理工學(xué)院 信息工程學(xué)院，貴州 貴陽550003)

基于CU變換的非齊次泊松過程的統(tǒng)計驗證模型

范朝霞1，趙 明2，楊劍鋒3*

(1.貴州大學(xué) 理學(xué)院，貴州 貴陽 550025;2.耶夫勒大學(xué) 技術(shù)與可持續(xù)發(fā)展學(xué)院，瑞典 80176;3.貴州理工學(xué)院 信息工程學(xué)院，貴州 貴陽550003)

開源軟件中故障輸入數(shù)據(jù)的非齊次泊松過程(NHPP)的檢驗和均值函數(shù)的估計是排隊論分析中重要的子問題。本文依據(jù)均值函數(shù)的連續(xù)性和緩慢變化性，把NHPP轉(zhuǎn)化為小區(qū)間泊松性的分段常數(shù)非齊次泊松過程(PC NHPP)進行研究，主要運用條件均勻(CU)變換法對數(shù)據(jù)調(diào)整變換后，驗證泊松性。本文還運用Bugzilla的錯誤報告數(shù)據(jù)對模型實踐分析，結(jié)果顯示CU變換后的驗證效果可以接受。最后運用R語言擬合參數(shù)，計算出故障輸入過程參數(shù)，為排隊論模型構(gòu)造和求解奠定基礎(chǔ)。

NHPP檢驗；CU變換；K-S檢驗； R語言

針對軟件排錯過程中故障輸入過程的統(tǒng)計分析方法[1]，國內(nèi)目前應(yīng)用廣泛的是假設(shè)錯誤輸入過程是NHPP，然后根據(jù)三類常見累積強度函數(shù)模型：線性模型、冪率模型、對數(shù)線性模型等，擬合比較最優(yōu)的累積強度函數(shù)模型，求解參數(shù)，運用到可靠性分析[2]中。但是這些忽略了一個重要的問題，即研究都是提前假設(shè)了NHPP，沒有對其合理性給出證明。因為NHPP的強度函數(shù)λ(t)形式變化的多樣性和不確定性，NHPP的驗證也成了一個難題。國外針對NHPP的假設(shè)檢驗有一定程度的研究，Brown[3]提出了用PC NHPP近似NHPP的驗證方法，根據(jù)λ(t)對于時間變化的連續(xù)性和緩慢性，在一定小的區(qū)間上可以認(rèn)為λ(t)不變，然后在每個充分小子區(qū)間上驗證泊松性。Kim[4]總結(jié)了四種NHPP的檢驗方法，對小區(qū)間的泊松過程做恰當(dāng)?shù)霓D(zhuǎn)換，變成我們易處理的分布類型，然后進行K-S檢驗，最終驗證數(shù)據(jù)的NHPP性。

參數(shù)估計是在NHPP得到驗證后，對擬合模型的估算參數(shù)[5-6]，求解均值函數(shù)和強度函數(shù)，進而為排隊論模型分析做準(zhǔn)備?？煽啃灾谐Ｓ玫氖亲訁^(qū)間上的分段函數(shù)擬合[7]，但是鑒于簡潔性，本文在整區(qū)間上運用最小二乘估計。

1 常見NHPP模型

1.1 Goel-Okumoto模型

G-O模型是NHPP類軟件可靠性最基本的模型[8]，該模型的適用環(huán)境較理想化，但是模型簡單，便于求解，故而應(yīng)用廣泛。該模型假設(shè)故障被檢測到的概率b是恒定不變的。

強度函數(shù)為：λ(t)=abe-bt，表示的是單位時間故障輸入率的變化趨勢。

均值函數(shù)為：m(t)=a(1-e-bt)，其表示的實際意義就是到時間t為止的累積頻數(shù)，故而本文中用均值函數(shù)來擬合原數(shù)據(jù)的累積頻數(shù)。

1.2 Yamada Delayed S-Shaped模型

Y-D模型是G-O模型的進一步推導(dǎo)，該模型適當(dāng)?shù)胤艑捔瞬糠掷碚摷僭O(shè)，更符合實際情況。模型假設(shè)故障被檢測到的概率b是與時間有關(guān)的函數(shù)，且檢測率隨著時間的增加、系統(tǒng)內(nèi)故障的減少而越來越大。

強度函數(shù)為：λ(t)=ab2te-bt，

均值函數(shù)為：m(t)=a[1-(1+bt)e-bt]，

a>0，b>0

1.3 Inflected S-Shaped模型

I-S模型也是假設(shè)故障檢測率b是一個與時間有關(guān)的函數(shù)，只是函數(shù)的形式發(fā)生了變化，自變量時間t與檢測率函數(shù)b(t)不再是多項式關(guān)系，變成了負(fù)指數(shù)形式的變化趨勢。

2 基于CU變換的NHPP統(tǒng)計檢驗?zāi)Ｐ?/h2>

在呼叫中心系統(tǒng)中，一天24小時的呼入率是不一樣的。因此，在以往的排隊論分析中，均假設(shè)每小時的呼入率是常數(shù)，在每個時間段上建立排隊論模型，綜合分析對應(yīng)指標(biāo)。類似于常用的分段常數(shù)的理念，直接把NHPP的強度函數(shù)轉(zhuǎn)化為PC NHPP進行檢驗分析。

2.1 CU(conditional uniform)轉(zhuǎn)換K-S檢驗

定理1 泊松分布事件發(fā)生時刻條件均勻原理：在[0，T]時間段內(nèi)，已知事件發(fā)生了n次的前提條件下，各事件發(fā)生的時刻{Ti，i∈1∶n}在不考慮順序依存關(guān)系發(fā)生的情況下，可看做相互獨立的U[0，T]隨機變量[9]。

本文的CU變換思想就源自泊松分布事件發(fā)生時刻條件均勻原理。當(dāng)不考慮Ti的順序時，{Ti，i∈1∶n}～U[0，T]，變形為{Ti/T，i∈1∶n}～U[0，1]，則經(jīng)驗累積分布函數(shù)Fn(x)：

(1)

同時，由均勻分布知對應(yīng)的理論累積分布函數(shù)F(x)：

F(x)=x，0≤x≤1

2.2 Log(logarithmic)變換K-S檢驗

定理2 泊松分布事件發(fā)生的時間間隔Xn，n=1，2，…服從參數(shù)為λ的指數(shù)分布，且相互獨立。

由定理2知，欲檢驗小區(qū)間上的泊松性質(zhì)，也可通過檢驗事件發(fā)生的時間間隔是不是服從指數(shù)分布來檢驗數(shù)據(jù)的泊松性。Brown構(gòu)造了率1指數(shù)分布隨機變量的Log變換：

1{Xjlog，n≤x}理論累積分布函數(shù)F(x)：

F(x)=1-e-λx，x≥0

2.3 K-S檢驗

1）掃描矢量化陜西省民政廳提供的1998年鄉(xiāng)級行政區(qū)域界線協(xié)議書附圖，獲取全省鄉(xiāng)級行政區(qū)域界線、界址點、界樁點和三交點等原始界線矢量數(shù)據(jù)，并對矢量數(shù)據(jù)進行坐標(biāo)轉(zhuǎn)換，即1954年北京坐標(biāo)系、1956黃海高程系轉(zhuǎn)換至2000國家大地坐標(biāo)系、1985國家高程基準(zhǔn)。

定義 Kolmogorov分布函數(shù)：

K-S檢驗是檢驗小樣本數(shù)據(jù)的泊松性。K-S檢驗分為單樣本K-S檢驗和兩樣本K-S檢驗，單樣本是檢驗一組樣本數(shù)據(jù)和已知的概率分布類型的擬合優(yōu)度的。

由定理3，統(tǒng)計量的漸近分布為：

3 實例驗證

3.1 數(shù)據(jù)來源與清洗

Bugzilla失效數(shù)據(jù)是由Mozilla公司開發(fā)的錯誤追蹤系統(tǒng)Bugzilla(http://www.bugzilla.org/)在版本升級和內(nèi)測期間由于系統(tǒng)失效而產(chǎn)生的一系列被詳細(xì)記錄的故障數(shù)據(jù)組成。本文數(shù)據(jù)即2010年11月Bugzilla的第四版本正式內(nèi)測上線開始記錄的4.0版本歷史故障數(shù)據(jù)。統(tǒng)計每月接收到的故障頻數(shù)。在Bugzilla V4.0中，故障輸入時間共延續(xù)45個單位，其中12個空數(shù)據(jù)，且都集中于序列尾部，可參照圖1。這就是統(tǒng)計中的拖尾現(xiàn)象。本文采取連續(xù)兩次遇零截尾，即Freq連續(xù)兩次為零時，截去后面的數(shù)據(jù)。

圖1 故障輸入過程頻數(shù)

3.2 NHPP的檢驗

由于K-S檢驗需要隨機變量是連續(xù)的，故而針對單個到達的離散變量，可以考慮其到達時刻、到達時間間隔等連續(xù)型變量，然后采用適合于小樣本檢驗的K-S檢驗準(zhǔn)確推導(dǎo)結(jié)論。

3.2.1 CU-KS檢驗

CU變換基于定理1對子區(qū)間樣本數(shù)據(jù)進行轉(zhuǎn)換。針對故障輸入時間數(shù)據(jù)。記到達時刻Opened_Day為Ti，區(qū)間長度T為31天，則根據(jù)CU變換理論，此處有Ti/T～U[0，1]。針對第四個子區(qū)間的CU變換數(shù)據(jù)整理如表1：

表1 第四個子區(qū)間CU變換數(shù)據(jù)

運用統(tǒng)計軟件中的R軟件結(jié)合K-S檢驗理論檢驗數(shù)據(jù)的泊松性質(zhì)，運行結(jié)論如下：

>ks.test(Temp，"punif")

>D = 0.12054， p-value = 0.9453

顯然，子區(qū)間上的泊松性是明顯的。但樣本量是19，在原數(shù)據(jù)中是相對較多的，故而檢驗效果相對較好，也證明了在足夠數(shù)據(jù)量的情況下，小區(qū)間的泊松性質(zhì)驗證很好。

但是在表2中有數(shù)據(jù)量Freq為1或0的情況，此時CU變換后小區(qū)間的泊松性質(zhì)并不樂觀。結(jié)果說明子區(qū)間上小樣本情況下的泊松性質(zhì)相對于多樣本的表現(xiàn)較差。同時需要注意，多樣本的小區(qū)間數(shù)據(jù)也出現(xiàn)有泊松檢驗效果較差的現(xiàn)象，分析發(fā)現(xiàn)，故障會有偶然聚集出現(xiàn)的特殊情況，這對統(tǒng)計量的均勻性檢驗很不利，故而出現(xiàn)了少有的特例情況。

表2 CU變換后各子區(qū)間上的p值和D值

但是整體來講，E(p)=0.522，E(D)=0.397，整體檢驗效果還是達到了的，已經(jīng)滿足了不能拒絕原假設(shè)PC NHPP的條件。

3.2.2 Log-KS檢驗

Log變換是基于定理2對子區(qū)間數(shù)據(jù)進行的對數(shù)變換。第四個子區(qū)間Log變換后數(shù)據(jù)如表3所示，根據(jù)Log變換理論可知，數(shù)據(jù)表中的序列0，是因為前后兩個時間點一致使得真數(shù)為1造成的。無窮小量(Inf)則是因為真數(shù)趨于0所致。

表3 第四個子區(qū)間Log變換數(shù)據(jù)

針對第四子區(qū)間變換數(shù)據(jù)進行K-S率1指數(shù)檢驗，R程序運行結(jié)果如下：

>ks.test(xlog，"pexp")

>D = 0.31579， p-value = 0.04521

顯然的，同一組預(yù)處理的數(shù)據(jù)，Log變換后檢驗效果沒有CU檢驗效果明顯。這也就暗示我們Log變換的數(shù)據(jù)預(yù)處理更復(fù)雜，需要更多的細(xì)節(jié)調(diào)整。

Log變換后的整體結(jié)論是E(p)=0.378，E(D)=0.327，雖然通過檢驗，但Log變換后的數(shù)據(jù)的率1指數(shù)檢驗效果沒有CU變換好。這與數(shù)據(jù)選取的時間不是絕對連續(xù)，使得對數(shù)變換對結(jié)果的影響較大等因素有關(guān)。Log變換時，應(yīng)該注意盡可能使得每個樣本數(shù)據(jù)都不一樣，避免0值和無窮小值出現(xiàn)，影響檢驗結(jié)果。鑒于驗證NHPP性已經(jīng)達到，此處Log變換的數(shù)據(jù)處理方法不做更深研究。

3.3 參數(shù)估計

鑒于故障輸入過程NHPP性質(zhì)的檢驗成立，此處用常用的NHPP類軟件可靠性模型對數(shù)據(jù)進行擬合分析。用Yamada Delayed S-Shaped(Y-D)模型擬合分析：

>nls(X$CFreq ～ a*(1-(1+b*X$time)*exp(-b*X$time))，data=X[，c(1，3)]，start = list(a=200，b=0.5))

用R里的nls函數(shù)對均值函數(shù)的參數(shù)進行估計，得a=252，b=0.2044，且兩個參數(shù)是顯著性不等于0。此時的標(biāo)準(zhǔn)殘差為3.811，迭代次數(shù)為6，容差為7個分?jǐn)?shù)位。

均值函數(shù)：

m(t)=252·[1-(1+0.2044·t)·e-0.2044t]

圖2 Y-D模型擬合

由圖2可以看出，Y-D模型擬合了故障輸入過程累積頻數(shù)的趨勢。分析K-S擬合檢驗結(jié)果，p值為0.9794，接近1，且D值相對較小，接受Y-D擬合。

>ks.test(X$CFreq，f(X$time))

>D = 0.11111， p-value = 0.9794

4 結(jié)論

本文運用CU變換和Log變換對NHPP檢驗，效果顯著，但是針對數(shù)據(jù)的預(yù)處理問題，本文根據(jù)數(shù)據(jù)類型和R編程難易程度選擇了大區(qū)間以月為時間單位，小區(qū)間以天為時間單位的處理方法，故而出現(xiàn)了小區(qū)間中有0值出現(xiàn)的情況，所幸K-S檢驗過關(guān)。

CU變換和Log變化是用了泊松分布的兩個特性延伸推導(dǎo)出的NHPP驗證方法，其實還有直接驗證法(頻數(shù)統(tǒng)計法)，也稱為標(biāo)準(zhǔn)泊松檢驗。Lewis[11]也曾提出了比較復(fù)雜的Lewis變換對PP過程進行驗證，并給出了推導(dǎo)，所以選擇合適的、有效的數(shù)據(jù)變換方法針對NHPP的檢驗都會有一定的效益。

[1] Dohi T， Matsuoka T， Osaki S. An Infinite Server Queuing Model for Assessment of the Software Reliability[J]. Electronics and Communications in Japan， 2002， 85(3): 43-51.

[2] Huang C Y， Hung T Y. Software reliability analysis and assessment using queueing models with multiple change-points[J]. Computers & Mathematics with Applications， 2010， 60(7): 2015-2030.

[3] Brown L， Zhao L. Statistical Analysis of a Telephone Call Center: A Queueing-Science Perspective[J]. Journal of the American Statistical Association， 2005， 100(March):36-50.

[4] Kim S H， Whitt W. Choosing arrival process models for service systems: Tests of a nonhomogeneous Poisson process[J]. Naval Research Logistics， 2014， 61(1):66-90.

[5] Massey W A， Parker G A， Whitt W. Estimating the parameters of a nonhomogeneous Poisson process with linear rate[J]. Telecommunication Systems， 1996， 5(4):361-388.

[6] 茆詩松. 高等數(shù)理統(tǒng)計[M]. 北京：高等教育出版社， 1998.

[7] 徐仁佐，劉蓮君，潘志宏，等. NHPP模型擬合質(zhì)量的改進[J]. 自然科學(xué)進展:國家重點實驗室通訊， 1991(6):535-542.

[8] 楊劍鋒. 復(fù)雜數(shù)據(jù)下的軟件可靠性分析方法[D].貴陽：貴州大學(xué)，2014.

[9] 張波，張景肖.應(yīng)用隨機過程[M].北京：清華大學(xué)出版社，2007.

[10] 朱力行. Kolmogorov統(tǒng)計量的精確分布及其在Bootstrap逼近中的應(yīng)用(英文)[J]. Journal of Mathematical Research with Applications， 1991， 11(2):163-164.

[11] Lewis P A W. Some results on tests for Poisson processes[J]. Biometrika， 1965， 36(52):67-77.

(責(zé)任編輯：曾 晶)

A Test Model of NHPP: Based on CU Transformation

FAN Chaoxia1，ZHAO Ming2，YANG Jianfeng3*

(1.College of Science， Guizhou University， Guiyang 550025，China; 2.Faculty of Technology and Sustainable Development，University of Gavle， Sweden 80176;3. College of Information Engineering， Guizhou Institute of Technology， Guiyang 550003，China)

The test of NHPP and the estimation of the mean function about the BUG input data in opened source software are important sub- problems in queuing theory analysis. In this paper， the test of NHPP was converted to the test of Piecewise- Constant NHPP as the continuity and slowness of the mean function， and finally problem converted to the test of the PP on subintervals firstly， and then combined to verify the NHPP on interval. The verification of PP on subinterval is mainly through the CU transformation of the data reasonablely. Also Bugzilla 's error report data was used to analyze the model， and the result shows acceptable. Finally， the parameters of the BUG input process are calculated using the R language fitting parameters， which lays the foundation for the construction and solution of the following queuing theory model.

NHPP test; CU transformation; K- S test; R language

1000-5269(2016)06-0010-04

10.15958/j.cnki.gdxbzrb.2016.06.03

2016-10-17

貴州省科學(xué)技術(shù)基金計劃(黔科合J字[2015]2064號)；高層次人才科研啟動經(jīng)費項目(XJGC20150106)

范朝霞(1990-)，女，在讀碩士，研究方向：應(yīng)用統(tǒng)計，Email：917855385@qq.com.

*通訊作者： 楊劍鋒，Email：jfyang1@163.com.

O211.6

A