Contributions to Industrial-Minerals Research

U.S. Department of the Interior

U.S. Geological Survey

Bulletin 2209–J

Issues and Challenges in the Application of 

Geostatistics and Spatial-Data Analysis to the 

Characterization of Sand-and-Gravel Resources

Chapter J of
Contributions to Industrial-Minerals Research

James D. Bliss, Phillip R. Moyle, and Keith R. Long, Editors

Bulletin 2209–J

U.S. Department of the Interior

U.S. Geological Survey 

Issues and Challenges in the Application of 

Geostatistics and Spatial-Data Analysis to the 

Characterization of Sand-and-Gravel Resources

By Daniel R. Hack

Chapter J of
U.S. Department of the Interior

Gale A. Norton, Secretary

U.S. Geological Survey

Charles G. Groat, Director

Version 1.0, 2005

This publication is available only online at:

http://pubs.usgs.gov/bul/b2209-j/

Manuscript approved for publication, July 24, 2003

Text edited by George A. Havach

Layout by Stephen L. Scott

Any use of trade, product, or firm names in this publication

is for descriptive purposes only and does not imply 

endorsement by the U.S. Government 
CONTENTS

Abstract--------------------------------------------------------------------------------------------------- 1

Introduction----------------------------------------------------------------------------------------------- 1

Background----------------------------------------------------------------------------------------------- 2

Current Approaches to Aggregate-Resource Exploration and Evaluation------------------------ 2

Concepts of Exploratory-Data Analysis and Spatial-Data Analysis------------------------------- 4

Basic Concepts of Geostatistics----------------------------------------------------------------------- 5

Existing Geostatistical Investigation------------------------------------------------------------------- 8

 Further Research Based on Existing Geostatistical Investigations--------------------------- 9

Limitations of Geostatistics----------------------------------------------------------------------------- 11

 Nonadditive Variables------------------------------------------------------------------------------ 11

 Uniform Sample Size------------------------------------------------------------------------------- 11

 Support and Selective-Mining-Unit Issues------------------------------------------------------ 12

 Other Considerations------------------------------------------------------------------------------- 12

How to Practically Confront These Limitations------------------------------------------------------- 12

 General Concepts----------------------------------------------------------------------------------- 12

 Drilling Programs----------------------------------------------------------------------------------- 13

Proposals for Future Work------------------------------------------------------------------------------ 13

References Cited----------------------------------------------------------------------------------------- 13

Figures

1. Basic features of the semivariogram------------------------------------------------------------ 6

2. Spherical, exponential, and Gaussian semivariogram model forms------------------------- 7

3. Semivariogram models demonstrating anisotropy and quasi-stationarity------------------ 7

4. Log-transformed semivariogram of density----------------------------------------------------- 10

Table

1. Correlation coefficients for the lower Tenmile Terrace, Idaho------------------------------- 11

III
1

Abstract

Sand-and-gravel (aggregate) resources are a critical 

component of the Nation’s infrastructure, yet aggregatemining technologies lag far behind those of metalliferous 

mining and other sectors. Deposit-evaluation and sitecharacterization methodologies are antiquated, and few 

serious studies of the potential applications of spatial-data 

analysis and geostatistics have been published. However, 

because of commodity usage and the necessary proximity 

of a mine to end use, aggregate-resource exploration and 

evaluation differ fundamentally from comparable activities for metalliferous ores. Acceptable practices, therefore, 

can reflect this cruder scale. The increasing use of computer 

technologies is colliding with the need for sand-and-gravel 

mines to modernize and improve their overall efficiency 

of exploration, mine planning, scheduling, automation, 

and other operations. The emergence of megaquarries 

in the 21st century will also be a contributing factor.

Preliminary research into the practical applications 

of exploratory-data analysis (EDA) have been promising. 

For example, EDA was used to develop a linear-regression 

equation to forecast freeze-thaw durability from absorption values for Lower Paleozoic carbonate rocks mined for 

crushed aggregate from quarries in Oklahoma. Applications 

of EDA within a spatial context, a method of spatial-data 

analysis, have also been promising, as with the investigation of undeveloped sand-and-gravel resources in the 

sedimentary deposits of Pleistocene Lake Bonneville, Utah.

Formal geostatistical investigations of sand-and-gravel 

deposits are quite rare, and the primary focus of those studies 

that have been completed is on the spatial characterization of 

deposit thickness and its subsequent effect on ore reserves. 

A thorough investigation of a gravel deposit in an active 

aggregate-mining area in central Essex, U.K., emphasized 

the problems inherent in the geostatistical characterization of 

particle-size-analysis data. Beyond such factors as common 

drilling methods jeopardizing the accuracy of the sizedistribution curve, the application of formal geostatistical 

principles has other limitations. Many of the variables used 

in evaluating gravel deposits, including such sedimentologic 

parameters as sorting and such United Soil Classification 

System parameters as gradation coefficient, are nonadditive. Also, uniform sampling methods, such as drilling, are 

relatively uncommon, and sampling is generally accomplished by a combination of boreholes, water-well logs, test 

pits, trenches, stratigraphic columns from exposures, and, 

possibly, some geophysical cross sections. When evaluated in consideration of the fact that uniform mining blocks 

are also uncommon in practice, subsequent complexities in 

establishment of the volume/variance relation are inevitable.

Several approaches exist to confront the limitations of 

geostatistical methods in evaluating sand-and-gravel deposits. Initially, we must acknowledge the practical requirements of the aggregate industry, as well as the limitations 

of the data collected by that industry, as a function of what 

the industry requires at the practical level, and consider 

that broader acceptance of formal geostatistics may require 

modifications of typical exploration and sampling protocols.

Future investigations should utilize data from the full 

spectrum of sand-and-gravel deposits (flood plain, glacial, 

catastrophic flood, and marine), integrate such other disciplines as sedimentology and geophysics into the research, 

develop commodity-specific approaches to nonadditive 

variables, and include the results of comparative drilling.

Introduction

Sand-and-gravel mines, as sources for concrete 

aggregate, exert a tremendous influence on the Nation’s 

infrastructure, yet few mining researchers or professionals have taken aggregate-mining technology seriously. 

Deposit-evaluation and site-characterization methodologies are at least 2 decades behind those in other mining 

sectors, and applications of geostatistics and spatial-data 

analysis have been little studied. This lack of commodity-specific engineering research consequently affects mine 

planning, scheduling, and other unit operations, compromising the overall efficiency of gravel-mining operations.

This chapter examines current approaches to aggregate-resource exploration and evaluation, describes the 

current status of research into applications of geostatistics and spatial-data analysis in the aggregate sector, 

reviews the methodologies developed and implemented (or 

Issues and Challenges in the Application of 

Geostatistics and Spatial-Data Analysis to the 

Characterization of Sand-and-Gravel Resources

By Daniel R. Hack1

Introduction

_____________________

1

Missouri Department of Transportation, Jefferson City, Mo.
2 Contributions to Industrial-Minerals Research

attempted) by previous workers, suggests the limitations 

inherent in using geostatistics in gravel-deposit evaluation, 

and proposes how exploration methodologies, sampling 

protocols, and data evaluation may need to be modified 

to characterize gravel deposits within a spatial context.

In this chapter, I define the basic concepts of exploratory-data analysis (EDA), spatial-data analysis (SDA), and 

formal geostatistics with respect to applications to aggregateresource characterization. I discuss the limitations of these 

applications, such as issues with the volume/variance relation, nonadditive variables, and the use of size-distribution 

curves. I propose potential solutions for practically confronting these limitations at the level of resource exploration and 

analysis, as well as concise suggestions for the design of 

sampling programs. I conclude by proposing measures to 

improve future results, such as comparative drilling exercises 

(that is, shell and auger versus reversed-circulation drilling) 

before beginning serious exploration, and assessing data from 

the full spectrum of sand-and-gravel deposits (flood plain, 

glacial, catastrophic flood, and marine) in order to better 

understand deposit-specific data-collection requirements.

Background

Arthur (1994), who presented the most comprehensive 

examination to date of the potential application of geostatistical principles in economic sand-and-gravel evaluation, 

remains the only serious such attempt, with a few minor 

exceptions. He proposed sampling methodologies that have 

never been implemented, or even taken seriously, by the 

aggregate sector. His summary comments were undeniably 

cynical (Arthur, 1994, p. 1):

* * * the lack of research on inplace resource estimation in the 

sand-and-gravel industry is highlighted by the lack of papers on 

the subject in the scientific literature, especially when compared to the rest of the mining industry. * * * there appeared to 

be several reasons for the lack of interest, the main ones being 

the lack of in-house knowledge of the techniques involved, 

the high cost of installing computerized systems [in 1994], 

scepticism concerning the advantages over the cost and time 

involved in obtaining data, and a general lack of understanding of the problems that can occur. * * * very little work has 

been carried out in the past to assess the precision of data 

obtained from the drilling and sampling of natural aggregates.

In the subsequent decade, little has changed, and the stumbling 

blocks that Arthur encountered still remain. I discuss his methods and results in the section below entitled “Existing Geostatistical Investigations.”

The increasing saturation of computer technologies in 

every facet of daily life, however, will inevitably influence 

the acceptance of computerized methods in the aggregate 

sector, particularly when this trend is taken into account 

with the emergence of megaquarries. Sand-and-gravel mines 

are becoming large independent operations that must plan 

and operate similarly to open-pit metal mines as the realities of government regulation, automation, and economies 

of scale catch up with them. Computerization in the aggregate sector, however, is may be 2 decades behind that in the 

metal-mining industry—easily several generations of computer technology. Given imminent modernization, however, 

mine-planning-software companies eventually must seriously 

address the requirements of this sector (Hack, 2002b).

The demands of 21st century infrastructure have compelled 

the emergence of megaquarries (Bliss and others, 2002), 

which are large mines accessible by rail or marine transport that produce more than 5 million tons of aggregate per 

year, with an expected lifetime of 50 years. Megaquarries 

are an inevitable result of two factors: the tremendous cost 

of permitting numerous small minesites, and an accelerating trend toward consolidation in the aggregate and concrete 

industry that began in earnest in the 1990s. This creation of 

aggregate companies with capacities for output scattered 

across the globe, and subsequent planning and operational 

requirements for these new economies of scale, ensure that 

the aggregate sector will require deposit evaluation (and the 

corresponding computer technologies) far advanced from 

what is presently in use. It is reasonable to propose that in the 

next decade, megaquarries identified as “type III” (consisting 

entirely of sand and gravel) will constitute most worldwide 

aggregate production (Bliss and others, 2002). The proximity of most type III megaquarries to maritime shipping 

lanes will broaden the scope of influence of these mines 

considerably. In addition, several large marine aggregate 

deposits are being developed to supplement onshore production, diminished because of environmental constraints.

Another development in the aggregate sector that will 

require deposit-evaluation methodologies comparable to those 

in other mining-industry sectors is the emergence of “type 

II” megaquarries. These mines produce rock not only for 

aggregate but also for various products suitable for such uses 

as chemicals, agricultural lime, filler, and cement. According to Bliss and others (2002, p. 306–307), “these quarries 

require considerable management to maximize return from 

the complex mixture of commodities they produce.” Aggregate miners, to date, have few (if any) industry-specific 

tools at their disposal to plan and manage such operations, 

whereas other mining-industry sectors have benefited from 

geostatistics and mine computing for decades. One of the 

first steps to the successful application of such techniques 

for sand-and-gravel mines is the evaluation of SDA methodologies, which will aid the evolution of aggregate-mine 

engineering from a haphazardly planned accessory to concrete production, ultimately to sophisticated three-dimensional planning, scheduling, and optimization routines.

Current Approaches to AggregateResource Exploration and Evaluation

Aggregate-resource exploration and evaluation differ 

fundamentally from comparable activities for metalliferous ores. Barksdale (1991) noted that the goal of aggregate 

exploration is to supply a given market at a competitive price. 

Although the goal in metal exploration can be described 
3

in comparable terms (that is, locating an economic-mineral resource), aggregate exploration is much more directly 

linked to the end use of the mined commodity. At the 

most basic scale, a single basalt bench can be located and 

mined for a few weeks to supply a paving project immediately adjacent. The factors in selecting of such a site are 

simple: material suitability, proximity to end use, landowner 

approval, and any applicable municipal land-use issues.

Likewise, deposit exploration and evaluation methodologies reflect this cruder scale. The various common aggregate 

references and guidebooks describe exploration and evaluation methodologies, although experience has demonstrated 

that the exploration procedure can easily be backward: the 

deposit finds the miner, rather than the other way around. For 

example, large gravel companies search out old mines and 

prospects, such as State department of transportation (DOT) 

pits, or promoters (“bird doggers”) call the sites to the companies’ attention. Also, a large producer can consolidate several 

small gravel pits in the same deposit into a single operation. 

With this methodology, at least a preliminary indication exists 

that the material meets some level of quality; and in terms 

of land-use restrictions, the value of the site for aggregateresource extraction has already been established. This final 

point implies an overwhelming factor in aggregate-resource 

development in the United States today: municipal land-use 

and zoning constraints. As a component of its Front Range 

Infrastructure Resources Project (FRIRP), the U.S. Geological Survey has developed geographic-information-system 

(GIS) techniques for the evaluation and protection of aggregate resources at the urban-planning stage, as well as for use 

as an exploration base (Melick and others, 1998). Lindsey 

(1997), in another report affiliated with FRIRP, provides an 

excellent introduction to the basic characteristics of deposit 

types, how models are developed, and the types of data used.

The “exploration guidance” of Barksdale (1991), which 

may be the most down to earth in terms of what is realistic 

for many companies, explains resource exploration in the 

most practical terms. His description of site selection begins 

with a discussion of the “principle of the weakest point” (and 

rightly so), pointing out the senselessness of spending money 

to develop a resource if one critical component in the exploration process (zoning, material quality, or cooperation of landowners) is missing. The importance of transportation routes is 

another factor that differentiates aggregate- from metalliferous-resource evaluation: aggregate must be mined near its end 

use, and haulage and handling costs compose a predominant 

percentage of final unit costs. Beyond these initial factors, 

Barksdale (1991) continues with further considerations, seemingly simplistic in opposition to something as sophisticated as 

the discovery and examination of porphyry copper deposits.

The exploration guidance provided by other workers, 

such as Smith and Collis (1993) and Primel and Tourenq 

(2000), is relevant, but the methodologies they outline are 

standard geologic procedures: literature review, reconnaissance-level field investigations, local sampling and mapping, data analysis, deposit modeling, and presentation of 

results. Both of these references, however, place particular emphasis on “inventories.” Smith and Collis describe 

British Geological Survey geologic- and geophysicalresource appraisals for various areas in the United Kingdom (including offshore resources). Primel and Tourenq 

describe small (1:100,000)- and medium (1:25,000)-scale 

inventories produced by the French Geological Survey 

and regional laboratories of the Department of Roads and 

Bridges. As does Barksdale (1991), Smith and Collis suggest the “three-legged stool” of municipal land use, transportation networks, and end use of the commodity.

In general, suitable information has been published on 

the traditional geologic aspects of aggregate-deposit exploration and evaluation, with notably progressive discussions of 

satellite imagery, remote sensing, and airborne geophysics. 

Many of these concepts, however, would be far beyond the 

understanding (or interest) of most regionally based concrete 

companies. Smith and Collis (1993, p. 55–58) skillfully 

summarize the various aspects of common drilling methods, with their general applications and characteristics:

• Hand augers: “the strata must be self-supporting, soft, 

 and with only very few pebbles.”

• Power augers: “numerous disadvantages when precise 

 sampling is required. Unconsolidated material, particu- 

 larly beneath the level of the water table, tends to col-

 lapse into the borehole *** tends to churn the strata, lead- 

 ing to contamination of one layer with another *** can 

 provide good guidance on the thickness of overburden 

 can quickly penetrate unstable strata without the need 

 for lining tubes.”

• Hollow-stem augers: “can penetrate unstable strata with-

 out the need for lining tubes *** sampling quality is likely 

 to be higher than that achieved by open-hole augering.”

• Light-cable percussion boreholes: “an adaptation of stand- 

 ard well-boring methods *** a limitation on the depth of 

 drilling is commonly imposed by the capacity of the rig to 

 extract the lining tubes from the ground strongly preferred 

 where sampling is the major objective because progress 

 cannot be maintained unless casing is inserted, and this 

 ensures samples remain uncontaminated.”

• Reverse circulation drills: “sampling is unconsolidated 

 and entirely representative *** the speed of the drilling 

 operation is such that refined sampling and logging may 

 not be achieved unless many hands are available.”

Arthur (1994) discusses some of the limitations highlighted 

for these methods, particularly with respect to geostatistical and 

other SDA approaches, in the course of his investigations. At 

the level of “prospecting,” the approach of Primel and Tourenq 

(2000) is not as firmly based on classical field geology as is that 

of Smith and Collis (1993). Primel and Tourenq particularly 

emphasize the use of topographic and isopach maps to properly 

evaluate and exploit variations in deposit and overburden thickness. They also astutely suggest that the most important data to 

extract from regional studies are mean grading within the zone, 

mean and extreme thicknesses of the overburden and gravel, 

mineralogic and petrographic compositions of various aggregate 

sizes, and principal heterogeneities.

Current Approaches to Aggregate-Resource Exploration and Evaluation
4 Contributions to Industrial-Minerals Research

Again, Barksdale (1991, p. 4.25)’s suggestions are 

generally candid and practical, most productive with 

respect to his target audience. Sampling methods are not 

confined only to drilling but also involve grab samples 

(“preliminary only and rarely can be considered typical of 

a deposit”) and pit samples (“relatively inexpensive”).

Principles of the sampling of particulate material must 

be regarded in tandem with these commodity- and industryspecific considerations. Resource data, as well as subsequent 

analysis and modeling, will be suspect unless workers have 

obtained representative samples in the field and retained 

their integrity throughout analytical procedures. Theories of 

the sampling of particulate solids were developed by Pierre 

Gy on the basis of seven types of sampling errors, along 

with concise techniques for their minimization. Examples of 

sampling errors include “increment-delimitation error” (error 

tied to inappropriate sampling design and the wrong choice 

of equipment) and “grouping and segregation error” (error 

due to nonrandom distribution of particles, usually by gravity). Gy’s theories became popular in practice because useful 

sampling and subsampling methods can be derived that are 

applicable for little or no extra expense. Pitard (1993)’s book 

is a fundamental reference on these theories and methods.

By the nature of the commodities in question, a fundamental stumbling block exists in the application of metalliferous-resource exploration and evaluation methods and 

protocols to the aggregate sector. “Gradation,” otherwise 

known as the characterization of particle-size distribution 

and also called granulometrics, is a key factor in identifying suitable aggregate prospects for a given market demand. 

In addition, the characterization within a spatial context of 

the gradation of a gravel deposit could enable gravel-mine 

planning to move closer to the level of planning expected at 

most other mines (metalliferous, coal, and most industrial 

minerals) in the present era. However, with a few exceptions that have met with mixed results, such investigations 

have never been taken seriously in the aggregate sector. 

The purpose of this chapter is to identify the analytical 

methods available, summarize the spatial-data investigations that have been completed, outline the limitations of 

these methods, and establish goals for further study.

Concepts of Exploratory-Data Analysis 

and Spatial-Data Analysis

The basic concepts of exploratory-data analysis (EDA) 

were summarized by Tukey (1977), in what remains the 

seminal EDA guidebook to this day, although many of these 

concepts had already been developed and used in Earth-science applications for some time (Krumbein and Graybill, 

1965). EDA is not a concise set of techniques or unifying 

principles (as with geostatistics) but a general analytical 

approach that usually employs graphical methods to allow the 

data to reveal their own model. The intent of EDA is to maximize insight into a dataset by (among other goals) uncovering underlying data structures and detecting outliers and 

anomalies. Typical techniques, as described by Tukey, (1977), 

are “stem-and-leaf displays,” “box-and-whisker plots,” and 

“parallel-schematic plots.” EDA has gained popularity in 

mining-operations research and simulations, because operating mines produce a tremendous amount of data that must 

be characterized if they are to be used efficiently. Kannan 

and others (1999) provided a good introductory demonstration of EDA applications in mine-operations research.

“Spatial-data analysis” (SDA) is a general term that 

encompasses methods as simple as rudimentary plotting 

and contouring, and as complex as lattices and geostatistics. In this chapter, the term includes methodologies that 

evaluate data in consideration of their position in two- or 

three-dimensional space. Although this definition includes 

the geostatistical spectrum, SDA here refers to techniques 

up to, but not including, geostatistics in level of complexity and sophistication. Utilization of both EDA and SDA 

requires a practical knowledge of such classical-statistical 

concepts as the four moments (mean, variance, skewness, 

and kurtosis), confidence intervals, the common statistical tests (t, z, f, and v2

), and calculation of correlation and 

covariance. Crow and others (1960), and Hammond and 

Householder (1962), and Krumbein and Graybill (1965) 

present a thorough introduction to background concepts.

EDA can be an approximate and simplistic method 

of data evaluation, and so it lends itself to gravel mining. 

Preliminary research on the practical applications of EDA 

in the aggregate sector has been promising. Bliss (1999) 

suggests that four types of model are used in regional 

aggregate assessment: descriptive models, distribution 

models, relational models, and spatial models. He derives 

a linear regression equation, by way of EDA, for forecasting freeze-thaw durability from the absorption values of 

Lower Paleozoic carbonate rocks mined for crushed aggregate from quarries in Oklahoma. He uses a “relational 

model,” whereby such statistical relations as regression 

are established among geologic variables to other variables, but with no consideration of spatial distribution. The 

results, however, are impressive because the regression 

equation explains 95 percent of the observed variation.

Bliss and Bolm (2001) provide a critical link between 

the use of EDA methodologies for practical data correlation (that is, predicting freeze-thaw durability to save on 

laboratory costs) and SDA modeling of aggregate-quality 

parameters within a spatial context. The two objectives of 

their study were to demonstrate how statistical methods can 

be used to examine the relations between geologic classification of sand-and-gravel deposits and material characteristics 

(such as Los Angeles wear and fines content). Bliss and Bolm 

(2001, p. 196) note that “the nature of undeveloped sand 

and gravel resources” can be forecasted by models prepared 

from statistical analysis. Their research was on the sedimentary deposits of Pleistocene Lake Bonneville, an important 

source of sand and gravel suitable for aggregate and construction in Utah. Lake Bonneville, which was comparable 

in size to Lake Michigan, straddles the boundary between 

the Basin and Range Province and the Wasatch Domain.
5

The result of Bliss and Bolm (2001)’s study was an 

application of EDA within a spatial context. The subsequent 

methodologies could be practically implemented, particularly in consideration of the interconnectedness of deposit 

exploration and land-use planning. Two classes of results 

were obtained. First, several trends were apparent, such as 

that sand and gravel in younger shorelines are slightly more 

durable and the deposits considerably larger in volume, and 

that the northern part of the Bonneville Basin contains slightly 

more durable sand and gravel than the southern part. Second, 

according to Bliss and Bolm (2001, p. 195), some results 

could have “immediate economic significance,” such as that 

the median sand-and-gravel deposit in the Wasatch Domain 

is 3 times larger than that in the Basin and Range Province.

As already mentioned, SDA encompasses a broad range 

of methodologies that evaluate data in consideration of their 

position in two- or three-dimensional space, such as trend-surface analysis (TSA), inverse-distance weighting, and geometric 

methods. Simply stated, TSA emphasizes the main features 

of a spatial distribution by smoothing over some of the local 

irregularities so that important trends can be isolated from the 

background more clearly. TSA mathematically separates spatial 

data into two components: regional and local (Davis, 1973).

Mineral-resource estimations performed manually on 

plan maps and cross sections typically belong to a class of 

“geometric methods” that are based on geometric weighting 

of assays (or other material properties). One of these methods is polygonal estimation, in which each assay grade for a 

given sampling unit is assigned its own polygon of influence 

(Noble, 1992). Polygons are drawn on maps on the basis of 

perpendicular bisectors of the line between each sample locality. Tonnage for each polygon is estimated by multiplying 

the area of each polygon by the thickness of the deposit, then 

dividing by the tonnage factor. “Nearest-neighbor estimation” is a computer approximation of polygonal estimation.

“Inverse-distance weighting” and “kriging” are methods of “moving-average estimation,” in which the basic 

procedure involves dividing the ore body into a matrix of 

rectangular blocks (taking into account any geologic controls that may affect grade assignment) and estimating 

the grade of each block by computing a weighted average 

based on the grades of surrounding blocks (Noble, 1992). 

Practical considerations involve determination of block 

size, anisotropies, and sample-selection criteria. Inversedistance weighting is based on the empirical observation 

that the weight of each sample from which the weightedaverage equation is derived is proportional to an inverse 

power of the distance from the location of the estimate 

to the sample. Kriging is described in the next section.

Cressie (1991) provides a useful survey of SDA 

methodologies, with an emphasis on the theoretical background rather than on concise procedures for field practitioners. He discusses such less common SDA topics 

as “lattice models,” “infill asymptotics,” and “spectral representations,” as well as the theoretical background behind such specialized geostatistical subjects as 

“Bayesian kriging” and “kriging-variance stability.”

Basic Concepts of Geostatistics

The term “geostatistics” refers to a collection of statistical and probabilistic methods that were developed in the 

South African mining industry in the 1950s (Krige, 1951) 

from forestry analysis in Scandinavia, although comparable 

mathematical concepts were being formulated by the Soviets with respect to the study of meteorologic fields (Gandin, 

1965). A theoretical framework for geostatistics was established by Matheron (1963) and his colleagues, termed by 

them “regionalized-variable theory.” The practice of geostatistics has been standardized in the mining of metalliferous ores since at least the 1970s, and has gained acceptance 

in such fields as hydrology, landscape ecology, and forest 

biometrics. Simply stated, geostatistics is used to describe 

and utilize correlations among spatially distributed data and 

to build models that can be adapted to practical situations. A 

geostatistical evaluation can drastically change the approach 

to the development of sampling methodology because spatial 

characterization is the intent, counter to traditional sampling, where the intend is to eliminate spatial dependencies. I 

avoid detailed mathematical and statistical discussions here, 

and the reader is referred to several excellent references, 

such as the reports by David (1977), Journel and Huijbregts 

(1978), Rendu (1978), Clark (1979), and Isaaks and Srivastava (1989). At this point, I am content to introduce general 

concepts and terminology derived from these references.

 A typical geostatistical analysis can be summarized in a 

few concise steps:

• Initial examination of the dataset, such as calculation 

 of basic statistics, production of simple contour maps, 

 and, most importantly, data validation

• Fitting of semivariogram models to the sample data, 

 allowing spatial and directional characterization of 

 variables

• On the basis of model fit, determination of such 

 semivariogram parameters as range, sill, and nugget

• Estimation

• Simulation, if desired

 I define a few necessary concepts and terms as follows. 

The term “semivariogram” is a fundamental expression of 

the spatiality in a dataset. Simply stated, a semivariogram 

aids in measuring or predicting how a given value (grain 

size, ore grade, and so on) will change in dependence on its 

location in two- or three-dimensional space. A semivariogram 

is a mathematical function that expresses the variation in 

a measured or predicted value among sample points separated by a distance h. Thus, a “directional” semivariogram 

is the function along a given azimuth in two or three dimensions, and an “omnidirectional” semivariogram does not 

take direction into account. If a dataset demonstrates different spatiality (expressed by the semivariogram parameters 

“sill” and “range”) in different directions, then the dataset 

is “anisotropic.” If no anisotropy is detected, then an omnidirectional semivariogram can be used. If the variable is 

nonadditive, then application of a suitable transformation 

Basic Concepts of Geostatistics
6 Contributions to Industrial-Minerals Research

before computing variograms, though not always possible, 

may be appropriate. Various transformations are available, of 

which the log transformation is the most commonly applied.

The term “nugget” (C0, fig. 1), which appears in the semivariogram as an apparent discontinuity at the origin, is an 

expression of either pure random variation or sampling error, 

or, more likely, both. The term “range” (a) identifies the distance h at which the semivariogram in a spherical model 

reaches a constant value, known as the “sill.” Intuitively, h is 

the maximum distance at which samples can be spatially correlated, beyond which no spatial correlation exists. The sill 

(C0+C1), which is the value of the semivariogram at a distance 

h beyond the range, represents the ordinary sample variance. 

These parameters are illustrated in figure 1.

Rather than trying to fit an equation to every possible 

form of semivariogram curve, several different mathematical 

models are used, most commonly “linear,” “spherical,” “exponential,” and “Gaussian.” (See the above-cited references for 

details.) The basic forms of these model types are illustrated in 

figure 2.

According to Clark (1979, p. 120), anisotropy is “perhaps the easiest ‘problem’ to tackle.” Basically, anisotropy 

represents variations in behavior or characteristics in different directions in two or three dimensions. In geostatistical 

investigations, the most common manifestation of anisotropy 

is a variation in the distance h in different directions with no 

change in the sill, which is a “geometric” anisotropy (fig. 3). 

A “zonal” anisotropy is one in which the sill changes but the 

variation in the distance h remains constant. Anisotropies are 

“corrected” in computations by means of various transforms, 

as discussed in detail throughout the above-cited references. In 

alluvial deposits, anisotropic characteristics of granulometric 

indices can suggest variations in depositional mechanics. For 

example, such indices as gradation coefficient evaluated for 

different strata (representing different flood events) in a catastrophic-flood deposit may exhibit varying degrees of anisotropy.

The issue concerning “stationarity” is another important consideration in variography and kriging. Without going 

into mathematical detail, if data are stationary, they have no 

definable trend. In the simplest terms, a sloping surface is not 

stationary. “Drift,” a simple polynomial function, is used to 

temporarily impose stationarity on data for the purposes of 

kriging. Even if data are nonstationary, they may be considered locally stationary (that is, to exhibit no definable trend) 

over a shorter distance, called the limit of quasi-stationarity 

(fig. 3). Kriging can be performed over this limit, rather than 

relying on a drift function.

The concept of “nonadditive variables” is an important 

consideration in the application of geostatistics to granulometric data and, potentially, with other aggregate parameters. 

An “additive model” is one in which the input variables have 

���������������

���

���

���

���

���

���

��� ��� ��� ���

�������������������

���������

����

�������������

������������

Figure 1. Basic features of the semivariogram.
7

�������� ��������

��������

������������

������������

������������

���������

��������

�����������

Figure 2. Spherical, exponential, and Gaussian semivariogram model forms.

��������

�������� ��������

������������

������������

������������

���������������������������

�������������������� ����������������

Figure 3. Semivariogram models demonstrating anisotropy and quasi-stationarity.

an additive effect on the output variables. As a simple example, if the input variable A has an effect of magnitude a on 

the output, and the variable B has an effect of magnitude b on 

the same output, then the combined effect of A and B would 

be a+b−a linear relation. Otherwise stated, the equation is 

additive because the contributions to the result are added 

together. This relation is not true, however, in a “nonadditive model,” where the combined effect of A and B cannot be 

expressed linearly, such as would be the case with, for example, 4A2

+2B3. Many variables used in the evaluation of metalliferous ores, such as troy ounces of Au per ton or weight 

percent Cu, are additive, whereas most granulometric indices, such as sorting or gradation coefficient, are nonadditive. 

(See section below entitled “Limitations of Geostatistics.”)

The concepts of “support” and the “volume/variance 

relation” are important in geostatistical analyses because 

the size (or support) of the dataset and the distribution of 

sample values are related. Intuitively, to use an example 

Basic Concepts of Geostatistics
8 Contributions to Industrial-Minerals Research

from Isaaks and Srivastava (1989), very small sample sizes, 

such as rock chips, would vary considerably if one nugget 

were pure gold and another were barren, whereas truckloadsize samples from the same deposit would vary much less. 

This concept is known as the volume/variance relation.

Evaluation and modeling of data under the assumption that the sample values are located at discrete points, 

rather than being part of a volume, such as a drill core, is 

termed a “point support.” On a point support, the sill ideally 

equals the variance. If the samples are not points, however, 

but drill cores of a given length, some of the point variation will be “smoothed out,” as in the previous example, 

where numerous individual values are effectively replaced 

with one average. Changing the semivariogram parameters 

in consideration of different support is termed “regularization” (Clark, 1979). During kriging, support volume is a 

consideration during estimation and “block modeling.”

The geostatistical procedure of estimating points or 

blocks (local versus global) is termed “kriging,” the designation of a set of linear-regression routines that minimize estimation of variance from a predefined covariance model (such as 

that represented by the semivariogram). “Kriging variance,” 

the variation in estimation error, is not data dependent but is 

determined by the semivariogram, the original sampling pattern, and the location of the point to be estimated relative to 

the sample locality. Estimates derived from kriging depend on 

the support of the estimated block, as well as on the sampling 

data.

The basic procedure of “ordinary kriging” in terms 

of “block modeling” can be summarized as follows:

• A search is made around the block to be estimated. 

 Sam ples located within the search window are used for 

 estimation of the block, and samples outside this window 

 are omitted.

• Samples within the search window are assigned weights 

 that reflect spatial characteristics, as manifested in the 

 semivariograms. In ordinary kriging, these weights sum 

 to unity.

• A weighted average is calculated to yield the block esti- 

 mate.

Numerous kriging techniques have been developed, 

depending on such factors as degree of stationarity, data density, and ultimate use of the estimates; common types are 

“ordinary kriging,” “simple kriging,” “universal kriging,” and 

“indicator kriging.” Simple kriging differs from ordinary kriging in that the sample weights do not sum to unity, using the 

average of the entire dataset rather than a local average and 

resulting in a less accurate but “smoother”-looking result. Both 

ordinary and simple kriging require stationarity or, at least, a 

practical limit of quasi-stationarity. Universal kriging can be 

used for nonstationary datasets by using the residuals (difference between the drift and the actual values of scatterpoints). 

Indicator kriging, which deals with threshold values, is basically a nonparametric counterpart to ordinary kriging. “Disjunctive kriging,” also known as full-indicator cokriging, is the 

cokriging of multiple indicators, beneficial over standard indicator kriging because an optimal estimation can be reached 

whereby each indicator profile is consistent with every other. 

“Multi-Gaussian kriging,” a generalization of lognormal kriging, is the application of kriging to the normal-score transformations of the input sample points. Thorough discussions of 

the numerous approaches to kriging are presented in the abovecited references.

Local and global estimates by way of kriging, which 

can be insufficient at the planning stage, can be beneficial to 

predict the variations in characteristics during mining and 

stockpiling. This estimation is the purpose of simulations, a 

field that has rapidly gained importance along with comparable increases in computing power. Simply stated, simulated 

deposits have the same values as those at the experimental 

locations and the same dispersion characteristics (up to order 

2) of the real deposit. The difference is that the objective of 

estimation is to provide, at each point, an estimate that is as 

close as possible to the true unknown grade. The simulation 

(or, better yet, “conditional simulation”) has the same two 

experimental moments (mean and covariance or variogram), 

as well as the same histogram, as the real grades—that is, 

it identifies the main dispersion characteristics of these true 

grades. In any single realization at each point, however, the 

simulated value is not the best possible estimate. “Conditioning” occurs when the simulations that are chosen are 

those that match the values at the experimental locations. 

(See Journel and Huijbregts, 1978, for additional details.)

Existing Geostatistical Investigations

Formal geostatistical investigations of sand-and-gravel 

deposits are quite rare, and the primary focus of completed 

studies has been on the spatial characterization of deposit 

thickness and its subsequent effect on ore reserves. Royle and 

Hosgit (1974) present a thorough treatment of the geostatistical evaluation of deposit thickness; their report generated 

some important preliminary results—in particular, that anisotropy and considerable nonstationarity should be expected in a 

deposit formed from multiple flood events.

Arthur (1994) investigated aggregate exploration by RMC 

Ltd. in an area of central Essex, U.K., that had been extensively worked for gravel. Objectives of his study included 

exploration and characterization of sand-and-gravel resources, 

covering such topics as optimum sampling pattern and density, 

errors inherent in different drilling and sampling techniques, 

geologic variation and its subsequent effect on plant feed, and 

the feasibility of formal geostatistical investigations. As stated 

above, such concepts are rudimentary in the evaluation of metalliferous resources, but their adoption in the aggregate sector 

has been slow to nonexistent.

Like Royle and Hosgit (1974), Arthur (1994) also examines the geostatistical characterization of deposit thickness, 

with applications toward computer-based extraction-volume 

calculation. A convincing example of infill drilling, using a 

25- versus a 100-m grid, is an important component of his 

examination. Arthur’s study also was somewhat successful 
9

in geostatistically analyzing the continuity of a deposit in the 

horizontal and vertical directions, to influence borehole spacing for the next round of drilling.

Arthur (1994) proceeds to examine the potential use of 

geos`Because aggregate-mining companies tend to report sizegrading data in terms of cumulative percent passing (owing 

to easier correlation to commercial screen sizes), most of the 

granulometric data used by Arthur (1994) were generally 

expressed in this format; however, the basic statistics of 

individual size fractions were calculated for percent-retained 

and percent-passing data. When percent-passing data were 

evaluated, the coarse fraction (20–75 mm) exhibited an inverse 

lognormal distribution, and only the medium (10, 14, and 20 

mm) fractions exhibited a normal population distribution, 

whereas when percent-retained data were evaluated, coarseand fine-fraction end members (28–75 mm and <75–300 

μm, respectively) exhibited lognormal distributions, and size 

fractions in the midrange (0.6–20 mm) exhibited normal 

distributions.

The most notable result of this analysis, in terms of 

further statistical and probabilistic examination, is the difference between the use of percent-passing and percent-retained 

data. For example, f and t tests comparing drilling data from 

the 25-m grid with those data from the 100-m grid showed 

a significant similarity for most size fractions when percentretained data were used; however, the same tests conducted by 

using percent-passing data showed that only the midrange size 

fraction <1.18 mm exhibited a significant similarity. Furthermore, correlations between individual size fractions can be 

derived only from percent-retained data. For example, significant negative correlations were measured between fractions in 

the gravel range (5, 10, 14, and 20 mm) and in the fine-sand 

range (75, 150, 212, and 300 μm), suggesting a subdivision 

of the study area into zones characterized by coarser and finer 

material and reinforcing results suggested by contouring.

Geostatistical analysis consisted of producing semivariograms on varying sample supports, borehole spacing, and 

locations within the study area. Drilling data from both the 25- 

and 100-m grids were used. Notable results are summarized as 

follows:

• Bench composites of 1 m produced from raw data did

 not lead to a significant increase in the quality of semivar- 

 iograms.

• Most semivariograms exhibited either a pure-nugget 

 effect or so much random behavior that no sill was defin- 

 able. Directional plots were mostly noise, owing to the 

 absence of usable sample pairs.

• Ranges modeled on the 100-m grid were from 400 to 800 

 m, whereas those modeled on the 25-m grid averaged 

 about 125 m.

• Semivariograms for size fractions above 14 mm exhibited 

 the pure-nugget effect for both datasets, and those for size 

 fractions between 75 and 300 μm required logarithmic 

 transformation before modeling.

The overall average correlation from a crossvalidation 

exercise performed on normally distributed midrange size 

fractions was only 0.52, suggesting that the model parameters do not accurately represent the original data.

Arthur (1994, p. 7–12) concludes that “it is pointless to carry on using this data to produce a geostatistical 

reserve estimate.” He emphasizes several considerations 

specific to the evaluation of sand-and-gravel deposits:

• Owing to the nature of deposition, most sand-and-

 gravel bodies will exhibit a proliferatio of small-scale 

 structures that coarse drilling grids may fail to detect, 

 including small channels, clay and silt pockets, and 

 organic inclusions. Changes in grain size from fine 

 sand to coarse gravel can occur over distances as small 

 as 0.3 m.

• A considerable drop in accuracy occurs for the 

 estimation of end-member (coarsand fine) fractions 

 from the shell-and-auger drilling method. In particular, 

 loss of fines can lead to underestimation by as much 

 as 75 percent. Therefore, grain-size distribution could 

 not be accurately modeled because actual cutpoints are 

 unknown. The shell-and-auger drilling method, however, 

 produced good depth estimates necessary for geostatisti 

 cal evaluation of horizon thickness.

Additional datasets based on reverse-circulation drilling were 

also suitable for quantifying vertical variations in the deposit, 

but much coarse material was lost.

Further Research Based on Existing 

Geostatistical Investigations

From the above discussions, it may be obvious that Arthur 

(1994) attempted several investigations on topics necessary 

needed for the further study of spatial-data applications in the 

aggregate sector, implying a twofold significance for his study: 

(1) although not all of his investigations were notably successful, the limitations of the methodology suggested, sometimes 

clearly, that common drilling methods jeopardize the accuracy 

of the size-distribution curve; and (2) no other researcher has 

seriously attempted that subject matter, at least in such detail.

A few more recent studies have investigated specific 

aspects of geostatistical analysis, with applications toward 

production accounting and land-use planning. Hack (2001) 

developed a correlation between δK" (cumulative percent 

passing through a K-in. mesh) and inplace density. His intent 

was to utilize changes in δK" (an index measure routinely 

measured by quality-control staff) to track changes within 

inplace density, necessary for proper accounting and verification of royalty payments. The logarithmically transformed 

semivariogram for density fitted a Gaussian model and was 

generally isotropic, except for a significant zonal anisotropy in the northwest/southeast (azm 135°±030°) direction. 

The semivariogram for δK" exhibited significant zonal and 

geometric anisotropies on azimuth 115°±030°. Following 

variography, ordinary kriging on a point support was used to 

estimate point values of δK" in several places where density 

samples had been taken. Kriging variance was approximately 

Existing Geostatistical Investigations
10 Contributions to Industrial-Minerals Research

1.2 percent squared, on a mean δK’’ of 36 percent. A regression equation of density on δK’’ was then derived, with the 

regression on δK’’ accounting for 80 percent of the observed 

variation in density (fig. 4). At the practical level, however, 

this relation should be applied with some degree of caution 

because of the small size of the datasets. Uncertainty in the 

final extracted quantity, however, is subject to other factors, such as vertical tolerances in topographic surveying.

SDA of aggregate-resource characteristics within the 

context of mine and development planning is the subject of 

Hack’s (2002a) report. This component of ongoing research 

involved Tenmile Terrace, the oldest major Pleistocene terrace 

in the lower Boise River Valley, Idaho. Data were collected 

by various parties, usually the Idaho Transportation Department, from 1953 to 1999. The terrace, whose largest outcrop 

is 19 km long by 5 km wide and generally about 10 m thick, 

is divided into two components: the main body designated as 

the “upper” terrace and the scattered remnants designated as 

the “lower” terrace. A total of 10 composite samples from the 

upper terrace, each representing an individual gravel pit or 

prospect, were derived from 46 individual samples, and a total 

of 13 composite samples for the lower terrace were derived 

from 39 individual samples. The composite samples were 

created to assess spatial variations across the entire Tenmile 

Terrace. Cumulative percent-passing data were entered from 

laboratory sheets into the software program SGMatrix, and 

then various sedimentologic, engineering, and concrete-industry factors, such as uniformity coefficient, gradation coefficient, sorting, and dust ratio, were calculated from cutpoints 

(D10, D25, and so forth). All the various factors were correlated with each other, as well as with (x, y) sample locations, 

to create a “matrix” of correlation coefficients (table 1).

The results were used in the preliminary identification 

of spatial structures. Long-range resource evaluation in the 

lower terrace indicated a strong negative correlation of easting 

and northing with cutpoints D50 and D75. Otherwise stated, the 

median and 75th percentile decrease as sampling progresses 

northeastward, reinforcing the previous geologic interpretation that deposition was from southwestward-flowing streams, 

because coarser material would have been deposited at the 

headwaters, with finer sediment settling out downstream. The 

intent of short-range resource evaluation was to appraise the 

spatiality of the individual sites (active pits, dormant pits, or 

prospects), with potential applications toward mine planning, 

as with Arthur (1994). Analysis at this scale suggested results 

of sedimentologic interest that were beyond the scope of 

Hack’s (2002a) study, such as consistent indications that sorting characteristics are related to the proportion of sand. Within 

a spatial context, cutpoints in the range D60−D75 exhibited 

strong correlations to northing and easting. Limited attempts at 

exploratory variography had mixed results, although midrange 

cutpoints generally exhibited linearity. Semivariograms for 

cutpoints D50 and D75 were spherical in form, with ranges typically 3 to 5 times those of sorting and gradation parameters.

Again using Tenmile Terrace as a test case, D.R. Hack 

(unpub. data, 2003) broadens the scope of long-range resource 

evaluations, emphasizing regional trends relevant to municipal land use and resource planning. The project area straddled 

two counties, Ada County and Canyon County, in southeastern 

Idaho. The terraces were generally split along county lines, with 

most of the lower terrace in Canyon County and most of the 

upper terrace in Ada County. In both counties (terraces), sampling tended to cluster in two general areas, about 20 km apart. 

Before the calculation of spatial correlations, sampling localities 

that potentially represented “anomalies” could be identified.

This approach highlighted significant differences 

between the two general areas. The samples from Canyon 

County (lower terrace) exhibited spatial correlations that 

could aid in deposit characterization. The above-noted correlation of cutpoint D75 with the southwest direction (≈0.70) 

decreased with a virtually linear trend through cutpoints down 

to D10, where the correlation was ≈0.20. Size fractions in 

the mesh-size range #10 (2,000 μm) to #30 (600 μm), as well 

as skewness and kurtosis, exhibited moderate correlations 

(≈0.60) with the northeast direction (180° from previously 

noted correlation). Therefore, skewness and kurtosis exhibit 

strong correlations with mesh sizes #10 and #30 (=0.60–

0.80), possibly providing insight into the sedimentologic 

environment. Also of potential sedimentologic interest, correlations of uniformity coefficient with mesh sizes #100 (150 

μm) and #200 (75 μm) were =0.68 and =0.85, respectively, 

and correlations of dust ratio with mesh sizes #50 (300 μm) 

and #100 were ≈−0.70.

The samples from Ada County (upper terrace) samples, 

evaluated as a whole, exhibited no notable spatial correlations 

(at least of >0.50): however, similarities in sedimentologic 

environment were noted between the two terraces. Correlations of uniformity coefficient with mesh sizes #100 and 

#200 were even stronger than in Canyon County. Dust ratio 

was also negatively correlated with mesh sizes #50 and #100, 

but not to so high a degree. In Ada County, however, sorting exhibited strong correlations with mesh sizes #10 and #30 

(=0.88 and =0.70, respectively).

Resource exploration or planning can benefit from such 

interpretations if granulometric data are available from prospects or operating pits in the area of interest. Using this result 

from Ada County, Idaho, an example, prospects that exhibit 

Figure 4. Log-transformed semivariogram of density (from Hack, 2001). higher proportions of mesh sizes #10 and #30 material will 

������������

�������������

���������������

�������������
11

statistically exhibit a higher level of sorting, a potentially 

important consideration when evaluating sites for concrete 

aggregate. Likewise, when evaluating the Canyon County sites 

within a spatial context, higher cutpoints (such as D75) are significantly correlated with the southwest direction, suggesting 

to a prospector that sites toward the southwest may be more 

likely to exhibit larger gravel (cobbles) at the coarse tail of the 

particle-size distribution.

Limitations of Geostatistics

A major obstacle in the adaptation of formal geostatistical principles to practical applications in aggregate 

exploration, evaluation, and mining is the rigorous theoretical standards to which the term “geostatistics” is normally 

applied. As noted above, however, data can be evaluated 

within a spatial context in the SDA spectrum, allowing 

acknowledgment of spatial-data structure without violating the academic constraints under which geostatistics 

operates, possibly providing the most practical route for 

workers in the aggregate sector. Nonetheless, geostatistical evaluation remains a reasonable goal to strive for, 

although its limitations must be confronted at the outset.

The most productive route to quantifying the role of 

geostatistical principles in aggregate-resource evaluation 

is to identify those aspects of geostatistical analysis that 

are taken for granted in the evaluation of metalliferous 

ores but would be challenging to apply to the sampling 

and modeling of a gravel pit. Most notable among these 

aspects are the additivity of the variables involved, the 

use of a uniform sample size (such as a drill-core barrel) 

in exploration, and the clear definition of a “selective 

mining unit.” Each of these terms is discussed below.

Nonadditive Variables

As suggested above, nonadditive variables are an 

important consideration in applications of geostatistics to 

granulometric data and other parameters used in the aggregate sector. Some of these variables are explicitly nonadditive, by function of their mathematical form—for example, 

the sedimentologic statistical moment identified as “skewness,” expressed as √[(D25×D75)/(D50)2

], or the United Soil 

Classification System soils-engineering parameter “gradation coefficient,” expressed as Cg=(D30)2

/(D60×D10). Other 

variables are not obviously nonadditive by their mathematical form but would be within the context of geostatistical evaluation and a subsequent block model or simulation—for example, such quality measures as Los Angeles 

wear and freeze-thaw durability. Not only could creation 

of suitable mathematical transforms be challenging, but 

also ensuring that the applications of these transforms are 

suitable with respect to other mathematical limitations.

Peroni and others (1999) provide a clever treatment of the geostatistical evaluation of nonadditive 

variables in a kaolinite deposit. Specifically, they use 

the Kubelka-Munk function to create a nonlinear transformation of the parameter brightness, a primary quality-control parameter in the evaluation of kaolin deposits. Its use avoids bias in block modeling and reserve 

estimation, and crossvalidation results are excellent.

Uniform Sample Size

Most practitioners of geostatistics in the metalliferous- or 

coal-mining industries take for granted the use of a uniform 

sample size, such as a drill-core barrel, during exploration and 

deposit evaluation. Likewise, most geostatistics textbooks, software, and ensuing methodologies are written with this aspect 

of geostatistical analysis in mind, because it applies to most 

conditions. In a gravel deposit, however, a uniform drill-core 

barrel is not practicable because the material is unconsolidated. 

Furthermore, borings are not the fundamental deposit-evaluation tool, for obvious reasons (see discussion above pertaining 

to Arthur, 1994). For many mine sites (if adequate sampling 

even exists), deposit evaluation will be a combination of 

boreholes, water-well logs, test pits, trenches, stratigraphic 

columns from exposures, and, possibly, some geophysical 

���� ����� ����� ��� ��� ��� ��� ��� ����

���� ������� ��� ����� ����� ������ ������ ������ ������ ������ ������

����������� ���� ��� ���� ����� ����� ����� ����� ����� ����

����� ������ ���� ���� ��� ����� ���� ����� ����� ���� �����

������������ ����� ����� ����� ��� ����� ���� ���� ����� ����

������������ ����� ����� ���� ����� ��� ���� ���� ���� �����

������������ ����� ����� ����� ���� ���� ��� ���� ���� �����

������������ ����� ����� ����� ���� ���� ���� ��� ���� ����

������������ ����� ����� ���� ����� ���� ���� ���� ��� �����

������������ ����� ���� ����� ���� ����� ����� ���� ����� ���

Table 1. Correlation coefficients for the lower Tenmile Terrace, Idaho.

[Bold numbers are diagonal elements. From Hack (2002a)].

Limitations of Geostatistics
12 Contributions to Industrial-Minerals Research

cross sections. A two-sided approach is needed to reconcile this 

issue: developing academically suitable methods to use such 

varied data (in widely varying three-dimensional geometric 

forms) in a proper geostatistical evaluation, and convincing 

those who prepare and utilize the end result that this diminished consistency in data types and methodology is acceptable, and that the results are still usable at a practical scale.

Support and Selective-Mining-Unit Issues

As implied by the preceding discussion, a nonuniform 

sample size (or, as in a trench or stratigraphic column, any 

definable sample size at all) will create subsequent complexities in establishing the volume/variance relation. This 

problem is not one sided, and it is complicated by the mode 

of operation in most aggregate mines. A uniform mining 

block, such as is known in the mining of metalliferous ores, 

is unheard of in a sand-and-gravel pit. The closest comparison would be to a standard bench height, based on equipment 

selection or regulatory requirements. Mine planning, at best, 

will be at the level of benches. However, the goal of geostatistical applications, or even SDA, is not necessarily the 

generation of a “block model” for a gravel pit that will not 

be mined in uniform blocks. Instead, the goal may to identify 

general trends, zones, or the sedimentologic environment at 

the inpit level, or to identify overall trends and environments 

at the level of exploration or resource-management plans.

Other Considerations

Beyond these obvious complexities in the use of “classical” geostatistics for aggregate-resource evaluation, the 

fundamental expression of data must be confronted. For 

example, will the use of granulometric data, expressed by 

way of a size-distribution curve, be theoretically sound in 

terms of its use in a geostatistical framework? Part of this 

problem has been dealt with previously, such as by the 

transformation of nonadditive variables. However, the basis 

of a size-distribution curve—finite points on a graph—presents a spectrum of academic questions beyond the obvious 

problems within common geostatistical practice. Initially, 

size fractions (expressed here as mesh sizes #4, #200, and 

so on) are only finite points on a distribution curve, and 

when they are evaluated as cumulative percent passing, 

they are seen to be interdependent—if one size fraction is 

removed, all the numbers change—and so they cannot be 

evaluated as independent variates. Cutpoints (expressed 

here as D25, D75, and so on) and the numerous service 

variables derived from them) are the converse side of this 

relation. The complexities inherent in the use of cutpoints 

are well examined by Arthur (1994). In the broader sense, 

the academic implications of SDA of particle-size data have 

been little studied, and skepticism could be considerable.

As suggested above, Arthur (1994) concludes that it is 

good to use percent retained when attempting to quantify 

the spatial distribution of individual size fractions, but it is 

better to use percent passing because it can provide information on the requirements of processing screens or classifiers. 

From a geostatistical point of view, percent-retained data 

can be more suitable because they generally conform to 

either a normal or lognormal distribution, whereas percent-passing data commonly exhibit an inverse lognormal 

distribution, particularly for the coarser size fractions.

How to Practically Confront These 

Limitations

General Concepts

Once the theoretical limitations of geostatistical methods 

in evaluating sand-and-gravel deposits have been identified, 

the next step is to clearly state how geostatistics in practical 

terms can be performed. With respect to the preceding discussion, the following three points must be acknowledged.

The practical requirements of the aggregate industry.—

Aggregate operators could benefit from deposit modeling, 

resource and mine planning, and computer-based scheduling. Expensive solutions at a level more complex than their 

operations will do them no good and will be met with immediate opposition.

The limitations of the data collected by the aggregate 

sector, as a function of what the mining industry requires 

at the practical level.—Aggregate operators do not mine 

their gravel pits on the basis of uniform mining blocks, 

and so they do not require intensive (and expensive) drilling campaigns that will result in block models which are 

irrelevant to their operating standards and likely unnecessary. However, aggregate operators have recently tended to 

collect more deposit data for requirements beyond deposit 

characterization, and useful data can be extracted from 

those examinations. For example, particle-size data obtained 

from a borehole drilled primarily to investigate aquifers, or 

from a test pit excavated to quantify overburden and topsoil 

thickness, can be used to aid in deposit characterization.

The possible requirement of modifications to the typical 

exploration and sampling protocol.—For formal geostatistics to be useful in aggregate-mining operations, sampling 

standards and protocols will need to be adapted toward this 

goal. Conversely, if geostatistical applications in aggregate 

are to be accepted by outside observers (skeptics), then sampling standards and protocols will need to be taken seriously.

If these requirements cannot be met—if the application 

of geostatistics to evaluating sand-and-gravel deposits is 

impractical—then one academic step backward in the SDA 

spectrum may need to be taken, such as earlier successful 

applications of EDA methodologies within a spatial context 

(for example, Bliss and Bolm, 2001). Indicator kriging is 

worth examining because the method is based on thresholds 

enabling evaluation of cutpoints. However, the successful 

demonstration of this methodology, as well as its ultimate 
13

application in the field, will require various robust datasets 

on which suitable crossvalidation examinations can be performed. Quantification of the variations of modeling requirements will require the compilation and analysis of data 

from the full spectrum of sand-and-gravel deposits (flood 

plain, glacial, catastrophic flood, and marine). If the chosen 

spatial methodology is initially successful, continues to 

validate the data, and yields beneficial results when implemented in the field, the concerns of skeptics can be reevaluated if shown to be irrelevant to actual field conditions.

Drilling Programs

With respect to drilling, Arthur (1994) proposes several 

approaches to improving the success of deposit exploration 

and characterization programs when geostatistical methods are 

being considered. Regarding phase I drilling: 

• Phase I drilling on a regular, wide-spaced coarse grid of 

 about 100 m will probably not define important character- 

 istics because gravel deposits are naturally highly diverse.

• Shell-and-auger drilling should be used in coarser, grav- 

 elly deposits, and reverse-circulation drilling in finer,

 sandy deposits.

• Drillholes should be sampled at 1-m intervals (or smaller 

 if suggested by lithology or type of sedimentary facies, or 

 both), and each sample analyzed for grading and petrogra- 

 phy.

• If the deposit is shallow, test pits should be used to direct- 

 ly assess the characteristics of the deposit, especially the 

 >75-mm size fraction. Additionally, a direct comparison 

 between drillhole samples and face samples could indi- 

 cate the type and extent of any sampling errors, such as 

 loss of fines by drilling. For example, holes could be 

 drilled behind the advancing face, and then adjacent sam- 

 ples collected when the face advances. This method, 

 however, could be prohibitively expensive.

• If the deposit is coarser, but the characteristics of the 

 fines are important to the final product, shell-and-auger 

 drilling should be supplemented with reverse-circulation 

 drilling. The opposite approach applies for sandier depos- 

 its containing coarse material of interest. 

 If these measures yield suitable geostatistical results for 

 phase I drilling, the following approaches should influ- 

 ence phase II drilling:

• Geostatistical ranges should influence drillhole spacing, 

 which would consist of at least three holes within the 

 range in each direction of the drilling grid. Small-scale 

 deposit-specific structures should be taken into account.

• Downhole semivariograms can aid in the identification of 

 horizontal structures, if used in conjunction with vertical- 

 grading profiles.

• Zones with different granulometric properties, as quan- 

 tified by SDA, can provide information on the final prod- 

 ucts that can be obtained from the deposit. This approach 

 is a first step toward mine planning and “block modeling” 

 of a sand-and-gravel operation.

Proposals for Further Work

To proceed beyond, and take best advantage of, the 

examinations that have been performed to date, I suggest 

the following steps:

1. Data from the full spectrum of sand-and-gravel deposits 

(flood plain, glacial, catastrophic flood, and marine) must be 

evaluated. Potentially, these data will present varying model 

types because they represent different depositional environments and mechanics.

2. Such scientists as sedimentologists must be integrated into 

this research. The geology and mechanics of gravel deposition differ from those of, for example, a porphyry copper ore 

body. Therefore, geostatistical evaluation of alluvial material will require far different academic considerations from 

those of metalliferous ore bodies. Pursuing rigorous SDA 

without understanding the science of orogenesis could lead to 

improper assumptions and ultimate failure.

3. Other professionals or scientists can assist mining engineers 

and geologists with insight into such conundrums as the 

question, because dense drilling grids are required to detect 

small-scale structures, how do we know where to expend 

more resources in looking for small-scale structures? From 

geophysics? From extrapolation of stratigraphic columns in 

available exposures?

4. With hindsight, Arthur (1994) proposed that serious exploration should have started with a comparative drilling 

exercise between shell-and-auger, reverse-circulation, and 

continuous-flight-auger methods. He then suggested that all 

1-m-interval samples should be sent for size-grading analysis, coupled with face sampling, to provide information on 

the precision and accuracy of each method and on the quality of size-grading data expectable from each method.

5. Evaluating data on a point-support basis can be suitable when 

the goal is to identify general spatial tendencies, trends, 

or basic zonation in terms of deposit characterization and 

reconnaissance-level resource planning. Beyond that goal, 

however, complexities with respect to the volume/variance 

relation will arise. Considerable theoretical and practical 

examinations will be needed to overcome this dilemma, and 

a mixture of theoretical and practical approaches is essential. 

A routine academic approach could be a formidable encounter. Suppose that the input data are from a few boreholes, 

test pits, and trenches; the output “spatial model” is, on 

the x, y plane, pit zones (sectors) and, in the z-direction, a 

bench height. The best approach may be to derive empirical 

methods and adapt them on the basis of validation results.

6. Even if support and volume/variance issues are reconciled, 

complexities from nonadditive variables may still arise. If 

the variables are explicitly nonadditive (such as skewness 

or gradation coefficient), mathematical transforms should 

suffice if given proper treatment. However, if the variables 

are quality based, such as Los Angeles wear, then such 

transforms as those used by Peroni and others (1999) will 

need to be created.

7. Additional sampling requires additional expenditure, which 

is likely the biggest stumbling block at present.

Proposals for Further Work
14 Contributions to Industrial-Minerals Research

References Cited

Arthur, John, 1994, Computerised sedimentological modelling 

and resource evaluation of sand and gravel deposits: Cardiff, 

U.K., University of Wales, Ph.D. thesis, 275 p.

Barksdale, R.D., ed., 1991, The aggregate handbook: Washington, D.C., National Stone Association, 670 p.

Bliss, J.D., 1999, Quantitative models for aggregate; some 

types and examples from Oklahoma carbonate rocks, in

Johnson, K.S., ed., Proceedings of the 34th Forum on the 

Geology of Industrial Minerals: Oklahoma Geological 

Survey Circular 102, p. 37–45.

Bliss, J.D., and Bolm, K.S., 2001, Statistical analysis of sand 

and gravel aggregate deposits of Late Pleistocene Lake 

Bonneville, Utah, in Bon, R.L., Riordan, R.F., and Krukowski, S.T., eds., Proceedings of the 35th Forum on the 

Geology of Industrial Minerals—the Intermountain West 

Forum, 1999: Utah Geological Survey Miscellaneous Publication 01–2, p. 195–214.

Bliss, J.D., Stanley, M.C., and Long, K.R., 2002, Role of megaquarries in future aggregate supply, in Proceedings of the 

53rd Annual Highway Geology Symposium: Sacramento, 

California Department of Transportation, p. 303–315.

Clark, Isobel, 1979, Practical geostatistics: London, Applied 

Science Publishers, 129 p.

Cressie, N.A.C., 1991, Statistics for spatial data: New York, 

John Wiley & Sons, 900 p.

Crow, E.L., Davis, F.A., and Maxfield, M.W., 1960, Statistics 

manual; with examples taken from ordnance development: 

New York, Dover Publications, 288 p.

David, Michel, 1977, Geostatistical ore reserve estimation 

(Developments in Geomathematics, v. 2): Amsterdam, Elsevier, 364 p.

Davis, J.C., 1973, Statistics and data analysis in geology, with

FORTRAN programs, by R.J. Sampson: New York, John 

Wiley and Sons, 550 p.

Gandin, L.S., 1965, Objective analysis of meteorological 

fields: Jerusalem, Israel Program for Scientific Translations, 

242 p.

Hack, D.R., 2001, Correlation of material quality and property 

variables in a gravel deposit, in Xie, Heping, Wang, Yuehan, 

and Jiang, Yaodong, eds., International Symposium on Computer Applications in the Mineral Industries, 29th, Beijing, 

2001, Proceedings: Rotterdam, Balkema, p. 31–34.

Hack, D.R., 2002a, A methodology for the geostatistical characterization and exploratory data analysis of sand and gravel 

deposits, in Bandopadhyay, Sukumar, ed., International 

Symposium on Application of Computers and Operations 

Research in the Mineral Industry, 30th, Phoenix, Ariz., 2002, 

Proceedings: Littleton, Colo., Society for Mining, Metallurgy, and Exploration, p. 49–56.

Hack, D.R., 2002b, Comparison of commercial mine planning 

packages available to the aggregate industry, in Bandopadhyay, Sukumar, ed., International Symposium on Application of Computers and Operations Research in the Mineral 

Industry, 30th, Phoenix, Ariz., 2002, Proceedings: Littleton, 

Colo., Society for Mining, Metallurgy, and Exploration, p. 

311–319.

Hammond, K.R., and Householder, J.E., 1962, Introduction to 

the statistical method; foundations and use in the behavioral 

sciences: New York, Alfred A. Knopf, 412 p.

Isaaks, E.H., and Srivastava, R.M., 1989, Applied geostatistics: 

New York, Oxford University Press, 561 p.

Journel, A.G., and Huijbregts, C.J., 1978, Mining geostatistics: 

New York, Academic Press, 600 p.

Kannan, G., Vorster, M., Martinez, J., and Rakes, T., 1999, 

Improving mine management through data mining, in

Dagdelen, Kadri, ed., International Symposium on Computer Applications in the Mineral Industries, 28th, Golden, 

Colo., Proceedings: Golden, Colorado School of Mines, p. 

637–645.

Krige, D.G., 1951, A statistical approach to some basic mine 

valuation problems on the Witwatersrand: Chemical Metallurgical and Mining Society of South Africa Journal, no. 52, 

p. 119–139.

Krumbein, W.C., and Graybill, F.A., 1965, Introduction to statistical models in geology: New York, McGraw-Hill, 475 p.

Lindsey, D.A., 1997, An introduction to sand and gravel 

deposit models, Front Range Urban Corridor: U.S. Geological Survey Open-File Report 97–81, 6 p.

Matheron, Georges, 1963, Principles of geostatistics: Economic 

Geology, v. 58, no. 8, p. 1246–1266.

Melick, Roger, Furman, T.L., Langer, W.H., and Van Brenk, 

W.H., 1998, Using geographic information systems with specific examples for the aggregate industry: U.S. Geological 

Survey Open-File Report 98–39, 12 p.

Noble, A.C., 1992, Ore reserve/resource estimation, in chap. 

5.6 of Hartman, H.L., ed., SME mining engineering handbook (2d ed.): Littleton, Colo., Society for Mining, Metallurgy, and Exploration, v. 1, p. 352–255.

Peroni, Rodrigo, Costa, J.F., Koppe, J.C., and Petter, C.O., 

1999, A novel methodology for modeling kaolian deposits, 

in Dagdelen, Kadri, ed., International Symposium on Computer Applications in the Minerals Industries, 28th, Golden, 

Colo., 1999, Proceedings: Golden, Colorado School of 

Mines, p. 545–553.

Pitard, F.F., 1993, Pierre Gy’s sampling theory and sampling 

practice: Boca Raton, Fla., CRC Press, 488 p.

Primel, Louis, and Tourenq, Claude, eds., 2000, Aggregates; 

geology, prospection, environment, testing, specifications, 

extraction, processing plants, equipments, quality control: 

Rotterdam, Balkema, 590 p.

Rendu, J.M., 1978, An introduction to geostatistical methods of 

mineral evaluation (Geostatistics Series, v. 2): Johannesburg, 

South African Institute of Mining and Metallurgy, 84 p.

Royle, A.G., and Hosgit, M.E., 1974, Local estimation of sand 

and gravel reserves by geostatistical methods: Institution 

of Mining and Metallurgy Transactions, sec. A, v. 83, p. 

A53–A61.

Smith, M.R., and Collis, L., eds., 1993, Aggregates; sand, 

gravel and crushed rock aggregates for construction purposes 

(2d ed.): London, Geological Society, 335 p.

Tukey, J.W., 1977, Exploratory data analysis: Reading, Mass., 

Addison-Wesley, 688 p.