CN103106332B - 一种测量不确定度的分析方法 - Google Patents

Abstract

本发明公开了一种测量不确定度的分析方法。利用参考测量系测量各标准样件，获取测量结果X，利用某一待分析的测量系统对标准样件进行测量，获取测量结果Y；分别计算测量结果X和Y的方差V(X)和V(Y)，并定义α＝V(X)/V(Y)；假设和是测量值X和Y的一个合理估计；并假设有关系存在；按照Mandel精度定义式计算求得计算定义的全局不确定度；再计算整体不确定度T，利用T²替换第5步公式α＝V(X)/V(Y)中的V(Y)，重新计算一个新的比例系数α¹；将α¹与α对比，如果满足迭代条件，则当前计算的T即为待分析系统的测量不确定度；否则继续迭代。该方法可以对测量系统的测量不确定度作更为合理的分析。

Description

一种测量不确定度的分析方法

技术领域

本发明属于测量技术领域，具体涉及一种新的分析测量不确定度的方法。相对于传统的单纯基于最小二乘法的测量不确定度分析方法，该方法能获得更加全面地考虑和衡量各种误差源对最终测量结果的影响。它适用于对测量系统的测量结果不确定度的分析。

背景技术

在用一个测量系统或测量仪器对标准样件进行测量时，测量的结果往往会在一定程度上偏离真实值。产生这些偏离的原因包括系统不确定度、环境随机噪声甚至样件的自身缺陷等等。实际中，往往将这些偏离在数学上描述为一个不确定度值，也就是该不确定度值在一定程度上确定了标准样件的真实值所在的范围。为了合理地计算该不确定度值，必须利用一定的数学工具来对各种不确定度或者误差来源作出一个合理的描述，并计算不确定度或误差的传递对最终测量结果的影响。目前已经使用的方法有最小二乘法，然而最小二乘法并没有完全考虑所有的测量不确定，相反地，最小二乘法武断地假定参考测量系的不确定度为零，这在一定程度上直接导致了最终测量结果不确定度的增加。因此，需要设计一种合理的测量不确定度分析方法，使得该方法能够尽可能全面地考虑输入不确定度对最终测量不确定度的影响。美国阿切尔等人(C.N.Archie et al.,US 7286247B2)提出了一种整体测量不确定度评估方法，该方法在Mandel精度分析方法基础上，提出了一种迭代计算全局不确定度TMU(Total MeasurementUncertainty)的思想。

发明内容

本发明的目的在于提供一种测量不确定度分析方法，该方法可以对测量系统的测量不确定度作更为合理的分析。

本发明提供的一种测量不确定度分析方法，包括下述步骤：

第1步利用参考测量系测量各标准样件，获取测量结果X＝(x₁,x₂,…,x_n)，测量分量x_k表示对标准样件组中的第k个标准样件进行测量的结果，n表示标准样件的个数；

参考测量系是独立于待分析的测量系统之外的一套已经比较可靠的测量系统；标准样件组含有一系列名义设计值各不相同的样件；

第2步利用某一待分析的测量系统对标准样件进行测量，获取测量结果Y＝(y₁,y₂,…,y_n)，每一个测量分量y_k表示对标准样件组中的第k个标准样件进行测量的结果；

第5步分别计算测量结果X＝(x₁,x₂,…,x_n)和Y＝(y₁,y₂,…,y_n)的方差V(X)和V(Y)，并定义一个比例参数α＝V(X)/V(Y)；

第6步由于样件的真实值无法获取，我们只能在最大限度上给出真实值的一个近似值。此处我们假设和是测量值X＝(x₁,x₂,…,x_n)和Y＝(y₁,y₂,…,y_n)的一个合理估计；并假设有关系存在，其中与是线性拟合方程系数的一个近似估计；

第7步按照Mandel精度定义式计算求得 $\widehat{\mathrm{\β}}=\frac{\mathrm{\α}w-u+\sqrt{{(u-\mathrm{\α}w)}^2+4{\mathrm{\αp}}^2}}{2\mathrm{\α}p},$ 进而得到 $\widehat{\mathrm{\α}}=\frac{\underset{i=1}{\overset{n}{\Σ}}{y}_i}{n}-\widehat{\mathrm{\β}}\frac{\underset{i=1}{\overset{n}{\Σ}}{x}_i}{n};$

第8步按照公式 ${\mathrm{\σ}}_{this}=\sqrt{\frac{\underset{i=1}{\overset{n}{\Σ}}{(y_i-{\hat{y}}_i)}^2}{n-2}}=\sqrt{\frac{\underset{i=1}{\overset{n}{\Σ}}{\left(\frac{y_i-\widehat{\mathrm{\α}}-\widehat{\mathrm{\β}}\·x_i}{1+\mathrm{\α}\·{\widehat{\mathrm{\β}}}^2}\right)}^2}{n-2}}$ 计算本发明中定义的全局不确定度σ_this；

第9步按照定义式计算整体不确定度T，其中σ_X为参考测量系的不确定度，定义为

第10步令V(Y)＝T²，利用，利用公式α＝V(X)/V(Y)，重新计算一个新的比例系数，记为α¹；

第11步将第10步中计算得到的新的比例系数α¹与最开始计算的α做对比，比较二者之间的差值是否小于某一个预先设定的阈值δ，如果有|α¹-α|＜δ，则迭代计算终止，输出当前计算的T即为待分析系统的测量不确定度；如果|α¹-α|＞δ，则继续迭代计算第5步～第10步，直到|α¹-α|＜δ。

其中，阈值δ的作用主要在于终止迭代计算，如δ＝10^-3α。一般来说，为了保证收敛的效果，δ应满足0<δ<＝10^-2α。

本发明是一种新的衡量某一种测量系统或仪器的测量不确定度的方法。该方法充分考虑了各种输入误差或不确定度的影响，并通过迭代的方式来确定系统的整体测量不确定度。该方法充分考虑了各种输入误差或不确定度的影响，并在一定程度上结合了最小二乘法的部分优点，从而实现了对测量系统的测量。本发明在Mandel精度和阿切尔等人提出的方法基础上，本发明提出了新的Mandel残差即σ_this的计算公式，从而更加合理地给出了待分析设备的整体测量不确定度。

附图说明

图1是本发明实例的实现流程图；

图2是本发明提出的全局不确定度与美国阿切尔等人提出的整体测量不确定度的图示对比。

具体实施方式

下面结合附图对本发明的具体实施方式作进一步说明。在此需要说明的是，对于这些实施方式的说明用于帮助理解本发明，但并不构成对本发明的限定。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。

下面将详细介绍本发明中涉及到的测量不确定度分析方法：

1.建立参考测量系；

参考测量系由一种或多种其它的测量仪器组成。建立参考测量系的目的是利用参考测量系对样件的测量结果与待分析的测量设备对样件的测量结果进行一个比对。由于参考测量系的目的是用作对比，因此需要保证组成参考测量系的测量仪器具有较高的测量精度和准确度。

2.标准样件组的准备；

标准样件组由一系列标准样件组成，这些样件应该具有不同大小的待测值。同时，在制作标准样件时，需要采用最成熟可靠的工艺，从而保证标准样件的实际值与名义设计值之间具有最小的偏差；

3.利用第1步中建立的参考测量系测量标准样件组中的每一个样件，获取测量结果；

设标准样件组中一共包含n个标准样件，将参考测量系测得的n个测量结果整理成向量形式X＝(x₁,x₂,…,x_n)，其中X的每一个分量为一个标准样件的测量结果。

其中，由于本方法是以统计学为基础的，因此从理论上来说n越大越好，但实际上n不可能取为无限大。此处，一般来说，n的取值最好大于30。

4.利用待分析的测量仪器对标准样件组中的每一个样件进行测量，获取测量结果：

同步骤3，将待分析测量仪器的n个测量结果整理成向量形式Y＝(y₁,y₂,…,y_n)，其中Y的每一个分量为一个标准样件的测量结果。

5.分别计算测量结果X＝(x₁,x₂,…,x_n)和Y＝(y₁,y₂,…,y_n)的方差V(X)和V(Y)，并定义一个比例参数α＝V(Y)/V(X)；

6.将X＝(x₁,x₂,…,x_n)和Y＝(y₁,y₂,…,y_n)组成n个测量数据点对{(x₁,y₁),(x₂,y₂),…,(x_n,y_n)}；

由于标准样件组中的每一个样件不可能获取准确值，我们只能给出准确值的一个近似，此处假设和是测量值X＝(x₁,x₂,…,x_n)和Y＝(y₁,y₂,…,y_n)的一个合理估计；并假设有关系存在，其中与是线性拟合方程系数的一个近似估计。对于两个完全理想的测量系统来说，其测量值组成的点对在笛卡尔坐标系中的曲线为Y＝X，即并且然而实际上由于测量过程中不可避免地存在各种误差的影响，导致我们不可能获取样件的准确值。

7.按照Mandel精度定义式计算求得 $\widehat{\mathrm{\&beta;}}=\frac{\mathrm{\&alpha;}w-u+\sqrt{{(u-\mathrm{\&alpha;}w)}^2+4{\mathrm{\&alpha;p}}^2}}{2\mathrm{\&alpha;}p},$ 进而得到 $\widehat{\mathrm{\&alpha;}}=\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{y}_i}{n}-\widehat{\mathrm{\&beta;}}\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{x}_i}{n};$

根据Mandel的理论，Mandel精度定义式为：

$WSS=\underset{i=1}{\overset{n}{\&Sigma;}}{\{\mathrm{\&alpha;}{(y_i-{\hat{y}}_i)}^2+(x_i-{\hat{x}}_i)\}}^2---\left(1\right)$

问题的关键在于最小化式(1)。由于我们已经假定将此式代入式(1)中，得：

$WSS=\underset{i=1}{\overset{n}{\&Sigma;}}{\{\mathrm{\&alpha;}{(y_i-\widehat{\mathrm{\&beta;}}\&CenterDot;x_i-\widehat{\mathrm{\&alpha;}})}^2+(x_i-{\hat{x}}_i)\}}^2---\left(2\right)$

定义：

$u=n\underset{i=1}{\overset{n}{\&Sigma;}}{(x_i-\stackrel{\&OverBar;}{x})}^2---\left(3\right)$

$w=n\underset{i=1}{\overset{n}{\&Sigma;}}{(y_i-\stackrel{\&OverBar;}{y})}^2---\left(4\right)$

$p=n\underset{i=1}{\overset{n}{\&Sigma;}}(x_i-\stackrel{\&OverBar;}{x})(y_i-\stackrel{\&OverBar;}{y})---\left(5\right)$

其中和分别为向量X和Y中所有元素的平均值。

进一步地，式(2)可以转化为：

$\mathrm{\&alpha;}\&CenterDot;p\&CenterDot;{\widehat{\mathrm{\&beta;}}}^2+(u-\mathrm{\&alpha;}\&CenterDot;w)\&CenterDot;\widehat{\mathrm{\&beta;}}-p=0---\left(6\right)$

解式(6)所示的一元二次方程可得：

$\widehat{\mathrm{\&beta;}}=\frac{\mathrm{\&alpha;}w-u+\sqrt{{(u-\mathrm{\&alpha;}w)}^2+4{\mathrm{\&alpha;p}}^2}}{2\mathrm{\&alpha;}p}---\left(7\right)$

进而可以得到：

$\widehat{\mathrm{\&alpha;}}=\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{y}_i}{n}-\widehat{\mathrm{\&beta;}}\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{x}_i}{n}=\stackrel{\&OverBar;}{y}-\widehat{\mathrm{\&beta;}}\&CenterDot;\stackrel{\&OverBar;}{x}---\left(8\right)$

8.按照公式 ${\mathrm{\&sigma;}}_{this}=\sqrt{\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{(y_i-{\hat{y}}_i)}^2}{n-2}}=\sqrt{\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{\left(\frac{y_i-\widehat{\mathrm{\&alpha;}}-\widehat{\mathrm{\&beta;}}\&CenterDot;x_i}{1+\mathrm{\&alpha;}\&CenterDot;{\widehat{\mathrm{\&beta;}}}^2}\right)}^2}{n-2}}$ 计算本发明中定义的全局不确定度σ_this；

第7步中计算出了直线拟合参数和将其代入本发明中定义的全局不确定度计算公式(9)中：

${\mathrm{\&sigma;}}_{this}=\sqrt{\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{(y_i-{\hat{y}}_i)}^2}{n-2}}=\sqrt{\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{\left(\frac{y_i-\widehat{\mathrm{\&alpha;}}-\widehat{\mathrm{\&beta;}}\&CenterDot;x_i}{1+\mathrm{\&alpha;}\&CenterDot;{\widehat{\mathrm{\&beta;}}}^2}\right)}^2}{n-2}}---\left(9\right)$

式(9)中的所有参量都是已经知道或者已经计算出来，因此可以直接计算得到本发明中定义的全局不确定度σ_this。相比较US 7286247 B2文献中全局不确定度的计算式 ${\mathrm{\&sigma;}}_{Mandel}=\sqrt{\frac{\underset{i=1}{\overset{n}{\&Sigma;}}\{{(y_i-{\hat{y}}_i)}^2+{(x_i-{\hat{x}}_i)}^2\}}{n-2}},$ 本发明剔除了这一项。两者之间的区别如图2所示，①为拟合得到的线性，②代表US 7286247 B2文献中涉及的计算式的计算方向，③代表本发明中涉及的计算式 ${\mathrm{\&sigma;}}_{this}=\sqrt{\frac{\underset{i-1}{\overset{n}{\&Sigma;}}{(y_i-{\hat{y}}_i)}^2}{n-2}}=\sqrt{\frac{\underset{i-1}{\overset{n}{\&Sigma;}}{\left(\frac{y_i-\widehat{\mathrm{\&alpha;}}-\widehat{\mathrm{\&beta;}}\&CenterDot;x_i}{1+\mathrm{\&alpha;}.{\widehat{\mathrm{\&beta;}}}^2}\right)}^2}{n-2}}$ 的计算方向。

根据公式(9)我们可以知道，参考测量系对标准样件组的测量结果X已经包含在全局不确定度σ_this的计算公式中，进而在中已经包含了由参考测量系引入的不确定度。因此，这一用来衡量参考测量系不确定度的子式需要剔除。这种计算方式就能够更合理地对待分析设备的整体测量不确定度进行分析。

9.按照定义式计算整体不确定度T，其中σ_X为参考测量系的不确定度，定义为

10.令V(Y)＝T²，利用公式α＝V(X)/V(Y)重新计算一个新的比例系数，记为α¹；

α¹＝T²/V(X)          (10)

11.将第10步中计算得到的新的比例系数α¹与最开始计算的α做对比，比较二者之间的差值是否小于某一个预先设定的阈值δ(如δ＝10^-3α)，如果有|α¹-α|＜δ，则迭代计算终止，输出当前计算的T即为待分析系统的测量不确定度；如果|α¹-α|＞δ，则继续迭代计算第5～10步，直到|α¹-α|＜δ。

以上所述为本发明的较佳实施例而已，但本发明不应该局限于该实施例和附图所公开的内容。所以凡是不脱离本发明所公开的精神下完成的等效或修改，都落入本发明保护的范围。

Claims (1)

1.一种测量不确定度的分析方法，该方法包括下述步骤：

第1步利用参考测量系测量各标准样件，获取测量结果X＝(x₁,x₂,…,x_n)，测量分量x_k表示对标准样件组中的第k个标准样件进行测量的结果，n表示标准样件的个数；

第2步利用某一待分析的测量系统对标准样件进行测量，获取测量结果Y＝(y₁,y₂,…,y_n)，每一个测量分量y_k表示对标准样件组中的第k个标准样件进行测量的结果；

第3步分别计算测量结果X＝(x₁,x₂,…,x_n)和Y＝(y₁,y₂,…,y_n)的方差V(X)和V(Y)，并定义一个比例参数α＝V(X)/V(Y)；

第4步假设和是测量值X＝(x₁,x₂,…,x_n)和Y＝(y₁,y₂,…,y_n)的一个合理估计；并假设有关系存在，其中与是线性拟合方程系数的一个近似估计；

第5步按照Mandel精度定义式计算求得 $\widehat{\mathrm{\&beta;}}=\frac{\mathrm{\&alpha;}w-u+\sqrt{{(u-\mathrm{\&alpha;}w)}^2+4{\mathrm{\&alpha;p}}^2}}{2\mathrm{\&alpha;}p},$ 进而得到 $\widehat{\mathrm{\&alpha;}}=\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{y}_i}{n}-\widehat{\mathrm{\&beta;}}\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{x}_i}{n};$

其中， $w=n\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{(y_i-\stackrel{\&OverBar;}{y})}^2,u=n\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{(x_i-\stackrel{\&OverBar;}{x})}^2,p=n\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}(x_i-\stackrel{\&OverBar;}{x})(y_i-\stackrel{\&OverBar;}{y});$

第6步按照公式 ${\mathrm{\&sigma;}}_{this}=\sqrt{\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{(y_i-{\hat{y}}_i)}^2}{n-2}}=\sqrt{\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{\left(\frac{y_i-\widehat{\mathrm{\&alpha;}}-\widehat{\mathrm{\&beta;}}\&CenterDot;x_i}{1+\mathrm{\&alpha;}\&CenterDot;{\widehat{\mathrm{\&beta;}}}^2}\right)}^2}{n-2}}$ 计算本发明中定义的全局不确定度σ_this；

第7步按照定义式计算整体不确定度T，其中σ_X为参考测量系的不确定度，定义为

第8步令V(Y)＝T²，利用公式α＝V(X)/V(Y)重新计算一个新的比例系数，记为α¹；

第9步将第8步中计算得到的新的比例系数α¹与最开始计算的α做对比，比较二者之间的差值是否小于某一个预先设定的阈值δ，如果有|α¹-α|＜δ，则迭代计算终止，输出当前计算的T即为待分析系统的测量不确定度；如果|α¹-α|＞δ，则继续迭代计算第2步～第8步，直到|α¹-α|＜δ。