Cokriging for optimal mineral resource estimates in mining operations

 

 

R.C.A. Minnitt^I; C.V. Deutsche^II

^IUniversity of the Witwatersrand, Johannesburg
^IIUniversity of Alberta, Canada

 

 

---

SYNOPSIS

Cokriging uses a sparsely sampled, but accurate and precise primary data-set, together with a more abundant secondary data-set, for example grades in a polymetallic orebody, containing both error and bias, to provide improved results compared to estimation with the primary data alone, as well as filtering the error and mitigating the effects of conditional bias. The method described here may also be applied in polymetallic orebodies and in other cases where the primary and secondary data could be collocated, and one of the data-sets need not be biased, unreliable, etc. An artificially created reference data-set of 512 lognormally distributed precious metal grades sampled at 25x25 m intervals constitutes the primary data-set. A secondary data-set on a 10x10 m grid comprising 3200 samples drawn from the reference data-set includes 30 per cent error and 1.5 multiplicative bias on each measurement. The primary and secondary non-collocated data-sets are statistically described and compared to the reference data-set. Variograms based on the primary data-set are modelled and used in the kriging of 10x10 m blocks using the 25x25 m and 50x50 m data grids for comparison against the results of the cokriged estimation. A linear model of coregionalization (LMC) is established using the primary and secondary data-sets and cokriging using both data-sets is shown to be a significant improvement over kriging with the primary data-set alone. The effects of the error and bias are filtered and removed during the cokriging estimation procedure. Thus cokriging using the more abundant secondary data, even though it contains error and bias, significantly improves the estimation of recoverable reserves.

Keywords: Cokriging, primary data-set, secondary data-set, linear model of coregionalization (LMC), ordinary kriging, optimal resource estimates.

---

 

 

Introduction

Cokriging is a widely proclaimed, but infrequently applied estimation technique using two or more sets of collocated or non-collocated data of quite different type and support, over the same domain. The method is applicable in polymetallic orebodies where one metal is well-sampled and the other less so, but may also be applied in orebodies where the primary and secondary data could be collocated, and one of the data-sets need not be biased, unreliable, etc. The method produces better estimates of recoverable reserves than ordinary kriging of either of the two data-sets individually. The application of cokriging demonstrated in this paper involves a primary random function {Z(u), uєA} of sparsely sampled, good quality measurements z(u_α), with mean m(u_α), α=1,..., n evenly distributed across a domain of interest for which grade and recoverable reserves are to be estimated. Over the same domain an often more densely distributed, but poorer quality, non-collocated secondary random function {Y(u), uєA} whose measurements y(ui>_β),with mean m(u_β), β=1,. ..,r include sampling errors and bias, is also available for use in the estimation procedure. Secondary data is often easier, quicker, or cheaper to collect than primary data, but its information content is suspect because of the concern that including poorer quality data with sampling error and bias would compromise the estimate. Examples include open pit mining operations where reverse circulation drilling is augmented by blast-hole sampling data, projects where both RC and diamond drill-hole data are available, tabular orebodies where regional exploration diamond drilling is supported by numerous face or chip samples in adjacent underground workings, or in exploration targets where a primary set of gold grades from diamond drilling is augmented by an abundance of secondary geochemical or geophysical information such as ground magnetic readings. Examples of cokriging two highly correlated variables sampled at different locations within a domain have been described in various applications by Meyers (1982), Wackernagel (2010), Goovaerts (1997), Chiles and Delfiner (2012), and Isaaks and Srivastava (1989).

 

Cokriging

Generally, data-sets will contain more than one variable of interest, and these are usually spatially cross-correlated as expressed in the cross-variogram. Cokriging is a method that takes advantage of the information embedded in the cross-correlation of a second variable in order to minimize the variance of the estimation error (David, 1977; Journel and Huijbregts, 1978). In this paper cokriging is considered in the context of simulated gold grades that have been sampled in two different campaigns, which are then used to estimate the recoverable reserves. Cokriging is the preferred method of estimation where there is an undersampled primary data-set and a closely sampled secondary data-set containing sampling error and bias. The advantage of cokriging with equally sampled data is unclear, but will be investigated in the future.

Isaaks and Srivastava (1989) indicate that in some situations cokriging will not improve an ordinary kriging estimate. Such situations arise when the primary and secondary variables are collocated, meaning that there is not one variable that is undersampled with respect to another, and that the auto- and the cross-variograms are proportional to the same basic model. In such cases the cokriging estimates will be identical to the ordinary kriging estimates. Isaaks and Srivastava conclude that, 'Thus if all the variogram models are "quite similar" in shape and the primary variable is not noticeably undersampled, cokriging will not improve things very much' (p. 405). Cokriging therefore is meant specifically for different, but correlated, non-collocated variables of different quality that are sampled at different densities across a domain. Deutsch and Journel (1998) note that variograms and cross-variograms can be used in the cokriging systems only provided that the requisite constraints on cokriging weights are met. Provided this is so, sample cross-variograms can be transformed into a corresponding cross-covariance C_ZY(h). The cross-covariance can then be modelled with a linear model of coregionalization to be used in the cokriging system. Where a single secondary variable (Y) is considered, the ordinary cokriging (OCK) estimator of a primary variable is:

where λ_α represents weights applied to the n primary data (z), X_β represents weights applied to the r secondary data (y), and the bold u identifies a vector.

Deutsch and Journel (1998) show that while kriging requires a model for the Z covariance, cokriging requires a joint model for the matrix of covariance functions including the Z covariance CZ(h), the Y covariance CY(h), and the cross Z-Y covariance C_zY(h) = Cov{Z(u), Y(u+h)}. In traditional ordinary kriging the sum of the weights applied to the primary variable is unity, while the sum of weights applied to the secondary variable is zero. Provided these constraints are met (Equation [2]) the estimator is unbiased:

The second constraint is that the sum of the secondary weights must be zero, which implies that the secondary data makes no net contribution to the estimate, a severe constraint on ordinary cokriging according to Deutsch and Journel (1998). The error variance is minimized under the constraints in Equation [1] as follows:

where the Lagrangian multipliers µ_α and µ_β account for the two unbiasedness constraints. The covariances are classically obtained from variograms as follows:

C_ZZ(0) and C_YY(0) are the variances of the Z and Y variables, respectively. C_ZY(0) is the estimated cross-covariance of collocated Z and Y data. The direct and cross-variograms are fitted with the well-known linear model of coregionalization (LMC).

where nst is the number of structures and Γl(h), l = 0,...; l is the l th nested structure defined by a shape, three angles, and three ranges; and C^l is the contribution of the l th structure to the variogram model, ZZ, ZY, or YY as appropriate.

In order for the model to be positive definite the following determinant must also be positive definite for all values of l and the cross-product of these coefficients must be greater than zero.

where the 0 identifies a scalar. Deutsch and Journel (1998) suggest that cokriging has not been extensively used in practice firstly because of the tedious joint modelling required by K² variograms when K variables are used, and secondly because of the screen effect of the better correlation among z(u_α) data compared with the weaker correlation between the z(u_α) andy(u_β) data. Furthermore, they describe two other cokriging techniques that are applied, namely simple cokriging and standardized ordinary cokriging.

Simple cokriging (SCK) with the prior assumption of a known mean is mathematically the best estimation method in that there are no constraints on the weights, it minimizes the estimation variance, and it produces estimates that are unbiased. The simple cokriging estimator at location (u₀) of a primary variable z using a secondary variable y is:

where λ _α and χ_βare the weights attributed to the primary and secondary data, and mz and my are the primary and secondary means in a stationary environment. The error variance at (u₀) is given by:

This is minimized by solving the following system of linear equations:

As with simple kriging, this version of cokriging employs the residuals from the data, or requires the standardization of the means such that they are zero. Although it is free of major constraints affecting other methods, the major concern in regard to application of this technique is the strong assumption of stationarity across the domain.

The standardized ordinary cokriging estimator (SOCK) uses a secondary variable that has been standardized such that the means for the primary and secondary variables are the same, and all the weights are constrained to add up to one. In this case Equation [1] can be written:

where λ_α and λ_βare cokriging weights obtained by solving the ordinary cokriging system of equations (Equation [3]) expressed in terms of correlograms (Goovaerts, 1998), with the single condition that:

where m_z=E{Z(u)} and m_y=E{Y(u)} are the stationary means of Z and Y (Deutsch and Journel, 1998).

In this case both the SOCK and SCK weights are linearly related. Goovaerts (1998) has shown that estimation with standardized or non-standardized data results in the same cokriging estimate. The non-standardized estimate is a rearrangement of Equation [4]:

With these constraints in mind, the use of cokriging and the necessity of a linear model of coregionalization (LMC) are now examined.

 

Linear model of coregionalization (LMC)

The linear model of coregionalization (LMC) informs the cokriging process as to how much of the secondary data to use. The cross-variogram carries the information content for the calibration that filters the error from the secondary data. The cross-variogram is constructed using the cross-covariance since the cross-variogram requires collocated data, but the cross-covariance does not; unequally sampled data is used directly. For a cross-covariance that is flat at zero distance it is possible to model an LMC, but all the secondary data would be assigned a weight of zero during the cokriging procedure. If the cross-covariance looks exactly like the direct variogram, it is again possible to model the LMC, but the primary and secondary data would be treated in exactly the same way, i.e. they would be assigned equal weights during the cokriging procedure.

The covariance matrix for jointly distributed primary and secondary variables must be positive definite. A simple way of constructing a valid cross-covariance is through the use of an LMC. While each variable has its own direct variogram and each pair of variables have their own cross-variogram, the LMC uses the same variogram structures and the same ranges for the direct variograms of Z and Y (although the contributions can vary), and for the cross-variograms of Z with Y (see Equation [3]). This ensures that the variance of any possible linear combination of these variables is always positive (Isaaks and Srivastava, 1989).

The model for each of the single variograms may consist of one or more components of the basic models, giving rise to so-called nested structures, which is acceptable since 'any linear combination of positive definite variogram models with positive coefficients is also a positive definite model (Isaaks and Srivastava, 1989, p. 375). These matrices must be positive definite, therefore the following relations must be true:

where l refers to the number of structures (nst), the 0 identifies a scalar, and C^l is the contribution of the l th structure to the variogram model ZZ, ZY, or YY as appropriate.

Isaaks and Srivastava (1989) state it is sufficient that all the eigenvalues of the matrices of b coefficients are positive, but generally the conditions shown in Equation [5] are more widely applied. They also warn that these constraints make the modelling of coregionalization difficult and that one of the single or cross-variograms may not fit the sample variogram well. They suggest that each individual model can be considered as being a small part of a total model and the overall fit judged accordingly. They also provide the helpful suggestion that when compiling a model of coregionalization, a basic model in the single variogram does not necessarily have to appear in the cross-variogram model. However, any basic model in the cross-variogram must be included in all the single variogram models.

 

Practical application of cokriging

In the progressing of a mineral prospect to a going mining concern it is common to accumulate primary data such as metal grade, as well as other types of data from a variety of studies. The latter may have different support and are termed secondary data. Directly or indirectly, the goal of such studies would be related to efforts to derive the best estimate of grades and their distribution across the orebody, the primary data source for such estimates being a limited number of carefully drilled diamond drill-holes. Although generally abundant, secondary data may be suspected of being in error or biased. Rather than applying ad-hoc calibration or correction factors to improve the quality of the secondary data, or simply rejecting secondary data as not usable, cokriging is an attractive alternative that allows secondary data to be used in a theoretically acceptable manner. However, the challenges with respect to the practical implementation of cokriging can be formidable, and it is the aim of this study to provide an understanding of and to document these challenges through a study of synthetic primary and secondary data-sets.

In a real-life mining situation it is essential to assemble the required modelling parameters for the primary data alone, considering issues such as stationarity, managing outliers, declustering, and variography. Attention to such details ensures that the appropriate parameters are available for investigating volume-variance relationships and the kriging of the data. Considerations about stationarity would answer questions such as whether the data belongs together, or does it need an adjustment to account for a trend. The way in which outliers are dealt with is always a contentious point, but where it is clear that outliers need attention, 'following local practice' or conforming to the mining personnel standard is a prudent position to take. When required, declustering would also be considered in order to produce an appropriately weighted distribution that, having taken account of sample redundancy and clustering, is summarized by a mean and a variance. Understanding and extracting the underlying variogram is essential to testing the sensitivity of the deposit to SMU size and volume-variance relationships, and in this case the kriging will require a variogram in original units. One might use normal scores, median indicators, correlograms, relative variograms, or any other tools available to extract information about the underlying spatial structures and to overcome the challenge of limited or skewed data.

In this particular exercise, the use of a simulated data-set and regularly spaced data means that none of the issues described has to be considered. The primary data is therefore suitable for incorporation, can be considered to be stationary, with no outliers, and the distribution within a regular grid obviates the issues that would normally have to be addressed in regard to declustering.

 

Reference data-set

A 400 m x 800 m Gaussian point data-set was simulated at a 1x1 m resolution to give 320 000 points. A lognormal distribution of fictitious metal grades with a mean of 0.70 g/t and a standard deviation of 1.08 g/t was created by transforming the normal scores data. The simulated lognormal reference data-set was sampled such that three highly correlated, non-collocated data-sets all having units of grams per ton (g/t) were created on regular square grids at 10x10 m, 25x25 m and 50x50 m, hereafter referred to as 10 m, 25 m, and 50 m grids. The relative position of these three grids in the first 100x100 m block is shown in Figure 1a

At a point scale (1x1x1 m), assuming a reef thickness of 1 m and a density of 3 t/m³, the project area contains 960 000 t at a grade of 0.70 g/t with a total metal content of 672 kg. True values for tons, grade, and contained metal were averaged for panels at 10x10 m, 25x25 m, and 50x 50 m grid size from the simulation reference data-set. Recoverable reserves were calculated at a cut-off of 1 g/t. The tons, average grade, and recoverable metal at point support, 10 m, 25 m, and 50 m panels are listed in Table I.

 

 

In proceeding from 1 m points to 50 m panels the grade of the recoverable reserves above 1 g/t decreases rapidly with a marginal increase in tons, due to the inclusion of increasing amounts of waste dilution in the SMUs. As a consequence the recoverable metal decreases from 418 kg to 324 kg as the size of the SMUs increases, the balance of the total metal in the deposit being lost in ore sent to waste and waste rock itself.

 

Sampling the reference data-set

Two primary data-sets (Z) sampled on 25x25 m and 50x50 m square grid spacings represent the 512 and 128 good quality drill intersections of a flat tabular orebody, carefully drilled, logged, and assayed diamond borehole cores, respectively. The more abundant secondary data (Y) was sampled at much closer 10x10 m square grid spacing representing 3200 measurements of more easily collected and probably cheaper geochemical or geophysical data. A 30 per cent random error and multiplicative bias of 1.5 were introduced to this secondary data during the sampling procedure. Descriptive statistics for each data-set as well as experimental, direct, and cross-variograms were calculated and modelled (Table II). Data on the 25 m and 50 m grids are cokriged independently with data from the 10 m grid for comparison with the true recoverable reserves.

The effect of introducing sampling errors and bias, compared to the true grades for the10 m grid data is shown in a scattergram in Figure 2. The material lying below the 1 g/t cut-off in the 3rd and 4th quadrants represents 10.9 per cent dilution (waste sent to the mill) and 3.1 per cent lost ore (ore sent to the waste dump), a total of 14.0 per cent.

 

 

The recoverable reserves in the true tons, grade, and recoverable metal for different panel sizes listed in Table I are compared with those of the sampled values listed in Table III.

Recoverable reserves of the sampled 10 m data are estimated to contain nearly 1.4 times more tonnage, and 1.85 times more contained metal, than the true block values. The reason for this is the sampled 10 m panels with introduced error and bias have 35 per cent higher grade, and 36 per cent higher tonnage, than the true values. For the 25 m panels (Table III), the tonnage is lower, but the higher grade means a 16 per cent higher metal content relative to the true values. For 50 m panels the much lower sampled ton is compensated by a much higher sampled grade to give equivalent amounts of metal.

 

Statistics and suitability of the secondary data

Assessing the suitability of the secondary data relative to the primary data for use in cokriging is essential. Cokriging depends heavily on two assumptions about the primary and secondary data: firstly, that data within the domain is stationary, and secondly, that it is sufficiently correlated for the data to be used together. Valid application of cokriging means there should also be a relationship between pairs of the primary and nearby secondary data, as reflected in their correlation coefficient. Where secondary data displays no evidence of correlation with the primary data, it is unlikely to contain suitable information and would have to be excluded from the analysis.

Cokriging also depends on primary and secondary data being near to, but not collocated with, one another. Although there are 3200 points on the 10 m grid with error and bias, only 512 secondary paired data lying within a 10 m radius of the primary data were extracted as shown in Figure 1b. The descriptive statistics for these two sets of data are shown in Table IV.

 

 

The paired data were plotted in a scattergram (Figure 3) to determine the correlation coefficient, examining the effect of the error and bias on the relationship between the data. The higher mean and variance of the secondary data relative to the primary data indicate a 30.6 per cent error in the mean, and 1.27 times bias in the variance. The scatterplot of the primary and secondary data is shown in Figure 3 and indicates a correlation coefficient of 0.689, which is considered significant at the 95 per cent level of confidence (we would reject the null hypothesis that there is no correlation between the primary and secondary data at r(P₅%) = 0.1449, since r = 0.689, and accept the alternative hypothesis: there is a correlation between the variables at 95 per cent significance).

 

 

In this case the difference in the mean of the primary and secondary data will give the bias associated with the secondary data, while the difference in variances provides an indication of the error associated with the secondary data which is due to the sampling error. The significant correlation (0.7) suggests that the secondary data is indeed suitable for incorporation in the estimation of recoverable reserves as part of a cokriging exercise. Permissible limits on the correlation coefficient are unclear, except to say that at or above 0.7 the method suggested here works very well, below 0.2 it does not work at all, and between 0.2 and 0.7 the results may be questionable.

 

Ordinary kriging with the primary data (25 m and 50 m)

Compiling an ordinary kriging model of the primary data-sets, i.e. the drill-hole data on the 25 m and the 50 m grids, is the best evaluation available to miners and mine planners if the benefits of cokriging are not applied. A model produced from ordinary kriging of the primary 25 m and 50 m data to be compared with the results from the cokriging is the only way in which the results of this exercise and purpose of this paper can be demonstrated. The basis of this comparison is t confirm that cokriging using the secondary data results not only in a meaningful improvement in the estimation of recoverable reserves, but also removes the bias and reduces the error in the final estimates. In a real-life analysis using cokriging, the number of samples, n, would be reasonably large, and should be subject to all stationarity considerations

The variograms shown in Figure 4a and 4b for use in the ordinary kriging estimation were modelled on the experimental variogram of the sampled 25 m and 50 m grid data respectively (see Figure 1a), the best option in a real-life situation where the only data available is what has been sampled. A nugget effect of 0.4 is used in both variograms following the C0 value of 0.4 in the point support simulation model. The variogram for the true regularized 25 m data and 50 m data does not provide any information at distances shorter than 25 m and 50 m respectively, but this data must be used since it is the only information available for these data-sets.

The distribution of recoverable reserves shown in Figures 5a and 5b is the complete kriged model for the project area using 25 m and 50 m data kriged into 10 m blocks, respectively. The grade, tons, and metal content of the true, sampled, and kriged models are listed for comparison in Table V.

The recoverable reserves estimated by ordinary kriging in the 10 m panels shown in Figure 5b and 5c are only grades above 1 g/t cut-off. The recoverable reserves for the true 10 m panels shown in Figure 5a are quite erratic, but there is a reasonable correspondence in the location of highs and lows in the true and kriged recoverable reserves for the 25 m data. The differences between the 10 m true and 50 m kriged data are quite marked (Figure 5a and c).

As the size of the panels increases from 10 m to 50 m the effect of smoothing due to the kriging procedure is more evident; the tonnage and grade estimates change predictably given the cut-off grade and distribution.

 

LMC for the case study

The aim of this study is to compare the benefits of cokriging against ordinary kriging using models of coregionalization at different data spacing. The first LMC is for a combination of 3200 secondary data on a 10 m grid with the 512 primary data on the 25 m grid. The second LMC uses a combination of the 3200 secondary data on the 10 m grid with the 128 primary data on the 50 m grid, as shown in Table VI.

 

 

The variograms for the primary (25 m and 50 m) and secondary (10 m) data are required. The variogram for the primary data is not very informative because there are only 512 data for the 25 m grid, and only 128 for the 50 m grid, but especially because at short distances there are no pairs closer than 25 m or 50 m (Figures 8a and 9a). The variograms of the secondary data for the 25 m LMC and the 50 m LMC are identical (Figures 8c and 9c), and have better continuity because of the abundance of data at distances of 10 m. Although the direct variograms do not appear very informative, we know that any structure in the cross-variogram must also be present in the direct variograms. Thi means that the structure seen in the cross-variogram is also present in the direct variograms, but this was not revealed before the cross-variogram was compiled; thus, the secondary data improves estimation of the variograms as well as local estimates.

Plotting the direct variograms for the primary and secondary data and the covariance for the combined primary-secondary data provides a route to the LMCs. The cross-covariance between the primary and secondary variables is plotted as red dots in Figure 6a and 7a. At short lag distances the cross-covariance is positive and decreases to zero with increasing lag distance. Once this structure is inverted (flipped over) using an Excel® spreadsheet to manipulate the data, it takes on the form of a variogram. The nugget effect on the cross-covariance does not appear because there is no collocated data. The intersection point on the y-axis acts as the pivotal horizontal axis around which the covariance is inverted. However, there is no way that the value of the intersection point on the y-axis can be calculated; it can be obtained only by extrapolating the trend of the covariance back to the y-axis. It is then necessary to subtract the values of 1.25 (25 m grid) and 0.93 (50 m grid) so that the nugget effect for the cross-variogram becomes zero as shown in Figures 6b and 7b. The nugget effect on the cross-variogram is not needed since there are no collocated Z-Y values.

It is necessary to fit an LMC which provides the quantitative calibration of the primary to the secondary data and permits the filtering of the error during the cokriging. Fitting of the LMC is done manually, although it is possible to fit the LMC automatically using the appropriate software, but this is usually done only if the LMC is fitted to more than two variables. The first step is to fit the cross-variogram as simply as possible in all directions. The sill is equal to the constant used to invert the covariance, and the nugget effect is zero. By applying the structure used in the cross-variogram it is now necessary to fit the direct primary and direct secondary variograms. There may be nested structures or completely different models in each of the direct primary and direct secondary variograms that do not appear in the cross-variogram. However, the constraint for the method is that any structures that appear in the cross-variogram must also appear in the direct primary and secondary variograms. This is the reason for starting with the cross-variogram, followed by the direct variograms.

Examples of this procedure for compiling the LMCs for the 25 m and the 50 m data are illustrated in Figures 6 and 7. In these figures the experimental points of the covariance are extrapolated back to intersect the y-axis at 1.25 for the 25 m data, and 0.93 for the 50 m data, when the distance is zero. These then become the points about which the rest of the covariances are rotated to give lines that look like variograms (Figures 8 and 9).

First LMC: The first model, using the 25 m primary data and the 10 m secondary data, requires the variograms for the primary and secondary data, and the covariance for the primary-secondary combination. The model variogram is used to identify the nugget effect and the sill for the direct primary and secondary data, as well as for the cross-variogram between primary and secondary data.

The equations for the variogram models shown in Figure 8 are as follows:

For a valid LMC it is also essential that the matrix be positive definite, meaning that the difference between the products of the primary and secondary direct variograms and cross-variogram (Equation [6]) must be positive (Goulard, 1989). This implies that all the direct variograms must be positive. Also, if there is some spatial structure in the cross-variogram between two variables it must be evident in the corresponding direct variograms as well, which implies that there is correlation between the two variables. On the diagonal of the kriging matrix the nugget of the primary and the nugget of the secondary are used, but the nugget effect of the cross-variogram is never used. The cross-variogram is also positive definite because there could be nugget effects on both direct variograms, but there is no nugget effect on the cross-variogram; only the structured portion of the covariance has a nugget effect. The positive definiteness condition for the first LMC is satisfied as:

The important implication here is that the cross-variogram will usually lie between the two direct variograms, so that in actually fitting the variograms it is convenient to fit the cross-variogram first so that the cross-variogram structures are defined, and then fit the direct variograms, because evidence of spatial structure in the cross-variograms constrains the direct variograms and it must also be evident in the direct variograms (Goulard and Vlotz, 1992). The fitting process is iterative but it is important to keep the range of the individual structures constant for all variograms.

Second LMC: The second model, using the 25x25 m primary data and the 10 m secondary data, also requires the variograms for the primary and secondary data, and the covariance for the primary-secondary combination. The model variogram is used in the same way as the first LMC to identify the nugget effect and the sill for the direct primary and secondary data, as well as for the cross-variogram between primary and secondary data.

The models for the direct variograms for primary and secondary data, and the cross-variogram for the combined primary and secondary data shown in Figure 9, are very acceptable fits of the equally acceptable experimental variograms. The equations for the variogram models shown in Figure 9 are as follows:

The positive definiteness condition for the second model is also satisfied as:

Having established that the cross-variograms are positive definite, it is now possible to perform the cokriging.

 

Comparing cokriging and ordinary kriging

The final step in the procedure is to create a block model that correctly accounts for all of the available data and at the same time demonstrates improvements in local and global estimates, by filtering the error and without transferring a bias. Firstly a single data file containing the primary and secondary data for use in the cokriging software routine is compiled. The compilation of the final data file is generally software-dependent, and many packages simply omit the cokriging capability because of the difficulties associated with the user interface to calculate the LMC; a further reason that cokriging is not used as often as it might be. Apart from additional implementation details, the cokriging and kriging algorithms are the same, and all the parameters for the cokriging routine are virtually identical to the kriging parameters. Usually the abundance of secondary over primary data in the cokriging routine means that a decision must be made about how much of the secondary data to use without overwhelming the process and disadvantaging the primary data. The main reason for limiting the amount of secondary data used is to avoid smoothing. If the model is being developed for a long-term resource estimate, over-smoothing may mean that the recoverable reserves will be estimated too low. A second reason to limit secondary data is that secondary data may be closer to the point being estimated than primary data, but the latter is more relevant than the secondary data. For this reason the cokriging software allows one to specify how many of each data type, usually a maximum of 12 each, will be used in the cokriging process.

The following analysis of the recoverable reserves above 1 g/t gold, derived from the kriging and cokriging of the 25 m and 50 m data-sets into 10 m panels, is compared against the known true data within the project area. The tonnage, grade, and recoverable metal are compared against the true data for the 10 m panels, in order to demonstrate that cokriging improves local and global estimation, that the bias has not been transferred, and that the error has been filtered out of the final result. Generally, ordinary kriging with the primary data will be unbiased, but the different varieties of cokriging with the primary and secondary data are also going to be unbiased. The caveat is that the mean of the primary, and particularly the secondary, data must be known, which is easy enough through the data.

The bias is dealt with through the means for the primary and secondary data, whereas the error is managed by essentially depressing the covariance between the primary and secondary variables C_zY(h).

The secondary variable has error content, and for this reason we would not want the variogram of the secondary data to have too much weight (Figure 10). The control of the weight attributed to the secondary variable is by reducing the sill. The amount of error decreases as the sill moves upwards and increases as the sill moves downwards. At the two extremes, i.e. when C_ZY(h) equals the variance of the primary variable, there is no error at all, and when C_ZY(h) is flat at zero lag there is only error and the data is not weighted at all. In effect this is a continuum that must be calibrated. The difference of the mean grades, (indicating bias), and the difference between variances, (indicating error), for the cokriging results using data on the 25 m and 50 m grid are closer to the true values than the ordinary kriged results. Again, these outcomes indicate that the cokriging procedure is highly efficient at removing the bias and filtering the error that were shown to be embedded in the secondary data.

Kriging and cokriging data on 25 m grid

Sample data on the 25 m grid has been kriged (OK) and cokriged into 10 m blocks (Figures 11b and 11c), and is visually compared with the known true values of 10 m blocks (Figures 11a).

The tons, grade, and recoverable metal for the 25 m kriged and cokriged panels are summarized in Table V.

Kriging and cokriging data on 50 m grid

Sampled data on a 50 m grid has been kriged (OK) and cokriged into 10 m blocks (Figures 12b and 12c, respectively), and compared with the known true distribution of 10 m blocks (Figures 12a) in the project area.

Again, it is evident that the co-kriging output is visually more similar to the known true distribution than the ordinary kriging outcome, shown in Figure 12.

 

Summary

A summary of the cokriging and ordinary kriging estimates using the data-sets on 25 m and 50 m grids to give tons, grade, and recoverable metal in the 25 m and 50 m panels is presented in Table VII. The reference data-set for the points and the 10 m, 25 m, and 50 m panels is the true recoverable reserves above 1 g/t in the deposit and represents the base case against which the sampled, ordinary kriged, and cokriged data is compared. The sampled data-set for the 10 m grid size is the poorer quality secondary data, while that for the 25 m, and 50 m grid sizes represents the higher quality primary data for use in cokriging.

Differences between the reference and the sampled data-sets at 10 m, 25 m, and at 50 m are a reflection of random selection at that grid spacing. It is noteworthy at each grid size that the sampled data-set marginally overstates the grade and contained metal.

The results of the kriging and cokriging into 10 m panels using sample data on the 25 m grid indicate that the tons predicted from ordinary kriging (195 000 t) are 11.0 per cent less than the known true tonnage (219 000 t) while cokriging predicts (221400 t), 1.1 per cent more than the known tonnage. The metal predicted in the ordinary kriged panels (315 kg) is 15.1 per cent lower than the true metal content (371 kg), while that for the cokriged panels (425 kg) is 14.6 per cent higher. The difference between the true known grade (1.69 g/t) and the ordinary kriged grade (1.61 g/t) is 4.7 per cent, but for the cokriged grade (1.92 g/t) it is 13.6 per cent higher. Cokriging provides tonnages and grade that are increased relative to the true values, representing an increase over the true metal content of 371 kg, but this is 0.5 per cent closer to the real content than the result obtained from ordinary kriging.

The predicted tons and metal for ordinary kriging and cokriging into 10 m panels using data sampled on the 50 m grid indicates that the tonnage from ordinary kriging (155 000 t) is 28.9 per cent lower than the true known tonnage (218000 t), while for the cokriging (198 300 t) it is only 9.0 per cent less than the true tonnage. On the 50 m grid size the amount of recoverable metal for the ordinary kriged panels (255 kg) is 21.3 per cent less than that for the true panels (324 kg), while for the cokriged panels (329 kg) it is 1.5 per cent more. The mean grade for the ordinary kriged estimates (1.64 g/t) is 10.1 per cent higher than the true known grade (1.49 g/t) of 10 m panels, while for the cokriged estimates (1.69 g/t) the difference between the means is 11.4 per cent higher.

In tons and recoverable metal the cokriging estimates are closer, if only marginally in the case of recoverable metal, to the known true values than the ordinary kriged estimates. The much better estimate of tonnage using cokriging relative to ordinary kriging is offset by the higher value for cokriged grade, but there is an overall marginal improvement in estimation of contained metal. These benefits due to cokriging have been further resolved by comparing the tons, grade, and recoverable metal for designations referred to as lost ore, dilution, ore, and waste estimated by ordinary and cokriging, and are given in Table VIII for the 25 m grid samples and 50 m grid samples. These figures, while useful in their own right, are better displayed as percentages in the relevant areas of Figure 13. The benefits of cokriging over ordinary kriging, provided the underlying assumptions are valid, are illustrated in the way that mine production will be allocated as lost ore, dilution, ore, and waste rock in a cross-plot of the kriged and cokriged estimates against the known true grade of materials in the designated areas shown in Figure 13.

By inspection of Figure 13 it is evident that for each classification category (lost ore, dilution, ore, and waste), standardized cokriging (SOCK) improves the estimation outcomes over ordinary kriging (OK). There are considerably fewer tons of lost ore for both the 25 m (29 400 t as against 67 500 t) and 50 m (47 100 t as against 130 500 t) sampled data, as well as fewer lost kilograms of metal (23.3 kg versus 47.8 kg) in the 25 m grid and in the 50 m grid (61.2 kg versus 230.8 kg).The same is true for dilution, in that there is 18 600 t and 36 900 t less dilution, and 32.9 kg and 62.7 kg improved metal recovery in the 25 m and 50 m sampled grids respectively, due to cokriging. In the case of ore, cokriging yields more tons, better grade, and more metal compared to ordinary kriging at the 25 m and 50 m sampling size. For the 25 m grid size there are more tons, lower grade, and less metal sent to waste, and the same is true for the 50 m grid size except that 8 kg of metal are lost to waste. The improvements due to cokriging over ordinary in the four categories lost ore, dilution, ore, and waste are shown as percentages in Figure 13.

The benefits of cokriging are particularly noticeable in the area representing lost ore, where cokriging gives fewer lost tons (3.1 per cent and 4.9 per cent) and less recoverable metal (3.4 per cent and 8.7 per cent) compared to ordinary kriging (7.0 per cent and 13.6 per cent) for tons (6.7 per cent and 29.0 per cent) and for metal in the 25 m and 50 m sampling grids (Figure 13). Lost ore is estimated to be below the cut-off grade, but it is truly economic in that it could make a positive contribution to the mining operation. However the company does incur a lost opportunity cost. This ore will never be accounted for in the balance sheet, nor will it add any value to the mine, except perhaps towards the end of the mine life when plant superintendents attempt to feed the plant from the low-grade stockpile.

Cokriging also yields less dilution in terms of tons (4.7 per cent versus 6.6 per cent ) and recoverable metal (8.5 per cent versus 12.6 per cent ) compared to ordinary kriging for the 25 m grid. The same is true for the 50 m grid, where tons (3.9 per cent versus 7.8 per cent ) and recovered metal (6.5 per cent versus 13.7 per cent ) are improved by cokriging (Figure 12). Dilution is mine production that is estimated to be above the cut-off, but is truly uneconomic in that its contribution is far outweighed by the costs incurred during milling and processing. Such material does show up in the balance sheet in that it raises milling and processing costs, and constitutes an opportunity cost in that it occupies milling capacity that could be better utilized by truly high-grade ore.

 

The work flow in practice

Practical application of cokriging in mining operations would broadly follow a four-stage work flow involving the making of assumptions, exploratory data analysis, variography, and estimation by cokriging. In addition, the underlying assumptions that must be established before cokriging is applied are summarized here.

Assumptions

Some of the assumptions would preclude the use of this technique, but where there are opportunities to apply it the benefits over ordinary kriging are substantial. The chief weakness of the cokriging method is non-stationarity, for example if the bias changes across the domain. The Z variable may be accompanied by another variable, Y, that is biased high and not perfectly correlated, but the method works well here because the data is strongly correlated, and in particular the bias is stationary. The bias may be as much as 30 per cent but provided it has the same average (uniform) across the whole domain the method is still applicable. A difficult case occurs if the data is biased in a non-stationary way. It is possible that data in the domain may on average have a bias that is close to zero, but if the bias is non-stationary application of this method will be flawed. The main assumption therefore is that of stationarity across the domain under consideration.

Exploratory data analysis

Evaluation of the primary data: this would involve a detailed evaluation of the primary data in terms of understanding the distribution of the variables using diagrammatic representations of the data in histograms, cumulative and probability distribution functions (cdf and pdf), as well as compiling the basic descriptive statistics, mean, variance, and coefficient of variation. The data location and spacing will indicate the necessity of declustering of the data to ensure that representative statistics are obtained. Plotting the declustered mean against the cell size will give an indication of the optimal cell size to be used in the analysis. Alternatively, the cell size can be estimated from the drill-hole spacing in sparsely sampled areas.

Evaluation of the secondary data: secondary data is subject to the same analysis as primary data except that there should be an evaluation of the bias and error content of the secondary data.

Confirmation of the stationarity of the populations: The strongest and most important assumption underlying the use of cokriging is that of stationarity within the domain. Where there is a trend it would need to be modelled. Another assumption concerns the nature of the data in that there should be two sets of non-collocated data. This method will work where the data is collocated. The primary data-set generally consists of few, high quality data while the secondary data-set comprises abundant, but poorer quality data that may contain sampling errors and sampling bias. The method also assumes that there is sufficient primary data to give a meaningful variogram, and usually where there is abundant secondary data there is no restriction on obtaining a stable variogram.

Variography

Establishing the experimental variogram for the primary data: this is usually not as easy to achieve as might be expected. Generally the primary data will have been drilled on a widely spaced grid, which means that the first point on the variogram occurs at a long lag. In addition there is very little to inform the short-range variability of the data and in particular the nugget effect is hard to estimate; in reality it is usually less than might be expected from an inspection of the experimental variogram.

Establishing the experimental variogram for the secondary data: this is generally easier than for the primary data, simply because there is more data and consequently the experimental variogram for the secondary data is better-behaved than that for the primary data.

Establishing the experimental variogram for the cross-variogram: this is done by firstly calculating the covariance for the primary and secondary data. This is a decreasing function with increasing lag, but its importance lies in the smooth extrapolation through points on the covariance curve to their point of intersection on the y-axis. This is used to invert the covariance in order to produce the cross-variogram.

Fiting the linear model of coregionalization: the main requirement of the LMC is that it be positive definite. This can be achieved by following the advice presented above.

Kriging estimation

Ordinary kriging: the only reason that ordinary kriging is performed is to provide a basis for comparison with the cokriging results. It might even be necessary to cross-validate the output to confirm that the ordinary kriging produces acceptable results.

Standardized cokriging: the means for the primary and secondary data must be used in the cokriging procedure. In addition it is important that the primary and secondary data be combined in a specific way (depending on the software) during the cokriging process. An important step in the analysis of the output is a comparison between the efficiencies of ordinary and cokriging. This is probably best done by comparing the improvements in a scatterplot of the four categories of production materials: lost ore, dilution, ore, and waste.

 

Conclusions

This paper demonstrates a little-used geostatistical technique that combines different sets of data, from different sources, with different qualities, in a cokriging procedure. The data needs only to be weakly correlated and not necessarily be collocated either. Cokriging provides better estimates of recoverable resources at local and global scales, compared to ordinary kriging of the best available data, and is marked by the highly significant benefit that errors and bias due to poor sampling in the secondary set of data are not transferred to the resulting estimates.

 

Acknowledgements

This paper arose from the Citation Course in Geostatistics presented by the second author in Johannesburg, South Africa in August 2011. The comments, discussion, and inputs arising from the other 15 attendees, some of which are incorporated here, are acknowledged. The authors would like to acknowledge Dr Chris Prins for editing various drafts of this paper.

 

References

Chiles, J. and Delfiner, P. 2012. Geostatistics: Modelling Spatial Uncertainty, 2nd edition. John Wiley and Sons, New York. 576 pp.         [ Links ]

David, M. 1977. Geostatistical Ore Reserve Estimation. Elsevier, Amsterdam. 364 pp.         [ Links ]

Deutsch, C.V. and Journel, A. 1998. GSLIB: Geostatistical Software Library and User's Guide, 2nd Edition. Oxford University Press, New York. 369 pp.         [ Links ]

Goovaerts, P. 1997. Geostatistics for Natural Resources Evaluation. Oxford University Press. 483 pp.         [ Links ]

Goovaerts, P. 1998. Ordinary cokriging revisited. Mathematical Geology, vol. 30, no. 1. pp. 21-42.         [ Links ]

Goulard, M. 1989. Inference in a coregionalisation model. Geostatistics, vol. 1. Armstrong, M. (ed.). Kluwer, Dordrecht. pp. 397-408.         [ Links ]

Goulard, M. and Vlotz, M. 1992. Linear coregionalisation model: tools for estimation and choice of cross-variogram matrix. Mathematical Geology, vol. 24, no. 3. pp. 269-286.         [ Links ]

Isaaks, E.H. and Srivastava, R.M. 1989. An Introduction to Applied Geostatistics. Oxford University Press, New York. 561 pp.         [ Links ]

Journel, A.G. and Huijbregts, J. C. 1978, Mining Geostatistics. Academic Press, New York. 600 pp.         [ Links ]

Meyers, D.E. 1982.Matrix formulation of co-kriging. Mathematical Geology, vol 13, no. 3. pp. 249-257.         [ Links ]

Wackernagel, H. 2010. Multivariate Statistics. Springer-Verlag, Berlin. 387 pp.         [ Links ] ♦