Contents
- Function implementing the GUM uncertainty framework
- Step 1: Summarize the probability distributions for the input quantities
- Step 2: Apply the law of propagation of uncertainty
- Step 3: Characterize the output quantity by a probability distribution
- Step 4: Evaluate the coverage factor
- Step 5: Evaluate the expanded uncertainty and coverage interval
- Function for summarizing the probability distribution for an input quantity
- Function for evaluating a coverage factor
- Function for evaluating a coverage factor for the t-distribution
- Program identifier
Function implementing the GUM uncertainty framework
This function implements the GUM uncertainty framework for uncertainty evaluation for a measurement model
with a scalar output quantity and a general number of uncorrelated input quantities.
The inputs to the function are a specification model of the measurement model, a specification pdfin of the probability distributions assigned to the input quantities in the measurement model, and a coverage probability p. The function returns an estimate yGUM of the output quantity, the standard uncertainty uyGUM associated with the estimate, the endpoints IyGUM of the probabilistically symmetric coverage interval for the output quantity corresponding to coverage probability p, and a specification pdfGUM of the probability distribution used to characterize the output quantity.
The estimate is obtained by evaluating the model at the expectations of the input quantities and the associated standard uncertainty by applying the law of propagation of uncertainty. The output quantity is characterized by a Gaussian distribution or a scaled and shifted t-distribution, and the distribution is used as the basis of evaluating the coverage interval for the output quantity.
model is a string containing the name of a function for evaluating the measurement model for values of the input quantities and calculating sensitivity coefficients. For an example of such a function, see model_additive.
Each input quantity in the measurement model is assigned a probability distribution, described by a probability density function, which may be (1) a Gaussian (or normal) distribution, (2) a scaled and shifted t-distribution, (3) a rectangular (or uniform) distribution, (4) a curvilinear trapezoidal distribution, or (5) a U-shaped distribution. (This list includes all the univariate probability distributions encountered in the examples of GUM Supplement 1. There is no consideration of joint probability distributions, such as the multivariate Gaussian distribution.)
pdfin is a matrix of N rows and 4 columns, the ith row of which defines the probability distribution for the ith input quantity. For a measurement model with a single input quantity, pdfin can take the following forms:
pdfin = [1, mu, sigma, inf] defines a Gaussian distribution with expectation mu and standard deviation sigma (the fourth element of pdfin is not used);
pdfin = [2, mu, sigma, nu] defines a t-distribution with shift parameter mu, scale parameter sigma and degrees of freedom nu;
pdfin = [3, a, b, inf] defines a rectangular distribution with lower limit a and upper limit b (the fourth element of pdfin is not used);
pdfin = [4, a, b, r] defines a curvilinear trapezoidal distribution with lower limit a, upper limit b and (fractional) reliability r for the semi-width;
pdfin = [5, a, b, inf] defines a U-shaped distribution with lower limit a and upper limit b (the fourth element of pdfin is not used).
Guidance on the assignment, in some common circumstances, of probability density functions for the input quantities in a measurement model is given in GUM Supplement 1 [3, section 6].
The coverage probability p takes one of the values 0.90, 0.95 or 0.99.
Section 4.1 of M G Cox, P M Harris and I M Smith, Software specifications for uncertainty evaluation, 2010, describes the procedure implemented in the GUM uncertainty framework.
Section 3.1.1 of M G Cox, P M Harris and I M Smith, Software for GUM Supplement 1: User Manual, 2011, describes the specifications of the different probability distributions. (Note that this document describes other probability distributions that are not considered here.)
function [yGUM, uyGUM, IyGUM, pdfGUM] = GUM_calculation(model, pdfin, p)
Input
model string Function defining the measurement model pdfin N x 4 Probability density functions for the input quantities in the measurement model p Coverage probability (0.90, 0.95 or 0.99)
Output
yGUM Estimate of the output quantity uyGUM Standard uncertainty associated with the estimate of the output quantity IyGUM 1 x 2 Lower and upper endpoints of a coverage interval for the output quantity pdfGUM 1 x 4 Probability density function for the output quantity: [1 or 2, yGUM, uyGUM, veff]
Step 1: Summarize the probability distributions for the input quantities
Each input quantity is summarized by an estimate (expectation) x, associated standard uncertainty (standard deviation) ux, and associated degrees of freedom nu, obtained from the probability density function assigned to the quantity.
N = size(pdfin,1); for i = 1:N [x(i,1), ux(i,1), nu(i,1)] = convert_pdf(pdfin(i,:)); end
Step 2: Apply the law of propagation of uncertainty
The model is evaluated for the estimates of the input quantities to obtain an estimate yGUM of the output quantity. The first order partial derivatives of the model with respect to the input quantities are evaluated for the estimates of the input quantities to obtain the sensitivity coefficients c. The law of propagation of uncertainty is applied to obtain the standard uncertainty uyGUM associated with yGUM.
See sections 4.2 and 4.3 of M G Cox, P M Harris and I M Smith, Software specifications for uncertainty evaluation, 2010.
[yGUM, c] = eval([model, '(x)']);
uyGUM = sqrt(c*diag(ux.^2)*c');
Step 3: Characterize the output quantity by a probability distribution
The output quantity is characterized by a Gaussian distribution (for an effective degrees of freedom exceeding 100) or a scaled and shifted t-distribution otherwise. The Welch-Satterthwaite formula is used to evaluate the effective degrees of freedom veff.
See section 4.4 of M G Cox, P M Harris and I M Smith, Software specifications for uncertainty evaluation, 2010.
denom_sum = sum(((c'.*ux).^4)./nu); if denom_sum == 0 veff = inf; pdfGUM = [1, yGUM, uyGUM, veff]; else veff = fix((uyGUM^4)/denom_sum); if veff > 100 veff = inf; pdfGUM = [1, yGUM, uyGUM, veff]; else pdfGUM = [2, yGUM, uyGUM, veff]; end end
Step 4: Evaluate the coverage factor
The coverage factor kp is obtained as a percentage point of the standard Gaussian distribution (for an effective degrees of freedom exceeding 100) or a t-distribution otherwise.
See section 4.5 of M G Cox, P M Harris and I M Smith, Software specifications for uncertainty evaluation, 2010.
kp = coverage_factor(p, veff);
Step 5: Evaluate the expanded uncertainty and coverage interval
The expanded uncertainty Uy is the product of the coverage factor kp and the standard uncertainty uyGUM. The probabilistically symmetric coverage interval (also the shortest coverage interval) is an interval centred on the estimate yGUM with half-width equal to the expanded uncertainty.
See section 4.5 of M G Cox, P M Harris and I M Smith, Software specifications for uncertainty evaluation, 2010.
Uy = kp*uyGUM; IyGUM = [(yGUM - Uy) (yGUM + Uy)];
Function for summarizing the probability distribution for an input quantity
This function summarizes an input quantity in a measurement model by an estimate (expectation), associated standard uncertainty (standard deviation), and associated degrees of freedom, obtained from the probability distribution assigned to the quantity.
function [x, ux, nu] = convert_pdf(pdf)
Input
pdf 1 x 4 Probability density function for an input quantity in a measurement model
Output
x Estimate of the input quantity ux Standard uncertainty associated with the estimate of the input quantity nu Degrees of freedom associated with the standard uncertainty
The calculation depends on the nature of the probability distribution assigned to the input quantity.
switch pdf(1)
Gaussian distribution.
case {1} x = pdf(2); ux = pdf(3); nu = inf; %%% % t-distribution. case {2} x = pdf(2); ux = pdf(3); nu = pdf(4); %%% % Rectangular distribution. case {3} a = pdf(2); b = pdf(3); x = (b + a)/2; ux = (b - a)/(2*sqrt(3)); nu = inf; %%% % Curvilinear trapezoidal distribution. case {4} a = pdf(2); b = pdf(3); r = pdf(4); x = (b + a)/2; ux = (b - a)/(2*sqrt(3)); nu = (1/2)*(1/r)^2; %%% % U-shaped distribution. case {5} a = pdf(2); b = pdf(3); x = (b + a)/2; ux = (b - a)/(2*sqrt(2)); nu = inf; end
Function for evaluating a coverage factor
In the GUM uncertainty framework, the random variable
is characterized by the standard Gaussian distribution N(0, 1) in the case that the effective degrees of freedom veff associated with u(y) is infinite, or by a t-distribution with veff degrees of freedom otherwise. This function evaluates the coverage factor kp corresponding to coverage probability p such that
The coverage factor is evaluated as a percentage point of the distribution (standard Gaussian or t-distribution) used to characterize Z.
See section 4.5 of M G Cox, P M Harris and I M Smith, Software specifications for uncertainty evaluation, 2010.
function kp = coverage_factor(p, veff)
Input
p Coverage probability (0.90, 0.95 or 0.99) veff Effective degrees of freedom (inf for a Gaussian distribution, finite for a t-distribution)
Output
kp Coverage factor
Distinguish between a degrees of freedom that is infinite or finite.
if isinf(veff) if p == 0.90 kp = 1.645; elseif p == 0.95 kp = 1.960; elseif p == 0.99 kp = 2.576; end else kp = t_dist_coverage(p, veff); end
Function for evaluating a coverage factor for the t-distribution
This function evaluates the coverage factor kp corresponding to coverage probability p such that
where T is characterized by a t-distribution with nu degrees of freedom.
See Table G.2 of BIPM, IEC, IFCC, ISO, IUPAC, IUPAP and OIML, Guide to the expression of uncertainty in measurement (ISBN 92-67-10188-9, corrected and reprinted, 1995).
function kp = t_dist_coverage(p, nu)
Input
p Coverage probability (0.90, 0.95 or 0.99) nu Degrees of freedom
Output
kp Coverage factor for the t-distribution
Look-up table TABLE for the coverage factor.
TABLE = [1 6.31 12.71 63.66 2 2.92 4.30 9.92 3 2.35 3.18 5.84 4 2.13 2.78 4.60 5 2.02 2.57 4.03 6 1.94 2.45 3.71 7 1.89 2.36 3.50 8 1.86 2.31 3.36 9 1.83 2.26 3.25 10 1.81 2.23 3.17 11 1.80 2.20 3.11 12 1.78 2.18 3.05 13 1.77 2.16 3.01 14 1.76 2.14 2.98 15 1.75 2.13 2.95 16 1.75 2.12 2.92 17 1.74 2.11 2.90 18 1.73 2.10 2.88 19 1.73 2.09 2.86 20 1.72 2.09 2.85 25 1.71 2.06 2.79 30 1.70 2.04 2.75 35 1.70 2.03 2.72 40 1.68 2.02 2.70 45 1.68 2.01 2.69 50 1.68 2.01 2.68 100 1.660 1.984 2.626];
Identify which column of TABLE to use (column 2 for 90 %, column 3 for 95 %, column 4 for 99 %).
if p == 0.90 col = 2; elseif p == 0.95 col = 3; elseif p == 0.99 col = 4; end
Use MATLAB's interp1 function to linearly interpolate the values in the identified column of TABLE. (For degrees of freedom greater than 100, take nu = 100.)
if nu > 100 kp = TABLE(27, col); else kp = interp1(TABLE(:, 1), TABLE(:, col), nu); end
Program identifier
NPLUnc_101 Open-Source Software (Release 1.1): Software for measurement uncertainty evaluation
M G Cox, P M Harris and I M Smith, Mathematics and Scientific Computing Group, National Physical Laboratory, UK
Crown copyright 2011