标准不确定度A类评定中极差法的深入讨论

标准不确定度A类评定中极差法的深入讨论
陈凌峰
【摘要】JJF 1059.1-2012《测量不确定度评定与表示》与GUM的区别之一是在标准不确定度的A类评定中引入了极差法.假设总体分别服从正态分布和均匀分布,则总体标准差的极差估计量,以及用于实际计算的极差系数可以从样本极差的分布函数导出.理论分析表明:虽然用极差法估计的总体标准差是无偏的,但是估计的总体方差偏大,这将导致最终测量结果的合成标准不确定度偏大.同时JJF 1059.1中仅提供了总体接近正态分布时的极差系数,并不适用于所有情况.作为比较,不论总体分布如何,使用贝塞尔公式估计的总体方差总是无偏的,不会给测量结果的合成标准不确定度带来原理性误差.由于极差法存在概率统计学上的原理性误差以及适用性限制,建议在标准不确定度A类评定中应审慎使用极差法.
【期刊名称】《计量学报》
【年(卷),期】2019(040)002
【总页数】6页(P347-352)
【关键词】计量学;标准不确定度;极差法;无偏估计;正态分布;均匀分布
【作者】陈凌峰
【作者单位】北京理工大学光电学院,北京100081
【正文语种】中文
【中图分类】TB9
1 引言
国家计量技术规范JJF 1059.1—2012《测量不确定度评定与表示》中推荐了两种基本的标准不确定度A类评定方法，即贝塞尔公式法和极差法。

其中贝塞尔公式法对输入量X的分布没有限制，但极差法的应用前提是输入量X接近服从正态分布[1]。

在重复性或复现性条件下，对被测量X进行n次独立重复观测，测得值分别为
x1,x2,…,xn，n个观测值的算术平均值为则单次测量结果的实验标准差sn可用贝塞尔公式计算：
其中在重复性或复现性条件下，对被测量X进行n次独立重复观测，若n个测得值x1,x2,…,xn中的最大值与最小值之差为Dn，在被测量X接近正态分布的前提下，单次测量结果的实验标准差s可近似表示为：
(1)
式中：系数C称为极差系数，其与测量次数n有关。

表1列出了JJF 1059.1—2012给出的不同测量次数n时的极差系数C。

JJF 1059.1是等效采用了2008版GUM[2]，但是在GUM原文中我们却找不到任何关于极差法的描述。

追寻极差法的由来，在JJF 1059—1999[3]中已经有了。

其主要起草人李慎安介绍极差法来自于JJG 1027—1991《测量误差及数据处理》，并建议一般测量次数n不小于4为好[4]。

从估计标准偏差的相对误差考虑，倪育才建议使用极差法时的测量次数n为4～9次为好[5]。

表1 极差系数Tab.1 The coefficient of the range
methodn23456789C1.131.692.062.332.532.702.852.97
2 极差法原理及极差系数的计算
设随机变量ξ 的分布函数为F(x)，概率密度函数为f(x)，ξ1,ξ2,…,ξn为其样本，将n个样本值从小到大排序：其中令则称Dn为样本极差。

记Dn的分布函数为FDn(r)，概率密度函数为fDn(r)，显然r≥0。

用事件A表示“样本ξ1,ξ2,…,ξn中有一个落在区间(x,x+dx)，而另n-1个同时落在区间(x,x+r)”，则事件A发生的概率为[6]：
其中(x,x+r)为样本观测值的范围，r为极差Dn观测值的上限。

于是有：
FDn(r)=P{Dn<r}= n[F(x+r)-F(x)]n-1f(x)dx
(2)
F(x)]n-2f(x+r)f(x)dx.
(3)
于是样本极差Dn的数学期望为：
E(Dn)=rfDn(r)dr=n(n-1)·
r[F(x+r)-F(x)]n-2f(x+r)f(x)dxdr
(4)
若ξ服从标准正态分布N(0,1)，即式(4)成为：
cn=E(Dn)=rfDn(r)dr
(5)
式(5)中的Φ(x)为标准正态分布的分布积分函数，
式(5)不易得到解析解，但在n确定的情况下，cn的值可用数值积分计算得出。

以
n=5为例，其MATLAB计算程序如下：
n=5; %设定样本数n
syms x r %设定积分变量
s=int(int(exp((-(x+r)^2-x^2)/2)*r*(normcdf(x+r)-normcdf(x))^(n-2)*(n*(n-1))/(2*pi),x,-inf,inf),r,0,12);
%设定积分表达式
Cn=vpa(s) %取高精度计算结果
上述计算程序对无穷积分区间做了截断处理，即将极差r的上限设为12。

此时的截断误差Δn可表示为：
例如，当n=5时，利用数值积分可算得c5的截断误差Δ5<2.6×10-14。

在忽略截断误差的情况下，n≤9的cn计算结果(精确到10-6)如表2所示。

对照表1可见，实际JJF 1059.1—2012中的极差系数C是将cn修正至0.01的结果。

若利用式(5)算得样本极差Dn的期望值cn，则样本极差Dn的方差计算如下：
(r-cn)2[Φ(x+r)-
(6)
对不同的样本数的数值积分结果同样列于表2中。

表2 X～N(0,1)时，样本极差的期望与方差Tab.2 The expectation and variance of the sample range when X～N(0,1)结果n23456789样本极差的期望
cn1.128 3791.692 5682.058 7502.325 9282.534 4122.704 3562.847 2002.970 026样本极差的方差d2n0.726 7600.789 1980.774 0620.746 6380.719 1710.694 2310.672 1240.652 596
3 有关极差法的深入讨论
3.1 极差法对测量不确定度评定结果的影响
若ξ服从正态分布N(μ,σ)，样本ξ1,ξ2,…,ξn的极差为Dn，做标准化变换则ξ′服从标准正态分布的极差为
即：
于是：
(7)
(8)
(9)
由式(7)、式(8)可知，若总体服从正态分布N(μ,σ)，则统计量是总体标准差σ的无偏估计量，但统计量(Dn/cn)2并非总体方差σ2的无偏估计量。

由于式(8)中则E(Dn/cn)2>σ2，统计量(Dn/cn)2实际上是总体方差σ2的偏大估计量。

利用表2的数据，对不同的样本数的计算结果如表3所示。

可见，当n=4时，用极差法估计的总体方差约偏大18.3%;当n=9时，约偏大7.4%。

表3 统计量(Dn/cn)2对总体方差σ2的偏大倍数Tab.3 The multiplier factor of
σ2 relative to (Dn/cn)2样本数n23456789d2n/c2n+11.570 81.275 51.182 61.138 01.112 01.094 91.082 91.074 0
不论总体ξ服从何种分布，只要期望E(ξ)=μ，方差D(ξ)=σ2皆有限，用贝塞尔公式计算的样本标准差sn并非总体标准差σ的无偏估计，但却是总体方差σ2的无
偏估计[6,7]，即这种情况与极差法刚好相反：若总体服从正态分布N(μ,σ)，则
E(Dn/cn)=σ，E(Dn/cn)2≠σ2。

由于在大部分测量不确定度评定中，标准不确定度A类评定结果仅是其中的一个
或几个分量。

通常还需要与B类评定的分量合成，才能得到测量结果的合成标准
不确定度。

GUM的第5.1.2节明确指出，合成标准不确定度是合成方差的正平方
根[2]。

JJF 1059—1999的6.1节对此也有明确描述，但JJF 1059.1—2012的
4.4节“合成标准不确定度的计算”却没有对这一关键点进行明确。

由于即
(Dn/cn)2是总体方差σ2的偏大估计，在标准不确定度的A类评定中采用极差法，将导致测量结果的合成标准不确定度偏大，并进一步导致测量结果的扩展不确定度偏大。

而使用贝塞尔公式则不会给合成标准不确定度的计算带来原理性误差。

需要说明的是：通常输入量的最佳估计值取n次独立观测值的算术平均值，在标准不
确定度中引入倍乘因子并不影响上述结论。

极差法对合成标准不确定度的影响程度取决于A类分量和B类分量的相对大小。

考虑合成标准不确定度uc仅包含一个A类分量，一个B类分量的简单情形，即
由表3数据可以算出，若A类分量不占优势(uA≤uB)，极差法对测量结果的合成
标准不确定度的影响有限：当n=4时，用极差法计算uA导致合成标准不确定度
uc的最大相对误差约+4.4%。

当n=5时，导致uc的最大相对误差约+3.4%。

当
n=9时，导致uc的最大相对误差约+1.8%。

若A类分量占优势(uA≥uB)，当n=4时，用极差法计算uA将导致合成标准不确
定度uc偏大4.4%～8.7%。

当n=5时，用极差法计算uA将导致合成标准不确定度uc偏大3.4%～6.7%。

当n=9时，用极差法计算uA将导致合成标准不确定度uc偏大1.8%～3.6%。

3.2 极差分布参数的蒙特卡洛模拟
在输入量总体服从正态分布N(μ,σ)的情况下，(Dn/cn)2是总体方差σ2的偏大估计。

虽然表3中高精度的数值积分计算结果已经说明了这一点，但我们还是希望
使用其他方法来验证该结论。

可选的验证方法是蒙特卡洛模拟[8,9]，为简化计算，假设输入量的总体服从标准正态分布N(0,1)。

由式(7)，
(10)
由式(10)，我们只需要验证不同样本数n时，样本极差统计量的数学期望为即可。

以n=5为例，其蒙特卡洛模拟的MATLAB程序如下：
n=5; %设定样本数n
cycle=2000000; %设定蒙特卡洛模拟次数为200万次
data=zeros(cycle,1); %数组初始化
for i=1:cycle
m=normrnd(0,1,1,n); %生成n个服从标准正态分布N(0,1)的随机数
data(i)=(max(m)-min(m))^2; %计算n个样本的极差的平方
end
E_Dn_2=mean(data); %计算200万个数据的平均值
表4列出了200万次蒙特卡洛模拟的的数学期望值可见与样本极差分布参数的数
值积分计算结果极为接近。

所以蒙特
卡洛模拟结果不仅验证了表2中计算的样本极差分布参数的正确性，而且间接验
证了(Dn/cn)2是总体方差σ2的偏大估计量这一结论。

表4 X～N(μ,σ)时，蒙特卡洛模拟结果Tab.4 The Monte Carlo simulation result when X～N(μ,σ)结果23456789数值积分计算结果cn1.128 3791.692 5682.058 7502.325 9282.534 4122.704 3562.847 2002.970 026d2n0.726 7600.789 1980.774 0620.746 6380.719 1710.694 2310.672 1240.652
596c2n+d2n2.000 03.654 05.012 56.156 67.142 48.007 88.778 79.473 7蒙
特卡洛模拟结果E(D2n)2.000 63.654 65.007 86.157 97.147 38.003 58.778 79.466 4
3.3 非正态分布的极差法系数计算
由第2节可知，JJF 1059.1—2012中的极差系数C实质上是假定输入量服从标准正态分布的情况下导出的样本极差分布的期望值，因而在应用极差法时就需要对输入量的分布予以限制。

只有在接近正态分布的情况下才能使用表1中的极差系数，否则就会给总体标准差的估计引入误差。

实际上，公式(3)中的样本极差分布函数
并未限定总体的分布，这意味着其他分布类型同样可以用极差法来估计总体标准差。

设随机变量ξ服从区间[0,a]上的均匀分布，其概率密度函数为f(x)，有
(11)
容易求得均匀分布的数学期望为μ=a/2，方差为σ2=a2/12，标准差为
设样本ξ1,ξ2,…,ξn的极差为Dn，由式(3)，均匀分布下样本极差Dn的概率密度
函数为：
fDn(r) =n(n-1)[F(x+r)-F(x)]n-2·
f(x+r)f(x)dx
于是样本极差Dn的数学期望和方差分别为：
(12)
(13)
由式(12)得：
于是我们得到输入量的总体服从区间[0,a]上均匀分布的情况下，用样本极差Dn估计总体标准偏差的表达式为：
(14)
式中的值如表5所示。

表5 X～U[0,a]时，从样本极差估计总体标准差的系数Tab.5 The coefficient of the range method when X～U[0,a]样本数n23456789极差系数kn1.154 71.732 02.078 52.309 42.474 42.598 12.694 32.771 3
将kn精确到0.01，然后与表1中的极差系数C做比较，虽然数值不尽相同，但总体来看差别并不大。

进一步分析可知，由于
(15)
所以，由式(14)计算的是均匀分布总体标准差σ的无偏估计量，但却不是均匀分布总体方差σ2的无偏估计量。

式(15)中，因n(n+1)>(n-1)(n+2)，可知实际是均匀分布总体方差σ2的偏大估计量。

因此，若已知输入量的总体服从均匀分布，在标准不确定度的A类评定中也可用极差法来估计总体标准差(使用表5中的极差系数kn)，但同样会导致测量结果的合成标准不确定度偏大。

3.4 总体方差的无偏极差估计量构造
由式(8)可知，在输入量的总体服从正态分布N(μ,σ)的情况下，(Dn/cn)2是总体方差σ2的偏大估计，偏大倍数为
显然将统计量取为就可以实现对总体方差σ2的无偏估计。

以代替表1中的极差系数C，由表2的计算结果，我们得到方差修正后的极差系数C如表6所示：
表6 X～N(μ,σ)时，做方差修正的极差系数Tab.6 The variance-corrected coefficient of the range method when X～N(μ,σ)样本数n23456789方差修正的极差系数C1.414 21.911 52.238 92.481 22.672 52.829 82.926 93.077 9
取表6中的C值，虽然Dn/C并非总体标准差σ的无偏估计量(偏小)，但(Dn/C)2却是总体方差σ2的无偏估计量，从而在扩展不确定度的计算中不会带入原理性误差。

实际上我们不可能构造出总体标准差σ的一个无偏估计量Ψ(即E(Ψ)=σ)，使得估计量Ψ2同时为总体方差σ2的无偏估计(即E(Ψ2)=σ2)。

对随机变量Ψ来说，
D(Ψ)恒大于0，由于E(Ψ2)=D(Ψ)+[E(Ψ)]2，这意味着如果Ψ是总体标准差σ的
无偏估计量，则用Ψ2估计总体方差σ2必然偏大。

JJF 1059.1中给出的极差法就是这种情况。

反之，如果Ψ2是总体方差σ2的无偏估计量，则用Ψ估计总体标准差σ必然偏小。

贝塞尔公式就符合这种情况。

或者更一般的说，我们只能针对总
体标准差和总体方差分别构建无偏估计量。

4 结论
综上所述，本文的结论是：
(1)若输入量的总体服从或接近服从正态分布，由于用极差法估计的总体方差偏大，在标准不确定度的A类评定中采用极差法必然导致测量结果的合成标准不确定度
偏大。

A类分量所占比重越大，测量次数n越少，极差法对测量结果的合成标准
不确定度的影响就越大。

(2)不同的总体分布将导致不同的样本极差分布，因而从样本极差估计总体标准差
的极差系数也不相同。

JJF 1059.1—2012中仅列出了输入量总体服从正态分布时
的极差系数，并不适用于所有情况。

(3)即便从计算简单性考虑、并且在被测量接近正态分布的情况下选用了极差法，
只要合成标准不确定度时不止一个分量，那么极差系数C应优先选用表6中经过
方差修正的系数，以避免合成标准不确定度偏大。

测量不确定度的准确评估对合格评定至关重要[10]。

JCGM 106:2012《测量不确
定度在合格评定中的应用指南》引入了审慎区间(Guard band)以限制做出错误评
定的概率[11]，而将该国际指南等效采用为国家计量技术规范(JJF1059.3)也已经
列入日程。

由于存在概率统计学上的原理性误差，在标准不确定度A类评定中使
用极差法得到的测量不确定度偏大，将导致合格评定中的审慎区间扩大，从而加大做出错误评定的概率，因此建议在标准不确定度A类评定中应审慎使用极差法。

[参考文献]
【相关文献】
[1] 国家质量监督检验检疫总局.JJF 1059.1—2012 测量不确定度评定与表示[S].2012.
[2] ISO/IEC GUIDE 98-3:Uncertainty of measurement-Part 3:Guide to the expression of uncertainty in measurement[S].2008.
[3] 国家质量监督检验检疫总局.JJF 1059—1999测量不确定度评定与表示[S].1999.
[4] 李慎安. JJF 1059—1999《测量不确定度评定与表示》讨论之三十五不确定度评定中极差法的应用问题[J].工业计量，2011，21(4):55.
[5] 倪育才. 测量不确定度理解与应用(二) 极差法和贝塞尔法之问的比较[J].中国计量，2004，(8):78-79.
[6] 梁之舜. 概率论及数理统计(下册)[M].北京：高等教育出版社，1980.
[7] 沙定国.误差分析与测量不确定度评定[M].北京：中国计量出版社，2003.
[8] 刘园园,杨健,赵希勇,等.GUM法和MCM法评定测量不确定度对比分析[J].计量学
报,2018,39(1):135-139.
Liu Y Y, Yang J, Zhao X Y, et al. Comparative analysis of uncertainty measurement evaluation with GUM and MCM[J]. Acta Metrologica Sinica, 2018,39(1): 135-139.
[9] 方兴华,宋明顺,顾龙芳,等.基于自适应蒙特卡罗方法的测量不确定度评定[J].计量学
报,2016,37(4):452-456.
Fang X H, Song M S, Gu L F, et al. Application of adaptive Monte Carlo method on measurement uncertainty evaluation[J]. Acta Metrologica Sinica, 2016,37(4):452-456. [10] 王宏涛,王汉斌,陈晓怀,等.测量不确定度引起的批量产品检验误判风险评估[J].计量学
报,2017,38(5):559-562.
Wang H T, Wang H B, Chen X H, et al. Evaluation of misjudgment risks caused by measurement uncertainty for batch product inspection[J]. Acta Metrologica
Sinica,2017,38(5):559-562.
[11] JCGM 106:Evaluation of measurement data-The role of measurement uncertainty in conformity assessment[S].2012.。