CN100442270C - 一种用蒙特卡罗统计模拟计算合成不确定度的方法 - Google Patents

Abstract

本发明涉及计量测试技术领域，尤其是涉及一种用蒙特卡罗统计模拟计算合成不确定度的方法。它包括如下步骤：a、引入测量输入量和输出量函数关系的数学模型Y＝f(x₁，x₂，…，x_n)，上式中Y为输出量、x₁，x₂，…，x_n为输入量；b、引入输入量x₁，x₂，…，x_n不确定度的概率分布和参数；c、根据输入量x₁，x₂，…，x_n不确定度的概率分布和参数，选择随机数ξ₁，ξ₂，…，ξ_n；d、选择对随机数ξ₁，ξ₂，…，ξ_n的检验方法；e、在计算机上进行模拟计算输出量Y、输出量Y的标准偏差和与置信水平精确相关的合成不确定度。本发明可使最终的扩展不确定度与置信水平精确关联，特别是高置信水平的扩展不确定度(置信水平95%或99%等)。这将在许多领域，如工业、商业以及卫生和安全等领域发挥很重要的作用。

Description

一种用蒙特卡罗统计模拟计算合成不确定度的方法

技术领域

本发明涉及计量测试技术领域，尤其是涉及一种用蒙特卡罗统计模拟计算合成不确定度的方法。

背景技术

不确定度的概念在计量测试技术领域占有着重要的地位。任一有效的测量必须带有不确定度的评定，以便于人们可以对这一测量过程的优劣进行评估、比较和复现这一测量过程。它在科技领域成为不可缺少的部分，所以国际上于1993年由国际标准化组织(ISO)等7个国际组织的名义联合发布了《测量不确定度表示指南》(Guide to the Expression of Uncertainty in Measurement简称GUM)。我国也于1999年发布了JJF1059-1999《测量不确定度评定与表示》，以便于与国际同步，可对重要的测量过程进行统一的评估，也保证了计量领域的量值传递的可行性和统一性。

但是，不确定度伴随着任一测量过程。许多领域，如工业、商业以及卫生和安全等领域，需要伴随着这个测量值的不确定度有一个确切的置信区间，以期望将其作为依据，预测项目的可行性。人们期望得到的数据(不确定度)既不夸大可预知性，即在概率要求确定的情况下，缩小不确定度；也不缩小预知性(如果发生相反的情况)，而期望获得正确的预测。由于大多数测量结果的概率分布是很难用数学方法得出的，所以上述标准根据不确定度传播率导出的合成标准不确定度u_c乘以复盖因子k_p得到的扩展不确定度U_p，并不确切地与置信水平的概率相关联(如p＝95％，，表示大约应有95％的测量结果落入此区间)。

发明内容

本发明所要解决的技术问题是提供一种可使最终的合成不确定度与置信水平精确关联的用蒙特卡罗统计模拟计算合成不确定度的方法，它包括如下步骤：a、引入测量输入量和输出量函数关系的数学模型Y＝f(x₁，x₂，…x_n)，上式中Y为输出量、x₁，x₂，…x_n为输入量；b、引入输入量x₁，x₂，…x_n不确定度的概率分布和参数；c、根据输入量x₁，x₂，…x_n不确定度的概率分布和参数，选择随机数ξ₁，ξ₂，…ξ_n；d、选择对随机数ξ₁，ξ₂，…ξ_n的检验方法；e、在计算机上进行模拟计算输出量Y、输出量Y的标准偏差和与置信水平精确相关的合成不确定度，在计算机上进行模拟计算是将各随机数ξ₁，ξ₂，…ξ_n代入Y＝f(x₁，x₂，…x_n)得到y_i＝f(ξ_1i，ξ_2i，…ξ_ni)，各随机数的取值范围在各不确定度分量分布的全范围，相应得到的y_i即为输出量的可能值，i＝1，2，…n，在计算机上产生n个模拟值y₁，y₂，…y_n，输出量 $Y=\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i$ ，输出量Y的标准偏差 $s={[\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{(y_i-Y)}^2/(n-1)]}^{1/2},$ 并根据需要评定合成不确定度置信区间的概率p，使得p＝m/n，式中n为模拟总次数，m为|y_i-Y|＜δ的模拟次数，δ为需要最终评定的输出量Y的与置信水平精确相关的合成不确定度。

本发明可使最终的合成不确定度与置信水平精确关联，特别是高置信水平的合成不确定度(置信水平95％或99％等)。如果根据不确定度传播率导出的合成标准不确定度u_c乘以复盖因子k_p得到的扩展不确定度U_p，由于输出量的分布往往不可知，而不确切地与置信水平的概率相关联。本发明将在许多领域，如工业、商业、医药卫生等领域发挥很大的作用。

附图说明

图1为模拟计算的计算机程序流程。

图2为模拟计算Y＝f(D₀，M_x，D_x)的分布。

具体实施方式

表1为模拟计算的10次重复计算值(不同模拟次数)。

实施例：本实施例是一个校准扫描电子显微镜记录图像上显微标尺的实例，其中测量数学模型即为计算校准值的计算公式：

M＝D₀×M_x/D_x

式中：

M——显微标尺的校准值；

D₀——标准物质的长度值；

M_x——记录图像上显微标尺的长度测量值；

D_x——标准物质的长度测量值。

各输入量不确定度的确定和产生随机数的方法

在测量的数学模型f(D₀，M_x，D_x)＝D₀×M_x/D_x已建立的情况下，其中各输入量D₀，M_x，D_x的不确定度假设确定如下：

标准物质的长度值D₀＝4.6μm±0.05μm，不确定度服从正态分布N(a，σ)，这里a＝4.6μm；σ＝0.05μm。标准物质的长度值用随机数ξ_1i表示，选择通常蒙特卡罗方法介绍的正态分布N(a，σ)的随机数产生的方法，可以得到：

${\mathrm{\ξ}}_{1i}=\mathrm{\σ}(\underset{k=1}{\overset{48}{\mathrm{\Σ}}}{r}_k-24)/2+a$

$=0.05(\underset{k=1}{\overset{48}{\mathrm{\Σ}}}{r}_k-24)/2+4.6$ (其中r_k为(0，1)上均匀分布的随机数)

记录图像上显微标尺的长度测量值M_x＝15mm±0.5mm，不确定度服从(a，b)区间的均匀分布，这里a＝15mm-0.5mm＝14.5mm；b＝15mm+0.5mm＝15.5mm。记录图像上显微标尺的长度测量值用随机数ξ_2i表示，选择通常蒙特卡罗方法介绍的对任意区域(a，b)上均匀分布的随机数产生的方法产生，可以得到：

ξ_2i＝a+(b-a)r_2i＝14.5+r_2i(其中r_2i为(0，1)上均匀分布的随机数)

标准物质的长度测量值D_x＝36mm±0.5mm，不确定度服从在(a，b)区间的均匀分布，这里a＝36mm-0.5mm＝35.5mm；b＝36mm+0.5mm＝36.5mm。标准物质长度测量值用随机数ξ_3i表示，选择通常蒙特卡罗方法介绍的对任意区域(a，b)上均匀分布的随机数产生的方法产生，可以得到：

ξ_3i＝a+(b-a)r_3i＝35.5+r_3i(其中r_3i为(0，1)上均匀分布的随机数)随机数的检验

由于本实例是采用美国Microsoft公司的产品Visual Basic 6.0软件编制的计算机程序来进行模拟计算，随机数的产生是该软件所附的随机数发生器产生的。可以认为这样获得的在(0，1)上均匀分布的随机数系列是通过各项参数检验、均匀性检验和独立性检验的，所以在计算程序中，我们仅对在该随机数基础上构造起来的在任意区域(a，b)上均匀分布的随机数和正态分布N(a，σ)的随机数的平均值 $\stackrel{\‾}{y}=\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i$ 进行参数检验，其方法如下：

对于任意区域(a，b)上均匀分布随机数的数学期望和方差表示为：

E(y)＝(a+b)/2

σ²＝E(y²)-[E(y)]²＝(b³-a³)/3-(a+b)²/4

由中心极限定理知统计量：

$u=[\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i-\mathrm{nE}\left(y\right)]/(\mathrm{\σ}\×\sqrt{n})=\sqrt{n}[\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i-\frac{a+b}{2}]/{[\frac{b^3-a^3}{3}-\frac{{(a+b)}^2}{4}]}^{1/2}$

当N充分大时，渐近地服从N(0，1)分布。

对于正态分布N(a，σ)的随机数的数学期望和方差直接就是a和σ，所以由中心极限定理得到统计量：

$u=[\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i-\mathrm{nE}\left(y\right)]/(\mathrm{\σ}\×\sqrt{n})=\frac{\sqrt{n}}{\mathrm{\σ}}(\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i-a)$

当N充分大时，渐近地服从N(0，1)分布。

取显著性水平α＝0.1，则当|u|＞1.645时，计算机程序将重新产生随机数，一直到产生的随机数通过参数检验然后继续进行模拟计算。

模拟计算的计算机程序流程参见图1。

模拟计算结果：

模拟计算Y＝f(D₀，M_x，D_x)的分布参见图2。

Y的模拟计算平均值：E(y)＝1.91670410155901μm；模拟计算次数：n＝100000；Y的模拟计算标准偏差：S＝4.62415096429435E-02μm；模拟计算值的置信概率为95％的合成不确定度：U₉₅＝0.086317873244385μm。

本发明还可以采用多次重复计算取得平均值，以便获得更高的精度。模拟计算的10次重复计算值(不同模拟次数)参见表1。

如果取n＝100000的10次重复计算平均值为最终结果，可得校准值及不确定度为：

图像上显微标尺校准值M＝1.917μm，标准偏差s＝0.046μm，U₉₅(M)＝0.085μm。

根据《测量不确定度表示指南》(Guide to the Expression of Uncertainty inMeasurement简称GUM)和中国国家标准JJF1059-1999《测量不确定度评定与表示》中规定的不确定度的传播率 $u_c^2\left(y\right)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\left[\frac{\∂f}{\∂x_i}\right]}^2{u}^2\left(x_i\right)$ 计算可得校准值及不确定度为：图像上显微标尺校准值M＝1.917μm，合成标准不确定度u_c(M)＝0.045μm，扩展不确定度U(M)＝0.090μm；k＝2。

模拟次数n	模拟计算平均值E(y)	标准偏差S	合成不确定度(95％)
				1000	1.91999960561	0.04343572944	0.07719564026
1000	1.91431399520	0.04693235069	0.08898753965
				1000	1.91650754874	0.04502600662	0.08558322496
1000	1.91461375245	0.04853157404	0.09428388535
				1000	1.92130136760	0.04419618587	0.08387538353
1000	1.92078692475	0.04687640923	0.08692830841
				1000	1.91278062577	0.04706553280	0.08933426860
1000	1.92021422396	0.04675727676	0.08659899333
				1000	1.92007606678	0.04442139537	0.08172215466
1000	1.91660246479	0.04489642036	0.08437755711
				平均值	1.91771965757	0.04581388812	0.08588869559
标准偏差	0.00311924152	0.00162946694	0.00462960147
10000	1.91327348144	0.04485437348	0.08474359067
				10000	1.91471483038	0.04687660217	0.08655682969
10000	1.91336180236	0.04237230110	0.08064296981
				10000	1.91254033916	0.04510784008	0.08615734388
10000	1.91822420335	0.04547925248	0.08737239349
				10000	1.91401236881	0.04731933820	0.08725464655
10000	1.91526978118	0.04781703640	0.09169709697
				10000	1.91719845107	0.04706752749	0.08899017738
10000	1.92087012419	0.04491498246	0.08518137817
				10000	1.92042430766	0.04588298945	0.08492132702
平均值	1.91598896896	0.04576922433	0.08635177536
				标准偏差	0.00302165441	0.00160376252	0.00290713773

(续下页)

(接上页)

50000	1.91445048858	0.04601576090	0.08670459993
				50000	1.91490416155	0.04514804391	0.08322765646
50000	1.91764740887	0.04391093842	0.08278820786
				50000	1.91899204445	0.04468988012	0.08517379323
50000	1.92114890430	0.04589111150	0.08504693505
				50000	1.91788347790	0.04502908213	0.08452275663
50000	1.91936685187	0.04553810318	0.08506173059
				50000	1.91788666864	0.04758123999	0.08746756429
50000	1.91578996295	0.04386265770	0.08274164036
				50000	1.91312546807	0.04560360155	0.08617017878
平均值	1.91711954372	0.04532704194	0.08489050632
				标准偏差	0.00249433394	0.00108890950	0.00162074597
				100000	1.91645610604	0.04647044109	0.08786854215
100000	1.91847265710	0.04554862898	0.08565389052
				100000	1.91560375003	0.04577514149	0.08593928581
100000	1.91709709493	0.04626752522	0.08563484604
				100000	1.91789391940	0.04447459245	0.08307888582
100000	1.91695612417	0.04530520963	0.08518167415
				100000	1.91808268402	0.04508823893	0.08446613380
100000	1.91783382781	0.04513480701	0.08416785468
				100000	1.91522932069	0.04712540374	0.08666412406
100000	1.91642910731	0.04621093718	0.08563994081
				平均值	1.91700545915	0.04574009257	0.08542951779
标准偏差	0.00108284077	0.00078684744	0.00133235532

表1

Claims (1)

1.一种用蒙特卡罗统计模拟计算合成不确定度的方法，其特征是它包括如下步骤：

a、引入测量输入量和输出量函数关系的数学模型Y＝f(x₁，x₂，…x_n)，上式中Y为输出量、x₁，x₂，…x_n为输入量；

b、引入输入量x₁，x₂，…x_n不确定度的概率分布和参数；

c、根据输入量x₁，x₂，…x_n不确定度的概率分布和参数，选择随机数ξ₁，ξ₂，…ξ_n；

d、选择对随机数ξ₁，ξ₂，…ξ_n的检验方法；

e、在计算机上进行模拟计算输出量Y、输出量Y的标准偏差和与置信水平精确相关的合成不确定度，在计算机上进行模拟计算是将各随机数ξ₁，ξ₂，…ξ_n代入Y＝f(x₁，x₂，…x_n)得到y_i＝f(ξ_1i，ξ_2i，…ξ_ni)，各随机数的取值范围在各不确定度分量分布的全范围，相应得到的y_i即为输出量的可能值，i＝1，2，…n，在计算机上产生n个模拟值y₁，y₂，…y_n，输出量 $Y=\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i,$ 输出量Y的标准偏差为 $s={[\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{(y_i-Y)}^2/(n-1)]}^{1/2},$ 并根据需要评定合成不确定度置信区间的概率p，使得p＝m/n，式中n为模拟总次数，m为|y_i-Y|＜δ的模拟次数，δ为需要最终评定的输出量Y的与置信水平精确相关的合成不确定度。

Images