–nitrogen in the Regional hybrid geospatial modeling of soil nitrate

+ MODEL
ARTICLE IN PRESS
Geoderma xx (2006) xxx – xxx
www.elsevier.com/locate/geoderma
Regional hybrid geospatial modeling of soil nitrate–nitrogen in the
Santa Fe River Watershed
S. Lamsal, S. Grunwald ⁎, G.L. Bruland, C.M. Bliss, N.B. Comerford
Soil and Water Science Department, IFAS, University of Florida, 2169 McCarty Hall; PO Box 110290, Gainesville FL, 32611, USA
Received 1 July 2005; received in revised form 29 November 2005; accepted 27 December 2005
Abstract
Typically, regional assessment of the spatial variability and distribution of environmental properties are constrained by sparse
field observations that are costly and labor intensive. We adopted a hybrid geospatial modeling approach that combined sparsely
measured soil NO3–N observations collected in three seasons (Sept. 2003, Jan. and May 2004) with dense auxiliary environmental
datasets to predict NO3–N within the Santa Fe River Watershed (3585 km2) in north-east Florida. Elevated nitrate–nitrogen
concentrations have been found in this watershed in spring, surface and ground water. We collected soil samples at four depths (0–
30, 30–60, 60–120, 120–180cm) based on a random-stratified sampling design. Classification and regression trees were used to
develop trend models for soil NO3–N predictions based on environmental correlation and predict values at the watershed scale.
Residuals were spatially autocorrelated only for the Jan. 2004 sampling and regression kriging was used to combine the kriged
residuals with tree-based trend estimates for this event. At each step of the upscaling process, error assessment provided important
information about the uncertainty of predictions, which was lowest for the Jan. sampling event. Sites that showed consistently high
NO3–N values throughout the cropping season (Jan–May 2004) with values ≤5 μg g− 1 covered 95.7% (3363.9 km2) of the
watershed. Values in the 5–10 μg g− 1 range covered 4.3% (150.7km2), while values exceeding 10μg g− 1 covered only 0.59%
(20.7 km2) of the watershed. Elevated soil NO3–N on karst, unconfined areas with sand-rich soils, or in close proximity to streams
and water bodies pose the greatest risk for accelerated nitrate leaching contributing to elevated nitrogen found in spring, surface and
ground water in the watershed. This approach is transferable to other land resource problems that require the upscaling of sparse
site-specific data to large watersheds.
© 2006 Elsevier B.V. All rights reserved.
Keywords: Nitrate–nitrogen; Geospatial; Classification and regression trees; Regression kriging
1. Introduction
Typically, regional assessment of the spatial variability and distribution of soil properties are constrained
by sparse field observations that are costly and labor
intensive. Therefore, hybrid modeling techniques have
⁎ Corresponding author. Tel.: +1 352 392 1951x204; fax: +1 352 392
3902.
E-mail address: [email protected] (S. Grunwald).
been proposed to combine sparsely sampled field data
with dense auxiliary environmental datasets to predict
soil properties at regional scale (Goovaerts, 1999;
McBratney et al., 2000; van Meirvenne and van
Cleemput, 2005). Commonly, regression modeling is
used to characterize relationships between environmental variables whereas geospatial analysis is focused on
the interpretation of soil-landscape variability. Least
squares linear regression and multiple regressions are
among the most commonly used analytical techniques
0016-7061/$ - see front matter © 2006 Elsevier B.V. All rights reserved.
doi:10.1016/j.geoderma.2005.12.009
GEODER-02541; No of Pages 15
ARTICLE IN PRESS
2
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
of ecologists and soil scientists (Trexler and Travis,
1993). Linear models are associated with assumptions
about the statistical distribution of a response variable,
the form of variance structure and the independence of
observations which are often difficult to meet with
environmental data (James et al., 2004). Furthermore,
environmental datasets are often complex and unbalanced with non-linear relationships that vary in space
and time. Therefore, linear regression models often fail
to effectively describe landscape pattern to address
environmental issues. Numerous studies employed
multivariate non-linear methods to interpret environmental behavior and model the spatial distribution of
environmental properties. Odeh et al. (1995), McBratney et al. (2000) and Hengl et al. (2004) used a suite of
different environmental variables such as soil, land use,
topographic and other landscape properties to predict
target soil variables.
These studies were rooted in the conceptual framework called SCORPAN described by McBratney et al.
(2003). Grunwald and Lamsal (2005) provided an
overview over mapping techniques to collect environmental data that can be included in a SCORPAN model.
McBratney et al. (2000) compared 24 different
univariate and multivariate statistical, geostatistical and
hybrid models at field, sub-catchment and regional
scales and found that hybrid techniques outperformed all
other methods at all scales. Hybrid models have the
ability to mix deterministic and stochastic components,
link data collected at different spatial scales, and
integrate sparse point datasets (e.g. soil samples) with
exhaustively measured environmental datasets (e.g.
digital elevation model—DEM, satellite imagery).
Bishop and McBratney (2001) argued that models that
incorporate secondary data (e.g. regression kriging—
RK) are superior in terms of their prediction performance
to generic geostatistical models such as ordinary kriging.
The performance of regression kriging was found to be
superior when compared to other prediction methods by
Odeh et al. (1995), Knotters et al. (1995), Odeh and
McBratney (2000), Hengl et al. (2004) and Triantafilis et
al. (2001) to predict a variety of different soil properties.
To address elevated nitrogen values found in the
Santa Fe River Watershed (SFRW) we present a case
study that uses land uses, soils, and landscape
characteristics in a multi-tier geospatial hybrid approach to model seasonal and regional patterns of soil
NO3–N. Increasing concentrations of NO3–N were
found in the ground, spring and surface water of the
Suwannee Basin in north-east Florida (Hornsby et al.,
2001). The SFRW is one of the tributaries of the
Suwannee River and has greatly contributed to the
NO3–N loads delivered by the Suwannee River to the
Gulf of Mexico (Hornsby et al., 2002). The watershed
covers 3585 km2 (13.8 %) of the Suwannee Basin but
contributed 20% of the total NO3–N loads that
drained from the Suwannee Basin into the Gulf of
Mexico accounting for about 2900 tons NO3–N
(Suwannee River Water Management District, 2003).
Our objectives were to upscale site-specific NO3–N
measurements to the watershed scale using a mixed
modeling approach based on environmental correlation
and geostatistical modeling.
2. Materials and methods
2.1. Study area
The SFRW, a tributary of the Suwannee Basin that
drains into the Gulf of Mexico, spans over an area of
approximately 3585 km2 across eight counties in northeast Florida (Fig. 1). The soils of the SFRW are
predominantly sandy in texture with loamy to clayey
deposits, organics and sites with sand hill karst terrain
with many solution basins. Ultisols with 36.7% aerial
coverage, Spodosols (25.8%), and Entisols (14.7%) are
the dominant soil orders in the watershed. Less
prominent are Histosols (2.0%), Inceptisols (1.1%),
and Alfisols (1.0%). Land use consists of pine plantation
(32.2%), wetlands (16.2%), upland forest (14.7%),
improved pasture (14.0%), urban (8.8%), forest regeneration (6.0%), crops (5.0%), rangeland (3.7%) and a
variety of high intensity land uses such as tree groves,
dairies, and feeding operations (SRWMD and SJWMD,
1995). Agricultural land uses in the watershed are
diverse, ranging from corn, peanuts, tobacco, vegetables, watermelons, strawberries, blueberries, and
pecans. The elevation ranges from around 3 m to over
91 m above mean sea level. Generally, the land is level
(0–2% slopes) to gently sloping and undulating (0–5%
slopes), with the major exception to this pattern being
the moderately and strongly sloping land (5–12%
slopes) along the Cody Scarp. Two main physiographic
regions in the watershed are the Gulf Coastal Lowlands
and the Northern Highlands, which are separated from
one another by the Cody Scarp (Schmidt, 1997).
Underlying geologic units include Eocene limestone
(which occurs near the ground surface in the highrecharge, strongly karst-influenced Gulf Coastal Lowlands), capped by Miocene sediments which tend to be
rather clayey and phosphatic (occurring at or near the
surface along the Cody Scarp), in turn capped by
Pliocene and Pleistocene–Holocene sediments which
tend to be sandy at the surface but having loamy subsoils
ARTICLE IN PRESS
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
3
Fig. 1. The Santa Fe River Watershed in north-east Florida. The enlarged watershed (right) shows the boundary of the eight counties with sampling
locations (May 2004 event).
or substrate at varying depths (Randazzo and Jones,
1997; Brown et al., 1990). Based on the National
Oceanic and Atmospheric Administration (NOAA)
records at seven monitoring stations across the watershed, the total precipitation was 680.7 mm in Sept. 2003,
274.3 mm in Jan. 2004 and 231.4 mm in May 2004
representing “wet/end of cropping season” (Sept.), “dry
winter season” (Jan.) and “dry spring season” (May),
respectively. The mean annual precipitation based on
30 years of records from 1971 to 2000 across the
watershed was 1334 mm and the mean annual temperature 20.4 °C (National Climatic Data Center, NOAA).
selected privately owned land. Therefore, the number
of samples collected varied slightly by season and depth
(Table 1). The samples were extracted with 2 M KCl
(Keeney and Nelson, 1982) and the extraction solution
was analyzed for NO3–N content and expressed in μg
g− 1 of dry soil. Detection limits were 0.05 ppm NO3–N.
Samples below the detection limit were assigned the
mean value derived from averaging the lowest observed
value in the respective sampling event and Null.
Profile average NO3–N values (za) were derived
according to Eq. (1):
4
X
zðxi Þ⁎di
2.2. Field data collection
za ¼
We selected 151 sampling locations using a stratified-random sampling design targeted to areas with land
use–soil combinations proportional to their aerial extent
and other land uses (e.g. tree groves) that were expected
to show high NO3–N in soils according to expertknowledge from experienced local extension specialists
(personal communication Dr. M.W. Clark). Soil samples
were collected during Sept. 2003, Jan. 2004 and May
2004 from four depth increments (0–30, 30–60, 60–120
and 120–180 cm) in composites proportional to the
depth of sampling to ensure a constant sampling support
for each layer increment. Composite soil samples from
each site lumped the local variability of soil properties
providing a local signature of average site conditions.
Not all sampling depths could be sampled due to field
conditions (e.g. high water table during rainy season) or
due to lack in receiving sampling permissions on
i¼1
4
X
ð1Þ
di
i¼1
where, z(xi) is the measured NO3–N in μg g− 1 soil at
sampling sites (xi) for the ith layer (i = 1, 2, 3 and 4) and
di is the thickness in cm of the ith layer.
Table 1
Number of samples for three sampling events
Sampling
time
Soil layers
Layer 1
(0–30cm)
Layer 2
(30–60cm)
Sept.,
2003
Jan.,
2004
May,
2004
101
101
89
59
123
122
116
103
128
128
125
122
Layer 3
(60–120cm)
Layer 4
(120–180cm)
ARTICLE IN PRESS
4
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
2.3. Spatial datasets
A comprehensive suite of soil-landscape and environmental predictor variables were assembled using
Table 2
Summary of predictor variables used for CART modeling
Soil-landscape and environmental attributes
Soil data: a
• Taxonomic data: Soil order, suborder, great group; subgroup,
component name
• Soil properties: Available water capacity low (AWCL) (cm/cm),
available water capacity high (AWCH) (cm/cm), surface texture,
water table depth (cm), drainage class, hydric/non-hydric, b
clay class, c depth to clay (cm), d bulk density (g cm− 3), e
organic matter (%),e spodic horizonb
Climate data: f
• Precipitation amount sampling month (mm)
• Precipitation amount previous 3 sampling months (mm)
• Precipitation amount previous 12 sampling months (mm)
Land cover/land use/landscape data:
• Landsat ETM+ 7 spectral bands—local pixel value g
• Landsat ETM+ 7 spectral bands—mean within 3 × 3 and
10 × 10 window g
• Landsat ETM+ 7 spectral bands—std. dev. within 3 × 3 and
10 × 10 window g
• Normalized Difference Vegetation Index (NDVI) g h
• Land use class (2003) i
• Majority land use within 3 × 3 windowi
• Land use change class j [agricultural/non-agricultural]
(period: 1990–2003)
• Population density k
• Physiographic division l
• Ecoregion m
Topographic data:
• Elevation (m)—local pixel value n
• Elevation (m)—3 × 3 windown
• Slope (percentage)—local pixel value n, o
• Slope (percentage)—3 × 3 window n, o
• Upslope drainage area value (m2) n, o
• Upslope drainage area class n, o, p
• Compound topographic index (CTI)—local pixel value n, o, q
• CTI within 3 × 3 window n, o, q
• Distance to stream value (m) r
• Distance to stream classr
Geologic data: s
• Environmental geology
• Floridian aquifer drastic index
• Intermediate aquifer drastic index
• Surficial aquifer drastic index
Geographic data: t
• x coordinates (m)
• y coordinates (m)
• xy coordinates
• Squared xy coordinates
local and focal geographic information system (GIS)
methods (Table 2). We derived soil attributes from the
Soil Survey Geographic Database (SSURGO) using
ArcGIS version 9.0 (Environmental System Research
Institute, Redlands, CA) including taxonomic (e.g. soil
order, suborder) and soil characteristics such as
available water capacity, bulk density and organic
matter content. We gathered precipitation data for the
sampling month (e.g. Sept. 2003), the past three months
of sampling (e.g. June–August for the September 2003
sampling event) and average annual precipitation for the
watershed from NOAA. Spectral reflectance data
(digital numbers—DN) were derived from Landsat
Enhanced Thematic Mapper (ETM+) satellite images. In
order to account for the dominant signature and
variability in the proximity of sampling locations, the
mean and standard deviation of the band reflectance
values were calculated using moving window neighborhoods of 3 × 3 pixels representing a neighborhood of
8100 m2 and 10 × 10 pixels corresponding to a neighborhood of size 90,000 m2, respectively. Land use data
were derived from a supervised classification of Landsat
ETM+ 2003 imagery with an overall classification
a
Soil Survey Geographic Database (SSURGO), Natural Resources
Conservation Service–U.S. Dept. of Agriculture.
b
Presence (1) or absence (0) of soil property.
c
Three clay classes corresponding to the clay content in the soil
profile: 1 = sand, sandy loam, loamy sand, fine sand; 2 = sandy clay
loam; and 3 = clay.
d
Depth to clay—four categories corresponding to depth at which
clay layer was encountered. 1 = clay layer within 50cm of the soil
surface; 2 = clay layer within 100cm of the soil surface; 3 = clay layer
within 200cm of the soil surface; and 4 = no clay layer in the soil
profile.
e
Depth averaged value within soil profile.
f
National Climatic Data Center.
g
Landsat Enhanced Thematic Mapper (ETM+) satellite image
(Febr. 13th, 2003).
h
Normalized Difference Vegetation Index (NDVI) (Jensen, 1996).
i
Seven land use classes (Sabesan, 2004).
j
Nine land use change classes based on agriculture and nonagriculture land uses from 1990, 2000 to 2003 (Sabesan, 2004).
k
U.S. Census Bureau (http://www.census.gov/main/www/cen2000.
html).
l
St. Johns River Water Management District.
m
World Wildlife Fund.
n
National Elevation Dataset (NED) U.S. Geological Survey.
o
Derived topographic attributed from NED.
p
Five classes based on the size of upslope drainage area: A ≤ 1800,
B 1800–3600, C 3600–5400, D 5400–9000, E >9000m2.
q
Compound topographic index (Wilson and Gallant, 2000).
r
Euclidean distance from sampling site to closest stream.
s
Florida Dept. of Environment Protection (1998).
t
Geographic coordinates of sampling sites in meters (Albers Equal
Area Conic map projection).
ARTICLE IN PRESS
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
accuracy of 82% (Sabesan, 2004). Demographic data
(e.g. population density) that served as a proxy
representing anthropogenic induced impact in the
SFRW were derived from U.S. Census data (http://
www.census.gov/main/www/cen2000.html). We used
the National Elevation Dataset (NED) with a pixel
resolution of 30 m to delineate multiple secondary
topographic attributes including slope, upslope drainage
area, and compound topographic index (CTI) using the
ArcHydro extension of ArcGIS. We used geologic
datasets developed by the Florida Department of
Environment Protection such as the Drastic Index that
was derived as a weighting average using the following
variables: (i) Depth to water, (ii) Net recharge, (iii)
Aquifer media, (iv) Soil media, (v) Topography, (vi)
Impact of the vadose zone, and (vii) Hydraulic
conductivity (Florida Dept. of Environmental Protection, 1998). Other landscape metrics included the
distance of each pixel in the watershed to the closest
streams and the geographic coordinates (x and y
coordinates).
2.4. Upscaling methods
Regional predictions of soil NO3–N across the
SFRW were conducted separately for each sampling
event using field observations and a set of environmental spatial variables summarized in Table 2. The
upscaling procedure followed the following basic
assumptions that the prediction of a variable Ẑ (xi) at
an unvisited location (xi) is made by summing three
components that describe the spatial variation of an
environmental variable (Matheron, 1965) (Eq. (2)).
Z ̂ðxi Þ ¼ mðxi Þ þ eðxi Þ þ e V
ð2Þ
where, Ẑ(xi) is the predicted value at location xi, m(xi) is
the trend/drift for the region (or the structural component); ε is a stochastic, locally varying but spatially
dependent residual from m(xi) known as the variation of
the regionalized variable, and ε′ is a residual error term
(spatially independent noise term).
To describe m(xi) we used Classification and
Regression Trees (CART) relating za to environmental
datasets summarized in Table 2. Classification and
Regression Trees are statistical data mining procedures
that were introduced by Breiman et al. (1984). Treebased modeling is suited for the analysis of complex
environmental datasets that require flexible and robust
analytical methods to deal with nonlinear relationships,
high order interactions and missing values (DeAth and
Fabricius, 2000). The CART methodology is based on
5
binary recursive partitioning; binary because parent
nodes are always split into exactly two child nodes and
recursive because the process can be repeated by treating
each child node as a parent node. Classification and
Regression Trees use either categorical or continuous
data types which predict the data class (classification
tree) or the data values (regression tree). Optimized
splitting rules are identified at each level of the tree. The
goal of regression tree models is to partition the data into
relatively homogenous (low deviation) terminal nodes,
and the mean of the values in each node is the predicted
value for that node. The errors from the classification
tree model are expressed with two measures of relative
cost (Breiman et al., 1984): (i) V-fold validation relative
cost which sets aside independent cases that are used to
evaluate tree predictions, and (ii) resubstitution relative
cost, which is computed using the same data used to
train/develop the tree.
We used the CART software package (Salford
Systems, San Diego, CA) to test a variety of tree
models. Trees were pruned to optimize the tree models
for each soil sampling event separately. Our goal was
to identify parsimonious trees that included environmental variables with large power to predict soil
NO3–N. We adopted the concept of minimal costcomplexity measure that is derived using the resubstitution misclassification rate and a penalty imposed
per each additional terminal node to optimize trees
(Breiman et al., 1984). We calculated the resubstitution relative cost and used V-fold validation with a
90% to 10% split between model and validation
datasets. For each of the three sampling events, a tree
was built and the environmental predictor variables
and splitting rules identified which were subsequently
used to predict soil NO3–N across the SFRW
representing the deterministic trend component of
Eq. (2). At each sampling location the residuals were
derived by subtracting the tree-derived trend model
from the soil NO3–N observations according to the
procedure outlined in Fig. 2.
The Moran's Index and semivariogram models were
used to test the residuals for spatial autocorrelation
(Moran, 1948). The Moran's index (I) is calculated
according to Eq. (3) (compare Goodchild, 1986):
XX
XX
I¼
wij ½ðzi −z̄ Þðzj −z̄ Þ=s
wij
ð3Þ
i
j
i
j
where zi and zj refer to attribute values measured at
different locations i and j, z̄ denotes the mean of the
attribute value, s is the sample variance, and the wij
term represents the spatial proximity (or similarity) of
i's and j's locations, with wii = 0 for all i and wjj = 0
ARTICLE IN PRESS
6
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
Fig. 2. Flow diagram summarizing our adopted methodology for regression-kriging.
for all j. The Moran index values that are positive
spatial autocorrelation, whereas Moran index values
smaller than 1 indicate independent, dissimilar attribute
values.
We found that the residuals of only one sampling
event (January 2004) were spatially autocorrelated. For
these spatially autocorrelated residuals a semivariogram
was generated and ordinary kriging used to generate soil
NO3–N residual predictions. According to Eq. (2) the
residuals were added to the trend component to create
the final output map showing the predicted spatial
distribution of soil NO3–N. This procedure resembles
regression kriging as outlined by Odeh et al. (1995)
(called “model C—RK”) and Hengl et al. (2004). The
advantage of adopting a hybrid modeling approach is
that kriging of regression residuals may improve the
prediction performance compared to the performance
when kriging or regression is done solely. While the
regression model predicts the trend component across
the area for which the explanatory variables are
available, residuals are restricted to the region for
which stationarity of the residuals can be assumed.
Regression kriging is a hybrid technique that combines a
multi-linear regression model (or Regression Tree or
Generalized Linear Model) with kriging of the regression residuals (Odeh et al., 1995; Goovaerts, 1997). This
method is based on the assumption that the deterministic
component of the target (soil) variable is accounted for
by the regression model, while the model residuals
represent the spatially dependent component.
Spatial modeling and upscaling were conducted
using the ISATIS software (Geovariances Americas
Inc., Houston, TX) and visualizations in ArcGIS 9.0
(ESRI, Redlands, CA). The prediction quality consists of
two components: map precision that measures the
residual variability in prediction, and map accuracy
that measures the closeness of the predictions to true
conditions (Mueller et al., 2001). We used the following
error measures (Mueller et al., 2001; Schloeder et al.,
2001; Webster and Oliver, 2001; Erxleben et al., 2002):
(i) Mean Prediction Error (MPE) to assess the systematic
error of the model; Mean Absolute Error (MAE) that is
indicative of the most accurate global or large-scale
estimates; Mean Square Error (MSE) that is more
sensitive to outliers than the MPE; Precision which is the
inverse of the variance; and the G value which measures
how effective a prediction might be relative to that which
could have been derived by using the sample mean.
3. Results and discussion
3.1. The NO3–N dataset—exploratory data analysis
Statistical properties for the depth averaged soil
NO3–N content at the different sampling events are
summarized in Table 3. Variations in soil NO3–N in
different seasons were due to land use management
practices (e.g. fertilizer applications), crop growth,
climatic conditions, nitrogen cycling and leaching. The
mean NO3–N was highest in Jan. 2004 with 4.1 μg g− 1
ARTICLE IN PRESS
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
Table 3
Statistical properties of profile averaged soil NO3–N (μg g− 1) for
different seasons
Statistics
Sept. 2003
Jan. 2004
May 2004
Number of observations
Mean
Standard error of the mean
95 % confidence of mean
— Lower
— Upper
Median
Min.
Max.
Standard deviation
Skewness
Kurtosis
101
0.70
0.13
0.45
123
4.09
1.45
1.22
128
1.17
0.24
0.69
0.97
0.23
0.05
6.54
1.31
2.77
7.59
6.96
0.11
0.06
103.71
16.09
4.91
25.36
1.64
0.29
0.05
19.92
2.74
4.37
22.30
followed by May 2004 with 1.2 μg g− 1 and Sept. 2003
with 0.7μg g− 1. The low NO3–N content during Sept.
2003 could be attributed to the influence of precipitation. High rainfall during the rainy season, with 600 mm
from June–Aug. 2003, favored leaching losses of NO3–
N from soils. The highest NO3–N values with a
maximum of 103.7 μg g− 1 encountered in Jan. 2004
were attributed to fertilization of crops and pastures, low
plant uptake as well as low microbial immobilization
and nitrification during the winter period. Similar results
were found in two other studies that documented the
impact of fertilization on soil NO3–N (Pérez et al.,
2003) and NO3–N leaching due to rainfall–runoff
events (Guo et al., 2001).
7
All events showed positively skewed NO3–N
distributions with many measurements below the
detection limit (BDL). It is a common phenomenon
that environmental datasets show many more low
observed values and fewer high values. For example,
Johnson et al. (2003) observed positively skewed
distribution of NO3–N content in 19 ground water
monitoring wells with values ranging from BDL of 0.10
up to 350.0 mg L− 1 with a median concentration of
2.0 mg L− 1. The soil NO3–N content by land use was
diverse showing high mean values in crops with 2.0μg
g− 1 (Sept. 2003), 26.9μg g− 1 (Jan. 2004) and 3.5 μg g− 1
(May 2004). Improved pasture showed mean NO3–N of
1.2 (Sept. 2003), 11.4 (Jan. 2004) and 2.9 (May
2004) μg g− 1 . Soil nitrate–nitrogen values in this
study were similar to other studies conducted on
comparable soils in Florida (Woodard et al., 2002).
Throughout all sampling events, pine plantations
showed low NO3–N measurements with means of
0.2 μg g− 1. These low values are typical for a Florida
forest ecosystem with tight cycling of nitrogen and low
fertilizer input that may only apply 50kg N ha− 1 once in
the first 5 years (forest regeneration phase) and as much
as 200kg N ha− 1 once during the following 15 years (N.
B. Comerford personal. communication, 2005). Low
soil NO3–N values were also found in wetlands, with
means of 0.3μg g− 1 (Sept. 2003), 0.1 μg g− 1 (Jan. 2004)
and 0.5μg g− 1 (May 2004), respectively. The low NO3–
N concentrations in wetlands are assumed to be due to
denitrification as well as mixing of water from the river
Fig. 3. Classification tree model of soil nitrate–nitrogen (NO3–N) for the Sept. 2003 event. The non-terminal nodes contain the node number and
number of predicted sites. The terminal nodes (TN) are shown in bold and contain terminal node number, the number of data (N) and the predicted
class.
ARTICLE IN PRESS
8
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
Table 4
Variable splitting criteria for the classification tree to predict NO3–N
classes (Sept. 2003)
Node
Splitting
number variable
1
2
3
Splitting criteria
Routing to left
child node
Land use
Agriculture,
rangeland
Physiographic High flatwoods,
divisions
Lake city karst,
lower, mcalpin
plain
Soil order
Alfisols,
Inceptisols,
Ultisols
Routing to right child
node
Pine plantation, upland
forest, urban, wetland
Bell sand hills and
willford flats, Haile
limestone plain, lake
city ridge, newberry
sand hills, perched
lakes and prairies, pose
creek scarp, san felasco
hammock, trail ridge
Entisols, Histosols,
Spodosols
riverbed and terrace deposits was accounted for by
denitrification and the rest by mixing of river and
floodplain water.
The chronological assessment using age dating
techniques made by Katz et al. (2001) have shown
that the NO3–N concentrations of springs of the
Suwannee River Basin have responded to increased
nitrogen loads from various sources in the basin in the
order of decades. Long term trends of NO3–N in the
basin showed that the increasing NO3–N concentrations in spring waters followed the steady increase in
fertilizer use (Katz et al., 1999). Katz (2004) found that
inorganic fertilizers (agricultural land use) were the
dominant source of nitrogen in spring waters in the
Suwannee River Basin based on δ15N isotope tracers.
Our study confirmed that land use is closely related to
soil NO3–N.
3.2. Regional predictions of soil NO3–N
and groundwater by water level fluctuations. In a
comparative study, McMahon and Böhlke (1996)
showed that the median NO3–N concentration in
groundwater from adjacent floodplain deposits
(468 μmol L− 1 ) and riverbed sediments (461μmol
L− 1) were lower than the median concentration in the
terrace deposits (1857 μmol L− 1). Authors estimated
that 15–30% of the difference between floodplain/
We calculated Moran's indices and semivariograms
using observed soil NO3–N values for each sampling
event, but no spatial autocorrelation was found. Because
all three seasonal soil datasets showed spatially
independent observations, we proceeded with relating
environmental datasets to soil NO3–N measurements to
generate trend models using CART.
Fig. 4. Regression tree model of NO3–N for Jan. (a) and May (b) 2004 sampling events.
ARTICLE IN PRESS
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
9
Table 5
Variable splitting criteria for regression tree model for Jan. 2004
Table 7
Relative cost of classification tree model for Sept. 2003 sampling event
Node
Class
No. of
observations
Error method
Misclassified
counts
Cost
measure
A
30
B
20
C
26
D
25
Resubstitution
Validation
Resubstitution
Validation
Resubstitution
Validation
Resubstitution
Validation
17
20
3
11
4
11
6
14
0.567
0.667
0.150
0.550
0.154
0.423
0.240
0.560
1
2
3
Splitting variable
Band 4 reflectance
Land use
Available water
capacity low
(AWCL)
Splitting criteria
Routing to the
left node
Routing to
the right node
≤113.5
Pine plantation,
upland forest, urban,
wetland
≤0.025
>113.5
Agriculture,
rangeland
>0.025
3.2.1. Classification tree model—Sept. 2003 sampling
event
The regression tree model of Sept. 2003 produced
only 2 nodes, suggesting that (i) either the NO3–N
distribution could not be predicted based on the selected
predictor variables, (ii) the size of the dataset was too
small to build a predictive tree, or (iii) the range of the
target soil variable values was too small as soil NO3–N
was relatively homogenous throughout the watershed.
To optimize tree-based predictions, a classification tree
Table 6
Variable splitting criteria for regression tree model for the May 2004
event
Node Splitting variable
Splitting criteria
Routing to left
child node
Routing to
right child
node
Pine plantation,
rangeland, urban,
wetland
≤14.0
Agriculture
1
Land use
2
Mean compound
topographic index
(CTI) 3 × 3 window
Component name
Arredondo, Bigbee,
Bivans, Blanton,
Blichton, Bonneau,
Chiefland, Dorovan,
Fluvaquents, Kendrick,
Kershaw, Lake, Lakeland,
Leon, Lochloosa, Lynn
Haven, Mascotte,
Millhopper, Monteocha,
Myakka, Neilhurst,
Newnan, Ocilla, Oleno,
Penney, Plummer,
Pomona, Ridgewood,
Sapelo, Sparr, Tavares,
Wauchula, Zolfo
Band 7 avg.
≤113.5
reflectance in
10 × 10 window
3
4
>14.0
Albany,
Chipley,
Pelham
was built to predict the four quartile classes of NO3–N:
class A (≤ 0.01μg g− 1), B (0.01–0.23 μg g− 1), C (0.23–
0.56 μg g− 1) and D (≥ 0.56 μg g− 1). The optimum tree
(Fig. 3) consisted of four nodes using land use,
physiographic divisions and soil orders as predictor
variables with different splitting criteria (Table 4).
3.2.2. Regression tree model—Jan. and May 2004
sampling events
The optimal regression tree models for the Jan. and
May soil sampling events are shown in Fig. 4(a) and (b)
and their splitting criteria are summarized in Tables 5
and 6, respectively. The Jan. model had 4 predicted
values of NO3–N, while the May model had 5 values.
The seasonal regression tree models used different sets
of environmental variables to predict NO3–N for
sampling events reflecting different seasonal patterns.
In Jan., spectral reflectance (band 4), land use and
available water capacity (low) were used as predictor
variables (Fig. 4). According to Lillesand and Kiefer
(2000) band 4 (0.76–0.90 μm near-infrared) absorbs
reflectance from healthy green vegetation and is
commonly used to estimate biomass. This band
emphasizes soil–crop and land–water contrasts (Jensen,
1996). In May 2004, predictor variables included land
use, mean compound topographic index within a local
Table 8
Errors associated with different sized regression tree models for the
Jan. 2004 sampling campaign
No. of nodes Validation relative error Resubstitution relative error
7
6
5
4a
2
1
>113.5
0.89 ± 0.26
0.89 ± 0.26
0.89 ± 0.26
0.86 ± 0.25
1.13 ± 0.16
1.01 ± 0.00
0.46
0.46
0.46
0.47
0.71
1.00
The mean and variance associated with one node tree were 1.47 and
24.56, respectively.
a
The optimal tree model.
ARTICLE IN PRESS
10
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
Table 9
Errors associated with different sized regression tree models for the
May 2004 sampling campaign
No. of nodes
Validation
relative error
Resubstitution
relative error
7
6
5a
4
2
1
0.59 ± 0.21
0.55 ± 0.19
0.51 ± 0.18
0.52 ± 0.18
0.85 ± 0.07
1.01 ± 0.00
0.08
0.10
0.14
0.20
0.34
1.00
The mean and variance associated with one node tree were 1.67 and
7.43, respectively.
a
The optimal tree model.
neighborhood of 3 × 3 pixels, the main component soils
and mean reflectance band 7 within a 10 × 10 pixel
neighborhood (Fig. 4b). Band 7 has a spectral resolution
from 2.08 to 2.35 μm representing the mid-infrared
range and is important to discriminate geologic rock
formations (Jensen, 1996). Though band 7 showed a
statistical relationship to soil NO3–N, it has less
meaning to describe nitrogen cycling or leaching
behavior. Land use was consistently used as a predictor
variable in all three tree models (Sept. 2003, Jan. and
May 2004), illustrating the importance of land use
management (e.g. fertilizer input), density and structure
of vegetation that are likely impacting nitrate leaching.
3.2.3. CART model error
The costs associated with the classification tree
model are summarized in Table 7. The validation cost
measure is considered more rigorous and conservative
than the resubstitution cost measure. Class C showed the
lowest costs with 0.15 based on resubstitution and 0.42
based on validation indicating the lowest misclassification rate. From an environmental perspective, it is more
important to predict high soil NO3–N values (classes C
and D) more accurately when compared to lower values
(classes A and B) that are closer to “natural” background
concentrations.
The errors associated with multiple regression tree
models for the Jan. 2004 event are summarized in Table
8. We selected the optimum regression tree with 4 nodes
(i.e., 3 splitting variables) based on validation results
with a minimum variance of 0.86. The simplest tree
consisted of a single node, with mean value of 1.47 and
variance of 24.56. The errors associated with regression
tree models for the May 2004 sampling event are
summarized in Table 9. The one node tree model had a
mean and variance of 1.17 and 7.43, respectively, and
the optimum tree with 5 nodes had a validation relative
error of 0.51 and a resubstitution relative error of 0.14.
V-fold validation is considered more rigorous than
resubstitution to evaluate modeling success. The NO3–
N predictions made for the May event are more accurate
compared to Jan. predictions according to the validation
relative errors (Tables 8 and 9).
3.2.4. Upscaling of the tree model predictions
The predicted NO3–N classes for the Sept. 2003
event (Fig. 5) were upscaled to unsampled locations
(pixels) across the watershed based on predictor
variables in the respective tree model. This model
shows the deterministic trend of NO3–N across the
watershed. Class D of the Sept. sampling event was
solely based on land use (agriculture and rangeland).
Predictions of class C were derived from land use and
physiographic division whereas classes A and B were
Fig. 5. Classification tree predictions of NO3–N for Sept. 2003 upscaled to the watershed scale. The unmapped areas (white) are water bodies.
ARTICLE IN PRESS
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
distinguished based on soil orders. Class B extended
along a NW–SE axis corresponding to the Cody Scarp.
Overall, in the western part of the watershed higher
NO3–N values dominated (classes C and D) whereas
east of the Cody Scarp classes A and C were prominent.
Since soil NO3–N observations were relatively low in
the Fall season (Sept. 2003), generated patterns of soil
NO3–N were relatively homogenous throughout the
mixed-use watershed.
Predicted NO3–N values for the Jan. and May 2004
events (Fig. 6) were also upscaled to the unsampled
locations across the watershed. Both the Jan. and May
trend models showed mixed distributions of high and
11
low NO3–N predictions across the watershed. Land use
was an important predictor variable in both seasons.
Speckled patterns of soil NO3–N were predicted based
on sets of environmental variables including land use,
spectral data, soil and topographic attributes. Although
the predicted values differed for the Jan. and May
events, the spatial distributions of high and low values
were similar. In general, the predicted higher soil NO3–
N values in the western part of the watershed pose a high
risk to be leached into the aquifer due to predominantly
sand-rich soils that are unconfined coinciding with karst
terrain. On the Lowland side of the watershed (west of
the Cody Scarp) the Floridian aquifer is unconfined
Fig. 6. (a) Regression tree predictions of NO3–N for Jan. and (b) May 2004 upscaled to the watershed scale.
ARTICLE IN PRESS
12
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
(Fernald and Purdum, 1998) posing a high risk for soil
NO3–N to be leached into the aquifer. Because the
Lowland side is also dominated by karst terrain and
sand-rich soils potential risks of leaching is higher in the
western region than the eastern Highland part of the
watershed where the Floridian aquifer system is
confined. The Cody Scarp represents the area containing
the most mature karst features including artesian
springs, steep-sided sinkholes, lakes that periodically
drain downward, disappearing and resurgences of
streams (Fernald and Purdum, 1998). High soil NO3–
N found along the Cody Scarp possibly poses the
greatest risk of accelerated nitrate leaching through the
vadose zone into the aquifer.
The residuals derived from the regression tree model
predictions at sampled sites account for the error of the
tree model and their error statistics are summarized in
Table 10. Error statistics showed smaller MPE, MAE,
MSE, and Precision for the May NO3–N trend model
compared to the Jan. 2004 model. The G values for both
models were high. The Moran's I analysis of the NO3–
N residuals for the Jan. 2004 event indicated that there
was significant positive autocorrelation with a positive I
value of 0.24 at the first distance interval (0–4500m).
For all other distance intervals, Moran's I values showed
neither significant positive nor negative autocorrelation.
Moran's I analysis for the May 2004 event revealed that
there was no significant autocorrelation of residuals at
any distance intervals indicating that the CART-based
trend map provided the best estimate.
For the Jan. dataset, a semivariogram model using the
residuals was derived by subtracting the trend values
from the observed values. A spherical semivariogram
model was fitted with the following parameters: range
5000 m, partial sill 9.2 and nugget of 1.8 (Fig. 7).
According to Camberdella et al. (1994), there is a strong
spatial dependence of the observations if the ratio of
nugget to total semivariance is <25%, which was the
case for the Jan. soil NO3–N residual dataset but not for
the residual May dataset. The interpolated residual
Table 10
Error statistics for the regression tree and regression kriging
predictions in different seasons
Error
measures
Jan. 2004
Trend model
(CART)
RK model
Trend model
(CART)
May 2004
MPE
MAE
MSE
Precision
G value
−0.034
1.129
11.321
11.414
53.900
0.230
0.738
3.043
3.090
82.570
−0.057
0.537
1.053
1.061
85.833
Fig. 7. Experimental (thin line) and fitted (bold line) semivariogram of
the NO3–N residuals for Jan sampling event.
model was added to the trend model to produce the final
NO3–N distribution map for the Jan. sampling campaign (Fig. 8). The addition of kriged residual data
corrected for the over- and under-predictions made by
the trend model. Cross-validation of soil NO3–N
predictions derived via regression kriging showed a
MPE of 0.23, indicating that predictions slightly overestimated observed values. Overall, the RK model
performed slightly better than the trend model for the
Jan. sampling event which was confirmed by a lower
MAE, MSE and Precision (Table 10). The G value was
higher for the RK model (82.57) when compared to the
trend model (53.90), suggesting that the former
provided better regional predictions of soil NO3–N
across the watershed (Table 10). The G value for the RK
model was slightly lower than that for the May
predictions which could be attributed to the higher
variance in the Jan. NO3–N data compared to the May
data. With a given set of predictor variables, it becomes
increasingly difficult to predict the target variables when
the target variables become more variable.
Spatial patterns of soil NO3–N differed among
seasons. In Sept., class A covered an area of 11.3%,
class B 6.2%, class C 40.4% and class D 42.1%. For the
Jan. sampling event, 90.7% of the watershed was
predicted to have ≤ 5μg g− 1 soil NO3–N, while 4.2%
was predicted to be within 5–10μg g− 1 soil NO3–N,
and 5.1% had predicted values ≥ 10 μg g− 1. During the
May 2004 sampling event, 81.9% showed soil NO3–N
predictions ≤5 μg g− 1 while 14.0% were in the range
between 5 and 10μg g− 1 and 4.1% of sites ≥ 10 μg g− 1.
Identifying sites that showed consistently high soil
NO3–N predictions in different cropping seasons (Jan.
and May) may help to provide a stronger signature of
“high risk areas” that potentially contribute to excessive
nitrate leaching. Sites that showed consistently higher
values throughout the cropping season, as represented
ARTICLE IN PRESS
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
13
Fig. 8. Predictions derived from regression kriging of NO3–N distribution across the watershed for the Jan. 2004 event (search neighborhood:
1500m).
by the Jan and May sampling events, in the predicted
5–10μg g− 1 range covered 4.3% (150.7 km2) of the
watershed, while sites exceeding predicted values of
10 μg g− 1 occurred on only 0.59% (20.7 km2) of the
watershed. Over 95.7% of the watershed expanding
over an area of 3363.9 km2 predicted soil NO3–N
were ≤5 μg g− 1 (Jan. and May, 2004).
In the western part of the watershed that is
characterized by the presence of sinkholes and karst
terrain, 93.5% of the area had soil NO3–N predictions
≤ 5μg g− 1, 2.8% within predicted 5–10μg g− 1 and
3.7% exceeding 10 μg g− 1 in Jan. 2004. Katz et al.
(2004) pointed out that karst areas in Florida are
susceptible to accelerated leaching of NO3–N, allowing
for ground water pollution. High soil NO3–N values in
close proximity to streams and water bodies pose a
greater risk to increase NO3–N in surface waters. We
found that only a small area (2.5%) showed soil NO3–N
exceeding 10 μg g− 1 within a 500 m buffer of streams
and water bodies in Jan 2004 (area of the buffer = 110.8km2). Within the same buffer, 0.2% showed
predicted soil NO3–N in the 5–10 μg g− 1 range while
97.3% showed predicted values ≤ 5 μg g− 1. Similarly,
the May prediction model showed 2.5%, 8.3% and
89.1% of the buffer zone accounting for > 10, 5–10 and
< 5μg g− 1, respectively.
4. Summary and conclusions
We adopted a hybrid geospatial modeling approach
that combined sparsely measured soil NO3–N observations collected in three seasons (Sept. 2003, Jan.
and May 2004) with dense auxiliary environmental
datasets to predict NO3–N within the SFRW in northeast Florida. For each season, different sets of
environmental variables had predictive power. Based
on the relationships between environmental predictor
variables and soil NO3–N, we selected season-specific
upscaling strategies. The weakest relationships were
found for the Sept. 2003 sampling event for which a
classification tree provided the best estimate. Stronger
relationships between environmental variables, especially land use and soil NO3–N, were found for the
Jan. and May 2004 events using regression trees. For
the Jan. event we found that soil NO3–N residuals
were spatially autocorrelated and therefore proceeded
with kriging of the residuals. Overall, NO3–N
prediction quality was highest in Jan. and May
2004. At each step of the upscaling process, error
assessment provided important information about the
uncertainty of predictions. We found that no universal
methodology (e.g. CART, RK) was applicable to
predict seasonal and regional patterns of soil NO3–N
within the SFRW. This might be due to the fact that
NO3–N is a highly dynamic soil property that leaches
readily and is highly dependent on climate, soil,
topography and season-specific land use management.
Therefore, a mixed approach that combines exhaustive
spatial environmental data (e.g. spectral data derived
from remote sensing) and sparse, season-specific field
observations of NO3–N provides a compromise
between labor, costs, and accuracy for regional
space–time predictions of soil NO3–N across a large
watershed. We were able to identify sites that showed
persistent high values in soil NO3–N in different
seasons that may be targeted first to implement best
ARTICLE IN PRESS
14
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
management practices. Other sites that show high soil
NO3–N and occur in proximity to streams, water
bodies or karst features pose a potential high risk for
accelerated nitrate leaching.
Acknowledgements
This project was supported by the USDA grant #200200501 funded through the “Nutrient Science for Improved
Watershed Management” program. This research was
supported by the Florida Agricultural Experiment Station
and approved for publication as Journal Series No. R10792. We would like to thank the support of Dr. M.W.
Clark and Dr. D.A. Graetz, and Dr. R.B. Brown.
References
ArcGIS software (Vers. 9.0). Environmental Systems Research
Institute (ESRI), Redlands, CA.
Bishop, T.F.A., McBratney, A.B., 2001. A comparison of prediction
methods for the creation of field-extent soil property maps.
Geoderma 103, 149–160.
Breiman, L., Friedman, J.H., Olshen, R.A., Stone, C.J., 1984.
Classification and Regression Trees. Wadsworth International
Group, CA.
Brown, R.B., Stone, E.L., Carlisle, V.W., 1990. Soils. In: Myers, R.L.,
Ewel, J.J. (Eds.), Ecosystems of Florida. University of Central
Florida Press, Orlando.
Camberdella, C.A., Moorman, T.B., Novak, J.M., Parking, T.B.,
Karlen, D.L., Turco, R.F., Konopka, E.E., 1994. Field-scale
variability of soil properties in central Iowa soils. Soil. Sci. Soc.
Am. J. 58, 1501–1511.
DeAth, G., Fabricius, K.E., 2000. Classification and regression trees: a
powerful yet simple technique for ecological data analysis.
Ecology 81 (11), 3178–3192.
Erxleben, J., Elder, K., Davis, R., 2002. Comparison of spatial interpolation methods for estimating snow distribution
in the Colorado Rocky Mountains. Hydrol. Process. 16,
3627–3649.
Fernald, E.A., Purdum, E.D. (Eds.), 1998. Water Resources Atlas of
Florida. Institute of Science and Public Affairs, Florida Sate
University, FL.
Florida Department of Environmental Protection, 1998. State wide
geologic data for Florida. Available at: http://www.fgdl.org.
Goodchild, M.F., 1986. Spatial autocorrelation. Concepts and
Techniques in Modern Geography, vol. 47. Geo Books, Norwich.
Goovaerts, P., 1997. Geostatistics for Natural Resources Evaluation.
Oxford University Press, NY.
Goovaerts, P., 1999. Using elevation to aid the geostatistical mapping
of rainfall erosivity. Catena 34 (3–4), 227–242.
Guo, L.P., Zhang, F.S., Wang, X.R., Mao, D.R., Chen, X.P., 2001.
Effect of long-term fertilization on soil nitrate distribution.
J. Environ. Sci. 13 (1), 58–63.
Grunwald, S., Lamsal, S., 2005. Emerging geographic information
technologies and soil information systems. In: Grunwald, S. (Ed.),
Environmental Soil-Landscape Modeling—Geographic Information Technologies and Pedometrics. CRC Press, New York,
pp. 127–154.
Hengl, T., Heuvelink, G.B.M., Stein, A., 2004. A generic framework
for spatial prediction of soil variables based on regression-kriging.
Geoderma 120, 75–93.
Hornsby, D., Mattson, R., Mirti, T., 2001. Surface Water Quality and
Biological Annual Report. Suwannee River Water Management
District, Florida.
Hornsby, D., Mirti, T., Zwanka, W., 2002. Nitrate Profile of the Lower
Santa Fe River. Department of Water Resources, Suwannee River
Water Management District, Live Oak Florida.
ISATIS software (Vers. 4.0.1). Geovariances, Houston, TX.
James, M.R., Turner, M.G., Smithwick, E.A.H., Dent, C.L., Stanley,
E.H., 2004. Spatial extrapolation: the science of predicting
ecological patterns and processes. Bioscience 54 (4), 310–320.
Jensen, J.R., 1996. Introductory Digital Image Processing, a Remote
Sensing Perspective. Prentice-Hall, Upper Saddle River, New
Jersey.
Johnson, D.M., Osiensky, J.L., Miller, S.M., 2003. Geostatistical
ground water monitoring of a point source NO3–N plume
entering a restored riparian zone. J. Environ. Sci. Health 38 (5),
719–735.
Katz, B.G., 2004. Sources of nitrate contamination and age of water in
large karstic springs of Florida. Environ. Geol. 46, 689–706.
Katz, B.G., Hornsby, H.D., Bohlke, J.F., Mokray, M.F., 1999. Sources
and chronology of nitrate contamination in spring waters,
Suwannee River Basin, Florida. USGS Water-Resources Investigations Report 99-4252. Tallahassee, Florida.
Katz, B.G., Bohlke, J.K., Hornsby, H.D., 2001. Timescales for nitrate
contamination of spring waters, northern Florida, USA. Chem.
Geol. 179, 167–186.
Katz, B.G., Chelette, A.R., Pratt, T.R., 2004. Use of chemical and
isotopic tracers to assess nitrate contamination and ground water
age, Woodville Karst Plain, USA. J. Hydrol. 289, 36–61.
Keeney, D.R., Nelson, D.W., 1982. Nitrogen—inorganic forms. In:
Page, A.L., Miller, R.H., Keeney, D.R. (Eds.). Methods of Soil
Analysis: Part 2. American Society of Agronomy Inc., Soil Sci.
Soc. of Am. Inc., Madison, WI.
Knotters, M., Brus, D.J., Voshaar, J.H.O., 1995. A comparison of
kriging, co-kriging and kriging combined with regression for
spatial interpolation of horizon depth with censored observations.
Geoderma 67, 227–246.
Lillesand, T., Kiefer, R., 2000. Remote Sensing and Image
Interpretation. John Wiley and Sons, New York, NY.
Matheron, G., 1965. Les Variable Régionalisees et leur Estimation.
Masson, Paris, France.
McBratney, A.B., Odeh, I.O.A., Bishop, T.F.A., Dunbar, M.S., Shatar,
T.M., 2000. An overview of pedometric techniques for use in soil
survey. Geoderma 97, 293–327.
McBratney, A.B., Mendonça Santos, M.L., Minasny, B., 2003. On
digital soil mapping. Geoderma 117, 3–52.
McMahon, P.B., Böhlke, J.K., 1996. Denitrification and mixing in a
stream-aquifer system: effects on nitrate loading to surface water.
J. Hydrol. 186, 105–128.
Moran, P.A.P., 1948. The interpretation of statistical maps. J. R. Stat.
Soc., Ser. B Methodol. 10, 243–251.
Mueller, T.G., Pierce, F.J., Schabenberger, O., Warncke, D.D., 2001.
Map quality for site-specific fertility management. Soil Sci. Soc.
Am. J. 65, 1547–1558.
Odeh, I.O.A., McBratney, A.B., 2000. Using AVHRR imageries for
spatial prediction of clay content in the lower Namoi Valley of
eastern Australia. Geoderma 97, 237–254.
Odeh, I.O.A., McBratney, A.B., Chittleborough, D.J., 1995. Further
results on prediction of soil properties from terrain attributes:
ARTICLE IN PRESS
S. Lamsal et al. / Geoderma xx (2006) xxx–xxx
heterotopic cokriging and regression-kriging. Geoderma 67,
215–226.
Pérez, J.M.S., Antigüedad, I., Arate, I., Linares, C.G., Morell, I., 2003.
The influence of nitrate leaching through unsaturated soil on
groundwater pollution in an agricultural area of the Basque
country: a case study. J. Sci. Total Environ. 317, 173–187.
Randazzo, A.F., Jones, D.S., 1997. The Geology of Florida. University
Press of Florida, Gainesville, FL.
Sabesan, A., 2004. Geo-spatial assessment of the impact of land cover
dynamics and distribution of land resources on soil and water
quality in the Santa Fe River Watershed. M.S. thesis, University of
Florida, Gainesville, FL.
Schloeder, C.A., Zimmerman, N.E., Jacobs, M.J., 2001. Comparison
of methods for interpolating soil properties using limited data. Soil
Sci. Soc. Am. J. 65, 470–479.
Schmidt, W., 1997. Geomorphology and physiography of Florida. In:
Randazzo, A.F., Hones, D.S. (Eds.), The Geology of Florida.
University Press of Florida, Gainesville, FL.
SRWMD and SJRWMD (Suwannee River Water Management District
and St. Johns River Water Management District), 1995. Remote
Sensing Derived Land Use According to the Florida Land Use and
Cover Classification System (FLUCCS). Available at: http://www.
fgdl.org).
15
Suwannee River Water Management District, 2003. Surface Water
Quality and Biological Annual Report. WR-02/03-03.
Trexler, J.C., Travis, J., 1993. Nontraditional regression analyses.
Ecology 76 (6), 1629–1637.
Triantafilis, J., Odeh, I.O.A., McBratney, A.B., 2001. Five geostatistical models to predict soil salinity from electromagnetic induction
data across irrigated cotton. Soil Sci. Soc. Am. J. 65, 869–879.
Van Meirvenne, M., van Cleemput, I., 2005. Pedometrical techniques
for soil texture mapping at a regional scale. In: Grunwald, S. (Ed.),
Environmental Soil-Landscape Modeling—Geographic Information Technologies and Pedometrics. CRC Press, New York,
pp. 323–341.
Webster, R., Oliver, M.A., 2001. Geostatistics for Environmental
Scientists. John Wiley and Sons, Ltd., Chichester, England.
Wilson, J.P., Gallant, J.C., 2000. Terrain analysis: principles and
applications. John Wiley & Sons, Inc., NY.
Woodard, R.K., French, C.E., Sweat, A.L., Graetz, D.A., Sollenberger,
L.E., Macoon, B., Portier, M.P., Wade, L.B., Rymph, J.S., Prine, G.
M., van Horn, H.H., 2002. Nitrate removal and nitrate leaching
from forage systems receiving dairy effluent. J. Environ. Qual. 31,
1980–1992.