Approximated Uncertainty Propagation of Correlated Independent Variables Using the Ordinary Least Squares Estimator

Abstract

For chemical measurements, calibration is typically conducted by regression analysis. In many cases, generalized approaches are required to account for a complex-structured variance–covariance matrix of (in)dependent variables. However, in the particular case of highly correlated independent variables, the ordinary least squares (OLS) method can play a rational role with an approximated propagation of uncertainties of the correlated independent variables into that of a calibrated value for a particular case in which standard deviation of fit residuals are close to the uncertainties along the ordinate of calibration data. This proposed method aids in bypassing an iterative solver for the minimization of the implicit form of the squared residuals. This further allows us to derive the explicit expression of budgeted uncertainties corresponding to a regression uncertainty, the measurement uncertainty of the calibration target, and correlated independent variables. Explicit analytical expressions for the calibrated value and associated uncertainties are given for straight-line and second-order polynomial fit models for the highly correlated independent variables.

1. Introduction

For calibration purposes, regression analysis wherein (in)dependent variables with a complex uncertainty structure has been intensively studied to yield the best estimates of the regression coefficient and associated uncertainties [1,2,3,4,5,6,7,8,9,10,11,12,13]. In particular, the uncertainty propagation of measured scalar responses toward the final calibrated value is a critical issue for ensuring the compatibility of inter- and intra-laboratory measurement results [14]. In the case of non-zero covariances among responses (dependent variables) and stimuli (independent variables), generalized least squares (GLS) regression is required to solve for the matrix form of squared residuals (SR), which is weighted against the variance–covariance matrix of calibration data vectors [15]. The weighted SR matrix is then minimized by an iterative solver. When the uncertainties of independent variables and covariances among dependent and independent variables are nonexistent, the GLS regression can be reduced to the weighted least squares (WLS) regression [1,2,3,4,5,6,7,8,16]. As the off-diagonal elements of the SR matrix are zero, the analytical expression for the fitted parameters of the regression model can be analytically derived.

In the metrological field, several research papers [17,18,19,20,21,22,23,24] and standards [25,26,27] have been published on regression analysis, and relevant softwares have been developed [28,29,30,31,32]. Analytical formulations of GLS are available from NPL and INRiM [18,20]. INRiM’s software, CCC (v2.0), works for linear and nonlinear regression problems with a full covariance matrix, excluding covariance between independent and dependent variables [28]. However, NPL’s software, XLGENLINE (v1.1), implemented in EXCEL by calling FORTRAN DLL, can work only for linear regression problems with up to fourth-order polynomials and diagonal covariance matrices [29]. In contrast, XGENLINE (v8.1) [30], which is a MATLAB version of XLGENLINE, uses a representation in terms of the Chebyshev polynomials to obtain better numerical conditioning than the general polynomials [17]. ISO/TS 28037:2010(E) [26] and the relevant software B_LEAST [30] can address the case of correlation among dependent and independent variables, but only when considering a straight-line model. B_LEAST can work with different types of regression models, including power and exponential models, but it requires diagonal covariance matrices only [32].

Following the Guide to the Expression of Uncertainty in Measurement (GUM), a common text for the metrological determination of the measurement value [33], the uncertainty of the calibration target can be evaluated by applying the (linear) law of propagation of uncertainty (LPU) to the abscissa. Therefore, for a generalized approach with the non-diagonalized SR matrix, an analytical process to derive the variance–covariance matrix of regression coefficients has been demonstrated only for particular cases wherein the correlation between variables does not disturb the derivation of the explicit variance–covariance formulae [33]. For chemical or biological RMs prepared from the same mother material, one may reasonably assume a strong correlation between these uncertainties [34]. For instance, a batch of chemical mixtures in various amount-fractions are typically blended using the same raw materials or mother mixture in the same facility and perhaps even by the same worker. Therefore, their variances may be correlated to some extent. Here, we attempt to derive an explicit expression of the uncertainty of the calibrated value which was propagated from measurement uncertainties of references and calibration targets using the ordinary least squares (OLS) estimator according to LPU. The calibration result using the proposed method is compared to the Monte Carlo simulation, complementing the validity of the proposed method. In addition, depicted in the right panel of Figure 1, a dataset in which the uncertainties of y_i are equivalent to the sum of squared residuals (SSR) and the uncertainties of x_i are highly correlated (Figure 1) was considered using LPU.

2. Theory

2.1. Least Squares Methods

Figure 1 illustrates the regression problem of the present study, presenting the measured data vector d, d_i = (x_i, y_i), where y corresponds to the response value and x the stimuli (typically, the value of the reference material (RM)), and reference values x_i are fully correlated to each other. The minimization of regression variance for the data vector can be expressed as a generalized least squares (GLS) formulation as follows [13]:

(1)${({\mathit{\beta }}^T, {{\mathit{X}}^{*}}^T,{{\mathit{Y}}^{*}}^T)}^T=arg{min}_{{\mathit{x}}^{\mathit{*}},{\mathit{y}}^{\mathit{*}},\tilde{\mathit{\beta }}}{\left(\begin{array}{c}\mathit{x}-{\mathit{x}}^{\mathit{*}}\\ \mathit{y}-{\mathit{y}}^{\mathit{*}}\end{array}\right)}^{\mathrm{T}}{\mathit{C}}_{\mathit{d}}^{-1}\left(\begin{array}{c}\mathit{x}-{\mathit{x}}^{\mathit{*}}\\ \mathit{y}-{\mathit{y}}^{\mathit{*}}\end{array}\right),$

where the unobservable true vector describing the fitted line $\left({\mathit{x}}^{*},{\mathit{y}}^{*}\right)$ is a constraint upon $x_i^{*}=q(y_i^{*},\tilde{\mathit{\beta }}$), i = 1, ∙∙∙∙, n. $\tilde{\mathit{\beta }}$ of order 1 × (m + 1) denotes the input regression coefficient vector of the mth-order polynomial model function, and $({\mathit{\beta }}^T, {{\mathit{X}}^{*}}^T,{{\mathit{Y}}^{*}}^T)$ are auxiliary output quantities. C_d is the variance–covariance matrix of the measurement data vector d of order n × n, where C_d is assumed to have full rank. Equation (1) does not have a closed-form solution, and therefore, an iterative solving algorithm is required to derive the estimate $\widehat{\mathit{\beta }}$ that is a vector of maximum likelihood estimators under the normality assumption, that is, the Best Linear Unbiased Estimator (BLUE) [17]. Equation (1) can be reduced to yield the OLS estimator at the regression coefficient vector β in the case that the C_d is the identity matrix. The SSR at β is expressed as follows:

(2)$\mathit{S}\left(\mathit{\beta }\right)={(\mathit{y}-\mathit{X}\mathit{\beta })}^{\mathrm{T}}(\mathit{y}-\mathit{X}\mathit{\beta })$

where $\mathit{y}={\left[y_1 y_2 \cdots \cdots y_n\right]}^{\mathrm{T}}$, the design matrix X has n rows and m columns having elements of $x_n^m$, and $\mathit{y}=\mathit{X}\mathit{\beta }+\mathit{\epsilon }$. The Gauss–Markov theorem allows the unobservable random error vector ε to have same variance σ². The function S(b) is quadratic in $\mathit{\beta }$ with positive definite Hessian, and therefore, this function possesses a unique global minimum at $\mathit{\beta }$ = $\widehat{\mathit{\beta }}$, which can be given by the explicit expression as follows:

(3)$$\begin{gathered} \frac{\partial \mathit{S}}{\partial \mathit{\beta }}\left(\widehat{\mathit{\beta }}\right)={\left.\frac{d}{d\mathit{\beta }}\left({\mathit{y}}^{\mathrm{T}}\mathit{y}-{\mathit{\beta }}^{\mathrm{T}}{\mathit{X}}^{\mathrm{T}}\mathit{y}-{\mathit{y}}^{\mathrm{T}}\mathit{X}\mathit{\beta }+{\mathit{\beta }}^{\mathrm{T}}{\mathit{X}}^{\mathrm{T}}\mathit{X}\mathit{\beta }\right)\right|}_{\mathit{\beta }=\widehat{\mathit{\beta }}} \\ =-2{\mathit{X}}^{\mathrm{T}}\mathit{y}+2{\mathit{X}}^{\mathrm{T}}\mathit{X}\widehat{\mathit{\beta }}=0 \end{gathered}$$

Therefore, a BLUE estimator of regression coefficients in the OLS problem can be derived from Equation (3) by applying the inverse of the Gram matrix ${\left({\mathit{X}}^{\mathbf{T}}\mathit{X}\right)}^{-1}$ to both sides. The measurement uncertainty (errors-in-variable on abscissa of data vector d) can be considered identical to the standard deviation of fit residuals, and Equation (3) yields the OLS estimator that is the BLUE, satisfying the Gauss–Markov theorem as follows [35,36]:

(4)${\widehat{\mathit{\beta }}}_{OLS}={\left({\mathit{X}}^{\mathrm{T}}\mathit{X}\right)}^{-1}{\mathit{X}}^{\mathrm{T}}\mathit{y},$

Therefore, ${\widehat{\mathit{\beta }}}_{OLS}$ is an affine transformation of the response vector y onto a column space of the regression line by the orthogonal projection operator, ${\left({\mathit{X}}^{\mathrm{T}}\mathit{X}\right)}^{-1}{\mathit{X}}^{\mathrm{T}}$ [37]. The GLS method is required when there are errors in the independent variable (along the ordinate). However, the error propagation of the regression variance of the GLS method requires, in general, the linearization of the regression variance at the solution point of Equation (1). The following sections explore the explicit expressions for the uncertainty of fully correlated reference u(x_i) propagated via ${\widehat{\mathit{\beta }}}_{OLS}$.

2.2. Variance–Covariance Matrix of the Regression Coefficient

A full matrix expression of Equation (4) can be expanded as follows:

(5)${\left[\begin{array}{cc}\begin{array}{cc}n& {\sum }_{i=1}^n{x}_i\\ {\sum }_{i=1}^n{x}_i& {\sum }_{i=1}^n{x}_i^2\end{array}& \begin{array}{cc}\cdots \cdots & {\sum }_{i=1}^n{x}_i^m\\ \cdots \cdots & {\sum }_{i=1}^n{x}_i^{m+1}\end{array}\\ \begin{array}{cc}\begin{array}{c}⋮\\ ⋮\end{array}& \begin{array}{c}⋮\\ ⋮\end{array}\\ {\sum }_{i=1}^n{x}_i^m& {\sum }_{i=1}^n{x}_i^{m+1}\end{array}& \begin{array}{cc}\ddots & \begin{array}{c}⋮\\ ⋮\end{array}\\ \cdots \cdots & {\sum }_{i=1}^n{x}_i^{2m}\end{array}\end{array}\right]}^{-1}\left[\begin{array}{c}\begin{array}{c}\sum _{i=1}^n{y}_i\\ \sum _{i=1}^n{x}_i{y}_i\end{array}\\ \begin{array}{c}⋮\\ ⋮\end{array}\\ \sum _{i=1}^n{x}_i^m{y}_i\end{array}\right]-\left[\begin{array}{c}\begin{array}{c}{\beta }_0\\ {\beta }_1\end{array}\\ \begin{array}{c}⋮\\ ⋮\end{array}\\ {\beta }_m\end{array}\right]=\mathit{Q}$

Then, the cost function matrix Q is expressed as Q = (q₀, q₁, ……, q_m)^T = 0 with coefficient vector β taken by the OLS estimator ${\widehat{\mathit{\beta }}}_{OLS}$. Then, Equation (5) is reorganized to ${\mathit{X}}^{\mathrm{T}}\mathit{X}\mathit{\beta }-{\mathit{X}}^{\mathrm{T}}\mathit{y}=0$ to have the implicit function vector elements q_m, expressed as follows [38]:

(6a)$q_0={\beta }_0·n+{\beta }_1·{\sum }_{i=1}^n{x}_i+\cdots \cdots +{\beta }_m·{\sum }_{i=1}^n{x}_i^m-{\sum }_{i=1}^n{y}_i=0$

(6b)$q_1={\beta }_0·{\sum }_{i=1}^n{x}_i+{\beta }_1·{\sum }_{i=1}^n{x}_i^2+\cdots \cdots +{\beta }_m·{\sum }_{i=1}^n{x}_i^{m+1}-{\sum }_{i=1}^n{x}_i{y}_i=0$

(6c)$q_m={\beta }_0·{\sum }_{i=1}^n{x}_i^m+{\beta }_1·{\sum }_{i=1}^n{x}_i^{m+1}+\cdots \cdots +{\beta }_m·{\sum }_{i=1}^n{x}_i^{2m}-{\sum }_{i=1}^n{x_i^{m}y}_i=0$

Equation (6) does not require an estimation of the inverse matrix by organizing each element implicitly. The systematic equations of Equation (6) do not have an exact analytical solution. Instead, coefficient vector β fits the equation best in the sense of finding the minima of SSR as described above. According to DIN 1319-4, the uncertainty transformation matrix ${\mathit{T}}_{{\widehat{\mathit{\beta }}}_{OLS}}$ for accounting for variability in ${\widehat{\mathit{\beta }}}_{OLS}$ should be formed by the ordinate components of the regression coefficient vector β and the response vector y. To propagate the regression variance, an uncertainty propagation matrix requires the Jacobian matrices of the cost function matrix Q against coefficient vector β as an input quantity and response vector y as an output quantity. Jacobian matrices Q_i and Q_o, where the subscripts represent the input and output quantities, respectively, yield an uncertainty transformation matrix ${\mathit{Q}}_o^{-1}{\mathit{Q}}_i$. To propagate the variances–covariances of the regression coefficients, the uncertainty propagation matrix requires Jacobian matrices of the cost function matrix Q against coefficient vector β and response vector y as input and output quantities, respectively. To obtain a positive definite, the uncertainty transformation matrix for estimating the OLS regression variance is negatively given as follows [23,39]:

(7)${\mathit{T}}_{{\widehat{\mathit{\beta }}}_{OLS}}=-{\mathit{Q}}_{{\widehat{\mathit{\beta }}}_{OLS}}^{-1}{\mathit{Q}}_{\mathit{y}},$

where ${\mathit{Q}}_{\mathit{y}}$ and ${\mathit{Q}}_{{\widehat{\mathit{\beta }}}_{OLS}}$ are the Jacobian matrices of the Q matrix against the response vector y and the regression coefficient vector $\mathit{\beta }$, respectively. These two matrices are given as follows:

(8)${\mathit{Q}}_{\mathit{y}}=\left[\begin{array}{cc}\begin{array}{cc}\frac{\partial q_0}{\partial y_1}& \frac{\partial q_0}{\partial y_2}\\ \frac{\partial q_1}{\partial y_1}& \frac{\partial q_1}{\partial y_2}\end{array}& \begin{array}{cc}\cdots \cdots & \frac{\partial q_0}{\partial y_n}\\ \cdots \cdots & \frac{\partial q_1}{\partial y_n}\end{array}\\ \begin{array}{cc}\begin{array}{c}⋮\\ ⋮\end{array}& \begin{array}{c}⋮\\ ⋮\end{array}\\ \frac{\partial q_m}{\partial y_1}& \frac{\partial q_m}{\partial y_2}\end{array}& \begin{array}{cc}\ddots & \begin{array}{c}⋮\\ ⋮\end{array}\\ \cdots \cdots & \frac{\partial q_m}{\partial y_n}\end{array}\end{array}\right]=\left[\begin{array}{cc}\begin{array}{cc}-1& -1\\ {-x}_1& {-x}_2\end{array}& \begin{array}{cc}\cdots \cdots & -1\\ \cdots \cdots & {-x}_n\end{array}\\ \begin{array}{cc}\begin{array}{c}⋮\\ ⋮\end{array}& \begin{array}{c}⋮\\ ⋮\end{array}\\ {-x}_1^m& {-x}_2^m\end{array}& \begin{array}{cc}\ddots & \begin{array}{c}⋮\\ ⋮\end{array}\\ \cdots \cdots & {-x}_n^m\end{array}\end{array}\right]=-{\mathit{X}}^{\mathrm{T}}$

(9)${\mathit{Q}}_{{\widehat{\mathit{\beta }}}_{OLS}}=\left[\begin{array}{cc}\begin{array}{cc}\frac{\partial q_0}{\partial {\beta }_0}& \frac{\partial q_0}{\partial {\beta }_1}\\ \frac{\partial q_1}{\partial {\beta }_0}& \frac{\partial q_1}{\partial {\beta }_1}\end{array}& \begin{array}{cc}\cdots \cdots & \frac{\partial q_0}{\partial {\beta }_m}\\ \cdots \cdots & \frac{\partial q_1}{\partial {\beta }_m}\end{array}\\ \begin{array}{cc}\begin{array}{c}⋮\\ ⋮\end{array}& \begin{array}{c}⋮\\ ⋮\end{array}\\ \frac{\partial q_m}{\partial {\beta }_0}& \frac{\partial q_m}{\partial {\beta }_1}\end{array}& \begin{array}{cc}\ddots & \begin{array}{c}⋮\\ ⋮\end{array}\\ \cdots \cdots & \frac{\partial q_m}{\partial {\beta }_m}\end{array}\end{array}\right]=\left[\begin{array}{cc}\begin{array}{cc}n& {\displaystyle \sum _{i=1}^n}{x}_i\\ {\displaystyle \sum _{i=1}^n}{x}_i& {\displaystyle \sum _{i=1}^n}{x}_i^2\end{array}& \begin{array}{cc}\cdots \cdots & {\displaystyle \sum _{i=1}^n}{x}_i^m\\ \cdots \cdots & {\displaystyle \sum _{i=1}^n}{x}_i^{m+1}\end{array}\\ \begin{array}{cc}\begin{array}{c}⋮\\ ⋮\end{array}& \begin{array}{c}⋮\\ ⋮\end{array}\\ {\displaystyle \sum _{i=1}^n}{x}_i^m& {\displaystyle \sum _{i=1}^n}{x}_i^{m+1}\end{array}& \begin{array}{cc}\ddots & \begin{array}{c}⋮\\ ⋮\end{array}\\ \cdots \cdots & {\displaystyle \sum _{i=1}^n}{x}_i^{2m}\end{array}\end{array}\right]={\mathit{X}}^{\mathrm{T}}\mathit{X}$

The variance matrix of the regression coefficient ${\mathit{V}}_{{\widehat{\mathit{\beta }}}_{OLS}}$ is calculated as follows:

(10)${\mathit{V}}_{{\widehat{\mathit{\beta }}}_{OLS}}={\mathit{T}}_{{\widehat{\mathit{\beta }}}_{OLS}}\mathit{\Omega }{\mathit{T}}_{{\widehat{\mathit{\beta }}}_{OLS}}^{\mathrm{T}},$

where Ω is the diagonal residual matrix with identical entries of s² [2]. s² is a variance of fit residual of the sampled data from the unobservable ‘true’ regression line, of which the regression coefficient vector is b. Notably, the statistical parameter for describing the extent of dispersion of an unobservable error e_i ~ y_i − q(x_i,b) that has a constant variance σ², in which the error vectors are normally distributed $\mathit{\epsilon }~N(0, {\sigma }^2)$ (zero mean and standard deviation equaling σ). The standard deviation of the fit residuals of the sampled data, s, can be expressed by the following:

(11)$s=\sqrt{\frac{{\sum }_{i=1}^n{\left\{y_i-\left(b_0+b_1{x}_i+\cdots +b_m{x}_i^m\right)\right\}}^2}{n-m-1},}$

where the degrees-of-freedom correction in the denominator is restricted by the polynomial order of the fitting model m and the number of data vectors n. ${\widehat{\mathit{\beta }}}_{OLS}$ allows the residual to be an observable estimate of the unobservable error ε. The probability distribution of the residual error converges to a normal distribution following the Gauss–Markov assumption. The variance–covariance matrix of the OLS estimator ${\mathit{V}}_{{\widehat{\mathit{\beta }}}_{OLS}}$ is expressed as follows:

(12)${\mathit{V}}_{{\widehat{\mathit{\beta }}}_{OLS}}=s^2{\mathit{G}}^{-1},$

where G is the Gram matrix $\left({\mathit{X}}^{\mathrm{T}}\mathit{X}\right)$ of Equation (9). Each diagonal entry corresponds to the variance of the regression coefficient of each term, whereas the off-diagonal entries correspond to the covariance between regression coefficients.

2.3. Propagation of Variance–Covariance of the Regression Coefficients

For calibration, the measured response value of the target y₀ is inversely evaluated to obtain the corresponding stimulus value x₀, which can be numerically rooted by using Newton’s method. To propagate the measurement uncertainty of u(y₀) onto u(x₀), the inverse evaluation requires the estimated regression coefficient vector ${\widehat{\mathit{\beta }}}_{OLS}$ and estimated x₀ as input and output quantities, respectively. The Jacobian matrix against the input and output quantities (${\mathit{q}}_{{\widehat{\mathit{\beta }}}_{OLS}}$ and $q_x$, respectively) can then be expressed as follows:

(13)${\mathit{q}}_{{\widehat{\mathit{\beta }}}_{OLS}}={\left[\frac{\partial q_0}{\partial {\beta }_0} \frac{\partial q_0}{\partial {\beta }_1} \cdots \cdots \frac{\partial q_0}{\partial {\beta }_m}\right]}_{x_0}=\left[1 x_0 \cdots \cdots x_0^m\right]$

(14)$q_x={\left[\frac{\partial q_0}{\partial x}\right]}_{x_0}={\beta }_1+2·{\beta }_2·x_0+\cdots \cdots +{m·\beta }_m·x_0^{m-1}=q^{’}\left(x_0\right)$

Here, $q^{’}\left(x_0\right)$ is the slope of the tangent of the determined regression line in the case of the second-order polynomial fit model at the target coordinate (x₀, y₀). This implies that the regression variance will be weighted by the slope of the tangent of the non-first-order calibration curve. The uncertainty transformation matrix for the regression variance toward the abscissa of the output quantity ${\mathit{t}}_{{\widehat{\mathit{\beta }}}_{OLS}}$ is then given as follows:

(15)${\mathit{t}}_{{\widehat{\mathit{\beta }}}_{OLS}}=q_x^{-1}{\mathit{q}}_{{\widehat{\mathit{\beta }}}_{OLS}}$

To transform the variance–covariance matrix of the OLS estimator ${\mathit{V}}_{{\widehat{\mathit{\beta }}}_{\mathit{O}\mathit{L}\mathit{S}}}$ with the mth-order polynomial model function (Equation (12)) onto the abscissa corresponding to the output quantity x₀, the uncertainty propagation matrix ${\mathit{t}}_{{\widehat{\mathit{\beta }}}_{OLS}}$ (Equation (15)) is bracketed according to the LPU as follows:

(16)${\mathit{V}}_{\mathit{m}.{\widehat{\mathit{\beta }}}_{OLS}}={\mathit{t}}_{{\widehat{\mathit{\beta }}}_{OLS}}{\mathit{V}}_{{\widehat{\mathit{\beta }}}_{OLS}}{\mathit{t}}_{{\widehat{\mathit{\beta }}}_{OLS}}^{\mathbf{T}}$

(17)$={\left(\frac{s}{q^{’}\left(x_0\right)}\right)}^2\left[1 x x^2 \cdots \cdots x^m\right]\left[\begin{array}{cc}\begin{array}{cc}{g}_{00}& g_{01}\\ g_{10}& g_{11}\end{array}& \begin{array}{cc}\cdots \cdots & g_{0m}\\ \cdots \cdots & g_{1m}\end{array}\\ \begin{array}{cc}\begin{array}{c}⋮\\ ⋮\end{array}& \begin{array}{c}⋮\\ ⋮\end{array}\\ g_{m0}& g_{m1}\end{array}& \begin{array}{cc}\ddots & \begin{array}{c}⋮\\ ⋮\end{array}\\ \cdots \cdots & g_{mm}\end{array}\end{array}\right]{\left[1 x x^2 \cdots \cdots x^m\right]}^{\mathrm{T}}$

where g_ij is an element of the inverse Gram matrix G⁻¹. In this study, $q^{’}\left(x_0\right)$ has only positive values. The propagated uncertainty of regression coefficients into the calibrated value x₀ is given as follows:

(18)$u_m\left(x_0\right)=\frac{s}{q^{’}\left(x_0\right)}\sqrt{{\displaystyle \sum _{i=0}^m}{\displaystyle \sum _{j=0}^m}{g}_{ij}·x_0^{(i+j)}}$

Because the highest-order element of the inverse Gram matrix $g_{mm}$ is 2m order, for the second-order polynomial, u₂(x₀) is shown as a quartic curve, as shown in Figure 2. The variance of each component is estimated using the dataset in Table 1.

2.4. Propagation of the Measurement Uncertainty of the Calibration Target

To propagate the measurement uncertainty of the calibration target, its uncertainty value is treated as an input value of the reverse evaluation. The partial derivative of the cost function at y₀ can be expressed as follows:

(19)$q_y=\left[\frac{\partial q_0}{\partial y_0}\right]=-1$

Then, the uncertainty transformation matrix for the measurement uncertainty of the unknown sample can then be obtained to obtain the positive definite as follows:

(20)$t_t=-q_x^{-1}{q}_y=\frac{1}{q^{’}\left(x_0\right)}$

The propagation of the measurement uncertainty of the target u(y₀) yields the target variance V_t as follows:

(21)$V_t=t_t{u}^2\left(y_0\right)t_t^{\mathrm{T}}=\left[\frac{1}{q^{’}\left(x_0\right)}·\left(\frac{s^2}{p}\right)·\frac{1}{q^{’}\left(x_0\right)}\right],$

where s is the measurement uncertainty of the calibration target, which, in general, can be approximated to the same extent of sampling variance for the regression data vector, and p is the number of repetitions of the target measurement. The measurement variance s² can be weighted by the number of repetitions p according to the central limit theorem. The target uncertainty can then be expressed as follows:

(22)$u_t\left(x_0\right)=\frac{s}{q^{’}\left(x_0\right)}\sqrt{\frac{1}{p}}$

2.5. Propagation of Uncertainty of the Independent Variable

The uncertainty of the independent variable was propagated using same method for the estimation of the regression variance. To define an uncertainty transformation matrix toward the abscissa, if there are no correlations between x_i and y_i, the Jacobian matrix of the cost function against the stimulus x_i is defined as follows:

(23)${\mathit{Q}}_{\mathit{x}}=\left[\begin{array}{cc}\begin{array}{cc}\frac{\partial q_0}{\partial x_1}& \frac{\partial q_0}{\partial x_2}\\ \frac{\partial q_1}{\partial x_1}& \frac{\partial q_1}{\partial x_2}\end{array}& \begin{array}{cc}\cdots \cdots & \frac{\partial q_0}{\partial x_n}\\ \cdots \cdots & \frac{\partial q_1}{\partial x_n}\end{array}\\ \begin{array}{cc}\begin{array}{c}⋮\\ ⋮\end{array}& \begin{array}{c}⋮\\ ⋮\end{array}\\ \frac{\partial q_m}{\partial x_1}& \frac{\partial q_m}{\partial x_2}\end{array}& \begin{array}{cc}\ddots & \begin{array}{c}⋮\\ ⋮\end{array}\\ \cdots \cdots & \frac{\partial q_m}{\partial x_n}\end{array}\end{array}\right]$

The uncertainty transformation matrix of the reference uncertainty to the abscissa is given as follows:

(24)${\mathit{T}}_{\mathit{x}.{\widehat{\mathit{\beta }}}_{OLS}}={\mathit{Q}}_{{\widehat{\mathit{\beta }}}_{OLS}}^{-1}{\mathit{Q}}_{\mathit{x}},$

where ${\mathit{Q}}_{{\widehat{\mathit{\beta }}}_{\mathit{O}\mathit{L}\mathit{S}}}$ is given in Equation (9). As an analytical derivation of the above matrix is nontrivial. Then, the variances and covariances of the references are propagated back to the abscissa using following equation:

(25)${\mathit{V}}_{\mathit{x}.{\widehat{\mathit{\beta }}}_{OLS}}={\mathit{T}}_{\mathit{x}.{\widehat{\mathit{\beta }}}_{OLS}}{\mathit{V}}_{\mathit{r}\mathit{e}\mathit{f}}{\mathit{T}}_{\mathit{x}.{\widehat{\mathit{\beta }}}_{OLS}}^{\mathrm{T}},$

where the variance–covariance matrix of the reference ${\mathit{V}}_{\mathit{r}\mathit{e}\mathit{f}}$ is given as follows:

(26)${\mathit{V}}_{\mathit{r}\mathit{e}\mathit{f}}={\mathit{D}}_{\mathit{r}\mathit{e}\mathit{f}}{\mathit{R}}_{\mathit{r}\mathit{e}\mathit{f}}{\mathit{D}}_{\mathit{r}\mathit{e}\mathit{f}}$

where ${\mathit{R}}_{\mathit{r}\mathit{e}\mathit{f}}$ is the correlation matrix in which the elements are r_x, and ${\mathit{D}}_{\mathit{r}\mathit{e}\mathit{f}}$ is the diagonal matrix of the reference uncertainty. The uncertainty of the calibrated value of the target u_x(x₀) estimated from the OLS estimator is given as follows:

(27)$u_x\left(x_0\right)=\frac{1}{q^{’}\left(x_0\right)}\sqrt{{\displaystyle \sum _{i=0}^m}{\displaystyle \sum _{j=0}^m}{\mathit{V}}_{\mathit{r}\mathit{e}\mathit{f}}\left(i,j\right)·{x_0}^{\left(i+j\right)}}$

Unfortunately, the formulaic expression of the explicit form of u_x(x₀) is nontrivial for higher-order polynomials. In a condition of no correlation (r_x = 0), a curvature of u_x(x₀) is quartically variated for the quadratic calibration function because of the nature of the regression variance–covariance matrix ${\mathit{V}}_{\mathit{r}\mathit{e}\mathit{f}}$, of which the highest-order element is $x^{2m}$. (Figure 2).

As a quantitative investigation of the correlation coefficient is beyond the limitation of the realization in many cases, a conservative approach is required. In an extremely conservative case, for an all-ones correlation matrix (all $r_{x_{ij}}=1$), the measurement uncertainty of the target along the abscissa is coherently dispersed, as presented below:

(28)$u_x\left(x_0\right)~\frac{1}{2}\left|\begin{array}{c}q\left(x_1+u\left(x\right), x_2+u\left(x\right), \cdots \cdots , { x}_n+u\left(x\right), { y}_1, y_2, \cdots \cdots y_n, { y}_0\right)\\ -q\left(x_1+u\left(x\right), x_2+u\left(x\right), \cdots \cdots , { x}_n+u\left(x\right), { y}_1, y_2, \cdots \cdots y_n, { y}_0\right)\end{array}\right|,$

where q is the cost function of the fitted model of the OLS estimator and u(x) is the averaged uncertainty of u(x_i). It was assumed that the uncertainty values of references are close to each other, as shown in Table 1 and Figure 3. The fitted line spans out in the positive and negative directions by the extent of the width of the normal distribution, of which uncertainty is u(x) (Figure 1). As expected, u_x(x₀) is smoothly variated as a function of x₀ when all correlations are fully positive (i.e., r_x = 1) based on the assumption that u(x₀) is constant with respect to x₀ (Figure 2). The variation rate is lower at large values of x₀ because the slope of the tangent for the used dataset is higher at such values. In the case that x₀ is centered within the reference scale, Equation (26) can be simplified as follows:

(29)$u_x\left(x_0\right)\cong \overline{u_{ref}}$

Table 1

Datasets to test the proposed method in this study. Note that all u(y_i) values are close to the standard deviation of the fit residual of the corresponding dataset. The uncertainties of reference values were assumed to be fully correlated with each other.

Reference (xi)	u (xi)	Response	u(yi)
			Set A	Set B	Set C
317.35	0.231	0.983	0.00098	0.00098
321.68	0.206	0.992	0.00074	0.00074	0.00041
325.11	0.208	1.001	0.00080
328.54	0.209	1.009	0.00085	0.00085
331.92	0.233	1.019	0.00077	0.00077	0.00042
335.24	0.243	1.028	0.00085
338.56	0.252	1.037	0.00093	0.00093
338.93	0.231	1.037	0.00106	0.00106	0.00058
344.98	0.230	1.055	0.00089
351.03	0.229	1.072	0.00072	0.00072
351.66	0.246	1.072	0.00082	0.00082	0.00045
356.06	0.239	1.085	0.00090
360.45	0.232	1.098	0.00098	0.00098
360.97	0.253	1.099	0.00069	0.00069	0.00038
364.85	0.241	1.111	0.00085
Target	(p = 1)	1.003	n = 15	n = 10	n = 5

3. Discussions

Artificial datasets mimicking the measurements of ambient-level N₂O using a gas chromatograph with an electron capture detector (GC-ECD) were utilized to evaluate the uncertainty of calibrated value x₀ (Table 1). The artificial data are dispersed along the ordinate to have uncertainties, u(y_i), which are closely set to the standard deviation of fit residuals s = 0.00073, 0.00085, and 0.00041, respectively, for dataset A, B, and C. This aspect implies that u(y_i) are considered as the sampling uncertainties from the $\left({\mathit{x}}^{*},{\mathit{y}}^{*}\right)$ which was modeled by the OLS estimator. It should be noted that the OLS and the WLS methods yield similar results in the regression variance, where the variance–covariance matrix of the WLS estimator ${\mathit{V}}_{{\widehat{\mathit{\beta }}}_{WLS}}$ was computed by the equation of ${\left({\mathit{X}}^{\mathrm{T}}\mathit{W}\mathit{X}\right)}^{-1}$, where W is the weighting matrix and its elements are inversed u(y_i) [5]. The u_x(x₀) curves for r_x = 1 and r_x = 0 cross outside of the range of the reference scale, implying that a conservative approach can be achieved by considering the correlation effect of reference values within the reference scale (Figure 2). The post-evaluation of the correlation factors among reference values is hard to achieve; one may assume r_x = 1 to avoid underestimating the calibration uncertainty due to uncertainties of reference values.

To verify the derived explicit equation of each uncertainty value, the Monte Carlo (MC) simulation was conducted employing identical datasets. In the MC method, individual data points were randomly sampled from normally distributed reference values (x_i), which were fully correlated to each other (all $r_{x_{ij}}=1$), and normally distributed response values (y_i). The uncertainties of both x_i and y_i, denoted as u(x_i) and u(y_i), respectively, for each dataset are detailed in Table 1. A total of 10,500 d_i = (x_i, y_i) were randomly sampled from the population to compute the uncertainty value of the calibration target. Running the OLS estimator on each sampled data vector yielded the single regression coefficient vector β and calibrated value x₀. An average and standard deviation of each parameter resulted in the corresponding value and its uncertainty. Because it is a valid assumption that a result by the MC simulation will be in close proximity to the unbiased limit of variance of the resulting distribution [40], a comparison between the x-correlated OLS method (_xcOLS) and the MC method allowed us to test the reliability of the _xcOLS method (Figure 4). With LPU, it is hard to handle the GLS problem in which the SR matrix cannot be reduced to the implicit model equations [33]. Therefore, approximating the uncertainties according to the LPU is prevalent in many cases. The uncertainties u₂(x₀), u_t(x₀), and u_x(x₀) from both the _xcOLS and MC methods exhibit a high degree of agreement for various datasets in ranges of variance of fit residuals s². Covariances stemming from the coherent dispersion of reference values were observed to be minor, as seen in the comparison between the MC and _xcOLS methods in the datasets (Table 2 and Table 3 and Figure 4).

4. Conclusions

In the field of chemical and biological analysis, the OLS method is a prevalent choice for calibration tasks involving highly accurate RMs, particularly when the negligible uncertainty of references obviates the need for their consideration in the regression analysis. However, numerous scenarios, such as those encountered in graphite furnace atomic absorption spectroscopy, X-ray fluorescence spectroscopy, and spark emission spectroscopy, deviate from the simplicity of the OLS approach due to the presence of considerable variance–covariance in the independent variables, namely reference values. To address this challenge, an approximated propagation method for the variance–covariance of the data matrix becomes imperative. This study addresses a specific calibration scenario wherein the independent variables exhibit high correlation, and uncertainties along the ordinate are equivalent to the standard deviation of the fit residuals. If the uncertainties of the independent variables are deemed negligible, this scenario aligns with the traditional OLS case. Conceptually illustrated in Figure 1, the uncertainties of fully correlated independent variables were dispersed orthogonally along the abscissa. The validity of this concept is substantiated through the step-by-step derivation of explicit equations. This study meticulously budgets the contributions of variances arising from regression coefficients, target measurements, and references, offering an analytical framework for comprehending uncertainty propagation for the calibration scenario. The calibration task employing the _xcOLS and the MC methods was exemplified using artificial datasets, with standard deviations of fit residuals s being 0.00073, 0.00085, and 0.00041. The resulting uncertainty values, as a metric for comparing results, were estimated by the _xcOLS and the MC methods to show good agreement. This good agreement underscores the validity of the proposed approach in propagating uncertainties associated with highly correlated independent variables using the OLS estimator.

Author Contributions

J.S.L.: Formal Analysis, Investigation, Resources, Validation, Visualization, Writing—Original Draft, Review, and Editing, Funding Acquisition, Project Administration, Y.D.K.: Data Curation, Formal Analysis, J.-C.W.: Conceptualization, Methodology, Data Curation, Formal Analysis, Validation, Investigation, Supervision. All authors have read and agreed to the published version of the manuscript.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Data are contained within the article.

Acknowledgments

The authors thank Sangil Lee for the discussions in the initial stage of this study.

Conflicts of Interest

The authors declare no conflicts of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. (Left) Notation of uncertainties and errors addressed in this study. u(xi) and u(yi) denote input variances in the reference and corresponding response values, respectively. The residual error εi represents the deviation between its theoretical and input values along the ordinate. The measured response of the target is y0, which is calibrated to x0. (Right) Conceptual representation of linear dispersion of the calibrated value due to the fully correlated reference values. The tested dataset for the proposed calibration method is given in Table 1 in Section 2.5.

Figure 2. (Left) Propagated reference variance corresponds to the calibrated value when there is no correlation, rx = 0 (blue); and correlation is unity, rx = 1 (red). (Right) Uncertainties propagated from the uncertainties of reference ux(x0) (gray); target measurement, ut(x0) (magenta); and regression, u2(x0) (cyan). The combined standard uncertainty u(x0) is represented by the green line. For detailed discussions of ux(x0) and ut(x0), see Section 2.4 and Section 2.5.

Figure 3. Graphical representation of the artificial dataset A given in Table 1. Note that the uncertainties of reference values and responses are amplified three times to be seen. The uncertainties of reference values were assumed to be fully correlated with each other.

Figure 4. A comparison of the variance values corresponding to the reference, [Forumla omitted. See PDF.] (gray); the target measurement, [Forumla omitted. See PDF.] (magenta); and regression, [Forumla omitted. See PDF.] (cyan). Note that the number of target measurement repetitions p is 1. Calibrated values obtained through each method are indicated on top of the respective bars.

Table 2

Estimated regression coefficients for each dataset. Associated uncertainties are given in the parentheses.

Dataset	Method	Regression Coefficients (Uncertainty)
		β 0	β 1	β 2
A	MC	0.7643 (0.1441)	−0.001084 (0.000844)	0.0000056 (0.0000012)
	xcOLS	0.7642 (0.1445)	−0.001084 (0.000846)	0.0000056 (0.0000012)
B	MC	0.8084 (0.2368)	−0.001338 (0.001394)	0.0000059 (0.0000020)
	xcOLS	0.8080 (0.1724)	−0.001336 (0.001015)	0.0000059 (0.0000015)
C	MC	0.6767 (0.1866)	−0.000565 (0.001089)	0.0000048 (0.0000016)
	xcOLS	0.6785 (0.1355)	−0.000576 (0.000794)	0.0000048 (0.0000012)

Table 3

The combined standard uncertainty of the determined value of the calibration target under the assumption of unitive correlation of the reference values. β₁ is the slope of a fitted straight line, and $\overline{x}$ is the mean of the reference value. The measurement uncertainty of the target is assumed to be at a level similar to the measurement uncertainty of the reference. The unitive correlation among the reference variances is approximated to $u_x\left(x_0\right) \cong \overline{u_{ref}}$.

Polynomial Order	$\mathit{V}\mathit{a}\mathit{r}\left({\mathit{x}}_0\right)$
m > l	$\frac{s^2}{{\left(q^{’}\left(x_0\right)\right)}^2}\left(\frac{1}{p}+{\displaystyle \sum _{i=0}^m}{\displaystyle \sum _{j=0}^m}{g}_{ij}·x^{(i+j)}\right)+{\overline{u_{ref}}}^2$
m = 1	$\frac{s^2}{{\beta }_1^2}\left(\frac{1}{p}+\frac{1}{n}+\frac{{\left(x_0-\overline{x}\right)}^2}{\sum _{i=1}^n{\left(x_i-\overline{x}\right)}^2}\right)+{\overline{u_{ref}}}^2$

References

1. Tellinghuisen, J. Least squares methods for treating problems with uncertainty in x and y. Anal. Chem.; 2020; 92, pp. 10863-10871. [DOI: https://dx.doi.org/10.1021/acs.analchem.0c02178] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32678579]

2. Tellinghuisen, J. Calibration: Detection, quantification, and confidence limits are (almost) exact when the data variance function is known. Anal. Chem.; 2019; 91, pp. 8715-8722. [DOI: https://dx.doi.org/10.1021/acs.analchem.9b00119] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/31180654]

3. Riu, J.; Ruis, F.X. Assessing the accuracy of analytical methods using linear regression with errors in both axes. Anal. Chem.; 1996; 68, pp. 1851-1857. [DOI: https://dx.doi.org/10.1021/ac951217s]

4. Tellinghuisen, J. A simple, All-purpose nonlinear algorithm for univariate calibration. Analyst; 2000; 125, pp. 1045-1048. [DOI: https://dx.doi.org/10.1039/b002171g]

5. Tellinghuisen, J. Weighted least-squares in calibration: What difference does it make?. Analyst; 2007; 132, pp. 536-543. [DOI: https://dx.doi.org/10.1039/b701696d] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/17525810]

6. Tellinghuisen, J. Least squares in calibration: Dealing with uncertainty in x. Analyst; 2010; 135, pp. 1961-1969. [DOI: https://dx.doi.org/10.1039/c0an00192a]

7. Tellinghuisen, J. Least-squares analysis of data with uncertainty in x and y: A Monte Carlo methods comparison. Chem. Int. Lab. Syst.; 2010; 103, pp. 160-169. [DOI: https://dx.doi.org/10.1016/j.chemolab.2010.07.003]

8. Francios, N.; Govaerts, B.; Boulanger, B. Optimal designs for inverse prediction in univariate nonlinear calibration models. Chem. Int. Lab. Syst.; 2004; 74, pp. 283-292. [DOI: https://dx.doi.org/10.1016/j.chemolab.2004.05.005]

9. Tellinghuisen, J. Simple algorithms for nonlinear calibration by the classical and standard additions methods. Analyst.; 2005; 130, pp. 370-378. [DOI: https://dx.doi.org/10.1039/b411054d]

10. Mulholland, M.; Hibbert, D.B. Linearity and the limitations of least squares calibration. J. Chrom. A.; 1997; 762, pp. 73-82. [DOI: https://dx.doi.org/10.1016/S0021-9673(96)00976-4]

11. Rehman, Z.; Zhang, S. Three-dimensional elasto-plastic damage model for gravelly soil-structure interface considering the shear coupling effect. Comput. Geotech.; 2021; 129, 103868. [DOI: https://dx.doi.org/10.1016/j.compgeo.2020.103868]

12. Rehman, Z.; Khalid, U.; Ijaz, N.; Mujtaba, H.; Haider, A.; Farooq, K.; Ijaz, Z. Machine learning-based intelligent modeling of hydraulic conductivity of sandy soils considering a wide range of grain sizes. Eng. Geol.; 2022; 311, 106899. [DOI: https://dx.doi.org/10.1016/j.enggeo.2022.106899]

13. Ni, Y.; Wang, H.; Hu, S. Linear Fitting of Time-Varying Signals in Static Noble Gas Mass Spectrometry Should Be Avoided. Anal. Chem.; 2023; 95, pp. 3917-3921. [DOI: https://dx.doi.org/10.1021/acs.analchem.3c00128] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/36786555]

14. Bremser, W.; Viallon, J.; Wielgosz, R.I. Influence of correlation on the assessment of measurement result compatibility over a dynamic range. Metrologia; 2007; 44, pp. 495-504. [DOI: https://dx.doi.org/10.1088/0026-1394/44/6/009]

15. Cox, M.G.; Jones, H.M. An algorithm for least-squares circle fitting to data with specified uncertainty ellipses. IMA J. Numer. Anal.; 1981; 1, pp. 3-22. [DOI: https://dx.doi.org/10.1093/imanum/1.1.3]

16. Tellinghuisen, J. Goodness-of-Fit Tests in Calibration: Are They Any Good for Selecting Least-Squares weighting formula. Anal. Chem.; 2022; 94, pp. 15997-16005. [DOI: https://dx.doi.org/10.1021/acs.analchem.2c02904]

17. Milton, M.J.T.; Harris, P.M.; Smith, I.M.; Brown, A.S.; Goody, B.A. Implementation of a generalized least-squares method for determining calibration curves from data with general uncertainty structures. Metrologia; 2006; 43, pp. S291-S298. [DOI: https://dx.doi.org/10.1088/0026-1394/43/4/S17]

18. Markovsky, I.; Van Huffel, S. Overview of total least-squares methods. Signal Process; 2007; 87, pp. 2283-2302. [DOI: https://dx.doi.org/10.1016/j.sigpro.2007.04.004]

19. Malengo, A.; Pennecchi, F. A weighted total least-squares algorithm for any fitting model with correlated variables. Metrologia; 2013; 50, pp. 654-662. [DOI: https://dx.doi.org/10.1088/0026-1394/50/6/654]

20. Krystek, M.; Anton, M. A weighted total least-squares algorithm for fitting a straight line. Meas. Sci. Technol.; 2007; 18, pp. 3438-3442. [DOI: https://dx.doi.org/10.1088/0957-0233/18/11/025]

21. Amiri-Simkooei, A.R.; Zangeneh-Nejad, F.; Asgari, J. On the covariance matrix of weighted total least-squares estimates. J. Surv. Eng.; 2016; 142, 04015014. [DOI: https://dx.doi.org/10.1061/(ASCE)SU.1943-5428.0000153]

22. Krystek, M.; Anton, M. A least-squares algorithm for fitting data points with mutually correlated coordinates to a straight line. Meas. Sci. Technol.; 2011; 22, 035101. [DOI: https://dx.doi.org/10.1088/0957-0233/22/3/035101]

23. Bremser, W.; Hasselbarth, W. Controlling uncertainty in calibration. Anal. Chim. Acta; 1997; 348, pp. 61-69. [DOI: https://dx.doi.org/10.1016/S0003-2670(97)00205-5]

24. Forbes, A.B. Generalised regression problems in metrology. Numer. Algorithms; 1993; 5, pp. 523-534. [DOI: https://dx.doi.org/10.1007/BF02108667]

25. ISO 6143:2001(E) Gas Analysis—Comparison Methods for Determining and Checking the Composition of Calibration Gas Mixtures, 2nd ed.; ISO: Geneva, Switzerland, 2001.

26. ISO/TS 28037:2010(E) Determination and Use of Straight Line Calibration Functions, 1st ed.; ISO: Geneva, Switzerland, 2010.

27. DIN 1319-4 1985 Basic Concepts of Measurements; Treatment of Uncertainties in the Evaluation of Measurements; DIN: Berlin, Germany, 1985.

28. Malengo, A.; Pennecchi, F.; Spazzini, P.G. Calibration Curve Computing (CCC) software v2.0: A new release of the INRIM regression tool. Meas. Sci. Technol.; 2015; 31, 114004. [DOI: https://dx.doi.org/10.1088/1361-6501/ab7d6e]

29. Smith, I.M.; Onakunle, F.O. XLGENLINE, Software for Generalized Least Squares Fitting. 2007; Available online: https://www.npl.co.uk/resources/software/xlgenline-and-xgenline (accessed on 27 February 2024).

30. XGENLINE. Software for Generalized Least Squares Fitting. Available online: https://www.npl.co.uk/resources/software/xgenline (accessed on 27 February 2024).

31. Software to Support ISO/TS 28037:2010(E). Available online: https://www.npl.co.uk/resources/software/iso-ts-28037-2010e (accessed on 27 February 2024).

32. Bremser, W. Calibration Tool B LEAST Software Supporting Implementation of ISO Standard 6143. 2001; Available online: https://www.iso.org/obp/ui/es/#iso:std:iso:6143:dis:ed-3:v1:en:sec:C (accessed on 27 February 2024).

33. Klauenberg, K.; Martens, S.; Bošnjaković, A.; Cox, M.G.; van der Veen, A.M.H.; Elster, C. The GUM perspective on straight-line errors-in-variables regression. Measurement; 2022; 187, 110340. [DOI: https://dx.doi.org/10.1016/j.measurement.2021.110340]

34. Khalid, U.; Rehman, Z.; Liao, C.; Farooq, K.; Mujtaba, H. Compressibility of Compacted Clays Mixed with a Wide Range of Bentonite for Engineered Barriers. Arab. J. Sci. Eng.; 2019; 44, pp. 5027-5042. [DOI: https://dx.doi.org/10.1007/s13369-018-03693-7]

35. Johnson, R.A.; Wichern, D.W. Applied Multivariate Statistical Analysis; 6th ed. Pearson: London, UK, 2002; Chapter 7 pp. 360-538.

36. Hill, C.; Griffiths, W.E.; Lim, G.C. Principles of Econometrics; 4th ed. John Wiley & Sons: Hoboken, NJ, USA, 1971; pp. 119-124.

37. Rao, C.R.; Toutenburg, H.; Shalabh,; Heumann, C. Linear Models and Generalizations; 3rd ed. Springer: Berlin/Heidelberg, Germany, 2007; pp. 36-38.

38. Walpole, R.E.; Myers, R.H.; Myers, S.L.; Ye, K.E. Probability and Statistics for Engineers and Scientists; 9th ed. Pearson: London, UK, 2014; 501p.

39. DIN 1319-4:1999-02 Fundamentals of Metrology: Part 4. Evaluation of Measurements; Uncertainty of Measurement; DIN: Berlin, Germany, 1999.

40. ISO/IEC Guide 98-3:2008/Suppl.1:2008 Propagation of Distributions Using a Monte Carlo Method; ISO: Geneva, Switzerland, 2008.