宋明順， 張俊亮， 方興華

（中國(guó)計(jì)量學(xué)院，浙江杭州 310018）

基于分位數(shù)函數(shù)和Bootstrap分布的不確定度評(píng)定研究

宋明順， 張俊亮， 方興華

（中國(guó)計(jì)量學(xué)院，浙江杭州 310018）

最大信息熵方法是基于概率分布評(píng)定測(cè)量不確度的主要方法之一。其所依賴的高階矩需要較大樣本的測(cè)量數(shù)據(jù)，而校準(zhǔn)／檢測(cè)實(shí)驗(yàn)室的測(cè)量一般為小樣本，故用最大熵方法評(píng)定小樣本測(cè)量不確定度缺乏一定的可靠性。提出了基于分位數(shù)函數(shù)和概率權(quán)重矩作為約束條件的最大信息熵不確定度評(píng)定法，把矩的計(jì)算從高次降為一次，并結(jié)合遺傳算法求解概率分布，用Bootstrap分布估計(jì)擴(kuò)展不確定度和包含區(qū)間，解決了由分位數(shù)區(qū)間估計(jì)分布不對(duì)稱所致的復(fù)雜計(jì)算問(wèn)題。

計(jì)量學(xué)；不確定度評(píng)定；分位數(shù)；概率權(quán)重矩；最大熵原理；Bootstrap分布

1 引 言

最大信息熵方法（PME）是基于概率分布評(píng)定測(cè)量不確定度的常用方法［1］。利用PME確定概率密度函數(shù)（PDF）是在已知先驗(yàn)信息約束條件下通過(guò)使信息熵最大來(lái)得到的，即：

且式（1）有以下約束條件：

式中，H（x）為熵函數(shù)，f（x）為隨機(jī)變量的概率密度函數(shù)，mi為第i階原點(diǎn)矩，通常mi由測(cè)量樣本給出，即：

通過(guò)求條件極值可得到測(cè)量值x的概率密度函數(shù)f（x），即：

式中，λi為拉格朗日乘子。

λi的計(jì)算依賴于樣本矩，樣本矩的計(jì)算又依賴于給定的測(cè)量樣本。理論上，樣本矩階次越高，所包含的信息量越大，求解的PDF越接近于真實(shí)分布。但高樣本矩階次帶來(lái)兩個(gè)問(wèn)題：一是λi的計(jì)算問(wèn)題，在實(shí)際的計(jì)算過(guò)程中，樣本矩越高，求解方程的非線性化程度越嚴(yán)重，求解就越復(fù)雜，計(jì)算誤差就比較大，如傳統(tǒng)λi的計(jì)算方法L-M法［2］和擬牛頓法［3］，均不能很好地解決上述問(wèn)題；二是樣本量的問(wèn)題，樣本矩階次越高，所需要的測(cè)量樣本量越大，但在實(shí)際校準(zhǔn)／檢測(cè)工作中，樣本量通常都是較少的，屬于小樣本測(cè)量（樣本容量小于50），在小樣本測(cè)量條件下，用PME求解PDF的可靠性大大下降。M.D.Pandey指出，在計(jì)算分位數(shù)的樣本矩時(shí)，階次取2～4偏差接近于0，效果比較好。

本文提出用概率權(quán)重矩（PWMs）為約束條件的分位數(shù)最大熵方法來(lái)求解小樣本條件下的PDF，然后得到分位數(shù)的點(diǎn)估計(jì)，即可以由PDF求出最佳估計(jì)值及標(biāo)準(zhǔn)差，并用Bootstrap方法得到給定包含概率下的擴(kuò)展不確定度及包含區(qū)間。

2 概率權(quán)重矩及分位數(shù)函數(shù)

針對(duì)小樣本測(cè)量中PME評(píng)定測(cè)量不確定度存在的高階矩問(wèn)題，通過(guò)分位數(shù)函數(shù)及PWMs將矩的高階次轉(zhuǎn)化為低階次，使其適用于小樣本測(cè)量，降低了計(jì)算的復(fù)雜性，并且保證測(cè)量不確定度評(píng)定的可靠性。

2.1 概率權(quán)重矩（PWM s）

現(xiàn)在常用的PWMs的定義是由Greewood提出的［4］：

1

其中，以下兩種PWM形式是簡(jiǎn)單而常用的。第1種：

比較式（8）和式（7）可以知道，當(dāng)樣本容量為k時(shí)，αk，βk分別與最小、最大順序統(tǒng)計(jì)量的期望有關(guān)，其關(guān)系式可表示為：

對(duì)于一個(gè)有n個(gè)順序統(tǒng)計(jì)量的樣本，ak，bk分別為αk和βk的無(wú)偏估計(jì)：

比較式（10）和式（4）可以看出，作為樣本信息的提取方式，樣本原點(diǎn)矩（mi）是求樣本中每個(gè)數(shù)的冪的均值，即對(duì)樣本數(shù)據(jù)進(jìn)行非線性處理，它對(duì)樣本的取值范圍很敏感，特別是在小樣本條件下，高階矩不能夠從小樣本中準(zhǔn)確地得到；PWMs是通過(guò)對(duì)樣本信息作線性計(jì)算來(lái)提取樣本信息的，它的優(yōu)點(diǎn)在于對(duì)樣本的邊界值不敏感，高階PWMs能從樣本中比較準(zhǔn)確地得出。因此，在小樣本測(cè)量的條件下，高階PWMs約束下得到的概率分布更為客觀、合理。

2.2 分位數(shù)函數(shù)的矩

PWMs是數(shù)據(jù)順序統(tǒng)計(jì)量的期望值，也可以解釋為非負(fù)隨機(jī)變量分位數(shù)函數(shù)的矩，分位數(shù)函數(shù)是分布函數(shù)的反函數(shù)。在非負(fù)隨機(jī)變量的情況下，常用βk作為分位數(shù)函數(shù)的矩，這也是計(jì)算中常用的一種情況，由傳統(tǒng)統(tǒng)計(jì)矩的定義知，k（k≥1）階矩的定義為：

把式（3）與式（12）比較，x（p）類似于f（x），而概率p類似于隨機(jī)變量x。分位數(shù)函數(shù)為一個(gè)連續(xù)函數(shù)，而且是非負(fù)，由于0≤x（p）≤∞，0≤p≤1。設(shè)：

根據(jù)式（12）和式（14），（βk／β0）可看作分位數(shù)函數(shù)x（p）的k階矩。在實(shí)際的計(jì)算中，如果對(duì)測(cè)量數(shù)據(jù)樣本作歸一化處理，則β0＝α0＝1。對(duì)于這樣的樣本，βk便為分位數(shù)函數(shù)的階矩。

3 測(cè)量不確定度的評(píng)定

現(xiàn)在來(lái)研究基于概率權(quán)重矩最大熵原理下的分位數(shù)函數(shù)測(cè)量不確定度的評(píng)定方法。

隨機(jī)變量x的分位數(shù)函數(shù)x（p）的熵定義為：

式中，bk為βk的無(wú)偏估計(jì)，m為考慮的概率權(quán)重矩的最高階數(shù)。

考慮到其約束條件，分位數(shù)的最大熵函數(shù)可以寫為：

對(duì)式（17）求導(dǎo)后得到PWMs為約束條件下用PME確定的分位數(shù)函數(shù)的表達(dá)式，即：

式中，拉格朗日乘子可以通過(guò)式（16）的約束條件來(lái)確定，即把式（18）代入式（16），得

對(duì)式（19）變形，且令殘差rk為：

殘差rk表征著拉格朗日乘子估計(jì)的水平，殘差越小，表明估計(jì)越準(zhǔn)確。為求得拉格朗日乘子λk的值，需要對(duì)殘差優(yōu)化處理，故令

式中，ε無(wú)限接近于0，以此為優(yōu)化過(guò)程的極值。

式（21）的求解，本文采用遺傳算法予以解決［5］。通過(guò)分位數(shù)函數(shù)可以對(duì)測(cè)量結(jié)果進(jìn)行不確定度評(píng)定。這種方法在小樣本條件下不確定度評(píng)定比傳統(tǒng)樣本原點(diǎn)矩估計(jì)的評(píng)定效果更好。

分位數(shù)函數(shù)是分布函數(shù)的反函數(shù)。因此，通過(guò)對(duì)分位數(shù)函數(shù)求反函數(shù)，就可以得到隨機(jī)變量的PDF，即：F（x）＝x－1（p），也可以直接通過(guò)分位數(shù)函數(shù)確定測(cè)量結(jié)果的最佳估計(jì)值和標(biāo)準(zhǔn)不確定度u：

在給定包含概率p的條件下確定包含區(qū)間Up時(shí)，傳統(tǒng)方法是先判斷其分布是否對(duì)稱，若對(duì)稱，可根據(jù)式（24）來(lái)求解其指定概率下的包含區(qū)間：

Up＝［x（（1－p）／2），x（（1＋p）／2）］（24）

若分布函數(shù)不對(duì)稱，通過(guò)文獻(xiàn)［6，7］可發(fā)現(xiàn)，Bootstrap在區(qū)間估計(jì)方面應(yīng)用很廣，效果也佳，故本文采取Bootstrap估計(jì)來(lái)確定其包含區(qū)間。

4 Bootstrap區(qū)間估計(jì)

Bootstrap方法是一種再抽樣統(tǒng)計(jì)方法。其基本思想是用已知的經(jīng)驗(yàn)分布函數(shù)代替未知母體的分布函數(shù)，并通過(guò)再抽樣技術(shù)獲得統(tǒng)計(jì)量的Bootstrap分布，在此基礎(chǔ)上就可進(jìn)行統(tǒng)計(jì)推斷，如求解統(tǒng)計(jì)量的標(biāo)準(zhǔn)差和包含區(qū)間等。

Bootstrap再抽樣方法很多，本文采取的是非參數(shù)方法［8］。其主要思想是利用經(jīng)驗(yàn)分布函數(shù)代替總體分布函數(shù)，從經(jīng)驗(yàn)分布函數(shù)中抽取相同容量的樣本來(lái)估計(jì)參數(shù)的抽樣分布。從經(jīng)驗(yàn)分布函數(shù)中抽樣與對(duì)原始樣本作有放回的抽樣是等效的。

Bootstrap的區(qū)間估計(jì)基本過(guò)程是首先抽取一組Bootstrap樣本，用遺傳算法計(jì)算相應(yīng)的統(tǒng)計(jì)量＝θ（，，…，Xn），如分位數(shù)x（p），稱為Bootstrap估計(jì)量。重復(fù)再抽樣B次，就可以獲得B個(gè)Bootstrap估計(jì)量：

這B個(gè)Bootstrap估計(jì)量就提供了統(tǒng)計(jì)量的Bootstrap分布信息。有了統(tǒng)計(jì)量分布，就可以在此基礎(chǔ)上進(jìn)行推斷統(tǒng)計(jì)量的包含區(qū)間。本文采用經(jīng)驗(yàn)百分位數(shù)方法。Timmerman指出，經(jīng)驗(yàn)百分位數(shù)方法的估計(jì)誤差將以1的速度趨近于0。

將B個(gè)Bootstrap估計(jì)量)θ*i由小到大進(jìn)行排序，由經(jīng)驗(yàn)百分位數(shù)方法得到包含概率為1－α的包含區(qū)間估計(jì)為：

5 實(shí)例分析

已知大型工具顯微鏡3倍物鏡長(zhǎng)度重復(fù)測(cè)量樣本，參見(jiàn)表1［9］。

表1 顯微鏡物鏡長(zhǎng)度測(cè)量數(shù)據(jù)mm

根據(jù)樣本對(duì)測(cè)量結(jié)果進(jìn)行不確定度評(píng)定。

由于樣本容量為20，屬于小樣本情況，因此用PWMs為約束條件確定分位數(shù)函數(shù)，通過(guò)分位數(shù)函數(shù)來(lái)對(duì)樣本進(jìn)行不確定度評(píng)定，通過(guò)Bootstrap來(lái)對(duì)區(qū)間進(jìn)行估計(jì)。此實(shí)例的計(jì)算，均采用Matlab程序完成。

（1）分位數(shù)函數(shù)的確定

在確定分位數(shù)的矩時(shí)，階次越高，精度越高，但其復(fù)雜性也越大。故根據(jù)樣本，利用式（10），計(jì)算得到分位數(shù)函數(shù)的前4階PWMs，分別為：

然后利用PME確定分位數(shù)函數(shù)表達(dá)式，如式（27）所示，對(duì)分位數(shù)函數(shù)求反函數(shù)即得到樣本的分布函數(shù)。

（2）根據(jù)分位數(shù)函數(shù)進(jìn)行不確定度評(píng)定

在得到分位數(shù)函數(shù)后，通過(guò)式（22）和式（23）來(lái)確定樣本的最佳估計(jì)值x＾1和標(biāo)準(zhǔn)不確定度u1：

取Bootstrap再抽樣次數(shù)B＝200，即對(duì)20個(gè)樣本進(jìn)行再抽樣200次，得到200組樣本量20的樣本，分別對(duì)這200組樣本進(jìn)行基于概率權(quán)重矩的分位數(shù)函數(shù)計(jì)算，得到200個(gè)分位數(shù)函數(shù)xi（p），然后取對(duì)應(yīng)p＝1－α的分位數(shù)值xi（α），將這200個(gè)新統(tǒng)計(jì)量xi（α）排序，依據(jù)式（25）得到包含概率1－α的包含區(qū)間，即在給定包含概率95%的條件下，得包含區(qū)間：

（3）評(píng)定結(jié)果比較分析

在GUM評(píng)定方法中［10］，最佳估計(jì)值和標(biāo)準(zhǔn)不確定度u2分別用樣本均值和多次測(cè)量的標(biāo)準(zhǔn)差來(lái)表示。

在給定包含概率95%的條件下，GUM中的方法按正態(tài)分布對(duì)樣本進(jìn)行擴(kuò)展不確定度評(píng)定，取包含因子k＝2，則此時(shí)的包含區(qū)間為］，最后計(jì)算得到包含概率95%下的包含區(qū)間為：［69.666 4，70.341 6］mm。

表2給出了GUM方法和用PWMs確定分位數(shù)函數(shù)進(jìn)行不確定度評(píng)定的結(jié)果的對(duì)比。

表2 與GUM中不確定度評(píng)定結(jié)果比較mm

根據(jù)式（27）得到用GUM中正態(tài)分布假設(shè)的95%包含概率下包含區(qū)間兩端點(diǎn)的概率為：p（69.666 4）＝0.136 5，p（70.341 6）＝0.959 7，也就是說(shuō)，在此區(qū)間下包含概率p＝96%，超出了假設(shè)條件95%。因此在此例中，用正態(tài)分布假設(shè)來(lái)進(jìn)行不確定度評(píng)定是不合理的。GUM方法是假設(shè)測(cè)量值服從正態(tài)分布的，只有在測(cè)量值真正服從正態(tài)分布時(shí)測(cè)量不確定度的結(jié)果才可靠；而PWMs確定分位數(shù)函數(shù)的方法沒(méi)有任何概率分布的假設(shè)，無(wú)論測(cè)量值服從什么概率分布，都適用。

6 結(jié) 論

當(dāng)測(cè)量數(shù)據(jù)樣本為小樣本（樣本容量小于50）時(shí)，在PWMs約束條件下，用PME確定樣本分位數(shù)函數(shù)，通過(guò)分位數(shù)函數(shù)來(lái)確定測(cè)量結(jié)果的最佳估計(jì)值和標(biāo)準(zhǔn)不確定度，用Bootstrap區(qū)間估計(jì)方法確定擴(kuò)展不確定度及包含區(qū)間，不依賴數(shù)據(jù)的分布形式，有更廣泛的應(yīng)用范圍，比GUM評(píng)定測(cè)量不確定度的方法要好。

［1］ Jaynes E T.Information theory and statisticalmechanics［J］.Phys Rev，1957，104（6）：620－630.

［2］ D′Antona G.Expanded Uncertainty and Coverage Factor Computation by High Order Moment Analysis［C］／／IMTC2004－Instrumentation and Measurement Technology Conference.Italy：Institute of Electrical and Electronics Engineers Inc，2004，234－238.

［3］ Pandey M D.Directestimation of quantile functions using the maximum entropy principle［J］.Structural Safety，2000，22（1）：61－79.

［4］ Greenwood J A，Landwehr J M，Matalas N C，et al. Probability weighted moments：definition and relation to parameters of several distributions expressible in inverse form［J］.Water Resourc Res，1979，15（5）：1049－1054.

［5］ Fang X H，Song M S.Estimation of Maximum－entropy Distribution Based on Genetic Algorithms in Evaluation of the Measurement Uncertainty［C］／／2010 Second WRI Global Congress on Intelligent Systems.Wuhan，China：Computer Society，2010，292－297.

［6］ 方興華.在測(cè)量不確定度評(píng)定中用PME確定PDF的研究［D］.杭州：中國(guó)計(jì)量學(xué)院.2009.

［7］ Efron B.Bootstrap methods：another look at the jackknife［J］.The Annalsof Statistics，1979，7（1）：1－26.

［8］ Efron B.Nonparametric estimates of standard error and confidence intervals［J］.The Canadian Journal of Statistics，1981，9（2）：137－172.

［9］ Timmerman M E，TerBrak C J F.Bootstrap confidence intervals for principal response curves［J］.Computational Statistic＆Data Analysis，2008，52（4）：1837－1849.

［10］ 孫永厚，周洪彪，黃美發(fā).基于最大熵的測(cè)量不確定度的貝葉斯評(píng)估方法［J］.統(tǒng)計(jì)與決策，2008，（12）：143－145.

［11］ BIPM，IEC，IFCC，ISO，IUPAC，IUPAP，OIML.Guide to the expression of uncertainty in measurement［S］. 2008.

Uncertainty Evaluation Based on Quantile Function and Bootstrap

SONG Ming-shun， ZHANG Jun-liang， FANG Xing-hua
（China Jiliang University，Hangzhou，Zhejiang 310018，China）

Maximum entropy is one of the chief methods in the evaluation of measurement uncertainty based on probability distribution，and the higher moments it depended on a larger sample of measurement data.However，measurement in calibration and testing laboratories is generally a small sample survey，so evaluation of measurement uncertainty based on the maximum entropy method for small samples lacks a certain degree of reliability.A method for uncertainty evaluation based on quantile function and probabilityweightingmoment is put forward as the constraint condition ofmaximum information entropymethod.By thisway，the high ordermoment is fell down for amoment，the probability distribution is solved with the combination of genetic algorithm，and the complicated calculation problem arising from quantile interval estimation of asymmetric distribution ismanaged aboutwith Bootstrap distribution estimation for expanded uncertainty and coverage interval.

Metrology；Uncertainty evaluation；Quantile functions；PWMs；Maximum entropy principle；Bootstrap distribution

TB9

A

1000-1158（2014）03-0300-05

10．3969／j．issn．1000-1158．2014．03．22

2013-02-21；

2013-07-23

國(guó)家自然科學(xué)基金（71071147）

宋明順（1961－），男，山東蓬萊人，中國(guó)計(jì)量學(xué)院教授，主要從事計(jì)量基礎(chǔ)、標(biāo)準(zhǔn)化和質(zhì)量管理的研究。smsqm＠cjlu.edu.cn