Software to support ISO/TS 28037:2010(E)
Contents
Introduction
This document runs the numerical example of weighted least squares (WLS) with known unequal weights described in Clause 6 (Model for uncertainties associated with the y_i), and performs the prediction described in 11.1 EXAMPLE 2.
Users of MATLAB may run the code in the corresponding M-file directly to obtain the results given in ISO/TS 28037:2010(E) and may also modify the data and run the code on the new data.
For users who do not have access to MATLAB, the software may be used as the basis for preparing implementations in other programming languages.
The software is provided with a software licence agreement (REF: DS/L/24/001) and the use of the software is subject to the terms laid out in that agreement. By running the M-code, the user accepts the terms of the agreement.
close all clear format short
Assign measurement data
Assign x-values.
x = [1.0 2.0 3.0 4.0 5.0 6.0]'; m = length(x);
Assign y-values.
y = [3.2 4.3 7.6 8.6 11.7 12.8]';
Assign uncertainties associated with y-values.
uy = [0.5 0.5 0.5 1.0 1.0 1.0]';
Obtain estimates of the straight line calibration function parameters and associated standard uncertainties and covariance
Solve the weighted least squares problem to obtain best fit straight-line parameters.
Step 1.
w = ones(m, 1)./uy; F2 = sum(w.*w);
Step 2.
g0 = (sum(w.*w.*x))/F2; h0 = (sum(w.*w.*y))/F2;
Step 3.
g = w.*(x - g0); h = w.*(y - h0);
Step 4.
G2 = sum(g.*g);
Step 5.
b = (sum(g.*h))/G2; a = h0 - b*g0;
Step 6.
u2a = 1/F2 + g0^2/G2; u2b = 1/G2; uab = -g0/G2;
Step 7.
r = w.*(y - a - b*x);
Step 8.
chi_sq_obs = sum(r.*r); nu = m - 2;
Step 9.
chi_sq = calc_chi_sq_95_percent_quantile(nu);
Display information on screen and generate figures
Measurement model.
fprintf('\nMODEL FOR UNCERTAINTIES ASSOCIATED WITH THE YI \n\n') fprintf('ISO/TS 28037:2010(E) CLAUSE 6 \n') fprintf('EXAMPLE (UNEQUAL WEIGHTS) \n\n')
MODEL FOR UNCERTAINTIES ASSOCIATED WITH THE YI ISO/TS 28037:2010(E) CLAUSE 6 EXAMPLE (UNEQUAL WEIGHTS)
Measurement data.
fprintf('FITTING \n') fprintf(['Data representing ', num2str(m),' measurement points, unequal weights \n']) for i = 1:m fprintf('%8.1f %8.1f %8.1f \n', [x(i), y(i), uy(i)]) end fprintf('\n')
FITTING Data representing 6 measurement points, unequal weights 1.0 3.2 0.5 2.0 4.3 0.5 3.0 7.6 0.5 4.0 8.6 1.0 5.0 11.7 1.0 6.0 12.8 1.0
Calculation tableau: see write_wls_tableau.m.
write_wls_tableau(x, y, w, g0, h0, g, h, a, b, r, '%8.3f ');
Calculation tableau 0.000 0.000 0.000 0.000 2.600 6.233 0.000 0.000 0.885 0.000 2.000 4.000 4.000 12.800 -3.200 -6.067 10.240 19.413 0.516 0.266 2.000 4.000 8.000 17.200 -1.200 -3.867 1.440 4.640 -1.398 1.955 2.000 4.000 12.000 30.400 0.800 2.733 0.640 2.187 1.088 1.183 1.000 1.000 4.000 8.600 1.400 2.367 1.960 3.313 -0.513 0.263 1.000 1.000 5.000 11.700 2.400 5.467 5.760 13.120 0.530 0.281 1.000 1.000 6.000 12.800 3.400 6.567 11.560 22.327 -0.427 0.182 0.000 15.000 39.000 93.500 0.000 0.000 31.600 65.000 2.057 4.131
Solution estimates.
fprintf('Estimate of intercept \n'), fprintf('%8.3f \n\n', a) fprintf('Estimate of slope \n'), fprintf('%8.3f \n\n', b)
Estimate of intercept 0.885 Estimate of slope 2.057
Standard uncertainties associated with solution estimates.
fprintf('Standard uncertainty associated with estimate of intercept \n'), fprintf('%8.3f \n\n', sqrt(u2a)) fprintf('Standard uncertainty associated with estimate of slope \n'), fprintf('%8.3f \n\n', sqrt(u2b)) fprintf('Covariance associated with estimates of intercept and slope \n'), fprintf('%8.3f \n\n', uab)
Standard uncertainty associated with estimate of intercept 0.530 Standard uncertainty associated with estimate of slope 0.178 Covariance associated with estimates of intercept and slope -0.082
Validation of the model.
fprintf('VALIDATION \n') fprintf('Degrees of freedom \n'), fprintf('%4u \n\n', nu) fprintf('Observed chi-squared value \n'), fprintf('%9.3f \n\n', chi_sq_obs) fprintf('95 %% quantile of chi-squared distribution with %u degrees of freedom', nu), fprintf('\n%9.3f \n\n', chi_sq) if chi_sq_obs > chi_sq fprintf('CHI-SQUARED TEST FAILED - STRAIGHT-LINE MODEL IS REJECTED \n\n') else fprintf('CHI-SQUARED TEST PASSED - STRAIGHT-LINE MODEL IS ACCEPTED \n\n') end
VALIDATION Degrees of freedom 4 Observed chi-squared value 4.131 95 % quantile of chi-squared distribution with 4 degrees of freedom 9.488 CHI-SQUARED TEST PASSED - STRAIGHT-LINE MODEL IS ACCEPTED
Figures.
set(0, 'DefaultLineLineWidth', 2) set(0, 'DefaultAxesFontSize', 12) set(0, 'DefaultAxesFontWeight', 'bold') figure, hold on errorbar(x, y, uy, 'ko', 'MarkerFaceColor', 'k', 'MarkerSize', 6) plot(x, a + b*x, 'b-') xlabel('\it x', 'FontName', 'Times', 'FontSize', 14) ylabel('\it y', 'FontName', 'Times', 'FontSize', 14) axis1 = axis; figure, hold on for i = 1:m plot([x(i), x(i)], [0, r(i)], 'k-') plot(x(i), r(i), 'ko', 'MarkerFaceColor', 'w', 'MarkerSize', 6) end plot([axis1(1), axis1(2)], [0, 0], 'b--') xlabel('\it x', 'FontName', 'Times', 'FontSize', 14) ylabel('\it r', 'FontName', 'Times', 'FontSize', 14)


Prediction
yp = 10.5; uyp = 1.0; xp = (yp - a)/b; ca = -1/b; cb = -(yp - a)/(b^2); cy = 1/b; u2xp = ca^2*u2a + cb^2*u2b + 2*ca*cb*uab + cy^2*uyp^2; fprintf('ISO/TS 28037:2010(E) 11.1 \n') fprintf('EXAMPLE 2 \n\n') fprintf('PREDICTION \n') fprintf('Measured y-value \n'), fprintf('%8.3f \n\n', yp) fprintf('Standard uncertainty associated with measured y-value \n'), fprintf('%8.3f \n\n', uyp) fprintf('Estimate of x \n'), fprintf('%8.3f \n\n', xp) fprintf('Sensitivity coefficient wrt a \n'), fprintf('%8.3f \n\n', ca) fprintf('Sensitivity coefficient wrt b \n'), fprintf('%8.3f \n\n', cb) fprintf('Sensitivity coefficient wrt y \n'), fprintf('%8.3f \n\n', cy) fprintf('Standard uncertainty associated with estimate of x \n'), fprintf('%8.3f \n\n', sqrt(u2xp))
ISO/TS 28037:2010(E) 11.1 EXAMPLE 2 PREDICTION Measured y-value 10.500 Standard uncertainty associated with measured y-value 1.000 Estimate of x 4.674 Sensitivity coefficient wrt a -0.486 Sensitivity coefficient wrt b -2.272 Sensitivity coefficient wrt y 0.486 Standard uncertainty associated with estimate of x 0.533
Acknowledgements
The work described here was supported by the National Measurement System (NMS) delivered through the UK Government’s Department for Science, Innovation and Technology (DSIT).