基于蒙特卡洛法的砝碼校準(zhǔn)不確定度評(píng)定

韓志國(guó)，李鎖印，馮亞南，許曉青

(中國(guó)電子科技集團(tuán)公司第十三研究所，河北 石家莊050051)

0 引言

測(cè)量是以確定量值為目的的一組操作。測(cè)量結(jié)果受人員、測(cè)量系統(tǒng)、測(cè)量程序、環(huán)境以及其它的因素影響而存在測(cè)量不確定度。測(cè)量不確定度是與測(cè)量結(jié)果相關(guān)聯(lián)的、表征合理的、賦予被測(cè)量值分散性的參數(shù)。ISO/IEC 發(fā)布的Guide98 －3 －2008 《測(cè)量不確定度表示指南》 (GUM)給出了測(cè)量不確定度評(píng)定與表示的通用方法。2008 年，計(jì)量指南聯(lián)合委員會(huì)(JCGM)在GUM 方法的基礎(chǔ)上推出了補(bǔ)充件——采用蒙特卡洛法評(píng)定測(cè)量不確定度，將蒙特卡洛法作為對(duì)GUM 方法的重要補(bǔ)充。我國(guó)新頒布的JJF1059.2 －2012《采用蒙特卡洛法評(píng)定測(cè)量不確定度》就是依據(jù)該補(bǔ)充件進(jìn)行編制的。該規(guī)范規(guī)定了用蒙特卡洛法評(píng)定與表示測(cè)量不確定度的方法，并就GUM 中未涉及的概率分布傳播的問題提供指導(dǎo)，擴(kuò)大了測(cè)量不確定度評(píng)定與表示的適用范圍，同時(shí)該規(guī)范提供了檢查GUM 法是否適用的驗(yàn)證方法。

1 蒙特卡洛法介紹

在測(cè)量不確定度的評(píng)定中采用的蒙特卡洛法(MCM)是一種通過重復(fù)采樣實(shí)現(xiàn)分布傳播的數(shù)值方法。與GUM 法利用線性化模型傳播不確定度的解析方法不同，MCM 通過對(duì)輸入量Xi的概率密度函數(shù)(PDF)離散采樣，由測(cè)量模型傳播輸入量的分布，計(jì)算獲得輸出量Y 的PDF 的離散采樣值，進(jìn)而由輸出量的離散分布數(shù)值直接獲取輸出量的最佳估計(jì)值、標(biāo)準(zhǔn)不確定度和包含區(qū)間。該輸出量的最佳估計(jì)值、標(biāo)準(zhǔn)不確定度和包含區(qū)間等特性的計(jì)算質(zhì)量隨PDF 采樣數(shù)增加可得到改善。

圖1 描述的是，由輸入量Xi(i =1，…，N)的PDF，通過模型傳播，給出輸出量Y 的PDF 的一個(gè)過程示意。圖1 列出了分別為相互獨(dú)立的正態(tài)分布、三角分布和正態(tài)分布的3 個(gè)輸入量，而輸出量的PDF 顯示為分布不對(duì)稱的情形。

圖1 輸入量獨(dú)立時(shí)分布傳播的描述

當(dāng)測(cè)量模型是線性的，且輸出量概率分布是正態(tài)分布時(shí)，GUM 方法可以提供準(zhǔn)確的不確定度和包含區(qū)間。但是，如果輸出量的結(jié)果是圖1 所示的不對(duì)稱的情形，那么采用GUM 方法進(jìn)行不確定度評(píng)定時(shí)可能就會(huì)得出不切實(shí)際的包含區(qū)間，此時(shí)采用蒙特卡洛法就是有效的替代方法。由于使用蒙特卡洛法進(jìn)行不確定度評(píng)定時(shí)需要對(duì)概率分布進(jìn)行隨機(jī)采樣，因此進(jìn)行不確定度評(píng)定時(shí)，應(yīng)利用合適的軟件進(jìn)行計(jì)算 (例如MATLAB)。

2 用蒙特卡洛法評(píng)定砝碼校準(zhǔn)不確定度

雖然GUM 法在許多情況下被認(rèn)為是非常適用的，但是確定是否滿足其所有應(yīng)用條件并不總是一件容易的事。由于MCM 的適用范圍比GUM 法更廣泛，建議用MCM 及GUM 法兩種方法同時(shí)進(jìn)行不確定度評(píng)定，并對(duì)結(jié)果進(jìn)行比較，如果比較結(jié)果較好，則GUM 法適用于此場(chǎng)合及今后足夠類似的情形，否則，應(yīng)考慮采用MCM 或者其它合適的替代方法。

下面以本實(shí)驗(yàn)室2012 年在砝碼能力驗(yàn)證中校準(zhǔn)100g 砝碼為例，介紹一下如何使用MCM 方法進(jìn)行不確定度評(píng)定。本次校準(zhǔn)使用的標(biāo)準(zhǔn)砝碼為E2等級(jí)100g砝碼，天平為TG328A 型機(jī)械天平。

2.1 建立測(cè)量模型

本次測(cè)量使用的是雙次替代法，用密度為ρr的標(biāo)準(zhǔn)砝碼R 對(duì)密度為ρt的砝碼進(jìn)行校準(zhǔn)，兩個(gè)砝碼標(biāo)稱值相同，需要在密度為ρa(bǔ)的空氣中進(jìn)行平衡配重。通常ρt和ρr不同，因此需考慮浮力的影響，應(yīng)用阿基米德定律，得到模型如下

式中:ms是密度為ρr靈敏度小砝碼質(zhì)量，用來測(cè)量機(jī)械天平的分度值。ΔI 為被測(cè)砝碼與標(biāo)準(zhǔn)砝碼標(biāo)尺讀數(shù)偏差;ΔIs為放上靈敏度小砝碼后，標(biāo)尺偏離的讀數(shù)。以折算質(zhì)量mct，mcr，mcs為變量，模型(1)轉(zhuǎn)變?yōu)?/p>

將公式(2)整理得

設(shè)Δmc= mct－ mnom為mct與標(biāo)稱質(zhì)量mnom=100 g的偏差。因此本例中使用的模型為

2.2 確定輸入量的概率密度分布函數(shù)(PDF)

可根據(jù)一系列測(cè)量值的分析，或根據(jù)某些歷史數(shù)據(jù)、校準(zhǔn)數(shù)據(jù)和專家判斷之類的信息所得到的科學(xué)判斷，為各輸入量設(shè)定PDF。各輸入量服從的分布、最佳估計(jì)值以及標(biāo)準(zhǔn)不確定度均列入表1。其中ρa(bǔ)0=1.2 kg/m3，是約定空氣密度，不包含不確定度。

表1 關(guān)于質(zhì)量校準(zhǔn)模型的輸入量及其概率密度函數(shù)

表1 中，標(biāo)準(zhǔn)砝碼的相關(guān)信息是通過檢定證書和JJG99 －2006 《砝碼檢定規(guī)程》得到。ΔI 與ΔIs的相關(guān)信息是通過使用雙次替代法重復(fù)測(cè)量10 次得到。砝碼的密度信息是通過JJG99 －2006 《砝碼檢定規(guī)程》中給出的典型材料密度表中查表得到。空氣密度ρa(bǔ)的信息是由JJG99 －2006 《砝碼檢定規(guī)程》中得到。

2.3 蒙特卡洛試驗(yàn)次數(shù)

在規(guī)定的數(shù)值容差下MCM 所提供的結(jié)果所需的試驗(yàn)次數(shù)跟輸出量的PDF “形狀”及包含概率有關(guān)，應(yīng)合理選擇蒙特卡洛試驗(yàn)次數(shù)及樣本量的大小M，一般情況下取M=106，這樣通常會(huì)為輸出量提供95%包含區(qū)間，該包含區(qū)間長(zhǎng)度被修約為1 或2 位有效十進(jìn)制數(shù)字。

2.4 程序

本程序是使用MATLAB 進(jìn)行編制的，具體如下:

clc;clear; %清除MATLAB 中的變量。

M=1000000; %設(shè)定蒙特卡洛試驗(yàn)次數(shù)。

m0 =100; %砝碼標(biāo)稱質(zhì)量(單位:g)。

dr=0.04; %標(biāo)準(zhǔn)砝碼折算質(zhì)量修正值(單位:mg)。

ur=0.033; %標(biāo)準(zhǔn)砝碼不確定度(單位:mg)。

mr=normrnd ((m0* 10^3 +dr)，ur，1，M); %標(biāo)準(zhǔn)砝碼采樣值(單位:mg)。

ms=normrnd (10，0.002，1，M); %靈敏度小砝碼采樣值(單位:mg)。

dI= normrnd (4.2，0.2，1，M); % 讀數(shù)重復(fù)性 (單位:格)。

dIs = normrnd (99.6，0.3，1，M); % 靈敏度讀數(shù)重復(fù)性(單位:格)。

rou0 =normrnd (8，0.07，1，M); %標(biāo)準(zhǔn)砝碼材料密度(g/cm3)。

rou1 =normrnd (8，0.07，1，M); %被檢砝碼材料密度(g/cm3)。

ra=0.24* rand (1，M) +1.08; %空氣密度(mg/cm3)。

r0 =1.2; %標(biāo)準(zhǔn)空氣密度(mg/cm3)。

dm= (mr+dI.* ms. /dIs). * (rou0* 10^3 －ra). * (rou1*10^3 －r0). / ((rou1* 10^3 －ra).* (rou0* 10^3 －r0)) －m0

* 10^3;%雙次交換衡量法數(shù)學(xué)模型。

s=std (dm); %輸出量標(biāo)準(zhǔn)不確定度。

avg=mean (dm); %輸出量估計(jì)值。

y=sort (dm); %將試驗(yàn)結(jié)果從小到大進(jìn)行排列。

y_ low=prctile (y，2.5); %包含概率為95%時(shí)，包含區(qū)間的下限值。

y_ high=prctile (y，97.5); %包含概率為95%時(shí)，包含區(qū)間的上限值。

2.5 報(bào)告結(jié)果

MATLAB 程序運(yùn)行得到如下結(jié)果:

被校砝碼折算質(zhì)量偏差為0.46 mg;標(biāo)準(zhǔn)不確定度為0.04 mg;包含概率為95%;包含區(qū)間為［0.38 mg，0.54 mg］。圖2 給出了上述程序得到的輸出量的近似概率分布。圖中，柱狀圖是輸出量的近似概率分布，實(shí)線是正態(tài)分布曲線，可見輸出量近似服從正態(tài)分布。

圖2 輸出量近似概率分布

3 蒙特卡洛法與其它不確定度評(píng)定方法的比較

3.1 GUM 法評(píng)定砝碼校準(zhǔn)不確定度

由于使用偏導(dǎo)的解析計(jì)算對(duì)公式(4)求偏導(dǎo)過于復(fù)雜，因此將公式(4)作近似處理，得

對(duì)上述模型求偏導(dǎo)，得到各輸入量靈敏系數(shù)如下:

則合成標(biāo)準(zhǔn)不確定度為

將表1 中輸入量的相關(guān)信息帶入上式，計(jì)算得uc=0.04mg，取包含概率為95%，得到如下結(jié)果:

被校砝碼折算質(zhì)量偏差為0.46 mg;標(biāo)準(zhǔn)不確定度為0.04 mg;包含概率為95%;包含區(qū)間為［0.38 mg，0.54 mg］。

3.2 JJG99 －2006 《砝碼檢定規(guī)程》中的方法評(píng)定砝碼校準(zhǔn)不確定度

規(guī)程中認(rèn)為砝碼的校準(zhǔn)不確定度由以下幾個(gè)方面引入:衡量過程，標(biāo)準(zhǔn)砝碼，空氣浮力修正，天平。各不確定度分量列入表2。

表2 測(cè)量不確定度分量一覽表

取包含概率為95%，得到如下結(jié)果:

被校砝碼折算質(zhì)量偏差為0.46 mg;標(biāo)準(zhǔn)不確定度為0.05 mg;包含概率為95%;包含區(qū)間為［0.36 mg，0.56 mg］。

三種方法得到的結(jié)果見表3。

表3 砝碼校準(zhǔn)計(jì)算結(jié)果

由三種方法得到的比較結(jié)果可知，MCM 與GUM得到的結(jié)果一致，而規(guī)程給出的不確定度略顯保守。這主要是由于某些不確定度分量被重復(fù)計(jì)入和計(jì)量人員對(duì)不確定度分量的不同理解程度而造成的。按照J(rèn)JF1059.2 －2012 《用蒙特卡洛法評(píng)定測(cè)量不確定度》中使用MCM 驗(yàn)證GUM 的要求，GUM 法通過了驗(yàn)證，在以后進(jìn)行砝碼校準(zhǔn)不確定度評(píng)定時(shí)，可以使用GUM法進(jìn)行評(píng)定。

4 結(jié)論

MCM 尤其適用于以下兩種情況:不宜對(duì)測(cè)量模型進(jìn)行線性化等近似的場(chǎng)合，在這種情況下，按GUM 確定輸出量的估計(jì)值和標(biāo)準(zhǔn)不確定度可能會(huì)變得不可靠;輸出量的概率密度函數(shù)較大程度地偏離正態(tài)分布或t 分布，例如分布明顯不對(duì)稱的場(chǎng)合，在這種情況下，可能會(huì)導(dǎo)致對(duì)包含區(qū)間或擴(kuò)展不確定度的估計(jì)不切實(shí)際。另外，如果GUM 法通過了MCM 的驗(yàn)證，即GUM 法若明顯適用，則GUM 法依然是不確定度評(píng)定的主要方法。在實(shí)際工作中可以根據(jù)實(shí)際情況選用合適的方法進(jìn)行測(cè)量不確定度的評(píng)定。

［1］國(guó)家質(zhì)量監(jiān)督檢驗(yàn)檢疫總局.JJF1059.2 －2012 用蒙特卡洛法評(píng)定測(cè)量不確定度［S］. 北京:中國(guó)計(jì)量出版社，2012.

［2］Evaluation of measurement data — Supplement 1 to the“Guide to the expression of uncertainty in measurement”— Propagation of distributions using aMonte Carlo method ［S］. ISO/IEC，2008.

［3］國(guó)家質(zhì)量監(jiān)督檢驗(yàn)檢疫總局.JJG99 －2006 砝碼檢定規(guī)程［S］. 北京:中國(guó)計(jì)量出版社，2006.

［4］陳懷艷，曹蕓，韓潔. 基于蒙特卡羅法的測(cè)量不確定度評(píng)定［J］. 電子測(cè)量與儀器學(xué)報(bào)，2007，25.

［5］倪育才. 實(shí)用測(cè)量不確定度評(píng)定［M］. 北京:中國(guó)計(jì)量出版社，2010.