Uncertainties Associated with Arithmetic Map Operations in GIS

Abstract

Arithmetic map operations are very common procedures used in GIS to combine raster maps resulting in a new and improved raster map. It is essential that this new map be accompanied by an assessment of uncertainty. This paper shows how we can calculate the uncertainty of the resulting map after performing some arithmetic operation. Actually, the propagation of uncertainty depends on a reliable measurement of the local accuracy and local covariance, as well. In this sense, the use of the interpolation variance is proposed because it takes into account both data configuration and data values. Taylor series expansion is used to derive the mean and variance of the function defined by an arithmetic operation. We show exact results for means and variances for arithmetic operations involving addition, subtraction and multiplication and that it is possible to get approximate mean and variance for the quotient of raster maps.

Key words
GIS; Interpolation variance; Map Algebra; Propagation of Uncertainty; Taylor series

INTRODUCTION

A map resulting from interpolation of field data must have some assessment of uncertainty (^{e.g. Heuvelink et al. 1989, Crosetto et al. 2000, Atkinson and Foody 2002}). When estimates are derived from field data, they have associated uncertainties caused by spatial variation of continuous variables, interactions among them and the effect of neighbor data (^{Wang et al. 2005}). Generally, kriging is used to predict values of the variable of interest at unsampled locations, because it provides an assessment of uncertainty as given by the kriging variance. However, the kriging variance does not depend on real local data values (^{Atkinson and Foody 2002}) and therefore this is not a reliable measure of local accuracy. According to Wang et al. (2005), local estimates are strongly associated with neighbor data. Actually, the kriging variance is just a measure of the ranking of data configurations given by the variogram model (^{Journel and Rossi 1989}). A measure of local accuracy should take into account the data value. In this sense, Yamamoto (²⁰⁰⁰) proposed interpolation variance as a measure of the reliability of kriging estimates taking into account both data configuration and data values.

Moreover, interpolation variance can also be used for quantifying uncertainty associated with interpolation of categorical variable (^{Yamamoto et al. 2012}). Therefore, interpolation variance is valid for both continuous and discrete variables. Actually, the interpolation variance can be used with any interpolation method based on weighted average formula (for instance, inverse of distance weighting, multiquadric equations, etc.). Yamamoto et al. (2012) used interpolation variance as a measure of uncertainty associated with multiquadric interpolation. Because interpolation is based on limited information given by a sample, the resulting raster map is always subject to uncertainty. Usually, the main source of uncertainty is related to a lack of knowledge, but there are other factors affecting the final uncertainty in a raster map. Burrough (1986) lists 14 factors affecting the uncertainty. These factors were classified into three groups: I) obvious sources of error; II) errors resulting from natural variation of original measurements; III) errors arising through processing (^{Burrough 1986}). Wellmann et al. (²⁰¹⁰) classified sources of uncertainty affecting geological data into three categories: I) imprecision and measurement error; II) stochastic nature of the geological variable; III) imprecise knowledge. Considering the field data with negligible uncertainty, the measured variance at unsampled location takes into account uncertainty due to lack of knowledge because of insufficient sampling (three factors of Burrough (1986), affect the interpolation variance: density of observations; natural variation and interpolation).

Considering that we have an interpolated raster map and its associated uncertainty, we can perform arithmetic operations with raster maps. Arithmetic map operations are very common procedures used in geographic information system. In general, two raster maps are combined using arithmetic operators (addition, subtraction, multiplication and division) to derive a new raster map. Let X and Y be interpolated raster maps from the same data set and let us perform some arithmetic operation between them. For each variable we know not only the interpolated value but also the uncertainty given by the interpolation variance. Thus, we have to calculate the resulting raster map after an arithmetic operation and the uncertainty as well. Uncertainty is magnified when cartographic overlay operations involve more than two steps (Burrough 1986). We can define a function of variables X and Y: $f(x,y)$ that can be expanded in a Taylor series about the mean values of X and Y. From this expansion we can apply the mathematical expectation operator to find the mean value of the function and the definition of the variance as a measurement of the uncertainty. Addition and subtraction are very simple arithmetic operators, but multiplication and division are much more complex operations and they require correcting the combined estimate and the calculation of the variance is even more complicated.

ESTIMATES AND VARIANCES OF ARITHMETICALLY COMBINED VARIABLES

Given a random variable X, the mathematical expectation is:

$$E\left[X\right]=\frac{1}{n}\sum _{i=1}^n{X}_i$$

The variance is a measure of the uncertainty associated with the mathematical expectation:

$$𝑉𝑎𝑟\left[X\right]=\frac{1}{n}\sum _{i=1}^n(X_i-E\left[X\right]{)}^2$$

which can be written as:

$$𝑉𝑎𝑟\left[X\right]=E{[X-E\left[X\right]]}^2=E\left[X^2\right]-{\left(E\left[X\right]\right)}^2$$

Given two random variables (X and Y), we can compute the covariance as:

$$𝐶𝑜𝑣(X,Y)=E\left[𝑋𝑌\right]-E\left[X\right]E\left[Y\right]$$

The correlation coefficient as a measure of the mutual relationship between two random variables is calculated as:

$${\rho }_{X,Y}=\frac{𝐶𝑜𝑣(X,Y)}{\sqrt{𝑉𝑎𝑟\left[X\right]𝑉𝑎𝑟\left[Y\right]}}=\frac{𝐶𝑜𝑣(X,Y)}{S_X{S}_Y}$$

Where $S_x$ and $S_y$ are standard deviations.

Let us call $f(x,y)$ a function resulting from arithmetically combined random variables X and Y. We are interested in the mean value of the function $f(x,y)$ around the mean values of random variables X and Y. To find the mean and variance of this function we use Taylor expansion around $\theta$, which is the point near the mean values of X and Y: $\theta =({\mu }_x,{\mu }_y)$. For the following development let us consider a simplified notation in which: ${\mu }_x=E\left[X\right];{\mu }_y=E\left[Y\right];{\sigma }_x^2=𝑉𝑎𝑟\left[X\right];{\sigma }_y^2=𝑉𝑎𝑟\left[Y\right]$ and ${\sigma }_{𝑥𝑦}=𝐶𝑜𝑣(X,Y)$. Notice that for propagation of uncertainty we need to know the means, variances, the covariance, and the type of operation (^{Heuvelink et al. 1989}).

The second-order Taylor expansion is (^{Weir and Hass 2014}):

$$\begin{array}{cc}\hfill f(x,y)=& f\left(\theta \right)+f_x\left(\theta \right)(x-{\mu }_x)+f_y\left(\theta \right)(y-{\mu }_y)\hfill \\ \hfill & +\frac{1}{2}{f}_{𝑥𝑥}\left(\theta \right){(x-{\mu }_x)}^2+f_{𝑥𝑦}\left(\theta \right)(x-{\mu }_x)(y-{\mu }_y)+\frac{1}{2}{f}_{𝑦𝑦}\left(\theta \right){(y-{\mu }_y)}^2+𝑟𝑒𝑚𝑎𝑖𝑛𝑑𝑒𝑟\hfill \end{array}$$

Applying the expectation operator in (5) and neglecting the remainder, we have:

$$E\left[f(x,y)\right]=f\left(\theta \right)+\frac{1}{2}{f}_{𝑥𝑥}\left(\theta \right){\sigma }_x^2+f_{𝑥𝑦}\left(\theta \right){\sigma }_{𝑥𝑦}+\frac{1}{2}{f}_{𝑦𝑦}\left(\theta \right){\sigma }_y^2$$

This is the general expression for computing the mathematical expectation of the function $f(x,y)$ about $\theta =({\mu }_x,{\mu }_y)$.

The variance of the function $f(x,y)$ around $\theta =({\mu }_x,{\mu }_y)$ can be computed as:

$$𝑉𝑎𝑟\left(f(x,y)\right)=E\left[{(f(x,y)-E\left(f(x,y)\right))}^2\right]$$

Replacing definitions (6) and (7) in (3), we have:

$$\begin{array}{cc}\hfill 𝑉𝑎𝑟\left[f\right(x,y\left)\right]=& E\left[\right(f_x\left(\theta \right)(x-{\mu }_x)+f_y\left(\theta \right)(y-{\mu }_y)+\frac{1}{2}{f}_{𝑥𝑥}\left(\theta \right)(x-{\mu }_x{)}^2+f_{𝑥𝑦}(x-{\mu }_x\left)\right(y-{\mu }_y)\hfill \\ \hfill & +\frac{1}{2}{f}_{𝑦𝑦}\left(\theta \right)(y-{\mu }_y{)}^2-\frac{1}{2}{f}_{𝑥𝑥}(\theta ){\sigma }_x^2-f_{𝑥𝑦}(\theta ){\sigma }_{𝑥𝑦}-\frac{1}{2}{f}_{𝑦𝑦}(\theta \left){\sigma }_y^2{)}^2\right]\hfill \end{array}$$

This expression cannot be simplified because it depends on partial derivatives of the function $f(x,y)$.

The mean and variance of the function $f(x,y)$ are presented in Table I. All formulas result from second-order Taylor expansion (except for equations (13) and (14)) and are equal to those derived by Mood et al. (1974). Equation (13) results from third-order Taylor expansion. This equation should provide a better approach for the mean of quotient of random variables. The variance after equation (14) uses only first-order terms (^{Mood et al. 1974}).

TABLE I
Mean and variance for the function $f(x,y)$ around $\theta =({\mu }_x,{\mu }_y)$.

$f(x,y)$	Derivatives	Mean	Variance
$X\pm Y$	$f_x(x,y)=1$ $f_y(x,y)=\pm 1$ $f_{𝑥𝑥}(x,y)=0$ $f_{𝑦𝑦}(x,y)=0$ $f_{𝑥𝑦}(x,y)=0$	$E\left[f\right(x,y\left)\right]={\mu }_x\pm {\mu }_y\left(10\right)$	$𝑉𝑎𝑟\left[f\right(x,y\left)\right]={\sigma }_x^2+{\sigma }_y^2\pm 2{\sigma }_{𝑥𝑦}\left(11\right)$
$𝑋𝑌$	$f_x(x,y)=y$ $f_y(x,y)=x$ $f_{𝑥𝑥}(x,y)=0$ $f_{𝑦𝑦}(x,y)=0$ $f_{𝑥𝑦}(x,y)=1$	$E\left[f\right(x,y\left)\right]={\mu }_x{\mu }_y+{\sigma }_{𝑥𝑦}\left(12\right)$	$𝑉𝑎𝑟\left[f\right(x,y\left)\right]={\mu }_y^2{\sigma }_x^2+{\mu }_x^2{\sigma }_y^2+E\left[\right(x-{\mu }_x{)}^2(y-{\mu }_y{)}^2]-{\sigma }_{𝑥𝑦}^2+2{\mu }_x{\mu }_y{\sigma }_{𝑥𝑦}+2{\mu }_{y}E\left[\right(x-{\mu }_x{)}^2(y-{\mu }_y)]+2{\mu }_{y}E[(x-{\mu }_x)(y-{\mu }_y{)}^2]\left(13\right)$
$X/Y$	$f_x(x,y)=\frac{1}{y}$ $f_y(x,y)=\frac{-x}{y^2}$ $f_{𝑥𝑥}(x,y)=0$ $f_{𝑦𝑦}(x,y)=\frac{2x}{y^3}$ $f_{𝑥𝑦}(x,y)=\frac{-1}{y^2}$ $f_{𝑥𝑥𝑥}(x,y)=0$ $f_{𝑥𝑥𝑦}(x,y)=0$ $f_{𝑥𝑦𝑦}(x,y)=\frac{2}{y^3}$ $f_{𝑦𝑦𝑦}(x,y)=\frac{-6x}{y^4}$	$E\left[f\right(x,y\left)\right]\approx \frac{{\mu }_x}{{\mu }_y}-\frac{{\sigma }_{𝑥𝑦}}{{\mu }_y^2}+\frac{{\mu }_x{\sigma }_y^2}{{\mu }_y^3}\left(14\right)$ $E\left[f\right(x,y\left)\right]\approx \frac{{\mu }_x}{{\mu }_y}-\frac{{\sigma }_{𝑥𝑦}}{{\mu }_y^2}+\frac{{\mu }_x{\sigma }_y^2}{{\mu }_y^3}+\frac{E\left[\right(x-{\mu }_x\left)\right(y-{\mu }_y{)}^2]}{{\mu }_y^3}+\frac{{\mu }_{x}E\left[\right(y-{\mu }_y{)}^3]}{{\mu }_y^4}\left(15\right)$	$𝑉𝑎𝑟\left[f\right(x,y\left)\right]\approx \frac{{\sigma }_x^2}{{\mu }_y^2}+\frac{{\mu }_x^2}{{\mu }_y^4}{\sigma }_y^2-2\frac{{\mu }_x}{{\mu }_y^3}{\sigma }_{𝑥𝑦}\left(16\right)$ $𝑉𝑎𝑟\left[f\right(x,y\left)\right]\approx \frac{{\sigma }_x^2}{{\mu }_y^2}+\frac{{\mu }_x^2}{{\mu }_y^4}{\sigma }_y^2-2\frac{{\mu }_x}{{\mu }_y^3}{\sigma }_{𝑥𝑦}+\frac{1}{{\mu }_y^4}E\left[\right({\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y){)}^2]+\frac{{\mu }_x^2}{{\mu }_y^6}E[\left(\right(y-{\mu }_y{)}^2-{\sigma }_y^2{)}^2]+\frac{2}{{\mu }_y^3}E[(x-{\mu }_x)[{\sigma }_{𝑥𝑦}-(x-{\mu }_x\left)\right(y-{\mu }_y\left)\right]]+2\frac{{\mu }_x}{{\mu }_y^4}E[(x-{\mu }_x)\left[\right(y-{\mu }_y{)}^2-{\sigma }_y^2\left]\right]-2\frac{{\mu }_x}{{\mu }_y^4}E\left[\right(y-{\mu }_y\left)\right[{\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y)\left]\right]-2\frac{{\mu }_x^2}{{\mu }_y^5}E\left[\right(y-{\mu }_y\left)\right[(y-{\mu }_y{)}^2-{\sigma }_y^2]]+2\frac{{\mu }_x}{{\mu }_y^5}E[[{\sigma }_{𝑥𝑦}-(x-{\mu }_x\left)\right(y-{\mu }_y\left)\right]\left[\right(y-{\mu }_y{)}^2-{\sigma }_y^2\left]\right]\left(17\right)$

It can be noted that for the quotient of X and Y, the mean and variance are calculated by approximate formulas. For a better approximation of the variance of the quotient of X and Y, we considered second-order terms according to Maskey and Guinot (2003). The variance of the quotient after equation (15) is better than the variance computed after equation (14). Indeed, the algebraic complexity increased very much with the inclusion of second-order terms (^{Maskey and Guinot 2003}). The mean and variance for the quotient of X and Y will always result in approximate formulas because the function $f(x,y)=X/Y$ is infinitely differentiable with respect to Y. An exact formula for the calculation of the variance of the quotient of X and Y was proposed by Frishman (¹⁹⁷¹). However, this alternative uses the variance of the product $X*1/Y$, which is implied in estimating $1/Y$ instead of $Y$. Details of the mathematical development of Taylor expansion to find the mean and variance of the function $f(x,y)$ can be found in the Appendix. This development was based entirely on statistics (Mood et al. 1974) and calculus (Weir and Hass 2014).

Heuvelink et al. (1989) proposed a general formula based on second order Taylor expansion for n variables. The mean of the function $f\left(M\right)$ around $M=\{{\mu }_1,{\mu }_2,\dots ,{\mu }_n\}$ is computed as (Heuvelink et al. 1989):

Where ${\rho }_{𝑖𝑗}$ is the correlation coefficient between variables $i$ and $j$; ${\sigma }_i$ is the standard deviation for variable $i$ and ${\sigma }_j$ is the standard deviation for variable $j$.

The general formula for variance of the function $f\left(M\right)$ is (Heuvelink et al. 1989):

$$𝑉𝑎𝑟\left(\acute{x}\right)=E\left\{{[\sum _{k=1}^n\left[(x_k-{\mu }_k)f_k\left(M\right)\right]+\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n[(x_i-{\mu }_i)(x_j-{\mu }_j)-{\rho }_{𝑖𝑗}{\sigma }_i{\sigma }_j?f_{𝑖𝑗}(M\left)\right]]}^2\right\}$$

$$\begin{array}{cc}\hfill 𝑉𝑎𝑟\left(\acute{x}\right)=& \sum _{k=1}^n\sum _{l=1}^n{\rho }_{𝑘𝑙}{\sigma }_k{\sigma }_l{f}_k\left(M\right)f_l\left(M\right)\hfill \\ \hfill & +?\sum _{k=1}^n\sum _{i=1}^n\sum _{j=1}^n\{E\left[(x_k-{\mu }_k((x_i-{\mu }_i)(x_j-{\mu }_j)-{\rho }_{𝑘𝑙}{\sigma }_k{\sigma }_l))\right]\left(f_k\left(M\right)\right)*f_{𝑖𝑗}(M\left)\right\}\hfill \\ \hfill & +?\frac{1}{4}\sum _{i=1}^n\sum _{j=1}^n\sum _{k=1}^n\sum _{l=1}^n\left\{E[(x_i-{\mu }_i)(x_j-{\mu }_j)(x_k-{\mu }_k)(x_l-{\mu }_l)-{\rho }_{𝑖𝑗}{\sigma }_i{\sigma }_j{\rho }_{𝑘𝑙}{\sigma }_k{\sigma }_l]f_{𝑖𝑗}\right(M\left)f_{𝑘𝑙}\right(M\left)\right\}\hfill \end{array}$$

We proved that all formulas presented in Table I (except equations (12) and (14)) are equal to those developed from Heuvelink’s general formulas for mean and variance (Appendix).

LOCAL ESTIMATES AND UNCERTAINTIES USING MULTIQUADRIC EQUATIONS

This paper concerns building raster maps from field data, which means interpolation of regularly spaced cells based on neighbor data. Multiquadric equations (^{Hardy 1971}) were chosen for interpolation of the variable at unsampled locations. Usually, several variables are measured in each sample location. For instance, in geochemical exploration, major and trace elements are analyzed simultaneously for each sample. Once we have built raster maps with these variables, we can combine variables to get another raster map. This is a common procedure in GIS. We present a method that allows combination of random variables and also uncertainty quantification based on equations summarized in Table I. It is important to emphasize that it depends on a reliable measure of uncertainty. In this sense, we propose the use of the interpolation variance (Yamamoto 2000) as an approach for uncertainty quantification. Interpolation variance can be used with any estimate based on weighted average formula.

For each cell location within the raster map, we have to find the nearest neighbor data. The dual form of multiquadric interpolation can be written as (^{Yamamoto 2002}):

$$Z_{𝑀𝑄}^{*}\left(x_o\right)=\sum _{i=1}^n{w}_{i}Z\left(x_i\right)$$

The weights $\{w_i,i=1,\dots ,n\}$ come from the solution of a system of linear equations (Yamamoto 2002):

$$\sum _{j=1}^n{w}_j\varnothing (x_i-x_j)+\mu =\varnothing (x_i-x_o)𝑓𝑜𝑟i=1\sum _{j=1}^n{w}_j=1$$

Where $\mu$ is the Lagrange multiplier and $\varnothing \left(x\right)=\sqrt{{\left|x\right|}^2+C}$ is the multiquadric kernel as proposed by Hardy (1971).

Because equation (17) is a weighted average formula, we can use the interpolation variance for quantifying the uncertainty associated with multiquadric interpolation. The interpolation variance is simply (Yamamoto 2000):

$$S_o^2\left(x_o\right)=\sum _{i=1}^n{w}_i{[Z\left(x_i\right)-Z_{𝑀𝑄}^{*}(x_o\left)\right]}^2$$

Raster map operations consist of combining two variables at once. Let us call $Z_1\left(x\right)$ and $Z_2\left(x\right)$ the two variables to be arithmetically combined. Both variables are measured on same data location and the data set has the same number of data points. This is the essential condition to perform arithmetical operation between variables $Z_1\left(x\right)$ and $Z_2\left(x\right)$. Since we are using multiquadric equations with the same multiquadric kernel, then the weights are identical: $\{w_i^1=w_i^2,i=1,\dots ,n\}$. In light of this fact we will consider hereafter the set of weights $\{w_i,i=1,\dots ,n\}$ for both variables.

Multiquadric interpolations for variables $Z_1\left(x\right)$ and $Z_2\left(x\right)$ are given as:

$$Z_{{𝑀𝑄}_1}^{*}\left(x_o\right)=\sum _{i=1}^n{w}_i{Z}_1\left(x_i\right)$$

$$Z_{{𝑀𝑄}_2}^{*}\left(x_o\right)=\sum _{i=1}^n{w}_i{Z}_2\left(x_i\right)$$

The term $E\left[Z_1\left(x\right)Z_2\left(x\right)\right]$ as a local statistic can be calculated as:

$$E\left[Z_1\left(x\right)Z_2\left(x\right)\right]=\sum _{i=1}^n\left(w_i{Z}_1\right(x_i\left)Z_2\right(x_i\left)\right)$$

Thus, the local covariance can be computed as:

$$𝐶𝑜𝑣(Z_1\left(x\right),Z_2\left(x\right))=\sum _{i=1}^n{w}_i{Z}_1\left(x_i\right)Z_2\left(x_i\right)-\left(\sum _{i=1}^n{w}_i{Z}_1\left(x_i\right)\right)\left(\sum _{i=1}^n{w}_i{Z}_2\left(x_i\right)\right)$$

See that equation (27) allows computing cell-by-cell covariance. This statistic is fundamental for propagation of uncertainty. The interpolation covariance in (27) is calculated directly from neighbor data using the same weights for multiquadric interpolation. Therefore, this approach is different from the spatial covariance derived from the variogram model. Actually, for $Z_1\left(x\right)=Z_2\left(x\right)$, this formula is equivalent to the interpolation variance.

Interpolation variances are computed as:

$$S_{o_1}^2=\sum _{i=1}^n{w}_i{\left[Z_1\right(x_i)-Z_{𝑀𝑄}^{*}(x_o\left)\right]}^2$$

$$S_{o_2}^2=\sum _{i=1}^n{w}_i{\left[Z_2\right(x{?}_i)-Z_{𝑀𝑄}^{*}(x_o\left)\right]}^2$$

Means and variances in Table I concern global statistics. But, raster maps imply local statistics and the purpose is to show that local statistics can be combined to compute mean and variances after mathematical operations between two random variables.

Keeping the same notation as used in Table I, we consider the following equivalence: ${\mu }_x=Z_{{𝑀𝑄}_1}^{*}\left(x\right);{\mu }_y=Z_{{𝑀𝑄}_2}^{*}\left(x\right);{\sigma }_x^2=S_{o_1}^2\left(x_o\right);{\sigma }_y^2=S_{o_2}^2\left(x_o\right)$ and ${\sigma }_{𝑥𝑦}=𝐶𝑜𝑣(Z_1\left(x\right),Z_2(x\left)\right)$.

In addition, high-order moments in Table I are calculated as:

$$E\left[.\right]=\sum _{i=1}^n{w}_i\left(.\right)$$

MATERIALS AND METHODS

For this study we consider a stratified random sample drawn from an exhaustive data set composed of 2500 data points arranged in an array of 50 rows by 50 columns. This sample, with 100 data points (4% of the exhaustive) presents two measured variables (VP and VS2) in each location (Figure 1a). These variables present a positive correlation (Figure 1b). In addition to performing arithmetic operations between VP and VS2, we have to prove that mean and variances of combined variables are valid. Thus, we generated four new variables: VP+VS2; VP-VS2; VPxVS2 and VP/VS2. Sample statistics for all variables are presented in Table II.

Figure 1
Location map of sample data (a) and VP X VS2 scattergram (b) (see the colors in the online version).

TABLE II
Sample statistics for variables VP, VS2, VP+VS2, VP-VS2, VPxVS2 and VP/VS2.

Statistics	VP	VS2	VP+VS2	VP-VS2	VPxVS2	VP/VS2
Nb. of data	100	100	100	100	100	100
Mean	15.568	25.051	40.619	-9.483	403.047	0.615
Variance	21.130	14.236	61.486	9.245	28571.847	0.017
Coef. Var.	0.295	0.151	0.193	-0.321	0.419	0.215
Max.	27.863	37.620	64.377	-3.087	1006.605	0.893
UQ	18.421	27.203	45.624	-7.647	499.767	0.694
Median	15.774	24.666	40.951	-9.326	394.319	0.625
LQ	12.757	22.391	35.219	-11.679	288.545	0.543
Min.	4.049	13.386	21.134	-17.389	86.794	0.189

• UQ=upper quartile and LQ=lower quartile.

VP and VS2 will be used for building raster maps, from which we want to obtain new raster maps by combining them arithmetically. Mean and variance for new random variables are listed in Table III. Comparing these results with directly computed statistics in Table II, we verify that except for the quotient of random variables, all others are exactly the same.

The sample data with six variables are the materials of this study. Original variables VP and VS2 will be used to build raster maps based on multiquadric interpolation. Equations (24) and (28) will be used to compute cells of the VP raster maps and equations (25) and (29) for VS2 raster maps. The covariance between VP and VS2 will be calculated after equation (27). These statistics will be combined according to formulas in Table I to derive means and uncertainties of resulting arithmetically combined raster maps. Now, we can compare resulting raster maps with those computed directly from new variables based on equations (16) and (17). Next, these equations ((21) and (23)) will be checked against equations (10) and (11) for variables VP±VS2, equations (12) and (13) for variable VPxVS2 and equations (14), (15), (16) and (17) for variable VP/VS2.

TABLE III
Mean and variance for combined random variables.

$f(x,y)$	Mean	Variance
$X+Y$	$15.56775+25.05096=40,61871\left(10\right)$	$21.12972+14.23615+26.12052=61.48639\left(11\right)$
$X-Y$	$15.56775-25.05096=-9.48321\left(10\right)$	$21.12972+14.23615-26.12052=9.24535\left(11\right)$
$𝑋𝑌$	$15.56775*25.05096+13.06026=403.04734\left(12\right)$	$13259.97187+3450.19941+636.64517-170.57033+10186.66456+729.07689+479.85943=28571.847\left(13\right)$
$X/Y$	$0.62144-0.02081+0.01410=0.61473\left(14\right)$ $0.62144-0.02081+0.01410+0.0009256476-0.0004044825=0.61525\left(15\right)$	$0.03367+0.00876-0.02587=0.01656\left(16\right)$ $0.03367+0.00876-0.02587+0.00118+0.00060-0.00185+0.00122+0.00122-0.00050-0.00151=0.01692\left(17\right)$

RESULTS AND DISCUSSION

Raster maps for variables VP and VS2 are presented in Figure 2. For the purpose of displaying interpolation variances and covariance, we considered a maximum equal to 95% of their distributions. Now, raster maps for VP and VS2 can be combined by applying arithmetic operators. Before proceeding it is important to show raster maps (Figure 3) computed directly from new variables.

The first operation is the sum: VP+VS2. Results of this operation are presented in Figure 4. The result (Figure 4a) is simply the sum of raster maps (for VP and VS2). The variance of the sum (Figure 4b) is equal to the sum of variances for VP $\left({\sigma }_x^2\right)$ and VS2 $\left({\sigma }_y^2\right)$ plus twice the covariance between VP and VS2 $\left({\sigma }_{𝑥𝑦}\right)$. The variance of the sum has great influence of the variable VP and the covariance as well. Results illustrated in Figure 4 are exactly the same as shown in Figure 3a.

Results of the subtraction of VS2 from VP are presented in Figure 5. The mathematical expectation and variance follow the same raster operations as for the sum. However, the variance is equal to the sum of variances minus twice the covariance between VP and VS2. The variance of the difference between random variables VP and VS2 shows a reasonable correlation with the variance of VP, but almost no correlation with the variance of VS2 and the covariance. Mean and variances computed after equations (10) and (11) for addition and subtraction match those values computed directly from equations (21) and (23). Once again, results can be compared graphically with images presented in Figure 3b.

Figure 2
Raster maps for variables VP (a) and for variable VS2 (b). Interpolated variables in column I, uncertainties in column II and covariance between VP and VS2 in column III. Interpolation variances and covariance were displayed for a maximum equal to 95% of their distributions (for the interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Figure 3
Raster maps computed directly from random variables VP+VS2 (a); VP-VS2 (b); VPxVS2 (c) and VP/VS2 (d). I) Interpolated variables; II) variances (for the interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Figure 4
Results of arithmetic operation VP+VS2. The result (a) is equal to the sum of VP $\left({\mu }_x\right)$ and VS2 $\left({\mu }_y\right)$. The associated uncertainty (b) is equal to the variance of VP $\left({\sigma }_x^2\right)$ plus variance of VS2 $\left({\sigma }_y^2\right)$ plus twice the covariance between VP and VS2 $\left({\sigma }_{𝑥𝑦}\right)$.

Figure 5
Results of arithmetic operation VP-VS2. The result (a) is equal to the sum of VP $\left({\mu }_x\right)$ and VS2 $\left({\mu }_y\right)$. The associated uncertainty (b) is equal to the variance of VP $\left({\sigma }_x^2\right)$ plus variance of VS2 $\left({\sigma }_y^2\right)$ minus twice the covariance between VP and VS2 $\left({\sigma }_{𝑥𝑦}\right)$ (for the interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

The results of the multiplication of VP and VS2 are presented in Figure 6, which can be compared with results of direct multiquadric interpolation of the variable VPxVS2 (Figure 3c). The mathematical expectation of this operation $E\left[f(x,y)\right]$ is equal to the product of VP and VS2 plus the covariance between VP and VS2. Note that the mathematical expectation has the contribution of the covariance. $E\left[f(x,y)\right]$ is highly correlated with the means ${\mu }_x$ and ${\mu }_y$ and low correlation with the covariance. The uncertainty $𝑉𝑎𝑟\left[f(x,y)\right]$ is equal to term 1 $\left({\mu }_y^2{\sigma }_x^2\right)$ plus term 2 $\left({\mu }_x^2{\sigma }_y^2\right)$ plus term 3 $\left(E\left[{(x-{\mu }_x)}^2{(y-{\mu }_y)}^2\right]\right)$ minus term 4 $\left({\sigma }_{𝑥𝑦}^2\right)$ plus term 5 $\left(2{\mu }_x{\mu }_y{\sigma }_{𝑥𝑦}\right)$ plus term 6 $\left(2{\mu }_{y}E\left[{(x-{\mu }_x)}^2(y-{\mu }_y)\right]\right)$ plus term 7 $\left(2{\mu }_{x}E\left[\right(x-{\mu }_x\left){(y-{\mu }_y)}^2\right]\right)$. Next, the variance was calculated after second-order Taylor expansion (equation (13)) and therefore higher order moments have been included. Equations (12) and (13) provided the same results as applying equations (21) and (23) for the new random variable VPxVS2.

Figure 6
Results of arithmetic operation VPxVS2. The result $\left(E\left[f(x,y)\right]\right)$ is equal to the product of VP $\left({\mu }_x\right)$ and VS2 $\left({\mu }_y\right)$ plus the covariance between VP and VS2 $\left({\sigma }_{𝑥𝑦}\right)$. The associated uncertainty $\left(𝑉𝑎𝑟\left[f(x,y)\right]\right)$ is equal to the term 1 $\left({\mu }_y^2{\sigma }_x^2\right)$ plus the term 2 $\left({\mu }_x^2{\sigma }_y^2\right)$ plus the term 3 $\left(E\left[{(x-{\mu }_x)}^2{(y-{\mu }_y)}^2\right]\right)$ minus the term 4 $\left({\sigma }_{𝑥𝑦}^2\right)$ plus the term 5 $\left(2{\mu }_x{\mu }_y{\sigma }_{𝑥𝑦}\right)$ plus the term 6 $\left(2{\mu }_{y}E\right[(x-{\mu }_x{)}^2(y-{\mu }_y\left)\right])$ plus the term 7 $\left(2{\mu }_{y}E\right[(x-{\mu }_x)(y-{\mu }_y{)}^2])$ (for the interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

The calculation of mean and variance for the quotient of VP and VS2 is more complex since the function $X/Y$ is infinitely differentiable with respect to $Y$. The mathematical expectation is not only equal to the quotient between VP and VS2 $\left({\mu }_x/{\mu }_y\right)$, but includes two other terms: ${\sigma }_{𝑥𝑦}/{\mu }_y^2$ and ${\mu }_x{\sigma }_y^2/{\mu }_y^3$(Figure 7a). This is the result after equation (14), which results from second-order Taylor expansion. Results from application of the equation (15) are shown in Figure 7b. Although results look very similar to each other, the mean according to equation (15) is more accurate than that provided by equation (14). Figures 7a, b basically display the same pattern as shown in Figure 3d-I. For calculating the variance of the ratio of X to Y, we have equation (16) that results from first-order Taylor expansion and another equation (17) derived after second-order Taylor expansion (Figure 8). It is virtually impossible to see the difference between results of these approaches in Figure 8. Thus, means and variances can be compared by being displayed in a scatter plot (Figure 9). Although, mean values are very close to each other, we can see that equation (15) provides a better approach than the equation (14). On the other hand, the variance after equation (17) is much better than that provided by equation (16).

We can follow arithmetic operations taking as an example the cell with coordinates CX=29.50 and CY=23.5 (Table IV). With eight neighbor data we interpolated this cell using field data for VP and VS2. All other variables (VP+VS2, VP-VS2, VPxVS2 and VP/VS2) were derived from these measurements (VP and VS2). Note that we can compute mean and variance directly for derived variables and thus these values can be checked against statistics computed using formulas listed in Table I. In Table IV, the data point with coordinates (CX=26.5, CY=25.5) has a weight equal to zero. This is because we applied an algorithm for correcting negative weights (^{Rao and Journel 1997}), in which a constant equal to the modulus of the largest negative weight was added to all weights and then restandardized to a sum equal to one.

Figure 7
Mean values for arithmetic operation VP/VS2. (a - upper row) mean after equation (14); (b - lower row) mean according to equation (15) (for the interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Figure 8
Variances of the arithmetic operation VP/VS2. a) resulting variance from first order Taylor expansion (equation (15)); b) variance calculated using second order Taylor expansion (equation (17)) (for the interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Figure 9
Scatterplots for means of the ratio of X to Y: a = equation (14) and b = equation (15); variances of the ratio: c = equation (16) and d = equation (17) (for the interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

TABLE IV
Neighboring data around a cell with coordinates (CX=29.50, CY=23.5), mean and variances for variables VP, VS2, VP+VS2, VP-VS2, VPxVS2 and VP/VS2.

CX	CY	W	VP	VS2	VP+VS2	VP-VS2	VPxVS2	VP/VS2
30.50	26.50	0.227	13.795	25.726	39.521	-11.931	354.890	0.53623
37.50	26.50	0.018	16.141	24.986	41.127	-8.845	403.299	0.64600
26.50	25.50	0	17.309	26.100	43.409	-8.791	451.765	0.66318
28.50	23.50	0.573	15.470	24.568	40.038	-9.098	380.067	0.62968
26.50	16.50	0.030	4.049	21.438	25.487	-17.389	86.802	0.18887
24.50	21.50	0.021	12.915	23.497	36.412	-10.581	303.464	0.54964
32.50	18.50	0.043	13.347	21.240	34.587	-7.894	283.490	0.62839
33.50	20.50	0.088	15.054	21.954	37.008	-6.900	33.496	0.68571
$Z_{MQ}^{*}\left(x_o\right)$			14.577	24.349	38.926	-9.773	355.831	0.59867
$S_o^2\left(x_o\right)$			4.129	1.656	7.594	3.977	2901.194	0.00730

Mean and variances for the cell with coordinates (CX=29.50, CY=23.5) after combining random variables VP and VS2 are presented in Table V. As can be seen, statistics for VP+VS2, VP-VS2 and VPxVS2 match with statistics computed directly from derived variables listed in Table IV. However, for VP/VS2 we obtain only approximate results.

TABLE V
Mean and variance for combined random variables.

$f(x,y)$	Mean	Variance
$X+Y$	$14.57656+24.34917=38.92573\left(10\right)$	$4.12927+1.65640+2*0.90406=7.59379\left(11\right)$
$X-Y$	$14.57656-24.34917=-9.77261\left(10\right)$	$4.12927+1.65640-2*0.90406=3.97755\left(11\right)$
$𝑋𝑌$	$14.57656*24.34927+0.900406=355.83120\left(12\right)$	$2448.17033+351.94621+29.44879-0.81733+641.75335-473.39666-95.91113=2901.19356\left(13\right)$
$X/Y$	$\frac{14.57656}{24.34917}-\frac{0.90406}{(24.34917{)}^2}+\frac{14.57656*1.65640}{(24.34917{)}^3}=0.59879\left(14\right)$	$\frac{4.12927}{(24.34917{)}^2}+\frac{(14.57656{)}^2*1.65640}{(24.34917{)}^4}-\frac{(2*14.57656)*0.90406}{(24.34917{)}^3}=0.0061402737\left(16\right)$ $0.0061402737+0.0000814529+0.0000073085+0.0013467544-0.0002728552-0.0002728552+0.0001317723-0.0000235991=0.0071382523\left(17\right)$

CONCLUSIONS

Propagation of uncertainty coming from arithmetic operations between random variables is well known in statistics (Mood et al. 1974). Heuvelink et al. (1989) and Heuvelink (¹⁹⁹⁸) have established the method of calculating mean and variances of the output raster map from several input maps based on Taylor method. This paper showed that propagation of uncertainty depends on a reliable measure of local accuracy and local covariance. We proved that it is possible to get exact mean and variance of the output raster map for map operations involving addition, subtraction and multiplication. For the mean of quotient, we proposed to compute it from a third-order Taylor expansion that provides a slightly better result than the mean computed from second-order Taylor expansion. For the variance of the ratio of X to Y, we demonstrated that the formula based on second-order Taylor expansion provides much better results than the first order.

ACKNOWLEDGMENTS

The authors wish to thank Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq) for supporting this research.

• ¹
  ATKINSON PM AND FOODY GM. 2002. Uncertainty in remote sensing and GIS: fundamentals. In: Foody GM and Atkinson PM (Eds), Uncertainty in Remote Sensing and GIS, p. 1-18.
• ²
  BURROUGH PA. 1986. Principles of geographical information system. Oxford, Oxford University Press, 194 p.
• ³
  CROSETTO M, TARANTOLA S AND SALTELI A. 2000. Sensitivity and uncertainty analysis in spatial modelling based on GIS. Agr Ecosyst Environ 81: 71-79.
• ⁴
  FRISHMAN F. 1971. On the arithmetic means and variances of products and ratios of random variables. National Technical Information Service, AD-785623, p. 331-345.
• ⁵
  HARDY RL. 1971. Multiquadric equations of topography and other irregular surfaces. J Geophys Res 76: 1905-1915.
• ⁶
  HEUVELINK GB. 1998. Error propagation in environmental modelling with GIS. London, Taylor & Francis, 127 p.
• ⁷
  HEUVELINK GB, BURROUGH PA AND STEIN A. 1989. Propagation of errors in spatial modelling with GIS. Int J Geogr Inf Syst 3: 303-322.
• ⁸
  JOURNEL AG AND ROSSI M. 1989. When do we need a trend model in kriging? Math Geolo 21: 715-739.
• ⁹
  MASKEY S AND GUINOT V. 2003. Improved first-order second moment method for uncertainty estimation in flood forecasting. Hydrolog Sci J 48: 183-196.
• ¹⁰
  MCCALLUM WG. 1998. Multivariable calculus. New York, J Wiley & Sons, 270 p.
• ¹¹
  MOOD AM, GRAYBILL FA AND BOES DC. 1974. Introduction to the theory of statistics. Tokyo, McGraw Hill, 3rd ed., 564 p.
• ¹²
  RAO SE AND JOURNEL AG. 1997. Deriving conditional distributions from ordinary kriging. In: Baafi EY and Schofield NA (Eds), Geostatistics Wollongong 96: 92-102.
• ¹³
  WANG G, GERTNER GZ, FANG S AND ANDERSON AB. 2005. A methodology for spatial uncertainty analysis of remote sensing and GIS products Photogramm Eng Rem S 71: 1423-1432.
• ¹⁴
  WEIR MD AND HASS J. 2014. Thomas’ calculus – multivariable. Boston, Person Education Inc, 13th ed., p. 572-1044.
• ¹⁵
  WELLMANN FJ, HOROWITZ GF, SCHILL E AND REGENAUER-LIEB K. 2010. Towards incorporating uncertainty of structural data in 3D geological inversion. Tectonophysics 490: 141-152.
• ¹⁶
  YAMAMOTO JK. 2000. An alternative measure of the reliability of ordinary kriging estimates. Math Geol 32: 489-509.
• ¹⁷
  YAMAMOTO JK. 2002. Ore reserve estimation using radial basis functions. Rev Inst Geol 23: 25-38.
• ¹⁸
  YAMAMOTO JK, MAO XM, KOIKE K, CROSTA AP, LANDIM PMB, HU HZ, WANG CY AND YAO LQ. 2012. Mapping an uncertainty zone between interpolated types of a categorical variable. Comput Geosci 40: 146-152.

APPENDIX

Mean and variance of the function $f(x,y)$ using Taylor expansion

In this appendix, we develop equations for mean and variance of the function $f(x,y)$ as presented in Table I. We will use Taylor expansion of $f(x,y)$ to find the mean and variance of this function around mean values ${\mu }_x$ and ${\mu }_y$. Besides mean values, variances of random variables $({\sigma }_x^2$ and ${\sigma }_y^2)$ and the covariance between random variables X and Y $\left({\sigma }_{𝑥𝑦}\right)$ must be known. The development for the mean and variance of product and quotient of random variables was based on Taylor expansion, according to Mood et al. (1974).

According to Weir and Hass (2014), the general equation of Taylor’s Formula for two variables at the point $\theta =(a,b)$ is:

Where $h=(x-a)$, $k=(y-b)$, $0<c<1$ and the last term is the remainder.

Neglecting the remainder term and taking n=2, we have:

$$\begin{array}{cc}\hfill f(x,y)=& f(a,b)+\frac{\partial f(a,b)}{\partial x}(x-a)+\frac{\partial f(a,b)}{\partial y}(y-b)\hfill \\ \hfill & +\frac{1}{2!}(\frac{{\partial }^{2}f(a,b)}{\partial x^2}{(x-a)}^2+2\frac{{\partial }^{2}f(a,b)}{\partial x\partial y}(x-a)(y-b)+\frac{{\partial }^{2}f(a,b)}{\partial y^2}{(y-b)}^2)\hfill \end{array}$$

Letting $(a,b)=({\mu }_x,{\mu }_y)$ in equation A1.2 we have:

$$\begin{array}{cc}\hfill f(x,y)=& f({\mu }_x,{\mu }_y)+f_x({\mu }_x,{\mu }_y)(x-{\mu }_x)+f_y({\mu }_x,{\mu }_y)(y-{\mu }_y)\hfill \\ \hfill & +\frac{1}{2!}(f_{𝑥𝑥}({\mu }_x,{\mu }_y){(x-{\mu }_x)}^2+2f_{𝑥𝑦}({\mu }_x,{\mu }_y)(x-{\mu }_x)(y-{\mu }_y)+f_{𝑦𝑦}({\mu }_x,{\mu }_y){(y-{\mu }_y)}^2)\hfill \end{array}$$

Applying the expectation operator to expression A1.3:

$$\begin{array}{cc}\hfill E\left[f(x,y)\right]=& E\left[f({\mu }_x,{\mu }_y)\right]+E\left[f_x({\mu }_x,{\mu }_y)(x-{\mu }_x)\right]+E\left[f_y({\mu }_x,{\mu }_y)(y-{\mu }_y)\right]\hfill \\ \hfill & +\frac{1}{2!}(E\left[f_{𝑥𝑥}({\mu }_x,{\mu }_y){(x-{\mu }_x)}^2\right]+2E\left[f_{𝑥𝑦}({\mu }_x,{\mu }_y)(x-{\mu }_x)(y-{\mu }_y)\right]+E\left[f_{𝑦𝑦}({\mu }_x,{\mu }_y){(y-{\mu }_y)}^2\right])\hfill \end{array}$$

Finally:

$$E\left[f(x,y)\right]=f({\mu }_x,{\mu }_y)+\frac{1}{2}({\sigma }_x^2{f}_{𝑥𝑥}({\mu }_x,{\mu }_y)+2{\sigma }_{𝑥𝑦}{f}_{𝑥𝑦}({\mu }_x,{\mu }_y)+{\sigma }_y^2{f}_{𝑦𝑦}({\mu }_x,{\mu }_y))$$

This is the general expression for mathematical expectation of $f(x,y)$. The variance of the function $f(x,y)$ can be written as:

$$𝑉𝑎𝑟\left(f\right(x,y\left)\right)=E.$$

Replacing A1.2 and A1.3 in A1.4, we get:

$$𝑉𝑎𝑟\left(f\right(x,y\left)\right)=E.$$

Calculation of mean and variance for sum and subtraction of random variables X and Y

For $f(x,y)=x\pm y$, we have the following derivatives:

$$\begin{array}{c}\hfill f_x(x,y)=1\hfill \\ \hfill f_y(x,y)=\pm 1\hfill \\ \hfill f_{𝑥𝑥}(x,y)=0\hfill \\ \hfill f_{𝑥𝑦}(x,y)=0\hfill \\ \hfill f_{𝑦𝑦}(x,y)=0\hfill \end{array}$$

Applying them on equations A1.4 and A1.6, we have equations (10) and (11) of Table I:

$$E\left[f(x,y)\right]=f({\mu }_x,{\mu }_y)$$

$$E\left[f(x,y)\right]={\mu }_x\pm {\mu }_y$$

$$𝑉𝑎𝑟\left(f(x,y)\right)=E\left[{((x-{\mu }_x)\pm (y-{\mu }_y))}^2\right]$$

$$𝑉𝑎𝑟\left(f(x,y)\right)=E[{(x-{\mu }_x)}^2\pm 2(x-{\mu }_x)(y-{\mu }_y)+{(y-{\mu }_y)}^2]$$

$$𝑉𝑎𝑟\left(f(x,y)\right)={\sigma }_x^2+{\sigma }_y^2\pm 2{\sigma }_{𝑥𝑦}$$

Calculation of the mean and variance for the product of random variables X and Y

For $f(x,y)=𝑥𝑦$, we have the following derivatives:

$$\begin{array}{c}\hfill f_x(x,y)=y\hfill \\ \hfill f_y(x,y)=x\hfill \\ \hfill f_{𝑥𝑥}(x,y)=0\hfill \\ \hfill f_{𝑦𝑦}(x,y)=0\hfill \\ \hfill f_{𝑥𝑦}(x,y)=1\hfill \end{array}$$

Applying them on equation A1.4, we get equation ((12) - Table I):

$$E\left[f(x,y)\right]={\mu }_x{\mu }_y+{\sigma }_{𝑥𝑦}$$

Replacing derivatives in equation A1.6, we have:

$$𝑉𝑎𝑟\left(f(x,y)\right)=E\left\{{\left[\begin{array}{c}\hfill {\mu }_y(x-{\mu }_x)+{\mu }_x(y-{\mu }_y)+(x-{\mu }_x)(y-{\mu }_y)-{\sigma }_{𝑥𝑦}\hfill \end{array}\right]}^2\right\}$$

Notice that:

$${(a+b+c+d)}^2=a^2+b^2+c^2+d^2+2𝑎𝑏+2𝑎𝑐+2𝑎𝑑+2𝑏𝑐+2𝑏𝑑+2𝑐𝑑$$

Expanding A1.10 according to A1.11 we obtain:

$$𝑉𝑎𝑟\left(f(x,y)\right)=E.$$

Applying definitions of variance, covariance and the general equation of variance A1.6, we get equation ((13) - Table I):

$$\begin{array}{cc}\hfill 𝑉𝑎𝑟\left(f(x,y)\right)=& {\mu }_y^2{\sigma }_x^2+{\mu }_x^2{\sigma }_y^2+2{\mu }_x{\mu }_y{\sigma }_{𝑥𝑦}+E\left[{(x-{\mu }_x)}^2{(y-{\mu }_y)}^2\right]-{\left({\sigma }_{𝑥𝑦}\right)}^2\hfill \\ \hfill & +2{\mu }_{y}E\left[{(x-{\mu }_x)}^2(y-{\mu }_y)\right]+2{\mu }_{x}E\left[(x-{\mu }_x){(y-{\mu }_y)}^2\right]\hfill \end{array}$$

Calculation of the mean and variance for the ratio of X to Y

For calculation of the mean of the quotient of random variables X and Y we will consider second and third order Taylor expansion as follows.

For function $f(x,y)=x/y$, we have the following derivatives:

$$\begin{array}{c}\hfill f_x=\frac{1}{y}\hfill \\ \hfill f_y=\frac{-x}{y^2}\hfill \end{array}\to \begin{array}{c}\hfill f_{𝑥𝑥}=0\hfill \\ \hfill f_{𝑥𝑦}=\frac{-1}{y^2}\hfill \\ \hfill f_{𝑦𝑦}=\frac{2x}{y^3}\hfill \end{array}\to \begin{array}{c}\hfill f_{𝑥𝑥𝑥}=0\hfill \\ \hfill f_{𝑥𝑥𝑦}=0\hfill \\ \hfill f_{𝑥𝑦𝑦}=\frac{2}{y^3}\hfill \\ \hfill f_{𝑦𝑦𝑦}=\frac{-6x}{y^4}\hfill \end{array}$$

Second order Taylor expansion for calculation of the mean of function $f(x,y)=x/y$:

Applying the derivatives (until second order) on the equation A1.4 we get (equation (14) - Table I):

$$E\left[f(x,y)\right]=\frac{{\mu }_x}{{\mu }_y}-{\sigma }_{𝑥𝑦}^2\frac{1}{{\mu }_y^2}+{\sigma }_y^2\frac{{\mu }_x}{{\mu }_y^3}$$

Third order Taylor expansion for calculation of the mean of function $f(x,y)=x/y$:

Developing equation A1.1 for n=3, we get:

$$\begin{array}{cc}\hfill f(x,y)=& f(a,b)+\frac{\partial f(a,b)}{\partial x}(x-a)+\frac{\partial f(a,b)}{\partial y}(y-b)\hfill \\ \hfill & +\frac{1}{2!}(\frac{{\partial }^{2}f(a,b)}{\partial x^2}{(x-a)}^2+2\frac{{\partial }^{2}f(a,b)}{\partial x\partial y}(x-a)(y-b)+\frac{{\partial }^{2}f(a,b)}{\partial y^2}{(y-b)}^2)\hfill \\ \hfill & +\frac{1}{3!}(\frac{{\partial }^{3}f(a,b)}{\partial x^3}{(x-a)}^3+3\frac{{\partial }^{3}f(a,b)}{\partial x^2\partial y}{(x-a)}^2(y-b)+3\frac{{\partial }^{3}f(a,b)}{\partial x\partial y^2}(x-a){(y-b)}^2+\frac{{\partial }^{3}f(a,b)}{\partial y^3}{(y-b)}^3)\hfill \end{array}$$

Letting $(a,b)=({\mu }_x,{\mu }_y)$ in equation A1.14, we have:

$$\begin{array}{cc}\hfill f(x,y)=& f({\mu }_x,{\mu }_y)+f_x({\mu }_x,{\mu }_y)(x-{\mu }_x)+f_y({\mu }_x,{\mu }_y)(y-{\mu }_y)\hfill \\ \hfill & +\frac{1}{2!}\left(f_{𝑥𝑥}\right({\mu }_x,{\mu }_y\left)\right(x-{\mu }_x{)}^2+2f_{𝑥𝑦}({\mu }_x,{\mu }_y)(x-{\mu }_x)(y-{\mu }_y)+f_{𝑦𝑦}({\mu }_x,{\mu }_y)(y-{\mu }_y{)}^2)\hfill \\ \hfill & +\frac{1}{3!}(\begin{array}{c}\hfill f_{𝑥𝑥𝑥}({\mu }_x,{\mu }_y)(x-{\mu }_x{)}^3+3f_{𝑥𝑥𝑦}({\mu }_x,{\mu }_y\left)\right(x-{\mu }_x{)}^2(y-{\mu }_y)\hfill \end{array}+3f_{𝑥𝑦𝑦}({\mu }_x,{\mu }_y\left)\right(x-{\mu }_x\left)\right(y-{\mu }_y{)}^2\hfill \\ \hfill & +f_{𝑦𝑦𝑦}(y-{\mu }_y{)}^3)\hfill \end{array}$$

Replacing derivatives in equation A1.15 we have (equation (15) - Table I):

$$E\left[f(x,y)\right]=\frac{{\mu }_x}{{\mu }_y}-{\sigma }_{𝑥𝑦}\frac{1}{{\mu }_y^2}+{\sigma }_y^2\frac{{\mu }_x}{{\mu }_y^3}+\frac{E\left[(x-{\mu }_x){(y-{\mu }_y)}^2\right]}{{\mu }_y^3}-\frac{{\mu }_x}{{\mu }_y^4}E\left[{(y-{\mu }_y)}^3\right]$$

First order Taylor expansion for calculation of the variance of function $f(x,y)=x/y$:

To obtain the expression for variance using First order Taylor expansion, we will develop the equation A1.1 without the summatory term, then we have:

$$f(x,y)=f(a,b)+\frac{\partial f(a,b)}{\partial x}(x-a)+\frac{\partial f(a,b)}{\partial y}(y-b)$$

Letting $(a,b)=({\mu }_x,{\mu }_y)$, we have:

$$f(x,y)=f({\mu }_x,{\mu }_y)+f_x({\mu }_x,{\mu }_y)(x-{\mu }_x)+f_y({\mu }_x,{\mu }_y)(y-{\mu }_y)$$

Then, we will get for the expectation:

$$E\left[f(x,y)\right]=f({\mu }_x,{\mu }_y)$$

And for the variance:

$$𝑉𝑎𝑟\left(f\right(x,y\left)\right)=E\left[{(f_x({\mu }_x,{\mu }_y)(x-{\mu }_x)+f_y({\mu }_x,{\mu }_y)(y-{\mu }_y))}^2\right]$$

$$𝑉𝑎𝑟\left(f\right(x,y\left)\right)=E[f_x^2({\mu }_x,{\mu }_y){(x-{\mu }_x)}^2+2f_x({\mu }_x,{\mu }_y)f_y({\mu }_x,{\mu }_y)(x-{\mu }_x)(y-{\mu }_y)+f_y^2({\mu }_x,{\mu }_y){(y-{\mu }_y)}^2]$$

Applying the corresponding derivatives (until first order) we get (equation (16) - Table I):

$$𝑉𝑎𝑟\left(f\right(x,y\left)\right)=\frac{{\sigma }_x^2}{{\mu }_y^2}-2\frac{{\mu }_x}{{\mu }^3}{\sigma }_{𝑥𝑦}+\frac{{\mu }_x^2}{{\mu }_y^4}{\sigma }_y^2$$

Second order Taylor expansion for calculation of the variance of function $f(x,y)=x/y$ (equation A1.6):

$$𝑉𝑎𝑟\left(f(x,y)\right)=E\{[\frac{1}{{\mu }_y}(x-{\mu }_x)-\frac{{\mu }_x}{{\mu }_y^2}(y-{\mu }_y){\frac{-1}{{\mu }_y^2}[(x-{\mu }_x)(y-{\mu }_y)-{\sigma }_{𝑥𝑦}]+\frac{{\mu }_x}{{\mu }_y^3}[{(y-{\mu }_y)}^2-{\sigma }_y^2]]}^2\}$$

$$𝑉𝑎𝑟\left(f(x,y)\right)=E\{[\frac{(x-{\mu }_x)}{{\mu }_y}-\frac{{\mu }_x}{{\mu }_y^2}(y-{\mu }_y){\frac{+[{\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y)]}{{\mu }_y^2}+\frac{{\mu }_x}{{\mu }_y^3}[{(y-{\mu }_y)}^2-{\sigma }_y^2]]}^2\}$$

Expanding A1.19 according to A1.11 we obtain:

$$𝑉𝑎𝑟\left(f(x,y)\right)=E.$$

Developing and simplifying, we get equation ((17) - Table I):

$$\begin{array}{cc}\hfill 𝑉𝑎𝑟\left[f(x,y)\right]\approx & \frac{{\sigma }_x^2}{{\mu }_y^2}+\frac{{\mu }_x^2}{{\mu }_y^4}{\sigma }_y^2-2\frac{{\mu }_x}{{\mu }_y^3}{\sigma }_{𝑥𝑦}+\frac{1}{{\mu }_y^4}E\left[{({\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y))}^2\right]+\frac{{\mu }_x^2}{{\mu }_y^6}E\left[{({(y-{\mu }_y)}^2-{\sigma }_y^2)}^2\right]\hfill \\ \hfill & +\frac{2}{{\mu }_y^3}E\left[(x-{\mu }_x)[{\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y)]\right]+2\frac{{\mu }_x}{{\mu }_y^4}E\left[(x-{\mu }_x)[{(y-{\mu }_y)}^2-{\sigma }_y^2]\right]\hfill \\ \hfill & -2\frac{{\mu }_x}{{\mu }_y^4}E\left[(y-{\mu }_y)[{\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y)]\right]-2\frac{{\mu }_x^2}{{\mu }_y^5}E\left[(y-{\mu }_y)[{(y-{\mu }_y)}^2-{\sigma }_y^2]\right]\hfill \\ \hfill & +2\frac{{\mu }_x}{{\mu }_y^5}E\left[[{\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y)][{(y-{\mu }_y)}^2-{\sigma }_y^2]\right]\hfill \end{array}$$

Heuvelink’s approach for calculation of the mean and variance of function $f(x,y)$based on second order Taylor expansion

In this appendix we wanted to show that the traditional approach used in Appendix A1 gives the same results as obtained by Heuvelink et al. (1989) and Heuvelink (1998). Actually, these references presented a general formula of second order Taylor expansion for n random variables.

The second order Taylor expansion of a function with n variables at the point M=(µ₁,...,µ_n) and neglecting the remainder can be written as (Heuvelink et al. 1989):

$$f\left(\acute{x}\right)=f\left(M\right)+\sum _{i=1}^n\left[(x_i-{\mu }_i)f_{x_i}\right(M\left)\right]+\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n\left[(x_i-{\mu }_i)(x_j-{\mu }_j)f_{x_i{x}_j}\left(M\right)\right]$$

Then, the mean around M is computed as (Heuvelink 1989):

$$E\left(f\right(M\left)\right)=E\left[f\left(M\right)\right]+E\left[\sum _{i=1}^n\left[(x_i-{\mu }_i)f_i\right(M\left)\right]\right]+\frac{1}{2}E\left[\sum _{i=1}^n\sum _{j=1}^n\left[(x_i-{\mu }_i)(x_j-{\mu }_j)f_{𝑖𝑗}\left(M\right)\right]\right]$$

$$E\left(f\right(M\left)\right)=f\left(M\right)+\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n\left[{\rho }_{𝑖𝑗}{\sigma }_i{\sigma }_j{f}_{𝑖𝑗}\left(M\right)\right]$$

Where ${\rho }_{𝑖𝑗}$ is the correlation coefficient between variables $i$ and $j$; ${\sigma }_i$ is the standard deviation for variable $i$ and ${\sigma }_j$ is the standard deviation for variable $j$.

$$𝑉𝑎𝑟\left(\acute{x}\right)=E\left\{{[\sum _{k=1}^n\left[(x_k-{\mu }_k)f_k\right(M\left)\right]+\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n\left[((x_i-{\mu }_i)(x_j-{\mu }_j)-{\rho }_{𝑖𝑗}{\sigma }_i{\sigma }_j)f_{𝑖𝑗}\left(M\right)\right]]}^2\right\}$$

$$\begin{array}{cc}\hfill 𝑉𝑎𝑟\left(\acute{x}\right)=& \sum _{k=1}^n\sum _{l=1}^n{\rho }_{𝑘𝑙}{\sigma }_k{\sigma }_l{f}_k\left(M\right)f_l\left(M\right)\hfill \\ \hfill & +\sum _{k=1}^n\sum _{i=1}^n\sum _{j=1}^n\left\{E\left[(x_k-{\mu }_k)((x_i-{\mu }_i)(x_j-{\mu }_j)-{\rho }_{𝑘𝑙}{\sigma }_k{\sigma }_l)\right]f_k\left(M\right)f_{𝑖𝑗}\left(M\right)\right\}\hfill \\ \hfill & +?\frac{1}{4}\sum _{i=1}^n\sum _{j=1}^n\sum _{k=1}^n\sum _{l=1}^n\left\{\begin{array}{c}\hfill E[(x_i-{\mu }_i)(x_j-{\mu }_j)(x_k-{\mu }_k)(x_l-{\mu }_l)-{\rho }_{𝑖𝑗}{\sigma }_i{\sigma }_j{\rho }_{𝑘𝑙}{\sigma }_k{\sigma }_l]*?*f_{𝑖𝑗}\left(M\right)f_{𝑘𝑙}\left(M\right)\hfill \end{array}\right\}\hfill \end{array}$$

Notice that for n=2, and setting x₁=x and x₂=y, the expression A2.2 becomes:

$$E\left(f\right(x,y\left)\right)=f({\mu }_x,{\mu }_y)+\frac{1}{2}({\sigma }_x^2{f}_{𝑥𝑥}({\mu }_x,{\mu }_y)+2{\sigma }_{𝑥𝑦}{f}_{𝑥𝑦}({\mu }_x,{\mu }_y)+{\sigma }_y^2{f}_{𝑦𝑦}({\mu }_x,{\mu }_y))$$

Expressions A1.3 and A2.4 are the same, so the development of the expectation for an arithmetic operation is the same as shown in Appendix A1. Therefore, equations (10), (12) and (14) of Table I can be derived from A2.4.

For the variance, we will need to expand the expression A2.3 and consider the arithmetic operation regarding the function $f(x,y)$. Then at the point $\mu =({\mu }_x,{\mu }_y)$ we have:

$$\begin{array}{cc}\hfill 𝑉𝑎𝑟(x,y)=& {\sigma }_x^2{f}_x^2\left(\mu \right)+2{\sigma }_{𝑥𝑦}{f}_x\left(\mu \right)f_y\left(\mu \right)+{\sigma }_y^2{f}_y^2\left(\mu \right)\hfill \\ \hfill & +E[{(x-{\mu }_x)}^3-(x-{\mu }_x){\sigma }_x^2]f_x\left(\mu \right)f_{𝑥𝑥}\left(\mu \right)\hfill \\ \hfill & +2E[{(x-{\mu }_x)}^2(y-{\mu }_y)-(x-{\mu }_x){\sigma }_{𝑥𝑦}]f_x\left(\mu \right)f_{𝑥𝑦}\left(\mu \right)\hfill \\ \hfill & +E[(x-{\mu }_x){(y-{\mu }_y)}^2-(x-{\mu }_x){\sigma }_y^2]f_x\left(\mu \right)f_{𝑦𝑦}\left(\mu \right)\hfill \\ \hfill & +E[{(x-{\mu }_x)}^2(y-{\mu }_y)-(y-{\mu }_y){\sigma }_x^2]f_y\left(\mu \right)f_{𝑥𝑥}\left(\mu \right)\hfill \\ \hfill & +2E[(x-{\mu }_x){(y-{\mu }_y)}^2-(y-{\mu }_y){\sigma }_{𝑥𝑦}]f_y\left(\mu \right)f_{𝑥𝑦}\left(\mu \right)\hfill \\ \hfill & +E[{(y-{\mu }_y)}^3-(y-{\mu }_y){\sigma }_y^2]f_y\left(\mu \right)f_{𝑦𝑦}\left(\mu \right)+\frac{1}{4}.\hfill \end{array}$$

Calculation of the variance for sum and subtraction of random variables X and Y

Since for $f(x,y)=x\pm y$ all second order derivatives are equal to zero, we have (equation (11) - Table I):

$$𝑉𝑎𝑟(x,y)={\sigma }_x^2+{\sigma }_y^2\pm 2{\sigma }_{𝑥𝑦}$$

Because all second order derivatives are equal to zero.

Calculation of the variance for the product of random variables X and Y

$$\begin{array}{cc}\hfill 𝑉𝑎𝑟(x,y)=& {\sigma }_x^2{\mu }_y^2+{\sigma }_y^2{\mu }_x^2+2{\sigma }_{𝑥𝑦}{\mu }_x{\mu }_y+2{\mu }_{y}E\left[{(x-{\mu }_x)}^2(y-{\mu }_y)\right]+2{\mu }_{x}E\left[(x-{\mu }_x){(y-{\mu }_y)}^2\right]\hfill \\ \hfill & +E\left[{(x-{\mu }_x)}^2{(y-{\mu }_y)}^2\right]-{\left({\sigma }_{𝑥𝑦}\right)}^2\hfill \end{array}$$

This is the same as equation ((13) - Table I).

Calculation of the variance for the ratio of X to Y

$$\begin{array}{cc}\hfill 𝑉𝑎𝑟(x,y)=& \frac{{\sigma }_x^2}{{\mu }_y^2}-2\frac{{\mu }_x}{{\mu }_y^3}{\sigma }_{𝑥𝑦}+\frac{{\mu }_x^2}{{\mu }_y^4}{\sigma }_y^2-\frac{2}{{\mu }_y^3}E[{(x-{\mu }_x)}^2(y-{\mu }_y)-(x-{\mu }_x){\sigma }_{𝑥𝑦}]+\frac{2{\mu }_x}{{\mu }_y^4}E[(x-{\mu }_x){(y-{\mu }_y)}^2-(x-{\mu }_x){\sigma }_y^2]\hfill \\ \hfill & +2\frac{{\mu }_x}{{\mu }_y^4}E[(x-{\mu }_x){(y-{\mu }_y)}^2-(y-{\mu }_y){\sigma }_{𝑥𝑦}]-\frac{2{\mu }_x^2}{{\mu }_y^5}E[{(y-{\mu }_y)}^3-(y-{\mu }_y){\sigma }_y^2]\hfill \\ \hfill & +\frac{1}{4}\{-4\frac{2{\mu }_x}{{\mu }_y^5}E[(x-{\mu }_x){(y-{\mu }_y)}^3-{\sigma }_y^2{\sigma }_{𝑥𝑦}]+4\frac{1}{{\mu }_y^4}E[{(x-{\mu }_x)}^2{(y-{\mu }_y)}^2-{\left({\sigma }_{𝑥𝑦}\right)}^2]+\frac{4{\mu }_x^2}{{\mu }_y^6}E[{(y-{\mu }_y)}^4-{\left({\sigma }_y^2\right)}^2]\}\hfill \end{array}$$

$$\begin{array}{cc}\hfill 𝑉𝑎𝑟(x,y)=& \frac{{\sigma }_x^2}{{\mu }_y^2}+\frac{{\mu }_x^2}{{\mu }_y^4}{\sigma }_y^2-2\frac{{\mu }_x}{{\mu }_y^3}{\sigma }_{𝑥𝑦}+\frac{2}{{\mu }_y^3}E\left[(x-{\mu }_x)[{\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y)]\right]\hfill \\ \hfill & +2\frac{{\mu }_x}{{\mu }_y^4}E\left[(x-{\mu }_x)[{(y-{\mu }_y)}^2-{\sigma }_y^2]\right]-2\frac{{\mu }_x}{{\mu }_y^4}E\left[(y-{\mu }_y)[{\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y)]\right]\hfill \\ \hfill & -2\frac{{\mu }_x^2}{{\mu }_y^5}E\left[(y-{\mu }_y)[{(y-{\mu }_y)}^2-{\sigma }_y^2]\right]-\frac{2{\mu }_x}{{\mu }_y^5}E[(x-{\mu }_x){(y-{\mu }_y)}^3-{\sigma }_y^2{\sigma }_{𝑥𝑦}]\hfill \\ \hfill & +\frac{E\left[{(x-{\mu }_x)}^2{(y-{\mu }_y)}^2\right]-{\left({\sigma }_{𝑥𝑦}\right)}^2}{{\mu }_y^4}+\frac{{\mu }_x^2}{{\mu }_y^6}E[{(y-{\mu }_y)}^4-{\left({\sigma }_y^2\right)}^2]\hfill \end{array}$$

As we can see, some terms of equation A2.7 do not match A1.21. Actually, we can rewrite these terms as follows.

$$E\left[{[{\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y)]}^2\right]={\left({\sigma }_{𝑥𝑦}\right)}^2-2{\sigma }_{𝑥𝑦}E\left[(x-{\mu }_x)(y-{\mu }_y)\right]+E\left[{(x-{\mu }_x)}^2{(y-{\mu }_y)}^2\right]=E\left[{(x-{\mu }_x)}^2{(y-{\mu }_y)}^2\right]-{\left({\sigma }_{𝑥𝑦}\right)}^2$$

$$E\left[{[{(y-{\mu }_y)}^2-{\sigma }_y^2]}^2\right]=E\left[{(y-{\mu }_y)}^4\right]-2{\sigma }_y^{2}E\left[{(y-{\mu }_y)}^2\right]+{\left({\sigma }_y^2\right)}^2=E\left[{(y-{\mu }_y)}^4\right]-{\left({\sigma }_y^2\right)}^2=E[{(y-{\mu }_y)}^4-{\left({\sigma }_y^2\right)}^2]$$

$$\begin{array}{cc}\hfill E\left[[{\sigma }_{𝑥𝑦}-(x-{\mu }_x)(y-{\mu }_y)][{(y-{\mu }_y)}^2-{\sigma }_y^2]\right]=& E\left[[(x-{\mu }_x)(y-{\mu }_y)-{\sigma }_{𝑥𝑦}][{\sigma }_y^2-{(y-{\mu }_y)}^2]\right]\hfill \\ \hfill =& E\left[(x-{\mu }_x)(y-{\mu }_y){\sigma }_y^2\right]-E\left[(x-{\mu }_x){(y-{\mu }_y)}^3\right]-E\left[{\sigma }_{𝑥𝑦}{\sigma }_y^2\right]\hfill \\ \hfill & +E\left[{\sigma }_{𝑥𝑦}{(y-{\mu }_y)}^2\right]=E\left[{\sigma }_{𝑥𝑦}{\sigma }_y^2\right]-E\left[(x-{\mu }_x){(y-{\mu }_y)}^3\right]\hfill \\ \hfill =& E[{\sigma }_{𝑥𝑦}{\sigma }_y^2-(x-{\mu }_x){(y-{\mu }_y)}^3]\hfill \end{array}$$

Replacing these terms into equation A2.7 we get exactly the same as A1.21, that is the equation ((17) - Table I).

Publication Dates

• Publication in this collection
  Aug 2018

History

• Received
  03 Aug 2017
• Accepted
  09 Oct 2017