- Email: [email protected]
Theoretical Basis of the Approach: The Theory of Regionalized Variables
Chapter
5
THEORETICAL BASIS OF THE APPROACH: THE THEORY OF REGIONALIZED VARIABLES
S Y N O P S I S — In this chapter, the theory of regionalized variables, as formal ized by Matheron (1965) will be presented, together with an attempt to show where it fits in the standard statistical literature. Bearing in mind that what is important to us is applications, the theoretical references will be limited and only what is needed for our purpose will be discussed. Thus, after defining a regionalized variable, the different hypotheses which one may have to use in order to make statistical inference possible, will be introduced. Since the estimators which will be later considered will be linear, properties of linear combinations of regionalized variables will be examined, and as mining samples are not point samples, we will see how to take care of non-point samples. Then the theoretical expression of estima tion variances and block variances as a function of the variogram will be given. Examples of the variograms most commonly used, will be shown and the general class of admissible functions for a variogram will be briefly reviewed. The notion of chaotic component or nugget effect will be intro duced and discussed. 5.1
FOREWORD
Most of the material in this chapter follows the presentation of Matheron ( 1 9 6 5 ) . In fact, this was also presented b y Matern ( 1 9 4 7 , 1 9 6 0 ) concerning forestry problems, and o n e reference t o t h e variogram can be traced back t o 1 9 2 6 in a paper b y Langsaetter, again a b o u t forest surveys. The sampling properties of regionalized variables, or t h e study of topographic variation as it has often been referred t o , have been discussed by Mahalanobis ( 1 9 4 4 ) , Nair ( 1 9 4 4 ) , Quenouille ( 1 9 4 9 ) , J o w e t t ( 1 9 5 5 ) , Williams ( 1 9 5 2 , 1 9 5 6 ) , Zubrycki ( 1 9 5 7 ) , Whittle ( 1 9 5 4 ) , Finney ( 1 9 5 0 ) . However, n o n e of these papers deal w i t h geological observations, let alone mining. A l t h o u g h the theory m a y n o t be n e w t o professional statisticians, o n e should recognize that very specific problems occur in the mining industry and this leads t o specific theoretical d e v e l o p m e n t s ajid different t e r m i n o l o g y . A m o n g other particularities, grade values are m u c h less c o n t i n u o u s than m o s t of biological or ecological data considered b y biometricians. Beside this, mines are truly three-dimensional b o d i e s and samples are n o t " p o i n t s " , while w h a t appeared in earlier literature m o s t l y concerns two-dimensional p o i n t processes. It is n o t the intention of this t h e o r y t o change the existing statistical t e r m i n o l o g y , but since t o mining professionals b o t h languages are equally u n k n o w n , it
92
d o e s n o t appear unreasonable t o try and d e v e l o p a presentation o f the t h e o r y w h i c h is specially tailored t o mining problems. D o i n g this has t h e usual drawback of disappointing b o t h statisticians and mining p e o p l e ; w e feel, however, that it is a possible avenue. 5.2 DEFINITION OF A R E G I O N A L I Z E D V A R I A B L E
Let x be a p o i n t in space R and z(x) t h e value of t h e function w e are interested in at p o i n t x. S u c h a f u n c t i o n will be called a regionalized variable. T h e grade of equal-size samples is such a regionalized variable. S o is t h e thickness of a d e p o s i t at o n e p o i n t or t h e specific gravity of t h e ore. Such a function is usually highly variable and n o n - c o n t i n u o u s and c a n n o t be studied directly. The structure w h i c h m a y appear in its variations will be studied b y examining its increments. The basic idea of the t h e o r y is t o consider such a function z(x) where x is a p o i n t or a vector o f R as o n e realization of a random f u n c t i o n ( R F ) Z(x). Thus w e are turning a perfectly well defined unique numerical value i n t o a realization of random process. We have o n l y o n e realization of that random function and t h e problem w e are faced w i t h consists in finding the characteristics of t h e R F Z(x) in order t o make the estimation of u n k n o w n points possible. 3
n
This decision t o l o o k at t h e values of t h e process w e are interested in as a realization of a random function is just an epistemological decision w h i c h is taken w h e n recognizing that deterministic m e t h o d s w o u l d n o t help us; it is n o t possible t o consider any e x p e r i m e n t which w o u l d c o n c l u d e that a deposit is n o t a realization of a random function (Matheron, 1 9 7 6 a ) . Since this is a very permissive m o d e l , t h e counterpart is o b v i o u s , it d o e s n o t h e l p m u c h ! H y p o t h e s e s have t o be introduced a b o u t t h e t y p e of random func t i o n which will be considered. 5.3 T H R E E P L A U S I B L E H Y P O T H E S E S
Going from t h e m o s t restrictive t o t h e more general, o n e has seen in t h e geostatistical literature three sets o f possible assumptions. 5.3.1 The weak-stationarity
assumption
(second-order stationarity)
This assumption,seldom found in natural p h e n o m e n a , consists in t w o conditions: — the e x p e c t e d value of t h e regionalized variable Z(x) is t h e same all over t h e field of interest; — t h e spatial covariance of t h e regionalized variable Z(x) is t h e same all over t h e field of interest.
93
Hence t h e e x p e c t e d value is: E[Z(x)]
= m
(5.1)
and t h e covariance is: E{[Z(x)-m][Z(x
+ h)-m]}
where h is a vector in R . random f u n c t i o n Z(x) is: n
V A R [Z(x)]
=
= K(x,x
+ h) = K(h)
(5.2)
/i is t h e length of vector /*. The variance of t h e
E{[Z(x)-m] }
= K(0)
2
(5.3)
h e n c e t h e process has a covariance o n l y if V A R [ £ ( # ) ] is finite. 5 . 3 . 2 The intrinsic
assumption
In m a n y deposits, as s h o w n b y Krige ( 1 9 5 1 ) , such a thing as a finite variance d o e s n o t exist; if o n e considers t h e variations of grade rather than t h e grade itself, it has a finite variance. H e n c e , considering t h e increments of t h e f u n c t i o n , Z(x) — Z(x + /*), rather than t h e function itself, o n e can make t h e following assumptions: E[Z(x
+h)-Z(x)]
VAR[Z(x
= 0
+ h)-Z(x)]
(5.4)
= 2y(h)
(5.5)
This last definition is t h e definition of t h e variogram. One can see that this definition can be rewritten: VAR[Z(x
+ h)-Z(x)]
= E[Z(x
+
h)-Z(x)-E{Z(x+h)-Z(x)}]
2
and since: E[Z(x
4- h) - Z(x)]
VAR[Z(x
=
+ h)-Z(x)]
0 = E[Z(x
+ h)-Z(x)]
(5.6)
2
This is the intuitive definition of t h e variogram w h i c h was given in Chapter 4. Also if a R F is second-order stationary it is also intrinsic and its variogram is: V A R [ Z ( x + h) - Z(x)]
= V A R [ Z ( x + h)] + VAR[Z(jc)] -2COV[Z(x
+
h),Z(x)]
or using t h e previous definition and dividing b y 2 : 7(A) = K(0)~K(h)
(5.7)
Remark on notation. As m a n y p e o p l e have already e x p e r i e n c e d , it is n o t l o n g in geostatistics before o n e encounters n o t a t i o n problems. We have just written our first formulae. Here is t h e first p r o b l e m : One m a y w o n d e r w h y
94
t h e factor 2 in t h e definition of 2y(h). It is because t h e quantity w h i c h is m o s t l y used is y(h) rather than 2y(h). C o n s e q u e n t l y y(h) has been called t h e semi-variogram or half-variogram. Experience s h o w s that almost every b o d y calls y(h) t h e variogram. One has even seen papers w h e r e 2y(h) was defined as the d o u b l e variogram. We could n o t favour that word b u t suggest that for the sake of simplicity, o n e should take as definition of the vario gram: V A R I O G R A M = y(h)
= \ V A R [ Z ( * + h) - Z(x)]
(5.8)
Point variogram. The previously defined variogram applies t o variable Z(x) defined at p o i n t x. One will often call such a variable a p o i n t variable and its variogram a p o i n t variogram as o p p o s e d t o n o n - p o i n t variables (see Section 5.4.2).
5 . 3 . 3 The hypotheses
of universal
kriging
A third kind of h y p o t h e s i s , less restrictive than t h e previous o n e s , assumes that t h e second m o m e n t of t h e R F or its increments has s o m e properties o f stationarity within a vicinity of restricted size and t h a t t h e e x p e c t a t i o n w h i c h is n o longer stationary varies in a regular manner in such a vicinity. If x and y = x + h are taken in the same vicinity: E[Z(x)]
= m(x)
or
E[Z(x)
-Z(y)]
= m{x)-m(y)
(5.9)
with: k
m{x)
=
X d^ix)
(5.10)
1 = 0
the fi(x) being k 4- 1 i n d e p e n d e n t functions and t h e d numerical coefficients. COV[Z(x),Z(x
+ h)]
= K(h)
or
V A R [ Z ( * + h) - Z(x)]
=
z
unknown
2y(h)
In this third kind of h y p o t h e s i s , n o t o n l y t h e covariance or variogram f u n c t i o n has t o be defined from the experimental values but also t h e size of t h e vicinity where the h y p o t h e s e s remain valid, the nature of t h e functions fx(x) as well as their n u m b e r and t h e values of t h e coefficients d which are a function of t h e position of t h e vicinity in the field. The progression in t h e h y p o t h e s e s can be seen. One m a y w o n d e r w h y s t o p here. In fact the last d e v e l o p m e n t s of the t h e o r y (Matheron, 1 9 7 3 ) n o w con sider generalized increments and stationarity properties for these generalized increments. It can also be seen that each h y p o t h e s i s covers t h e previous ones. N o w the t h e o r y of generalized random intrinsic functions of order k covers all t h e previous cases. It is far b e y o n d t h e scope of this v o l u m e and will n o t be discussed. A review of it is given in Delfiner ( 1 9 7 6 ) . u
95 5.4 L I N E A R C O M B I N A T I O N S A N D A V E R A G E
5.4.1 Statistical
properties
VALUES
of the linear combinations
of random
variables
Linear c o m b i n a t i o n s of random variables are random variables (see for instance, Feller, 1 9 5 7 ) . The c o m b i n a t i o n s m a y be defined from c o m p o n e n t s of a R F o n a finite set o f points S[: n
Z*
= E
di Z(x ),
jc,-eS =
t
[x . . x
.x ] n
i =l
or an infinite set o f points which is a d o m a i n V of v o l u m e V. Z*
= 1/V
f
Z(x) dx
(5.11)
Y.
The variance o f Z* can be expressed for each of t h e previous sets o f hypotheses. If Z(x) is a second-order stationary RF, whatever t h e a which is expressed from t h e covariance of Z(x):
i9
VAR[Z*]
= t
Z* has a variance
t<*i
i=i
(5.12)
y=i
or o n an infinite set: V A R [ Z * ] = 1/V
f dx f K{x-y)dy Jy_ JV
2
(5.13)
If Z(x) is an intrinsic RF, it has b e e n s h o w n (Matheron, 1 9 7 1 , p . 5 6 ) that Z* has a finite a-priori variance if and o n l y if t h e c o n d i t i o n 2 " a = 0 is realized. In this case: = 1
VAR[Z*]
= - £
£
i = l
Z* = 1/V f Z(x) v
t o dxjV
atajjixi-xj)
(5.14)
j = l
dx is a linear c o m b i n a t i o n where t h e weights are all equal
and since 1/V f
v
5.4.2 Non-point
x
variables:
dx = 1, t h e n Z* d o e s n o t have a finite variance. smoothing
When w e are concerned with t h e average grade of a v o l u m e V centered at point x, w e can write the grade of that b l o c k as in the previous paragraph: Z[V(x)]
= Z
v
Z(x)dx
V being the v o l u m e o f set V, as m e n t i o n e d before.
96
Z[V(.x:)] is s o m e t i m e s called a s m o o t h e d variable and it might be interest ing t o express its variogram as a function o f t h e variogram of the p o i n t variable Z(x). Let k(x) designate a function w h i c h is equal t o o n e if x is included in V and zero if n o t . T h e n w e call t h e geometric covariogram of V t h e quantity: (5.15) It can easily be s h o w n that t h e variogram of variables of size V can be written: 7 v ( * ) = f 7(A +x)P(x)dx-$
(5.16)
P(x)y(x)dx
Note. The geometric covariogram has an interesting geometrical expres sion. Let us consider a set V and a vector h. N o w let us define V as the set obtained b y translating V b y a quantity h. N o w the intersection of V and V is such t h a t j c e V a n d * + heV. H e n c e P ( h ) = M e a s ( V n V )/V . More simply t h e geometric covariogram of a piece of core of length / will be \ l — h |// . h
h
2
h
2
5.5 T H E O R E T I C A L E X P R E S S I O N O F V A R I A N C E S
5 . 5 . 1 The extension
variance
Let V be a v o l u m e t h e grade of w h i c h Z(V) is u n k n o w n and t o be esti m a t e d b y the k n o w n grade Z(v) of v o l u m e v. Then Z(V) and Z(v) can be written as stochastic integrals, t h e variances and covariance of w h i c h can be s h o w n t o be functions o f t h e variogram:
The variance o f t h e error of estimation Z(V) — Z(v) w h i c h is m a d e w h e n assigning t o V t h e grade o f v exists since this is a linear c o m b i n a t i o n where t h e sum of weights is zero and it is s h o w n t o be (Matheron, 1 9 6 5 ) .
VAR[Z(V)-Z(*)]
=
-r
T
Ax
7 | x — y | dy
(5.17)
97
This relationship is t h e theoretical k e y t o any estimation problem since it allows a numerical characterization of t h e magnitude of t h e error involved in an estimation procedure as s o o n as y(h) is k n o w n . N o t e that since x and y are points in R , t h e integrals are s e x t u p l e o n e s . This m a y involve s o m e c o m p u t a t i o n a l problems. T h e y will be discussed later. This c o n c e p t and formulation m a y s e e m fairly simple and straightforward, but apparently this n o t i o n had never b e e n clearly defined until Matheron did it. 3
Estimation by a series of samples Usually, o n e wishes t o give t o a block V n o t t h e grade of a b l o c k y_ b u t rather t h e average grade of a series of punctual values, say: 1 Z*
N
I
= -
Z(x ) t
In this case t h e previous formulation is reduced t o : °N
= J^Z
l(x -x)dx~^\^
dx\^j(x-y)dy
i
lyixi-xj)
(5.18)
This variance is usually called an e s t i m a t i o n variance. The programming of its c o m p u t a t i o n will be seen in a later chapter, this being an essential t o o l for t h e appraisal of t h e quality of a sampling pattern. 5 . 5 . 2 The variance
of a block and the Krige's
relationship
We are n o w interested in t h e variance o f t h e grade of b l o c k s of a given size in t h e deposit. For punctual " b l o c k s " w e will use t h e q u a n t i t y a ( 0 | F ) where 0 stands for t h e p o i n t and V t h e d e p o s i t . It can be s h o w n that: 2
a (0\V)
= ±
2
J
z
6x\
y\x-y\dy
v
(5.19)
In standard analysis of variance t e r m i n o l o g y , o n e should call this function t h e within block variance. For n o n degenerate b l o c k s of size v t h e variance of t h e grade of blocks of size v, a (v \ V), within a larger v o l u m e V can be calculated from t h e follow ing relationship: 2
o (v\V) 2
=
a (0|y)-a (0iz;) 2
2
(5.20)
98
In short o n e says "the variance of blocks v". Coming back t o t h e analysis of variance t e r m i n o l o g y o n e should call it t h e b e t w e e n block variance. N o t e h o w e v e r that there is an infinite number of blocks v within V. Consequently this is a variance strictly speaking o n l y w h e n v and V are h o m o t h e t i c . When t h e y are n o t , w e shall define w h a t w e call t h e variance of a block by this relationship. N o w , considering three blocks, v v , v , o n e has: u
o (v \v ) x
= o {v \v ) 2
3
x
2
+
o (v \v ) 2
2
3
2
3
(5.21)
This had been n o t i c e d experimentally b y D.G. Krige and is k n o w n as Krige's relationship in t h e mining industry. We could also state it as: The total variance is equal t o t h e sum o f the within b l o c k variance and the b e t w e e n block variance. N o t e again that w e are dealing w i t h an infinite n u m b e r of terms in each sum of squares. 5 . 5 . 3 The covariance
of two
blocks
Covariances were introduced in S e c t i o n s 4 . 4 . 3 and 4 . 4 . 4 . Given t w o blocks, V and V , o n e of t h e m eventually being reduced t o a p o i n t — a sample — o n e easily s h o w s t h a t t h e covariance of t h e grades of the t w o blocks within a deposit D can be expressed as a function of t h e variogram, using t h e .following integrals: x
2
(5.22)
A program will be given later t o numerically estimate t h e s e x t u p l e integrals involved. One c o u l d say that this is in fact the o n l y program required in geostatistics and t h e presentation of the t h e o r y could s t o p here. All w h a t is left are mere application details. T h e y are of course essential for t h e practical use of geostatistics, b u t all major ideas are contained in t h e preceding few paragraphs. This is probably w h y geostatistics until recently never attracted the a t t e n t i o n of professional statisticians. 5.6 T H E N U G G E T E F F E C T C
5.6.1
0
Generalities
In m a n y instances, t h e variogram will n o t tend t o zero as h d o e s . This means that if o n e goes back t o sample the same site, t h e results will differ.
99
This m a y occur for a n u m b e r of reasons, such as p o o r analytical precision, p o o r sample preparation or more often highly erratic mineralization at l o w scale. This later p h e n o m e n o n , like t h e presence or absence of a gold nugget, is obviously responsible for t h e n a m e of that intercept at t h e origin. Usually statisticians simply call this t h e Chaotic C o m p o n e n t . It can be considered as t h e variance o f a totally random c o m p o n e n t superimposed o n t h e regionalized variable. As such, it is inversely proportional t o t h e v o l u m e of samples and thus t h e nugget effect of t h e variogram of 5 0 ' samples will be five times less than t h e o n e of 1 0 ' samples. It has n o effect o n covariances but will increase p o i n t variances b y C . One should consider it as an un certainty o n t h e value of t h e sample itself. 0
5.6.2 Theoretical
approach
A t h o r o u g h theoretical treatment of t h e nugget effect s h o w s that o n e c a n n o t keep a simple form like j(h) = C + y (h). This is d u e t o t h e fact that t h e d i s c o n t i n u i t y is n o t really occurring at the p o i n t level. What o n e has in fact is a very short range transition and t h e c o m p u t a t i o n s have t o be made taking i n t o a c c o u n t t w o variograms w h i c h represent variations at t w o different levels of observation, o n e being t h e chaotic c o m p o n e n t and t h e other the structured c o m p o n e n t . H e n c e t h e variance of any sample should be t h e sum of t w o terms c o m i n g from each o f t h e t w o variograms since all integrals are linear in 7 . One can, however, s h o w that in c o m m o n cases, w h e n t h e size o f samples is larger than t h e range of t h e small transition, t h e overall result is simplified, t h e additional term being simply as s h o w n before inversely proportional t o t h e size of t h e sample. 0
0
The nugget term for a block variance The rapid transition f u n c t i o n can be written y(h) = C ~ C(h) where C(h) is a positive f u n c t i o n w h i c h is zero for distances larger than t h e small range. It corresponds t o t h e chaotic c o m p o n e n t . T h e n t h e variance of a sample v inside a large o n e V is increased b y a term 0% w h i c h is, if w e write h— \x —y\\ 0
or:
0%
100
since:
Let us c o m p u t e o n e of the remaining integrals: C\x—y\dy = C(h)dh = A, A being a constant, Jv Jv since the v o l u m e o n w h i c h C(h) is > 0 is smaller than V. Now: \
f
[\
C\x-y\dy]
dx = ~
f
A dx =
-
Thus the variance is increased b y : a
-
= A
(^i)
(5
-
23)
Usually 1/V is zero compared t o 1/v and thus the additional term is simply A/v. In practice a nugget effect will have t o be taken i n t o a c c o u n t o n l y for samples, since mining b l o c k v o l u m e s are e n o r m o u s in comparison t o sample v o l u m e s and the constant AIV again b e c o m e s absolutely negligible. The nugget term for a covariance It simply d o e s n o t exist since t h e structured c o m p o n e n t and the chaotic c o m p o n e n t are i n d e p e n d e n t of each other and t h e chaotic c o m p o n e n t s are uncorrelated b y definition. 5 . 6 . 3 Remarks
on the origin of a nugget
effect
Splitting t h e total variances of samples i n t o a chaotic "nugget" c o m p o n e n t and a structured c o m p o n e n t means that o n e considers the observed regionalized variable as t h e sum of a structured random variable, plus a random c o m p o n e n t , uncorrelated w i t h the grade. This random c o m p o n e n t m a y have m a n y origins. It m a y be a real small-scale transition d u e t o the mineralization, it can be due t o p o o r sample preparation and analysis, or simply lack of care in sample collecting. This is s o m e t i m e s referred t o as a "human nugget e f f e c t " and w h e n identified it should be eliminated if possible. On Figs. 6 9 and 7 0 t w o cases of such a h u m a n nugget effect are s h o w n . The first e x a m p l e b y Serra ( 1 9 6 7 ) clearly s h o w s t h e structured c o m p o n e n t of the variogram, w h i c h remains unaffected through t i m e , while the nugget effect sharply decreases w h e n o n e o n l y considers recent carefully collected samples. The s e c o n d e x a m p l e is from a m o l y b d e n u m deposit where t h e mineralization occurs in half an inch veinlets. In this case t h e reduction of the sample is w h a t is essential; o n o n e h o l e the preparation was careful
101
F i g . 6 9 . " H u m a n " n u g g e t e f f e c t i l l u s t r a t e d . T w o v a r i o g r a m s in t h e s a m e a r e a o f a n i r o n m i n e c o m p u t e d o n o l d ( A ) and n e w ( B ) s a m p l e s . A f t e r Serra ( 1 9 6 7 ) .
#(h)
0-05 0.04 0.03 0.02 0.01
10
20
30
40
50
Fig. 7 0 . " H u m a n " nugget effect
60 generated by poor sample preparation (A)
disappears
w i t h c o r r e c t s a m p l e r e d u c t i o n ( B ) . T h e t w o v a r i o g r a m s are f r o m a m o l y b d e n u m d e p o s i t .
but following a w r o n g s c h e m e , thus giving perfectly random results, while for the s e c o n d , Gy's sampling principles were used, and t h e structure is clearly revealed. More practical considerations a b o u t t h e nugget effect will be given in t h e chapter a b o u t variogram interpretation. We will also spend s o m e time later t o discuss sampling errors.
102 5.7 T H E O R E T I C A L M O D E L S O F ISOTROPIC V A R I O G R A M S
In practice y(h) is n o t k n o w n and has t o be estimated from t h e available samples. One thus obtains for various values of h a series of values 7 * ( / * ) . One t h e n has t o fit an e q u a t i o n t o these experimental values. So far it has b e e n possible t o estimate m o r e than a hundred deposits w i t h o n l y t w o m o d e l s and t h e t e n d e n c y is n o w t o use o n l y o n e : the spherical m o d e l . This m o d e l is m e n t i o n e d in Matern ( 1 9 6 0 , p. 3 0 ) . T o the author's k n o w l e d g e , however, it was n o t used in practical applications b y biometricians. T h e y m o s t l y tend t o use t h e e x p o n e n t i a l m o d e l although it was recognized that it c a n n o t claim any divine right (Whittle, 1 9 5 4 ) . In mining problems it was f o u n d that t h e e x p o n e n t i a l m o d e l has t o o gentle a growth t o satisfactorily represent grade variations. Without claiming more divine rights for t h e spherical m o d e l , let us present it and t h e m a n y e x a m p l e s where it can be used. 5 . 7 . 1 The spherical
model
Its shape is s h o w n o n Fig. 7 1 ; its equation is:
(5.24)
Fig. 7 1 . T h e spherical m o d e l variogram.
103
500
1000
1500
2000
2500
h(m)
Fig. 7 2 . H o r i z o n t a l variogram o f D . D . H . i n t e r s e c t i o n in iron ore in Lorraine After S e r r a ( 1 9 6 7 ) .
(France).
F i g . 7 3 . H o r i z o n t a l v a r i o g r a m o f v e r t i c a l p a c k s a k D . D . H . in a s t r a t a - b o u n d C a n a d i a n u r a n i u m d e p o s i t , s h o w i n g a r a n g e o f 6 0 m ( s a m p l i n g i n t e r v a l is 1 5 m ) . A f t e r B l a i s a n d Carlier ( 1 9 6 8 ) .
a is called t h e range, C is t h e nugget effect and C + C t h e sill. This m o d e l occurs naturally in all t h e deposits where grades b e c o m e i n d e p e n d e n t of each other o n c e a given distance, a, is reached. It is t h e c o m m o n rule in m o s t sedimentary deposits and also in porphyry copper deposits. T h e variograms s h o w n in Figs. 72—80 have b e e n satisfactorily described by such a m o d e l . Deposits as different as iron ore deposits in France, porphyry copper in Chile and Canada, lead—zinc in Ireland and Canada, have been f o u n d t o have their grade variations adequately represented b y this m o d e l . S o were gold deposits in S o u t h Africa, bauxite in France, lateritic nickel in N e w Caledonia, uranium in Canada and phosphates in Africa. This m o d e l is also k n o w n as t h e Matheron m o d e l and is said t o describe transition p h e n o m e n a 0
0
100
300
500
700
h (m)
Fig. 7 4 . A 1 0 h o r i z o n t a l variograms in the bauxite deposit of Mazauges Aval, showing isotropy. After Marechal and Roullier ( 1 9 7 0 ) . 2
3
100
300
500
2
150 4
125
H
25
H
100
200
300
400
h(m)
Fig. 7 5 . S i 0 horizontal variograms in the bauxite deposit of Mazauges Aval, show ing small nugget effect and isotropy. After Marechal and Roullier (1970).
Y(h)
100
700
h(ft)
Fig. 7 6 . Vertical variogram o f m a g n e t i t e grade in Canadian iron ore d e p o s i t .
105 Y(h)
100
300
500
700
900
1100
h(ft)
Fig. 7 7 . V a r i o g r a m s a l o n g strike and d o w n d i p for g o l d grade in a Canadian quartz vein deposit showing over 6 0 0 of continuity.
2Y(h)f
300 H
50
h(m)
Fig. 7 8 . H o r i z o n t a l variogram in an African g o l d m i n e s h o w i n g very short z o n e o f influ e n c e . A f t e r C o r t e z e t al. ( 1 9 7 4 ) .
V(h) 1200
300 H 150
500
1000
1400
h(m)
Fig. 7 9 . H o r i z o n t a l variogram o f t h e Silica a c c u m u l a t i o n in garnieritic n i c k e l of N e w Caledonia. After Journel ( 1 9 7 4 ) .
deposits
106 v(h)
r
r1
—Ti
150
1 —
1000
500
Fig. 8 0 . Variogram
of
the
nickel accumulation
—
i
—
1400
h(m)
in garnieritic nickel d e p o s i t s o f
New
Caledonia. After Journel ( 1 9 7 4 ) .
as it is t h e o n e w h i c h occurs w h e n o n e has structures i n d e p e n d e n t of each other but within each, t h e grades are highly correlated. A digression on the origin of the 'spherical' name and genetical implication of variograms A process having a spherical variogram can be generated as f o l l o w s : first take a random (Poisson) distribution of masses in space. T h e n assign t o each p o i n t of space a grade given b y the n u m b e r of masses within a radius a of the point (in other words, take t h e number of masses within a "sphere" . . .). Also n o t e that t h e v o l u m e of t h e intersection of t w o spheres of diameter a, v o l u m e V, t h e centers of w h i c h are h apart is:
if h < a and 0 otherwise. It is interesting t o be able t o construct (for simu lation purposes) a process w i t h a given variogram; w e will n o t e x p a n d o n t h e subject however, e x c e p t t o warn a t e m p t e d reader n o t t o draw totally h y p o t h e t i c a l genetic conclusions from mathematical games. The only accept able process w o u l d be the reverse: given an alleged mineralizing process t o theoretically c o m p u t e its variogram and later experimentally c h e c k it. 5.7.2 TheDe
Wijsian model
(after Prof. H.J. de Wijs)
This m o d e l is e n c o u n t e r e d in cases where there is n o such thing as a range of d e p e n d e n c e . y(h) keeps increasing as h increases. In m o s t of t h e hydrothermal deposits, t h e variogram plots as a straight line w h e n y(h) is p l o t t e d against ln(h). The equation is thus: y(h)
= A ln(fc) + B
(5.25)
Such variograms can be seen in Figs. 8 1 — 8 3 . The vogue of this m o d e l was certainly due t o t h e " g o o d " properties of the logarithmic function w h i c h
107
I
i
i
I
2
3
i
4
i
5
i
6
7
i
i
8
9
i
i
•
10
h
F i g . 8 1 . D e Wijsian v a r i o g r a m f o u n d f r o m c h a n n e l s a m p l e s a l o n g u p p e r d r i f t i n a Z a m b i a sulfide deposit.
I
i
i
I
2
3
i
4
i
5
6
i
7
i
8
i
9
i
i
•
10
h
F i g . 8 2 . D e Wijsian v a r i o g r a m f o u n d f r o m c h a n n e l s a m p l e s a l o n g l o w e r drift i n a Z a m b i a sulfide deposit. X(h)4
I
2
3
4 5 6 7 8 9 1 0
15
20
h
F i g . 8 3 . G e n e r a l v e r t i c a l v a r i o g r a m ( s o l i d l i n e ) o f d i a m o n d drill h o l e c o r e s a m p l e s a n d horizontal variogram ( b r o k e n line) o f c h a n n e l s a m p l e s o f t h e mineralization in the Mounana uranium deposit, Gabonese Republic. After Matheron ( 1 9 6 2 ) .
108
Y(h)f
0.8
1.6
2.4
3.2
4.0
4.8
h
Fig. 8 4 . Several e x p e r i m e n t a l variograms correctly represented b y a linear m o d e l .
a m o n g other things is harmonic. The first v o l u m e s of Matheron ( 1 9 6 3 ) and Carlier ( 1 9 6 4 ) were entirely d e v o t e d t o it. T h e y have lost a l o t of their interest w i t h t h e appearance of c o m p u t e r s , b u t should be remembered specially in case where a c o m p u t e r is n o t available and quick results are n e e d e d (see Exercise 4 . 6 . 2 ) . 5 . 7 . 3 Other
models
5 . 7 . 3 . 1 The linear model The m o s t simple m o d e l w h i c h o n e tries t o use as m u c h as possible for obvious reasons is just a linear o n e . Quite often it is possible at least within a limited n e i g h b o u r h o o d t o use a linear equation: y(h)
=
Ah+B
Several e x a m p l e s are s h o w n in Fig. 8 4 w h i c h is a p l o t o f Fig. 8 9 for distances smaller than 5 m. N o t e that this variogram is the theoretical variogram of a process called t h e Brownian m o t i o n . (Feller, 1 9 6 6 , p. 9 7 ) . 5 . 7 . 3 . 2 The ah model We have already seen t h e linear m o d e l w h e n p l o t t e d o n a log scale, w h i c h is the D e Wijsian m o d e l ; w e can also in m a n y cases m a k e the variogram straight by plotting it on a log-log scale. Then the equation is y(h) = ah since ln y(h) = X ln h + constant. This is often m e t in elevation variograms (Figs. 8 5 and 8 6 ) or in the s t u d y of mill feed variability (Fig. 8 7 ) . K
x
5 . 7 . 3 . 3 The exponential model The e x p o n e n t i a l m o d e l is n o t e n c o u n t e r e d t o o often in mining practice
109
¥(h)
8.0
6.0
4.0
2.0
0.0
Fig. 8 5 . Plot o n arithmetic scale o f the variogram of the elevation o f a hanging wall in a gold d e p o s i t .
Y(h)
Fig. 8 6 . S a m e variogram plotted o n a log-log scale.
as
in
Fig.
85,
4 6 10 20 4060 100 200 400 h( minutes) Fig. 8 7 . P l o t o n a log-log scale o f t h e variogram o f mill feed grade for the mill o f t h e SaintS e b a s t i e n M i n e ( F r a n c e ) . A f t e r P. G y ( 1 9 6 7 ) .
110
since its infinite range is associated t o a t o o c o n t i n u o u s process. It is h o w ever very practical. The equation is: y(h)
=C {
+ C
0
5.1.3.j/The
exp
hole effect
-
AM
a
model
The "hole e f f e c t " m o d e l s h o w n in Fig. 8 8 has for an e q u a t i o n : tW
=
/sin ah \ c^i—
It can be used t o represent fairly c o n t i n u o u s processes (the tangent at t h e origin is horizontal and it s h o w s a periodic behaviour w h i c h is o f t e n e n c o u n tered w h e n there exists a succession of rich and p o o r z o n e s , for instance).
F i g . 88. P l o t o f t h e t h e o r e t i c a l v a r i o g r a m 7 ( h ) = C [ l — s'm(ah)/ah] and a periodic hole effect.
showing slow growth
The e x a m p l e in Fig. 8 9 is from an Australian uranium m i n e where rich and p o o r layers alternate. In Fig. 9 0 w e see an e x a m p l e from Serra ( 1 9 7 0 ) w h i c h corresponds t o t h e e x i s t e n c e of waste concretions in t h e Lorraine iron ore. Figure 9 1 s h o w s a more elaborated m o d e l b y Marechal and D i d y k ( 1 9 7 1 ) w h i c h was d e v e l o p e d for t h e El Salvador copper m i n e in Chile.
Ill
Fig. 8 9 . Vertical variograms c o m p u t e d for increasing s a m p l e l e n g t h s in an uranium deposit showing a marked hole effect.
0
0.5
1.0
1.5
Australian
2.0 2.5 (METERS)
Fig. 9 0 . H o r i z o n t a l variogram o f iron, c a l c i u m and silica grade s h o w i n g s h o r t scale, transition p h e n o m e n o n and h o l e effect ( e x i s t e n c e of c o n c r e t i o n s ) in t h e Lorraine iron ore deposit. After Serra ( 1 9 7 0 ) .
The equations w h i c h were used are a c o m b i n a t i o n o f e x p o n e n t i a l and sine model: # ) =
P
( l - e - )
+
d
-
p
)
(
^
)
112 Y(h>r 0.54
OA A 3
9
15
21
27
33
39
45
51
57
h(m)
F i g . 9 1 . V a r i o g r a m o f c o p p e r g r a d e at t h e El S a l v a d o r c o p p e r d e p o s i t ( C h i l e ) . A f t e r Marechal and D i d y k ( 1 9 7 1 ) .
5 . 7 . 4 Admissible
functions
for a
variogram
N o t all functions are suitable t o represent a variogram. One has t o bear in mind that a variogram represents a physical process n a m e l y t h e grade at o n e p o i n t Z{x). In the estimation procedure w e will have t o use linear combina t i o n s o f k n o w n values. We m u s t m a k e sure that t h e variance of these c o m binations exists and that it is positive. It has b e e n s h o w n that a linear c o m bination Z* has a finite variance if and o n l y if t h e sum of t h e weights a is equal t o 0. In this case: t
VAR(Z*) = -I
£
a.ajjix.-xj)
1=1 7 = 1
This defines t h e c o n d i t i o n w h i c h a 7 - f u n c t i o n should satisfy in order t o be admissible as a variogram. This c o n d i t i o n o n — 7 is that it m u s t be condition ally positive definite: for any set o f n p o i n t s x x , . . x and any set o f weights a a , . . ., a such that 2 a = 0 o n e should have: u
u
n
-
2
n
2
n
f
n
Z Z aityyixi /=i y=i
~xj)
>
0
This is a very strict c o n d i t i o n and c o n s e q u e n t l y o n e c a n n o t simply fit which ever p o l y n o m i a l m o d e l o n e w o u l d like. One should c h e c k for t h e positivity o f the m o d e l or use a m o d e l w h i c h is k n o w n t o be satisfactory like t h e ones w e have just m e n t i o n e d . Many other equations have b e e n proposed b y Matheron ( 1 9 6 8 ) . He has also c o m p u t e d all t h e p o l y n o m i a l s o f order less than 6 w h i c h are admissible.
113
However, for practical purposes w e will see that w e should be h a p p y w i t h the few m o d e l s listed above and w i t h c o m b i n a t i o n s of t h e m , since a linear c o m b i n a t i o n of t h e previous equations is an admissible e q u a t i o n . We will see m a n y e x a m p l e s in t h e chapter a b o u t m o d e l fitting. For readers interested in t h e theoretical p r o b l e m of t h e class o f admissible isotropic covariance functions in R , w e refer t h e m t o t h e very detailed treatment of this p r o b l e m w h i c h was given b y Matheron in t h e first part of his thesis ( 1 9 6 5 ) , or t o Matern ( 1 9 6 0 ) w h o sums up t h e w o r k of B o c h n e r ( 1 9 3 2 ) , S c h o n b e r g ( 1 9 3 8 ) , Hartman and Wintner ( 1 9 4 0 ) , Lord ( 1 9 5 4 ) and Yaglom ( 1 9 5 7 ) and gives an extensive list of admissible functions. n
5.8 HOW C L O S E IS T H E I N T R I N S I C H Y P O T H E S I S T O R E A L I T Y ?
5 . 8 . 1 Empirical
discussion
of the
problem
It might be important at this stage t o spend m o r e t i m e discussing w h a t h y p o t h e s e s and m o d e l s are, and also introducing t h e c o n c e p t of robustness. First of all, a m o d e l is n o t reality; it is t h e expression o f our inability t o fully understand and unravel t h e deterministic processes w h i c h generated the present state of t h e p h e n o m e n o n (in our case, t h e ore) w h i c h is observed. Geostatistics is a probabilistic m o d e l ; t h e law o f gravity is a deterministic m o d e l ; w e could use statistics instead t o find t h e p o s i t i o n at a given t i m e of a falling b o d y . U n f o r t u n a t e l y , t o c o m p u t e t h e grade of a b l o c k such a simple equation as e = \gt c a n n o t be used. What is d o n e instead is t o consider that t h e grade of our b l o c k is o n e realization of a random variable, just like t h e result o f a die t h r o w is considered as o n e realization o f a uniform distri b u t i o n b e t w e e n 1 and 6. In fact t h e o u t c o m e of t h e die t h r o w is purely determined b y t h e original p o s i t i o n of t h e die, t h e w a y it is t h r o w n and possibly o t h e r parameters w h i c h are simply t o o n u m e r o u s t o be c o u n t e d and measured. N o w t h e probabilistic m o d e l w h i c h is used t o describe the out c o m e of t h e t h r o w o f a die is very simple and it is easy t o c h e c k w h e t h e r t h e die is biased or n o t , in other w o r d s , if t h e m o d e l is correct for this die. When w e are dealing w i t h ore grade, w e are dealing w i t h a random variable w i t h an infinite n u m b e r o f d i m e n s i o n s (grade at each p o i n t ) and t h e equivalent of k n o w i n g t h e probability of t h e o u t c o m e of a 1 for a die n o w b e c o m e s k n o w i n g for any n and any set of p o i n t s x (i = 1, . . . , n) t h e probability that simultaneously: 2
t
{Z(x ) x
< Z
u
Z(x ) 2
<
Z , . . . , Z(x ) 2
n
<
Z) n
This is clearly impossible, at least at present. The intrinsic h y p o t h e s i s intro duces s o m e suppositions, such as t h e stationarity of increments, w h i c h can t h e n lead t o t h e definition and c o m p u t a t i o n o f a variogram and t h e n t o any variance c o m p u t a t i o n s . One of t h e first criticisms o f t h e intrinsic h y p o t h e s i s was that nature is
114
n o t even w e a k l y stationary. This is true, but w h a t should be considered, rather than rejecting t h e m o d e l , is t o try and estimate h o w w r o n g t h e m e t h o d is w h e n applied in non-stationary c o n d i t i o n s . This is the c o n c e p t of robust ness. Since these early criticisms, m a n y d e v e l o p m e n t s have occurred and the universal kriging h y p o t h e s i s of Matheron ( 1 9 7 1 ) has already b e e n introduced in order t o deal with non-stationary p h e n o m e n a . More recently, the t h e o r y of t h e random intrinsic functions of order k has firmly established t h e problem o f t h e estimation of a random f u n c t i o n under any h y p o t h e s i s (Matheron, 1 9 7 3 ) . Thus t h e theoretical objections have all been answered, b u t this d o e s n o t m e a n that these sophisticated t o o l s have t o be used as s o o n as non-stationarity is suspected. In fact, m a n y experiments have been m a d e t o s h o w that using a simple variogram is very little affected b y very strong departure from h y p o t h e t i c a l c o n d i t i o n s . A n a c c o u n t of this is given in David and Dagbert ( 1 9 7 6 ) . Besides this, the n o t i o n of stationarity can o n l y be con sidered w i t h respect t o a given scale at w h i c h o n e wishes t o work, and m o s t of the time there exists a n e i g h b o u r h o o d for an area of interest where stationarity can be assumed. With a given set of samples, it is always possible t o c o m p u t e a variogram; this d o e s n o t mean that this variogram is t o be used if it has b e e n c o m p u t e d o n a n o n - h o m o g e n e o u s region; it will exhibit features w h i c h should send us back t o t h e geological definition of h o m o g e n e o u s z o n e s in an ore deposit. One sees at this p o i n t that g e o l o g y and numerical c o m p u t a t i o n s should c o m p l e m e n t each other in order t o p o i n t o u t the different z o n e s w h i c h o n e should consider in estimating t h e deposit. 5 . 8 . 2 Theoretical
discussion
of the
problem
A theoretical discussion of these problems is definitely important but m u c h b e y o n d t h e s c o p e of this manual. Interested readers will find all necessary justifications in Matheron ( 1 9 7 1 ) under t h e heading "Statistical Inference" and in his m o s t recent paper ( 1 9 7 6 a ) a b o u t t h e c h o i c e of m o d e l s in geostatistics. A m o n g the points that are stressed, considering t h e relevance of a m o d e l , four c o n d i t i o n s w h i c h a m o d e l m u s t satisfy are p o i n t e d o u t : Statistical inference m u s t be possible, i.e., it should be possible t o estimate t h e parameters o f t h e m o d e l . Then t h e m o d e l m u s t be operational, i.e., it m u s t give an answer t o the question asked. This is obvious but should rather be considered in t h e reverse; depending o n the questions o n e wants t o answer, o n e m a y use different m o d e l s for t h e same deposit. The m o d e l has t o be compatible w i t h t h e data, and finally t h e predictions o f t h e m o d e l must be checked b y experience. This will be an a-posteriori proof. N o w a m o n g all m o d e l s m e e t i n g these c o n d i t i o n s , h o w is t h e c h o i c e t o be made? The answer is short and clear: Take t h e simplest.
5
THEORETICAL BASIS OF THE APPROACH: THE THEORY OF REGIONALIZED VARIABLES
S Y N O P S I S — In this chapter, the theory of regionalized variables, as formal ized by Matheron (1965) will be presented, together with an attempt to show where it fits in the standard statistical literature. Bearing in mind that what is important to us is applications, the theoretical references will be limited and only what is needed for our purpose will be discussed. Thus, after defining a regionalized variable, the different hypotheses which one may have to use in order to make statistical inference possible, will be introduced. Since the estimators which will be later considered will be linear, properties of linear combinations of regionalized variables will be examined, and as mining samples are not point samples, we will see how to take care of non-point samples. Then the theoretical expression of estima tion variances and block variances as a function of the variogram will be given. Examples of the variograms most commonly used, will be shown and the general class of admissible functions for a variogram will be briefly reviewed. The notion of chaotic component or nugget effect will be intro duced and discussed. 5.1
FOREWORD
Most of the material in this chapter follows the presentation of Matheron ( 1 9 6 5 ) . In fact, this was also presented b y Matern ( 1 9 4 7 , 1 9 6 0 ) concerning forestry problems, and o n e reference t o t h e variogram can be traced back t o 1 9 2 6 in a paper b y Langsaetter, again a b o u t forest surveys. The sampling properties of regionalized variables, or t h e study of topographic variation as it has often been referred t o , have been discussed by Mahalanobis ( 1 9 4 4 ) , Nair ( 1 9 4 4 ) , Quenouille ( 1 9 4 9 ) , J o w e t t ( 1 9 5 5 ) , Williams ( 1 9 5 2 , 1 9 5 6 ) , Zubrycki ( 1 9 5 7 ) , Whittle ( 1 9 5 4 ) , Finney ( 1 9 5 0 ) . However, n o n e of these papers deal w i t h geological observations, let alone mining. A l t h o u g h the theory m a y n o t be n e w t o professional statisticians, o n e should recognize that very specific problems occur in the mining industry and this leads t o specific theoretical d e v e l o p m e n t s ajid different t e r m i n o l o g y . A m o n g other particularities, grade values are m u c h less c o n t i n u o u s than m o s t of biological or ecological data considered b y biometricians. Beside this, mines are truly three-dimensional b o d i e s and samples are n o t " p o i n t s " , while w h a t appeared in earlier literature m o s t l y concerns two-dimensional p o i n t processes. It is n o t the intention of this t h e o r y t o change the existing statistical t e r m i n o l o g y , but since t o mining professionals b o t h languages are equally u n k n o w n , it
92
d o e s n o t appear unreasonable t o try and d e v e l o p a presentation o f the t h e o r y w h i c h is specially tailored t o mining problems. D o i n g this has t h e usual drawback of disappointing b o t h statisticians and mining p e o p l e ; w e feel, however, that it is a possible avenue. 5.2 DEFINITION OF A R E G I O N A L I Z E D V A R I A B L E
Let x be a p o i n t in space R and z(x) t h e value of t h e function w e are interested in at p o i n t x. S u c h a f u n c t i o n will be called a regionalized variable. T h e grade of equal-size samples is such a regionalized variable. S o is t h e thickness of a d e p o s i t at o n e p o i n t or t h e specific gravity of t h e ore. Such a function is usually highly variable and n o n - c o n t i n u o u s and c a n n o t be studied directly. The structure w h i c h m a y appear in its variations will be studied b y examining its increments. The basic idea of the t h e o r y is t o consider such a function z(x) where x is a p o i n t or a vector o f R as o n e realization of a random f u n c t i o n ( R F ) Z(x). Thus w e are turning a perfectly well defined unique numerical value i n t o a realization of random process. We have o n l y o n e realization of that random function and t h e problem w e are faced w i t h consists in finding the characteristics of t h e R F Z(x) in order t o make the estimation of u n k n o w n points possible. 3
n
This decision t o l o o k at t h e values of t h e process w e are interested in as a realization of a random function is just an epistemological decision w h i c h is taken w h e n recognizing that deterministic m e t h o d s w o u l d n o t help us; it is n o t possible t o consider any e x p e r i m e n t which w o u l d c o n c l u d e that a deposit is n o t a realization of a random function (Matheron, 1 9 7 6 a ) . Since this is a very permissive m o d e l , t h e counterpart is o b v i o u s , it d o e s n o t h e l p m u c h ! H y p o t h e s e s have t o be introduced a b o u t t h e t y p e of random func t i o n which will be considered. 5.3 T H R E E P L A U S I B L E H Y P O T H E S E S
Going from t h e m o s t restrictive t o t h e more general, o n e has seen in t h e geostatistical literature three sets o f possible assumptions. 5.3.1 The weak-stationarity
assumption
(second-order stationarity)
This assumption,seldom found in natural p h e n o m e n a , consists in t w o conditions: — the e x p e c t e d value of t h e regionalized variable Z(x) is t h e same all over t h e field of interest; — t h e spatial covariance of t h e regionalized variable Z(x) is t h e same all over t h e field of interest.
93
Hence t h e e x p e c t e d value is: E[Z(x)]
= m
(5.1)
and t h e covariance is: E{[Z(x)-m][Z(x
+ h)-m]}
where h is a vector in R . random f u n c t i o n Z(x) is: n
V A R [Z(x)]
=
= K(x,x
+ h) = K(h)
(5.2)
/i is t h e length of vector /*. The variance of t h e
E{[Z(x)-m] }
= K(0)
2
(5.3)
h e n c e t h e process has a covariance o n l y if V A R [ £ ( # ) ] is finite. 5 . 3 . 2 The intrinsic
assumption
In m a n y deposits, as s h o w n b y Krige ( 1 9 5 1 ) , such a thing as a finite variance d o e s n o t exist; if o n e considers t h e variations of grade rather than t h e grade itself, it has a finite variance. H e n c e , considering t h e increments of t h e f u n c t i o n , Z(x) — Z(x + /*), rather than t h e function itself, o n e can make t h e following assumptions: E[Z(x
+h)-Z(x)]
VAR[Z(x
= 0
+ h)-Z(x)]
(5.4)
= 2y(h)
(5.5)
This last definition is t h e definition of t h e variogram. One can see that this definition can be rewritten: VAR[Z(x
+ h)-Z(x)]
= E[Z(x
+
h)-Z(x)-E{Z(x+h)-Z(x)}]
2
and since: E[Z(x
4- h) - Z(x)]
VAR[Z(x
=
+ h)-Z(x)]
0 = E[Z(x
+ h)-Z(x)]
(5.6)
2
This is the intuitive definition of t h e variogram w h i c h was given in Chapter 4. Also if a R F is second-order stationary it is also intrinsic and its variogram is: V A R [ Z ( x + h) - Z(x)]
= V A R [ Z ( x + h)] + VAR[Z(jc)] -2COV[Z(x
+
h),Z(x)]
or using t h e previous definition and dividing b y 2 : 7(A) = K(0)~K(h)
(5.7)
Remark on notation. As m a n y p e o p l e have already e x p e r i e n c e d , it is n o t l o n g in geostatistics before o n e encounters n o t a t i o n problems. We have just written our first formulae. Here is t h e first p r o b l e m : One m a y w o n d e r w h y
94
t h e factor 2 in t h e definition of 2y(h). It is because t h e quantity w h i c h is m o s t l y used is y(h) rather than 2y(h). C o n s e q u e n t l y y(h) has been called t h e semi-variogram or half-variogram. Experience s h o w s that almost every b o d y calls y(h) t h e variogram. One has even seen papers w h e r e 2y(h) was defined as the d o u b l e variogram. We could n o t favour that word b u t suggest that for the sake of simplicity, o n e should take as definition of the vario gram: V A R I O G R A M = y(h)
= \ V A R [ Z ( * + h) - Z(x)]
(5.8)
Point variogram. The previously defined variogram applies t o variable Z(x) defined at p o i n t x. One will often call such a variable a p o i n t variable and its variogram a p o i n t variogram as o p p o s e d t o n o n - p o i n t variables (see Section 5.4.2).
5 . 3 . 3 The hypotheses
of universal
kriging
A third kind of h y p o t h e s i s , less restrictive than t h e previous o n e s , assumes that t h e second m o m e n t of t h e R F or its increments has s o m e properties o f stationarity within a vicinity of restricted size and t h a t t h e e x p e c t a t i o n w h i c h is n o longer stationary varies in a regular manner in such a vicinity. If x and y = x + h are taken in the same vicinity: E[Z(x)]
= m(x)
or
E[Z(x)
-Z(y)]
= m{x)-m(y)
(5.9)
with: k
m{x)
=
X d^ix)
(5.10)
1 = 0
the fi(x) being k 4- 1 i n d e p e n d e n t functions and t h e d numerical coefficients. COV[Z(x),Z(x
+ h)]
= K(h)
or
V A R [ Z ( * + h) - Z(x)]
=
z
unknown
2y(h)
In this third kind of h y p o t h e s i s , n o t o n l y t h e covariance or variogram f u n c t i o n has t o be defined from the experimental values but also t h e size of t h e vicinity where the h y p o t h e s e s remain valid, the nature of t h e functions fx(x) as well as their n u m b e r and t h e values of t h e coefficients d which are a function of t h e position of t h e vicinity in the field. The progression in t h e h y p o t h e s e s can be seen. One m a y w o n d e r w h y s t o p here. In fact the last d e v e l o p m e n t s of the t h e o r y (Matheron, 1 9 7 3 ) n o w con sider generalized increments and stationarity properties for these generalized increments. It can also be seen that each h y p o t h e s i s covers t h e previous ones. N o w the t h e o r y of generalized random intrinsic functions of order k covers all t h e previous cases. It is far b e y o n d t h e scope of this v o l u m e and will n o t be discussed. A review of it is given in Delfiner ( 1 9 7 6 ) . u
95 5.4 L I N E A R C O M B I N A T I O N S A N D A V E R A G E
5.4.1 Statistical
properties
VALUES
of the linear combinations
of random
variables
Linear c o m b i n a t i o n s of random variables are random variables (see for instance, Feller, 1 9 5 7 ) . The c o m b i n a t i o n s m a y be defined from c o m p o n e n t s of a R F o n a finite set o f points S[: n
Z*
= E
di Z(x ),
jc,-eS =
t
[x . . x
.x ] n
i =l
or an infinite set o f points which is a d o m a i n V of v o l u m e V. Z*
= 1/V
f
Z(x) dx
(5.11)
Y.
The variance o f Z* can be expressed for each of t h e previous sets o f hypotheses. If Z(x) is a second-order stationary RF, whatever t h e a which is expressed from t h e covariance of Z(x):
i9
VAR[Z*]
= t
Z* has a variance
t<*i
i=i
(5.12)
y=i
or o n an infinite set: V A R [ Z * ] = 1/V
f dx f K{x-y)dy Jy_ JV
2
(5.13)
If Z(x) is an intrinsic RF, it has b e e n s h o w n (Matheron, 1 9 7 1 , p . 5 6 ) that Z* has a finite a-priori variance if and o n l y if t h e c o n d i t i o n 2 " a = 0 is realized. In this case: = 1
VAR[Z*]
= - £
£
i = l
Z* = 1/V f Z(x) v
t o dxjV
atajjixi-xj)
(5.14)
j = l
dx is a linear c o m b i n a t i o n where t h e weights are all equal
and since 1/V f
v
5.4.2 Non-point
x
variables:
dx = 1, t h e n Z* d o e s n o t have a finite variance. smoothing
When w e are concerned with t h e average grade of a v o l u m e V centered at point x, w e can write the grade of that b l o c k as in the previous paragraph: Z[V(x)]
= Z
v
Z(x)dx
V being the v o l u m e o f set V, as m e n t i o n e d before.
96
Z[V(.x:)] is s o m e t i m e s called a s m o o t h e d variable and it might be interest ing t o express its variogram as a function o f t h e variogram of the p o i n t variable Z(x). Let k(x) designate a function w h i c h is equal t o o n e if x is included in V and zero if n o t . T h e n w e call t h e geometric covariogram of V t h e quantity: (5.15) It can easily be s h o w n that t h e variogram of variables of size V can be written: 7 v ( * ) = f 7(A +x)P(x)dx-$
(5.16)
P(x)y(x)dx
Note. The geometric covariogram has an interesting geometrical expres sion. Let us consider a set V and a vector h. N o w let us define V as the set obtained b y translating V b y a quantity h. N o w the intersection of V and V is such t h a t j c e V a n d * + heV. H e n c e P ( h ) = M e a s ( V n V )/V . More simply t h e geometric covariogram of a piece of core of length / will be \ l — h |// . h
h
2
h
2
5.5 T H E O R E T I C A L E X P R E S S I O N O F V A R I A N C E S
5 . 5 . 1 The extension
variance
Let V be a v o l u m e t h e grade of w h i c h Z(V) is u n k n o w n and t o be esti m a t e d b y the k n o w n grade Z(v) of v o l u m e v. Then Z(V) and Z(v) can be written as stochastic integrals, t h e variances and covariance of w h i c h can be s h o w n t o be functions o f t h e variogram:
The variance o f t h e error of estimation Z(V) — Z(v) w h i c h is m a d e w h e n assigning t o V t h e grade o f v exists since this is a linear c o m b i n a t i o n where t h e sum of weights is zero and it is s h o w n t o be (Matheron, 1 9 6 5 ) .
VAR[Z(V)-Z(*)]
=
-r
T
Ax
7 | x — y | dy
(5.17)
97
This relationship is t h e theoretical k e y t o any estimation problem since it allows a numerical characterization of t h e magnitude of t h e error involved in an estimation procedure as s o o n as y(h) is k n o w n . N o t e that since x and y are points in R , t h e integrals are s e x t u p l e o n e s . This m a y involve s o m e c o m p u t a t i o n a l problems. T h e y will be discussed later. This c o n c e p t and formulation m a y s e e m fairly simple and straightforward, but apparently this n o t i o n had never b e e n clearly defined until Matheron did it. 3
Estimation by a series of samples Usually, o n e wishes t o give t o a block V n o t t h e grade of a b l o c k y_ b u t rather t h e average grade of a series of punctual values, say: 1 Z*
N
I
= -
Z(x ) t
In this case t h e previous formulation is reduced t o : °N
= J^Z
l(x -x)dx~^\^
dx\^j(x-y)dy
i
lyixi-xj)
(5.18)
This variance is usually called an e s t i m a t i o n variance. The programming of its c o m p u t a t i o n will be seen in a later chapter, this being an essential t o o l for t h e appraisal of t h e quality of a sampling pattern. 5 . 5 . 2 The variance
of a block and the Krige's
relationship
We are n o w interested in t h e variance o f t h e grade of b l o c k s of a given size in t h e deposit. For punctual " b l o c k s " w e will use t h e q u a n t i t y a ( 0 | F ) where 0 stands for t h e p o i n t and V t h e d e p o s i t . It can be s h o w n that: 2
a (0\V)
= ±
2
J
z
6x\
y\x-y\dy
v
(5.19)
In standard analysis of variance t e r m i n o l o g y , o n e should call this function t h e within block variance. For n o n degenerate b l o c k s of size v t h e variance of t h e grade of blocks of size v, a (v \ V), within a larger v o l u m e V can be calculated from t h e follow ing relationship: 2
o (v\V) 2
=
a (0|y)-a (0iz;) 2
2
(5.20)
98
In short o n e says "the variance of blocks v". Coming back t o t h e analysis of variance t e r m i n o l o g y o n e should call it t h e b e t w e e n block variance. N o t e h o w e v e r that there is an infinite number of blocks v within V. Consequently this is a variance strictly speaking o n l y w h e n v and V are h o m o t h e t i c . When t h e y are n o t , w e shall define w h a t w e call t h e variance of a block by this relationship. N o w , considering three blocks, v v , v , o n e has: u
o (v \v ) x
= o {v \v ) 2
3
x
2
+
o (v \v ) 2
2
3
2
3
(5.21)
This had been n o t i c e d experimentally b y D.G. Krige and is k n o w n as Krige's relationship in t h e mining industry. We could also state it as: The total variance is equal t o t h e sum o f the within b l o c k variance and the b e t w e e n block variance. N o t e again that w e are dealing w i t h an infinite n u m b e r of terms in each sum of squares. 5 . 5 . 3 The covariance
of two
blocks
Covariances were introduced in S e c t i o n s 4 . 4 . 3 and 4 . 4 . 4 . Given t w o blocks, V and V , o n e of t h e m eventually being reduced t o a p o i n t — a sample — o n e easily s h o w s t h a t t h e covariance of t h e grades of the t w o blocks within a deposit D can be expressed as a function of t h e variogram, using t h e .following integrals: x
2
(5.22)
A program will be given later t o numerically estimate t h e s e x t u p l e integrals involved. One c o u l d say that this is in fact the o n l y program required in geostatistics and t h e presentation of the t h e o r y could s t o p here. All w h a t is left are mere application details. T h e y are of course essential for t h e practical use of geostatistics, b u t all major ideas are contained in t h e preceding few paragraphs. This is probably w h y geostatistics until recently never attracted the a t t e n t i o n of professional statisticians. 5.6 T H E N U G G E T E F F E C T C
5.6.1
0
Generalities
In m a n y instances, t h e variogram will n o t tend t o zero as h d o e s . This means that if o n e goes back t o sample the same site, t h e results will differ.
99
This m a y occur for a n u m b e r of reasons, such as p o o r analytical precision, p o o r sample preparation or more often highly erratic mineralization at l o w scale. This later p h e n o m e n o n , like t h e presence or absence of a gold nugget, is obviously responsible for t h e n a m e of that intercept at t h e origin. Usually statisticians simply call this t h e Chaotic C o m p o n e n t . It can be considered as t h e variance o f a totally random c o m p o n e n t superimposed o n t h e regionalized variable. As such, it is inversely proportional t o t h e v o l u m e of samples and thus t h e nugget effect of t h e variogram of 5 0 ' samples will be five times less than t h e o n e of 1 0 ' samples. It has n o effect o n covariances but will increase p o i n t variances b y C . One should consider it as an un certainty o n t h e value of t h e sample itself. 0
5.6.2 Theoretical
approach
A t h o r o u g h theoretical treatment of t h e nugget effect s h o w s that o n e c a n n o t keep a simple form like j(h) = C + y (h). This is d u e t o t h e fact that t h e d i s c o n t i n u i t y is n o t really occurring at the p o i n t level. What o n e has in fact is a very short range transition and t h e c o m p u t a t i o n s have t o be made taking i n t o a c c o u n t t w o variograms w h i c h represent variations at t w o different levels of observation, o n e being t h e chaotic c o m p o n e n t and t h e other the structured c o m p o n e n t . H e n c e t h e variance of any sample should be t h e sum of t w o terms c o m i n g from each o f t h e t w o variograms since all integrals are linear in 7 . One can, however, s h o w that in c o m m o n cases, w h e n t h e size o f samples is larger than t h e range of t h e small transition, t h e overall result is simplified, t h e additional term being simply as s h o w n before inversely proportional t o t h e size of t h e sample. 0
0
The nugget term for a block variance The rapid transition f u n c t i o n can be written y(h) = C ~ C(h) where C(h) is a positive f u n c t i o n w h i c h is zero for distances larger than t h e small range. It corresponds t o t h e chaotic c o m p o n e n t . T h e n t h e variance of a sample v inside a large o n e V is increased b y a term 0% w h i c h is, if w e write h— \x —y\\ 0
or:
0%
100
since:
Let us c o m p u t e o n e of the remaining integrals: C\x—y\dy = C(h)dh = A, A being a constant, Jv Jv since the v o l u m e o n w h i c h C(h) is > 0 is smaller than V. Now: \
f
[\
C\x-y\dy]
dx = ~
f
A dx =
-
Thus the variance is increased b y : a
-
= A
(^i)
(5
-
23)
Usually 1/V is zero compared t o 1/v and thus the additional term is simply A/v. In practice a nugget effect will have t o be taken i n t o a c c o u n t o n l y for samples, since mining b l o c k v o l u m e s are e n o r m o u s in comparison t o sample v o l u m e s and the constant AIV again b e c o m e s absolutely negligible. The nugget term for a covariance It simply d o e s n o t exist since t h e structured c o m p o n e n t and the chaotic c o m p o n e n t are i n d e p e n d e n t of each other and t h e chaotic c o m p o n e n t s are uncorrelated b y definition. 5 . 6 . 3 Remarks
on the origin of a nugget
effect
Splitting t h e total variances of samples i n t o a chaotic "nugget" c o m p o n e n t and a structured c o m p o n e n t means that o n e considers the observed regionalized variable as t h e sum of a structured random variable, plus a random c o m p o n e n t , uncorrelated w i t h the grade. This random c o m p o n e n t m a y have m a n y origins. It m a y be a real small-scale transition d u e t o the mineralization, it can be due t o p o o r sample preparation and analysis, or simply lack of care in sample collecting. This is s o m e t i m e s referred t o as a "human nugget e f f e c t " and w h e n identified it should be eliminated if possible. On Figs. 6 9 and 7 0 t w o cases of such a h u m a n nugget effect are s h o w n . The first e x a m p l e b y Serra ( 1 9 6 7 ) clearly s h o w s t h e structured c o m p o n e n t of the variogram, w h i c h remains unaffected through t i m e , while the nugget effect sharply decreases w h e n o n e o n l y considers recent carefully collected samples. The s e c o n d e x a m p l e is from a m o l y b d e n u m deposit where t h e mineralization occurs in half an inch veinlets. In this case t h e reduction of the sample is w h a t is essential; o n o n e h o l e the preparation was careful
101
F i g . 6 9 . " H u m a n " n u g g e t e f f e c t i l l u s t r a t e d . T w o v a r i o g r a m s in t h e s a m e a r e a o f a n i r o n m i n e c o m p u t e d o n o l d ( A ) and n e w ( B ) s a m p l e s . A f t e r Serra ( 1 9 6 7 ) .
#(h)
0-05 0.04 0.03 0.02 0.01
10
20
30
40
50
Fig. 7 0 . " H u m a n " nugget effect
60 generated by poor sample preparation (A)
disappears
w i t h c o r r e c t s a m p l e r e d u c t i o n ( B ) . T h e t w o v a r i o g r a m s are f r o m a m o l y b d e n u m d e p o s i t .
but following a w r o n g s c h e m e , thus giving perfectly random results, while for the s e c o n d , Gy's sampling principles were used, and t h e structure is clearly revealed. More practical considerations a b o u t t h e nugget effect will be given in t h e chapter a b o u t variogram interpretation. We will also spend s o m e time later t o discuss sampling errors.
102 5.7 T H E O R E T I C A L M O D E L S O F ISOTROPIC V A R I O G R A M S
In practice y(h) is n o t k n o w n and has t o be estimated from t h e available samples. One thus obtains for various values of h a series of values 7 * ( / * ) . One t h e n has t o fit an e q u a t i o n t o these experimental values. So far it has b e e n possible t o estimate m o r e than a hundred deposits w i t h o n l y t w o m o d e l s and t h e t e n d e n c y is n o w t o use o n l y o n e : the spherical m o d e l . This m o d e l is m e n t i o n e d in Matern ( 1 9 6 0 , p. 3 0 ) . T o the author's k n o w l e d g e , however, it was n o t used in practical applications b y biometricians. T h e y m o s t l y tend t o use t h e e x p o n e n t i a l m o d e l although it was recognized that it c a n n o t claim any divine right (Whittle, 1 9 5 4 ) . In mining problems it was f o u n d that t h e e x p o n e n t i a l m o d e l has t o o gentle a growth t o satisfactorily represent grade variations. Without claiming more divine rights for t h e spherical m o d e l , let us present it and t h e m a n y e x a m p l e s where it can be used. 5 . 7 . 1 The spherical
model
Its shape is s h o w n o n Fig. 7 1 ; its equation is:
(5.24)
Fig. 7 1 . T h e spherical m o d e l variogram.
103
500
1000
1500
2000
2500
h(m)
Fig. 7 2 . H o r i z o n t a l variogram o f D . D . H . i n t e r s e c t i o n in iron ore in Lorraine After S e r r a ( 1 9 6 7 ) .
(France).
F i g . 7 3 . H o r i z o n t a l v a r i o g r a m o f v e r t i c a l p a c k s a k D . D . H . in a s t r a t a - b o u n d C a n a d i a n u r a n i u m d e p o s i t , s h o w i n g a r a n g e o f 6 0 m ( s a m p l i n g i n t e r v a l is 1 5 m ) . A f t e r B l a i s a n d Carlier ( 1 9 6 8 ) .
a is called t h e range, C is t h e nugget effect and C + C t h e sill. This m o d e l occurs naturally in all t h e deposits where grades b e c o m e i n d e p e n d e n t of each other o n c e a given distance, a, is reached. It is t h e c o m m o n rule in m o s t sedimentary deposits and also in porphyry copper deposits. T h e variograms s h o w n in Figs. 72—80 have b e e n satisfactorily described by such a m o d e l . Deposits as different as iron ore deposits in France, porphyry copper in Chile and Canada, lead—zinc in Ireland and Canada, have been f o u n d t o have their grade variations adequately represented b y this m o d e l . S o were gold deposits in S o u t h Africa, bauxite in France, lateritic nickel in N e w Caledonia, uranium in Canada and phosphates in Africa. This m o d e l is also k n o w n as t h e Matheron m o d e l and is said t o describe transition p h e n o m e n a 0
0
100
300
500
700
h (m)
Fig. 7 4 . A 1 0 h o r i z o n t a l variograms in the bauxite deposit of Mazauges Aval, showing isotropy. After Marechal and Roullier ( 1 9 7 0 ) . 2
3
100
300
500
2
150 4
125
H
25
H
100
200
300
400
h(m)
Fig. 7 5 . S i 0 horizontal variograms in the bauxite deposit of Mazauges Aval, show ing small nugget effect and isotropy. After Marechal and Roullier (1970).
Y(h)
100
700
h(ft)
Fig. 7 6 . Vertical variogram o f m a g n e t i t e grade in Canadian iron ore d e p o s i t .
105 Y(h)
100
300
500
700
900
1100
h(ft)
Fig. 7 7 . V a r i o g r a m s a l o n g strike and d o w n d i p for g o l d grade in a Canadian quartz vein deposit showing over 6 0 0 of continuity.
2Y(h)f
300 H
50
h(m)
Fig. 7 8 . H o r i z o n t a l variogram in an African g o l d m i n e s h o w i n g very short z o n e o f influ e n c e . A f t e r C o r t e z e t al. ( 1 9 7 4 ) .
V(h) 1200
300 H 150
500
1000
1400
h(m)
Fig. 7 9 . H o r i z o n t a l variogram o f t h e Silica a c c u m u l a t i o n in garnieritic n i c k e l of N e w Caledonia. After Journel ( 1 9 7 4 ) .
deposits
106 v(h)
r
r1
—Ti
150
1 —
1000
500
Fig. 8 0 . Variogram
of
the
nickel accumulation
—
i
—
1400
h(m)
in garnieritic nickel d e p o s i t s o f
New
Caledonia. After Journel ( 1 9 7 4 ) .
as it is t h e o n e w h i c h occurs w h e n o n e has structures i n d e p e n d e n t of each other but within each, t h e grades are highly correlated. A digression on the origin of the 'spherical' name and genetical implication of variograms A process having a spherical variogram can be generated as f o l l o w s : first take a random (Poisson) distribution of masses in space. T h e n assign t o each p o i n t of space a grade given b y the n u m b e r of masses within a radius a of the point (in other words, take t h e number of masses within a "sphere" . . .). Also n o t e that t h e v o l u m e of t h e intersection of t w o spheres of diameter a, v o l u m e V, t h e centers of w h i c h are h apart is:
if h < a and 0 otherwise. It is interesting t o be able t o construct (for simu lation purposes) a process w i t h a given variogram; w e will n o t e x p a n d o n t h e subject however, e x c e p t t o warn a t e m p t e d reader n o t t o draw totally h y p o t h e t i c a l genetic conclusions from mathematical games. The only accept able process w o u l d be the reverse: given an alleged mineralizing process t o theoretically c o m p u t e its variogram and later experimentally c h e c k it. 5.7.2 TheDe
Wijsian model
(after Prof. H.J. de Wijs)
This m o d e l is e n c o u n t e r e d in cases where there is n o such thing as a range of d e p e n d e n c e . y(h) keeps increasing as h increases. In m o s t of t h e hydrothermal deposits, t h e variogram plots as a straight line w h e n y(h) is p l o t t e d against ln(h). The equation is thus: y(h)
= A ln(fc) + B
(5.25)
Such variograms can be seen in Figs. 8 1 — 8 3 . The vogue of this m o d e l was certainly due t o t h e " g o o d " properties of the logarithmic function w h i c h
107
I
i
i
I
2
3
i
4
i
5
i
6
7
i
i
8
9
i
i
•
10
h
F i g . 8 1 . D e Wijsian v a r i o g r a m f o u n d f r o m c h a n n e l s a m p l e s a l o n g u p p e r d r i f t i n a Z a m b i a sulfide deposit.
I
i
i
I
2
3
i
4
i
5
6
i
7
i
8
i
9
i
i
•
10
h
F i g . 8 2 . D e Wijsian v a r i o g r a m f o u n d f r o m c h a n n e l s a m p l e s a l o n g l o w e r drift i n a Z a m b i a sulfide deposit. X(h)4
I
2
3
4 5 6 7 8 9 1 0
15
20
h
F i g . 8 3 . G e n e r a l v e r t i c a l v a r i o g r a m ( s o l i d l i n e ) o f d i a m o n d drill h o l e c o r e s a m p l e s a n d horizontal variogram ( b r o k e n line) o f c h a n n e l s a m p l e s o f t h e mineralization in the Mounana uranium deposit, Gabonese Republic. After Matheron ( 1 9 6 2 ) .
108
Y(h)f
0.8
1.6
2.4
3.2
4.0
4.8
h
Fig. 8 4 . Several e x p e r i m e n t a l variograms correctly represented b y a linear m o d e l .
a m o n g other things is harmonic. The first v o l u m e s of Matheron ( 1 9 6 3 ) and Carlier ( 1 9 6 4 ) were entirely d e v o t e d t o it. T h e y have lost a l o t of their interest w i t h t h e appearance of c o m p u t e r s , b u t should be remembered specially in case where a c o m p u t e r is n o t available and quick results are n e e d e d (see Exercise 4 . 6 . 2 ) . 5 . 7 . 3 Other
models
5 . 7 . 3 . 1 The linear model The m o s t simple m o d e l w h i c h o n e tries t o use as m u c h as possible for obvious reasons is just a linear o n e . Quite often it is possible at least within a limited n e i g h b o u r h o o d t o use a linear equation: y(h)
=
Ah+B
Several e x a m p l e s are s h o w n in Fig. 8 4 w h i c h is a p l o t o f Fig. 8 9 for distances smaller than 5 m. N o t e that this variogram is the theoretical variogram of a process called t h e Brownian m o t i o n . (Feller, 1 9 6 6 , p. 9 7 ) . 5 . 7 . 3 . 2 The ah model We have already seen t h e linear m o d e l w h e n p l o t t e d o n a log scale, w h i c h is the D e Wijsian m o d e l ; w e can also in m a n y cases m a k e the variogram straight by plotting it on a log-log scale. Then the equation is y(h) = ah since ln y(h) = X ln h + constant. This is often m e t in elevation variograms (Figs. 8 5 and 8 6 ) or in the s t u d y of mill feed variability (Fig. 8 7 ) . K
x
5 . 7 . 3 . 3 The exponential model The e x p o n e n t i a l m o d e l is n o t e n c o u n t e r e d t o o often in mining practice
109
¥(h)
8.0
6.0
4.0
2.0
0.0
Fig. 8 5 . Plot o n arithmetic scale o f the variogram of the elevation o f a hanging wall in a gold d e p o s i t .
Y(h)
Fig. 8 6 . S a m e variogram plotted o n a log-log scale.
as
in
Fig.
85,
4 6 10 20 4060 100 200 400 h( minutes) Fig. 8 7 . P l o t o n a log-log scale o f t h e variogram o f mill feed grade for the mill o f t h e SaintS e b a s t i e n M i n e ( F r a n c e ) . A f t e r P. G y ( 1 9 6 7 ) .
110
since its infinite range is associated t o a t o o c o n t i n u o u s process. It is h o w ever very practical. The equation is: y(h)
=C {
+ C
0
5.1.3.j/The
exp
hole effect
-
AM
a
model
The "hole e f f e c t " m o d e l s h o w n in Fig. 8 8 has for an e q u a t i o n : tW
=
/sin ah \ c^i—
It can be used t o represent fairly c o n t i n u o u s processes (the tangent at t h e origin is horizontal and it s h o w s a periodic behaviour w h i c h is o f t e n e n c o u n tered w h e n there exists a succession of rich and p o o r z o n e s , for instance).
F i g . 88. P l o t o f t h e t h e o r e t i c a l v a r i o g r a m 7 ( h ) = C [ l — s'm(ah)/ah] and a periodic hole effect.
showing slow growth
The e x a m p l e in Fig. 8 9 is from an Australian uranium m i n e where rich and p o o r layers alternate. In Fig. 9 0 w e see an e x a m p l e from Serra ( 1 9 7 0 ) w h i c h corresponds t o t h e e x i s t e n c e of waste concretions in t h e Lorraine iron ore. Figure 9 1 s h o w s a more elaborated m o d e l b y Marechal and D i d y k ( 1 9 7 1 ) w h i c h was d e v e l o p e d for t h e El Salvador copper m i n e in Chile.
Ill
Fig. 8 9 . Vertical variograms c o m p u t e d for increasing s a m p l e l e n g t h s in an uranium deposit showing a marked hole effect.
0
0.5
1.0
1.5
Australian
2.0 2.5 (METERS)
Fig. 9 0 . H o r i z o n t a l variogram o f iron, c a l c i u m and silica grade s h o w i n g s h o r t scale, transition p h e n o m e n o n and h o l e effect ( e x i s t e n c e of c o n c r e t i o n s ) in t h e Lorraine iron ore deposit. After Serra ( 1 9 7 0 ) .
The equations w h i c h were used are a c o m b i n a t i o n o f e x p o n e n t i a l and sine model: # ) =
P
( l - e - )
+
d
-
p
)
(
^
)
112 Y(h>r 0.54
OA A 3
9
15
21
27
33
39
45
51
57
h(m)
F i g . 9 1 . V a r i o g r a m o f c o p p e r g r a d e at t h e El S a l v a d o r c o p p e r d e p o s i t ( C h i l e ) . A f t e r Marechal and D i d y k ( 1 9 7 1 ) .
5 . 7 . 4 Admissible
functions
for a
variogram
N o t all functions are suitable t o represent a variogram. One has t o bear in mind that a variogram represents a physical process n a m e l y t h e grade at o n e p o i n t Z{x). In the estimation procedure w e will have t o use linear combina t i o n s o f k n o w n values. We m u s t m a k e sure that t h e variance of these c o m binations exists and that it is positive. It has b e e n s h o w n that a linear c o m bination Z* has a finite variance if and o n l y if t h e sum of t h e weights a is equal t o 0. In this case: t
VAR(Z*) = -I
£
a.ajjix.-xj)
1=1 7 = 1
This defines t h e c o n d i t i o n w h i c h a 7 - f u n c t i o n should satisfy in order t o be admissible as a variogram. This c o n d i t i o n o n — 7 is that it m u s t be condition ally positive definite: for any set o f n p o i n t s x x , . . x and any set o f weights a a , . . ., a such that 2 a = 0 o n e should have: u
u
n
-
2
n
2
n
f
n
Z Z aityyixi /=i y=i
~xj)
>
0
This is a very strict c o n d i t i o n and c o n s e q u e n t l y o n e c a n n o t simply fit which ever p o l y n o m i a l m o d e l o n e w o u l d like. One should c h e c k for t h e positivity o f the m o d e l or use a m o d e l w h i c h is k n o w n t o be satisfactory like t h e ones w e have just m e n t i o n e d . Many other equations have b e e n proposed b y Matheron ( 1 9 6 8 ) . He has also c o m p u t e d all t h e p o l y n o m i a l s o f order less than 6 w h i c h are admissible.
113
However, for practical purposes w e will see that w e should be h a p p y w i t h the few m o d e l s listed above and w i t h c o m b i n a t i o n s of t h e m , since a linear c o m b i n a t i o n of t h e previous equations is an admissible e q u a t i o n . We will see m a n y e x a m p l e s in t h e chapter a b o u t m o d e l fitting. For readers interested in t h e theoretical p r o b l e m of t h e class o f admissible isotropic covariance functions in R , w e refer t h e m t o t h e very detailed treatment of this p r o b l e m w h i c h was given b y Matheron in t h e first part of his thesis ( 1 9 6 5 ) , or t o Matern ( 1 9 6 0 ) w h o sums up t h e w o r k of B o c h n e r ( 1 9 3 2 ) , S c h o n b e r g ( 1 9 3 8 ) , Hartman and Wintner ( 1 9 4 0 ) , Lord ( 1 9 5 4 ) and Yaglom ( 1 9 5 7 ) and gives an extensive list of admissible functions. n
5.8 HOW C L O S E IS T H E I N T R I N S I C H Y P O T H E S I S T O R E A L I T Y ?
5 . 8 . 1 Empirical
discussion
of the
problem
It might be important at this stage t o spend m o r e t i m e discussing w h a t h y p o t h e s e s and m o d e l s are, and also introducing t h e c o n c e p t of robustness. First of all, a m o d e l is n o t reality; it is t h e expression o f our inability t o fully understand and unravel t h e deterministic processes w h i c h generated the present state of t h e p h e n o m e n o n (in our case, t h e ore) w h i c h is observed. Geostatistics is a probabilistic m o d e l ; t h e law o f gravity is a deterministic m o d e l ; w e could use statistics instead t o find t h e p o s i t i o n at a given t i m e of a falling b o d y . U n f o r t u n a t e l y , t o c o m p u t e t h e grade of a b l o c k such a simple equation as e = \gt c a n n o t be used. What is d o n e instead is t o consider that t h e grade of our b l o c k is o n e realization of a random variable, just like t h e result o f a die t h r o w is considered as o n e realization o f a uniform distri b u t i o n b e t w e e n 1 and 6. In fact t h e o u t c o m e of t h e die t h r o w is purely determined b y t h e original p o s i t i o n of t h e die, t h e w a y it is t h r o w n and possibly o t h e r parameters w h i c h are simply t o o n u m e r o u s t o be c o u n t e d and measured. N o w t h e probabilistic m o d e l w h i c h is used t o describe the out c o m e of t h e t h r o w o f a die is very simple and it is easy t o c h e c k w h e t h e r t h e die is biased or n o t , in other w o r d s , if t h e m o d e l is correct for this die. When w e are dealing w i t h ore grade, w e are dealing w i t h a random variable w i t h an infinite n u m b e r o f d i m e n s i o n s (grade at each p o i n t ) and t h e equivalent of k n o w i n g t h e probability of t h e o u t c o m e of a 1 for a die n o w b e c o m e s k n o w i n g for any n and any set of p o i n t s x (i = 1, . . . , n) t h e probability that simultaneously: 2
t
{Z(x ) x
< Z
u
Z(x ) 2
<
Z , . . . , Z(x ) 2
n
<
Z) n
This is clearly impossible, at least at present. The intrinsic h y p o t h e s i s intro duces s o m e suppositions, such as t h e stationarity of increments, w h i c h can t h e n lead t o t h e definition and c o m p u t a t i o n o f a variogram and t h e n t o any variance c o m p u t a t i o n s . One of t h e first criticisms o f t h e intrinsic h y p o t h e s i s was that nature is
114
n o t even w e a k l y stationary. This is true, but w h a t should be considered, rather than rejecting t h e m o d e l , is t o try and estimate h o w w r o n g t h e m e t h o d is w h e n applied in non-stationary c o n d i t i o n s . This is the c o n c e p t of robust ness. Since these early criticisms, m a n y d e v e l o p m e n t s have occurred and the universal kriging h y p o t h e s i s of Matheron ( 1 9 7 1 ) has already b e e n introduced in order t o deal with non-stationary p h e n o m e n a . More recently, the t h e o r y of t h e random intrinsic functions of order k has firmly established t h e problem o f t h e estimation of a random f u n c t i o n under any h y p o t h e s i s (Matheron, 1 9 7 3 ) . Thus t h e theoretical objections have all been answered, b u t this d o e s n o t m e a n that these sophisticated t o o l s have t o be used as s o o n as non-stationarity is suspected. In fact, m a n y experiments have been m a d e t o s h o w that using a simple variogram is very little affected b y very strong departure from h y p o t h e t i c a l c o n d i t i o n s . A n a c c o u n t of this is given in David and Dagbert ( 1 9 7 6 ) . Besides this, the n o t i o n of stationarity can o n l y be con sidered w i t h respect t o a given scale at w h i c h o n e wishes t o work, and m o s t of the time there exists a n e i g h b o u r h o o d for an area of interest where stationarity can be assumed. With a given set of samples, it is always possible t o c o m p u t e a variogram; this d o e s n o t mean that this variogram is t o be used if it has b e e n c o m p u t e d o n a n o n - h o m o g e n e o u s region; it will exhibit features w h i c h should send us back t o t h e geological definition of h o m o g e n e o u s z o n e s in an ore deposit. One sees at this p o i n t that g e o l o g y and numerical c o m p u t a t i o n s should c o m p l e m e n t each other in order t o p o i n t o u t the different z o n e s w h i c h o n e should consider in estimating t h e deposit. 5 . 8 . 2 Theoretical
discussion
of the
problem
A theoretical discussion of these problems is definitely important but m u c h b e y o n d t h e s c o p e of this manual. Interested readers will find all necessary justifications in Matheron ( 1 9 7 1 ) under t h e heading "Statistical Inference" and in his m o s t recent paper ( 1 9 7 6 a ) a b o u t t h e c h o i c e of m o d e l s in geostatistics. A m o n g the points that are stressed, considering t h e relevance of a m o d e l , four c o n d i t i o n s w h i c h a m o d e l m u s t satisfy are p o i n t e d o u t : Statistical inference m u s t be possible, i.e., it should be possible t o estimate t h e parameters o f t h e m o d e l . Then t h e m o d e l m u s t be operational, i.e., it m u s t give an answer t o the question asked. This is obvious but should rather be considered in t h e reverse; depending o n the questions o n e wants t o answer, o n e m a y use different m o d e l s for t h e same deposit. The m o d e l has t o be compatible w i t h t h e data, and finally t h e predictions o f t h e m o d e l must be checked b y experience. This will be an a-posteriori proof. N o w a m o n g all m o d e l s m e e t i n g these c o n d i t i o n s , h o w is t h e c h o i c e t o be made? The answer is short and clear: Take t h e simplest.