JEQ Grow Your Career With ASA
HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS
 QUICK SEARCH:   [advanced]


     


Published online 1 May 2008
Published in J Environ Qual 37:1145-1157 (2008)
DOI: 10.2134/jeq2007.0245
© 2008 American Society of Agronomy, Crop Science Society of America, and Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
This Article
Right arrow Abstract Freely available
Right arrow Figures Only
Right arrow Full Text (PDF) Free
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Citing Articles
Right arrow Citing Articles via HighWire
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Webb, R. M. T.
Right arrow Articles by Linard, J.
Right arrow Search for Related Content
PubMed
Right arrow Articles by Webb, R. M. T.
Right arrow Articles by Linard, J.
Agricola
Right arrow Articles by Webb, R. M. T.
Right arrow Articles by Linard, J.
Related Collections
Right arrow Pedotransfer Functions
Right arrow Pesticides
Right arrow Ecological Risk Assessment
Right arrow Organic Compounds
Right arrow Vadose Zone Processes and Chemical Transport

Variations in Pesticide Leaching Related to Land Use, Pesticide Properties, and Unsaturated Zone Thickness

Richard M. T. Webba,*, Michael E. Wieczorekb, Bernard T. Nolanc, Tracy C. Hancockd, Mark W. Sandstroma, Jack E. Barbashe, E. Randall Baylessf, Richard W. Healya and Joshua Linarda

a U.S. Geological Survey, Denver, CO
b U.S. Geological Survey, Baltimore, MD
c U.S. Geological Survey, Reston, VA
d U.S. Geological Survey, Richmond, VA
e U.S. Geological Survey, Tacoma, WA
f U.S. Geological Survey, Indianapolis, IN

* Corresponding author (rmwebb{at}usgs.gov).

Received for publication May 16, 2007.

    ABSTRACT
 TOP
 NOTES
 ABSTRACT
 INTRODUCTION
 Materials and Methods
 Results and Discussion
 Conclusions, Implications, and...
 REFERENCES
 
Pesticide leaching through variably thick soils beneath agricultural fields in Morgan Creek, Maryland was simulated for water years 1995 to 2004 using LEACHM (Leaching Estimation and Chemistry Model). Fifteen individual models were constructed to simulate five depths and three crop rotations with associated pesticide applications. Unsaturated zone thickness averaged 4.7 m but reached a maximum of 18.7 m. Average annual recharge to ground water decreased from 15.9 to 11.1 cm as the unsaturated zone increased in thickness from 1 to 10 m. These point estimates of recharge are at the lower end of previously published values, which used methods that integrate over larger areas capturing focused recharge in the numerous detention ponds in the watershed. The total amount of applied and leached masses for five parent pesticide compounds and seven metabolites were estimated for the 32-km2 Morgan Creek watershed by associating each hectare to the closest one-dimensional model analog of model depth and crop rotation scenario as determined from land-use surveys. LEACHM parameters were set such that branched, serial, first-order decay of pesticides and metabolites was realistically simulated. Leaching is predicted to be greatest for shallow soils and for persistent compounds with low sorptivity. Based on simulation results, percent parent compounds leached within the watershed can be described by a regression model of the form e–depth (a ln t1/2–b ln KOC) where t1/2 is the degradation half-life in aerobic soils, KOC is the organic carbon normalized sorption coefficient, and a and b are fitted coefficients (R2 = 0.86, p value = 7 x 10–9).

Abbreviations: ATR, atrazine • HYA, hydroxyatrazine • DEA, deethylatrzine • DIA, deisopropylatrazine • DDA, didealkylatrazine • MET, metolachlor • ESA, ethane sulfonic acid • OXA, oxanillic acid • MES, metolachlor ESA • MOX, metolachlor OXA • GLY, glyphosate • AMPA or AMP, aminomethylphosphonic acid • SIM, simazine • PAQ, paraquat • LEACHM, leaching estimation and chemistry model • t1/2, half-life • KOC, organic carbon normalized sorption coefficient • DTK model, depth, t1/2, KOC model


    INTRODUCTION
 TOP
 NOTES
 ABSTRACT
 INTRODUCTION
 Materials and Methods
 Results and Discussion
 Conclusions, Implications, and...
 REFERENCES
 
THE goal of the U.S. Geological Survey (USGS) National Water Quality Assessment Program (NAWQA) is to assess past, current, and future water-quality conditions and trends in river basins and aquifers across the nation (Hirsch et al., 1988). Beginning in 1991, standard techniques were used to document variations in water quality in undeveloped, rural, agricultural, and urban areas (Fuhrer et al., 1999). Now in its second decade, NAWQA is placing more emphasis on understanding how the observed variations can be explained by natural and anthropogenic factors, with specific attention given to recommendations made by the National Research Council (2002): conduct mass balances on constituents of concern; ensure uniformity in the use of models; address uncertainty associated with all aspects of its water-quality modeling evaluations; and use process-based models that can simulate interactions between surface water and ground water. Agricultural catchments in California, Washington, Nebraska, Indiana, and Maryland were the focus of research to better understand the fate and transport of agricultural chemicals through the landscape (Capel et al., 2008). This paper describes simulated leaching of five pesticides and seven metabolites (Table 1 ) at the Maryland site. Degradation and transformation of atrazine and metolachor observed at the Nebraska site were included to constrain model simulations in Maryland because 2 yr had passed since the application of these compounds at the Maryland site.


View this table:
[in this window]
[in a new window]

 
Table 1. Compound names, abbreviations, and chemical names for the five pesticides and seven metabolites modeled in this report.

 
This paper presents estimates of the application, transformation, and leaching of pesticides to ground water for the 32-km2 Morgan Creek watershed in Kent County on the Delmarva Peninsula (Fig. 1a ). The conceptual and numerical models described here address three questions: (i) How do the quantity and quality of water recharging ground water vary with time and across areas with variable unsaturated zone thickness? (ii) What are the seasonal and interannual variations of sorbed and dissolved masses of a variety of pesticides and degradates, each with distinct application times, rates, and properties in the unsaturated silty loams beneath corn (Zea mays L.) and soybean (Glycine max L.) fields on the Delmarva peninsula? and (iii) Is a regression model based on numerical model results capable of describing the relative ground water contamination potential of pesticides and metabolites for the Morgan Creek watershed?


Figure 1
View larger version (92K):
[in this window]
[in a new window]

 
Fig. 1. (A) Location of the Morgan Creek watershed on the Delmarva peninsula. Also shown are the locations of the three weather stations in the Chesapeake Bay region that were used to fill in data or to provide observations back to October 1994. (B) Land use for growing season in 2004 determined from field surveys. The soybean field surrounding the intensive study site containing the flow-path wells is located near the outlet. The soils and instruments in the cross-section below X-X' are presented in Fig. 1D. (C) Unsaturated-zone thickness computed as the difference between the land surface elevation and the steady-state potentiometric head simulated with a ground-water model (MODFLOW). (D) Cross-section of Maryland intensive study site showing unsaturated-zone thickness, and the location of lysimeters and peizometers. Volumetric soil moisture was measured with time-domain reflectometers (TDRs) collocated with the pan lysimeters at M21.

 

    Materials and Methods
 TOP
 NOTES
 ABSTRACT
 INTRODUCTION
 Materials and Methods
 Results and Discussion
 Conclusions, Implications, and...
 REFERENCES
 
Fifteen one-dimensional Leaching Estimation and Chemistry Models (LEACHM; Hutson and Wagenet, 1992; Hutson, 2005) were constructed for different crop rotations and different unsaturated zone thicknesses. LEACHM was selected as an appropriate model for upscaling given its capabilities, ease of use, adequate documentation, and open architecture (Nolan et al., 2005).

A regression model was also constructed from LEACHM simulation results to predict the percent of leached mass as a function of unsaturated zone thickness, organic carbon normalized sorption coefficient, solubility, and half-life. The objective of the regression model was to provide a reduced form, or metamodel, of the process-based model. The metamodel represents key concepts in LEACHM but is easier to use and has computational advantages. The underlying conceptual model is (i) the mass of pesticides and metabolites leaching to ground water is inversely related to residence time in the unsaturated zone, and (ii) residence times are less, and leaching greater, for persistent compounds with low sorptivity that are applied during wet periods than they are for short-lived, sorptive compounds applied during dry periods.

The LEACHM model used in our simulations employs the Richard's equation algorithm for water flow. Before simulating reactive transport of the pesticides, the model parameters related to the fate of water and pesticides were calibrated. Parameters describing flow were calibrated using observations of soil moisture and bromide concentrations made at the intensively studied field in the Morgan Creek watershed whereas the half-lives for atrazine, metolachlor, and their degradates were adjusted so that simulations fit observations of these compounds made at the Nebraska site (Webb and Sandstrom, 2008). All other pesticide properties were culled or derived from existing literature (Capel et al., 2008). The following sections describe the validation process for estimating soil water flow properties and the pesticide degradation constants.

Climate and Land Cover
During 2003 and 2004, local weather conditions, soil moisture and tension, and water quality of the unsaturated and saturated zones were observed on and beneath an intensively studied field, referred to in this paper as the intensive study site, subject to corn–soybean rotations with winter wheat (Triticum aestivum L.) ground cover (Capel et al., 2008). Field by field surveys completed in 2004 showed 75% of the 32-km2 watershed was comprised of agricultural fields (predominantly corn–soybean crop rotations) and 14% was in pasture and hay (Fig. 1b). The remaining areas not in corn–soybean rotation or pasture/hay include roughly equal areas of poultry, dairy, forest, woody wetlands, and the residential areas in and around the town of Kennedyville.

Daily precipitation totals along with evapotranspiration from soils and vegetation provide boundary conditions for the upper layers of the LEACHM models. Data from nearby weather stations were used to extend the 2 yr of weather records (2003–2004) at the intensive study site back to October 1994 to include a wider range of interannual climate variations. The 10 yr of daily observations of precipitation along with minimum and maximum temperatures were used as driving variables to estimate potential evapotranspiration using the Water, Energy, and Biogeochemical model (WEBMOD) (Webb et al., 2006) built for the Morgan Creek watershed (Linard et al., 2008). Precipitation and WEBMOD estimates of potential evapotranspiration were, in turn, used as driving variables for the hydrology simulated in LEACHM. WEBMOD simulations for the 10-yr period predict that approximately 60%, or 72 cm yr–1, of precipitation returns to the atmosphere as a result of evapotranspiration. The 2 mo with highest precipitation and stream discharge are March and September, corresponding to periods of spring storms and the passing extratropical depressions.

Soils
Soils for all models were represented by three layers of silty and sandy loams observed when installing lysimeters, tensiometers, and heat dissipation probes at the intensive study site: a silty loam with less than 1% organic carbon from the surface to a depth of 1 m; a sandy loam also with low percentage of organic carbon from 1 to 2 m depth; and a loam with about 1.4% organic carbon from 2 m to the bottom of the modeled core (site M21; Hancock et al., 2008). Spatial variations in soil type are not explicitly included in our approach because the watershed soils are predominantly horizontally extensive marine deposits where vertical heterogeneity is greater than horizontal heterogeneity. Available soil-series data (STATSGO) had insufficient accuracy and resolution compared with measured data from the soil column, i.e., more than 90% of the watershed is described by only two soil series characterized as fine silts and loams similar to those included in the three layers of the model presented here. Other data (SSURGO) have finer resolution but, like STATSGO, describe only the plant root zone (:2 m deep) and not the deeper unsaturated zone; additional data would be needed to describe soil variations from the bottom of the root zone to the water table.

It should be noted that soils were not sampled from the soil surface, which is rich in organic matter. Woody debris and compacted sillage can provide an important substrate to occlude parent pesticides like atrazine and metolachlor from bacterially mitigated degradation (Lerch et al., 1997). Sorption was simulated using linear isotherms; algorithms to simulate occlusion with associated product persistence are not currently part of the model. The organic-rich layer at the surface not represented in the models would likely further retard the infiltration of pesticides into the soil, further reducing the estimates of leaching through the unsaturated zone simulated here. In practice, pesticides retained in surface soils would have a greater likelihood of being washed from fields into surface drainage. Transport of pesticides sorbed to sediments and organic debris is not considered in this study; all of the mass applied is assumed to remain on the surface until it dissolves and infiltrates during precipitation or irrigation events. Therefore, the thickness of the unsaturated zone, rather than soil type, is used as the dominant factor to represent unsaturated zone variability in the watershed. Fluxes for the entire Morgan Creek watershed were estimated by relating each hectare in the watershed to the model results that most closely matched the crop rotation and unsaturated zone thickness at that point.

Unsaturated Zone Thickness
The thickness of the unsaturated zone (Fig. 1c) was estimated by subtracting the steady-state potentiometric heads, estimated with a MODular three-dimensional ground-water FLOW model (MODFLOW; Capel et al., 2008), from the land-surface digital elevation model (DEM), derived from Light Detection and Ranging (LIDAR, M.R. Nardi, U.S. Geological Survey, written communication, 2005). The LIDAR DEM had an original resolution of 2 m and was resampled to produce average elevations for each 30- by 30-m cell to compare with the 30-m discretization used in the MODFLOW simulation. As derived, 15% of the watershed had no unsaturated zone thickness (where MODFLOW heads exceeded land-surface elevation, producing negative thickness estimates). Therefore the vertical datum for the MODFLOW heads was lowered by 1.65 m, resulting in 5.7% of the watershed with no unsaturated zone thickness. The resulting simulation more closely matched the area classified as wetlands and riparian forest in the National Land Cover Dataset (NLCD) and produced an average unsaturated zone thickness of 7.5 m for the intensive study site. The thickness of the unsaturated zone beneath the intensive study site thus derived is consistent with observations (Fig. 1d) and appreciably thicker than 4.9 m, the mean unsaturated-zone thickness in the watershed. Beneath the flow path at the intensive study site (cross-section X-X' in Fig. 1d), the unsaturated zone thins from 10 m thick at the eastern, upgradient end of the flow path to 0 m in the riparian areas of Morgan Creek.

Final LEACHM Models
Fifteen LEACHM models describing five unsaturated zone thicknesses and three crop rotations were constructed. The five soil columns depths were 1, 2.5, 5, 7.5, and 10 m each discretized into 10 cm thick horizons. The three crop rotations and associated agricultural chemical applications are (i) fields planted in corn in 2004 were assumed to have been planted in corn in even-numbered years and in soybean in odd-numbered years beginning with soybean in 1995; (ii) fields planted in soybean in 2004 were assumed to have been planted in soybean in even years and corn in odd years beginning with corn in 1995; and (iii) the crop rotation at the intensive study site had soybean in 1995 but skipped the planned corn crop in 2004 in favor of soybean. For lack of specific information throughout the rest of the watershed, it was assumed that the 34 ha including and surrounding the weather station was the only place and time that the alternating corn–soybean crop cycle was disrupted.

As a soil conservation measure typical of the area, winter wheat is commonly planted in the fall and turned under in the spring. For the years in corn, the herbicides atrazine, metolachlor, simazine, and paraquat are typically applied to the field in late April at rates of 2.08, 1.36, 1.16, and 0.84 kg ha–1 respectively (Hancock et al., 2008). For the years in soybean, 1.1 kg ha–1 of paraquat is applied in early May with a subsequent application in late July of glyphosate and paraquat at rates of 1.53 and 0.31 kg ha–1, respectively. The rates were assumed typical for the Morgan Creek watershed and were used as the application rates for the target crops biannually in the model. Because the intensive study site was planted in soybean in two consecutive years, it received 20% less pesticide loading than other corn/soybean fields in the watershed (four applications instead of five).

Reactive Transport of Pesticides
The potential for a pesticide to leach through the unsaturated zone to ground water depends on several factors including: (i) the application timing and loads; (ii) fluxes of water into and out of the unsaturated zone; (iii) the properties of the compound; and (iv) the properties of the soils. The pesticide properties used in the LEACHM models are listed in Tables 2 and 3 . The sorbtivity, degradation half-life in aerobic soils, and solubility of the pesticides and metabolites span a wide range of values (Table 3). The detection frequency of two pesticides with similar physical and chemical characteristics may differ as a result of different application rates or different analytical detection thresholds in the laboratory.


View this table:
[in this window]
[in a new window]

 
Table 2. Transformation rates used to model degradation and production of pesticides and metabolites. Pesticide abbreviations are defined in Table 1.

 

View this table:
[in this window]
[in a new window]

 
Table 3. Physical properties of five pesticides and seven metabolites. Pesticide abbreviations are defined in Table 1.

 
LEACHM simulates pesticide transport and fate on the basis of the physical and chemical properties of the compound and a first-order decay model. Linear sorption was assumed. Organic carbon sorption coefficients (KOC), solubilities, vapor densities, and transformation rates at standard temperature and pressure for various pesticides used in our simulations were modified from Gilliom et al. (2006) and cover a wide range of values. Half-lives for atrazine, metolachlor, and their metabolites reported for non-sterile soil were adjusted to be consistent with a branched serial first-order decay model that accounts for partial degradation rates in cases where two or more daughter compounds are produced from the same parent compound (Webb and Sandstrom, 2008).

LEACHM
Hydrology and Transport of Conservative Solutes
Fluxes of water through the unsaturated zone in LEACHM respond to differences of hydraulic potential gradient and unsaturated conductivities as described by Richards' equation. The transport of a conservative solute is described by the convection dispersion equation (CDE), with smaller dispersivity values resulting in a delayed first arrival of the solute at a given depth.

Water and solutes are redistributed at each time step in the simulation in response to infiltration and hydraulic potential gradient. After that, variations in partitioned fractions or transformation rates for any specific compound in a given layer are a result of the intrinsic properties of the compound, the organic carbon content, temperature, or water content of that layer.

Sensitivity Analysis and Calibration
Richard's equation requires knowledge of the soil moisture retention curves for the soils and how unsaturated conductivities vary with water potential. van Genuchten (1980) provided a framework for describing soil moisture retention.

Formula 1[1]
where

Se is effective saturation = Formula 1, 0 ≤ Se ≤ 1;
{theta} is volumetric water content, or volume of water per unit bulk volume of soil;
{theta}r is residual water content;
{theta}s is saturated water content;
h is soil water tension, h ≥ 0;
{alpha} is a curve fitting parameter related to air entry pressure;
n,m are curve fitting parameters related to pore size distribution; the relationship, m = 1–1/n is often assumed.

Unsaturated conductivity is estimated using the van Genuchten hydraulic conductivity relationship, and based on the hydraulic conductivity model of Mualem (1976):

Formula 2[2]
where

Ks is saturated hydraulic conductivity.

Porosity, analagous to saturated water content, can be determined from particle density and bulk density. Bulk densities of core samples recovered during well installation were compressed either during sampling or during shipping such that values exceeded 2.1 g cm–3. Porosities associated with these densities would be unrealistically low so the bulk densities became a variable to calibrate to match observed variations in soil moisture. Bulk density (bd), residual soil moisture ({theta}r), the two van Genuchten fitting coefficients, {alpha} and n, and solute dispersivity (dsp) were calibrated to fit variations in simulated soil moisture and bromide concentrations with observations. During simulations, the variables were allowed to vary within ranges described by Meyer et al. (1997) for each soil textural class.

Dispersivity was estimated using results from a bromide tracer experiment. In the spring of 2004, a solution of KBr (15 mg m–2) was applied to the land surface directly above the lysimeters. Bromide was detected through the summer of 2004 in a pan lysimeter at 55 cm (1.8 ft).

These concentrations along with water content distribution were used to simultaneously optimize {alpha}, n, bd, and dsp using the PEST routine (Doherty, 2001). The model failed to converge when the bulk density of the three layers was allowed to vary independently. Therefore, the bulk densities of the top (bd1) and bottom (bd3) layers were fixed at 1.3 g cm–3, the mean value for silty loams (Meyer et al., 1997), while the density of the middle layer (bd2) was allowed to vary freely. Residual soil moisture ({theta}r) was found to be the most sensitive variable for all layers. Because 52 soil moisture measurements were included in the calibration dataset versus six measurements of bromide concentrations, the dispersivity (dsp) parameter was relatively insensitive. Matching bromide concentrations was difficult in light of the fact that extremely dry conditions were prevalent during the entire summer of intensive study and sampling (Fig. 2 ). Bromide concentrations had only started to rise in the final samples collected in the lysimeters at 55 cm in September 2004. Figure 2 shows the final match between observed and simulated variations of soil moisture at depths of 45 and 115 cm and bromide concentrations at a depth of 55 cm.


Figure 2
View larger version (29K):
[in this window]
[in a new window]

 
Fig. 2. Observed versus predicted moisture content (A) and bromide concentration (B) for the calibrated LEACHM model.

 
Extrapolation to Morgan Creek Watershed
The calibrated LEACHM model was extrapolated throughout the watershed by assigning predicted values to locations with similar land use and unsaturated zone thickness. Both the unsaturated zone thickness grid and the 2004 land-cover map were queried in a geographic information system at a spacing of 100 m, or one sample per hectare, resulting in 3275 sample points. Pesticides were assumed to be applied only to the 2330 ha that were in corn–soybean crop rotation. Summary values compiled for each of the 15 model runs (five depths for three crop scenarios) included amount of pesticide applied, the gains and subsequent losses of metabolites, the residual amounts of the compounds in the unsaturated zone, and the amount leached to ground water. Where the unsaturated zone was between 1 and 10 m thick (87.7% of crop area), summary values were linearly interpolated between the results for the five depth models. Where the unsaturated zone was less than 1 m thick (7.7% of the crop area) summary values of the 1-m column were used. Areas with unsaturated zones greater than 10 m thick (4.6% of the crop area) were assigned the values of the 10-m column.

Relative Potential Ground Water Contamination
The detailed soil, crop, climate, and irrigation data required for LEACHM simulations limit the model's usefulness for making regulatory decisions. Therefore, a variety of simpler screening and index models have been developed to rank the relative potential of various pesticides to contaminate ground water. To simplify the numerical code and improve on earlier screening methods, a regression model was fit to 20 simulation results in which LEACHM predicted leaching to occur from the bottom of the 1-, 2.5-, or 5-m column; no compound was predicted to leach from the 7.5- or 10-m columns by the end of the 10-yr model run. The regression model predicts the mass of pesticides and metabolites leached over a 10-yr period under fields planted in corn on even years as a function of depth, half-life (t1/2), KOC, and solubility (Sol). On the basis of observed relations, the variables were transformed as e–depth, ln(t1/2), ln(KOC), and ln(Sol). The response variable was the amount leached as predicted by LEACHM, as a percent of parent application. Squares of the variables were included as were interaction terms between any two variables to arrive at 14 terms plus an intercept for beginning a stepwise linear regression. An automated stepwise linear regression (Insightful Corporation, 2002) removed only one term, (e–depth)2, to produce a model with a multiple R2 of 0.97 and a p value of 6 x 10–4. Fewer terms were desired to increase the degrees of freedom for the residuals and to arrive at an overall more robust model that is transferable to other pesticides in similar settings. Therefore, beginning with the 14 terms plus an intercept, the least significant term (highest p value for each coefficient) was sequentially eliminated until only two interaction terms remained (e–depth|ln(t1/2) and e–depth| ln(KOC)). The index model will be referred to as the DTK model for depth, half-life (t1/2), and organic carbon normalized sorption coefficient (KOC).

The DTK model results are compared with the Gustafson index (Gustafson, 1991) that predicts the relative mobility of pesticide compounds as a function of t1/2 and KOC and the California method of Specific Numerical Values (Johnson, 1991; Clayton, 2005) that identifies a compound as a potential leacher or non-leacher on the basis of its solubility, half-life (aerobic and anaerobic), and organic carbon normalized sorption coefficient.

A variety of other indices, which will not be considered here, have been developed with increasing levels of sophistication to estimate attenuation of pesticides based on climate, soils, and geologic variations. See Díaz-Díaz and Loague (2001) and Bernard et al. (2005) for examples.


    Results and Discussion
 TOP
 NOTES
 ABSTRACT
 INTRODUCTION
 Materials and Methods
 Results and Discussion
 Conclusions, Implications, and...
 REFERENCES
 
The soil matrix above the saturated zone is referred to as the unsaturated zone, or alternatively as the vadose zone in recognition that during intense rain or snowmelt some pore spaces may become saturated. Perhaps because of its inaccessibility or difficulty in sampling, flow and transport through the unsaturated zone has not received the same attention as flow and transport on the surface or in aquifers. This is unfortunate, because the unsaturated zone influences atmospheric, ecological, and biogeochemical transport processes more than any other hydrologic compartment (Bronstert et al., 2005).

Recharge
The time that rainwater will reside in the Morgan Creek watershed ranges from less than a minute for that landing in the flowing stream or on saturated riparian areas near the outlet to millennia for rainfall that recharges deep, slow moving ground water. Much of the rainwater that infiltrates the soil never reaches the outlet at all, evaporating instead back into the atmosphere. Infiltration is only a fraction of total precipitation, and recharge is only a fraction of infiltration.

For the water years 1995 to 2004 (1 Oct. 1994 through 30 Sept.2004), the Morgan Creek watershed received an average of 115 cm of precipitation per year; the driest year was 2002 with 73 cm and the wettest year was 1996 with 161 cm. The largest single rainfall event was 16–17 Sept. 1999, when the tropical depression that was Hurricane Floyd dumped more than 32 cm of rain on the Delmarva Peninsula. LEACHM predicted that the hydraulic conductivity of the soils would permit only 5 cm of the deluge to infiltrate the soil, the remaining 27 cm becoming surface runoff. This conservative estimate of infiltration during large events is one reason that overall recharge rates are somewhat lower than other independent estimates. For the 10-yr simulation, recharge was greatest in spring seasons (April through June) that followed the wet winters of 1996 through 1998 (Fig. 3 ).


Figure 3
View larger version (48K):
[in this window]
[in a new window]

 
Fig. 3. (A) Seasonal variations in runoff, infiltration, evapotranspiration, and recharge simulated for a 5-meter thick unsaturated zone. The sum of runoff and infiltration are equal to the total rain plus snowmelt for the same period. A difference between infiltration and the sum of evapotranspiration and recharge will result in a change of moisture content in the unsaturated zone. (B) Seasonal recharge for simulated for 1-m, 5-m, and 10-m thick unsaturated zones. Recharge that peaks in the winter at 1-m lags until the following fall before reaching a peak recharge at 10 m. Astronomical seasons are used: Fall (Oct.-Dec.); Winter (Jan.-Mar.); Spring (Apr.-June); and Summer (July-Sept.).

 
All other factors being equal, ground water will be most vulnerable to contamination by pesticides and metabolites where residence time is shortest. For a conservative solute in the unsaturated zone, the mean residence time is equal to the average water content of the unsaturated zone (area normalized to units of length) divided by the recharge, or drainage rate in units of length time–1. Using the average soil moisture and drainage rates for 1995 through 2004, the mean residence time for the top 50 cm of soil was 180 d, increasing to 15 yr for the 10-m column. The average annual recharge for the 5-m column, close to the average mean unsaturated zone thickness was 12.4 cm yr–1. Fisher and Healy (2008) estimated recharge using a mass-balance technique at 31.5 cm yr–1, close to the 29.4 cm yr–1 estimated using the WEBMOD surface-water model. These estimates are within the range of 10 to 40 cm yr–1 computed by Böhlke and Denver (1995) and 34 to 173 cm yr–1 computed by Nolan et al. (2007) by the saturated zone chloride tracer method, and 5.6 to 23.2 cm yr–1 by the unsaturated-zone chloride tracer method. It is not surprising that our estimate (12.4 cm yr–1) agrees more closely with Nolan et al. (2007), because both the LEACHM water flux and the unsaturated zone chloride method are point estimates of drainage through the unsaturated zone. The other methods integrate over larger areas and therefore include focused recharge in the numerous surface water detention structures in the watershed.

To compare the behavior of conservative compounds with that of reactive compounds, LEACHM models included bromide applications with the same loading (208 mg m–2 = 2.08 kg ha–1) and date (28 April) as the atrazine applications. Results of the models that simulated corn crops on even years were combined with the results from models that simulated corn crops on odd years. In this fashion, variations in concentration and loads for 10 2-yr periods could be examined, even though only five applications would occur on any given corn–soybean field in the watershed.

Simulations of bromide transport were consistent with conceptual models. Approximately half of the applied mass was predicted to leach from the top 50 cm in the first 180 d, consistent with a steady-state recharge model. Leaching rates and concentrations were greatest in the weeks following application followed by seasonal variations that responded to climate. Through summer and fall, bromide concentrations rise (Fig. 4 ) and bromide leaching rates fall (Fig. 5 ) as evapotranspiration exceeds infiltration. In fact, recharge at 50-cm depth (as indicated by bromide leaching rate) was predicted to slow or even cease during the summer (July-Sept.) for most years, thus the mass transfer of pesticides and metabolites from the unsaturated zone to ground water during these seasons are insignificant even though concentrations are high. Concentrations predicted to reach equilibrium conditions within a few application cycles for the shallower horizons took longer to equilibrate in deeper horizons. As expected, the 10 m thick unsaturated zone model, with a residence time of approximately 15 yr, was not predicted to reach equilibrium conditions by the end of the 10-yr simulation.


Figure 4
View larger version (38K):
[in this window]
[in a new window]

 
Fig. 4. Variations in mean daily bromide concentrations at 50 cm simulated for 2 yr following a hypothetical application of 208 mg m–2 on 28 April of each year starting in 1995. The dark gray envelope and the black line shows the minimum, maximum, and mean concentrations simulated for ten 2-yr periods following springtime applications from 1995 to 2004. The white lines show daily variations in concentrations during the wettest period (1996–97) and the driest period (2001–02). Astronomical seasons are used: Sp, Spring (Apr.-June); Su, Summer (July-Sept.); F, Fall (Oct.-Dec.); and W, Winter (Jan.-Mar.).

 

Figure 5
View larger version (23K):
[in this window]
[in a new window]

 
Fig. 5. Variations of leached bromide, atrazine, and metabolites simulated for 50 cm and 180 cm, for 2 yr following a hypothetical application of 208 mg m–2 on 28 April of each year starting in 1995. Values are 30-d moving averages. Astronomical seasons are used: Sp, Spring (Apr.-June); Su, Summer (July-Sept.); F, Fall (Oct.-Dec.); and W, Winter (Jan.-Mar.).

 
Chemical Transport
Reactive compounds differ from conservative ones in that their transport is determined not only by application rates and climate but by interactions between the compounds and the soil column. For this study, climate and soil variations with depth were assumed to be the same throughout the Morgan Creek watershed, so predictions differ according to pesticide properties.

The mass leaching and concentrations of a conservative solute such as bromide serve as an upper limit for the other, reactive compounds. Whereas bromide was predicted to begin leaching from the 10-m column during the 10-yr simulation, glyphosate and paraquat (the least mobile of the pesticides) were not predicted to leach from even the 1-m column. The remaining 12 compounds were predicted to leach from the 1-m column, but only the most mobile of the metabolites, DEA, DDA, metolachlor ESA, and metolachlor OXA, were predicted to leach from the 5-m column. No pesticide or metabolite was predicted to leach from the 7.5-m column above the minimum model reporting level of 0.001 kg ha–1. In total, 20 compound-depth combinations had leaching greater than the minimum model reporting level.

The concentrations of atrazine and its metabolites at a given depth reflected the differing properties of the parent, daughter, and granddaughter products. Because of the greater production rate assigned to HYA relative to DEA or DIA, HYA was predicted to be the dominant atrazine metabolite at depths less than a meter. DEA, with its slower degradation rate, was predicted to become the dominant metabolite at deeper depths (Fig. 5). The first-order decay model resulted in log-linear relations for leached mass versus depth (Fig. 6 ). The relation between predicted leaching and other physical and chemical properties is not so clear, but when plotted as a linear-log relation (Fig. 7 ), a general trend of greater leaching with longer half-lives, lower organic sorption, and greater solubility can be discerned.


Figure 6
View larger version (14K):
[in this window]
[in a new window]

 
Fig. 6. Total leached mass of three parent pesticides and seven metabolites simulated by LEACHM for unsaturated zone thickness of 1.0, 2.5, and 5.0 m. Paraquat and glyphosate are not shown since simulations predicted no leaching even for the 1.0-m-thick unsaturated zone.

 

Figure 7
View larger version (14K):
[in this window]
[in a new window]

 
Fig. 7. Mass of five parent pesticides and seven metabolites simulated to leach past 1 m versus compound half-life, organic sorption coefficient, and solubility. [Abbreviations: ATR, Atrazine; HYA, Hydroxyatrazine; DEA, Deethylatrzine; DIA, Deisopropylatrazine; DDA, Didealkylatrazine; MET, Metolachlor; MES, Metolachlor ESA; MOX, Metolachlor OXA; GLY, Glyphosate; AMP, AMPA; SIM, Simazine; and PAQ, Paraquat.]

 
Depth, t1/2, KOC Model (DTK)
The final regression model had a multiple R2 of 0.86 and a p value of 7 x 10–9 (Fig. 8 ):

Formula 3[3]
where

PL is total leached mass of pesticide or metabolite for water years 1995–2004, as percent of total applied parent product, beginning with no pesticides in the soil column;
depth, is thickness of the unsaturated zone, in meters;
t1/2 is half-life, in days; and
KOC is the organic carbon normalized sorption coefficient, in milliliters per gram of organic carbon.


Figure 8
View larger version (26K):
[in this window]
[in a new window]

 
Fig. 8. Percent of applied pesticides and metabolites to leach to ground water simulated by the distributed LEACHM model and by the two-term regression model. [Abbreviations: ATR, Atrazine; HYA, Hydroxyatrazine; DEA, Deethylatrzine; DIA, Deisopropylatrazine; DDA, Didealkylatrazine; MET, Metolachlor; MES, Metolachlor ESA; MOX, Metolachlor OXA; AMP, AMPA; and SIM, Simazine.]

 
This regression equation predicts the mass of pesticides and metabolites to leach from the unsaturated zone beneath corn and soybean crops in Morgan Creek after a 10-yr simulation. The combination of t1/2 and KOC for glyphosate and AMPA were such that no leaching (PL ≤ 0) was predicted for even the thinnest unsaturated zones.

The DTK model is similar to the Gustafson index or 'Ground water Ubiquity Score (GUS):

Formula 4[4]
A GUS value greater than 1.8 would present a leaching risk; reasonable classification accuracies (see discussion of Specific Numeric Values below) for the DTK model were obtained by using a reference depth of 1 m. The DTK model predicts the same rank and relative mobilities as the GUS for five herbicides applied to two different soil types in France (Bernard et al., 2005). The DTK model has the additional advantage of including depth so the results can be mapped where pesticide application and unsaturated zone thickness are known.

Specific Numeric Values, or SNVs (Johnson, 1991; Clayton, 2005), are used to classify potential leaching compounds by the State of California; no national SNV criteria have yet been established by the U.S. Environmental Protection Agency. The California SNV criteria are:

  1. Water solubility greater than 3 ppm.
  2. Organic carbon normalized soil sorption coefficient (KOC) less than 1900 cm3 g–1.
  3. Hydrolysis half-life greater than 14 d.
  4. Aerobic soil metabolism half-life greater than 610 d.
  5. Anaerobic soil metabolism half-life greater than 9 d.

A pesticide is considered a potential leacher if the Boolean expression "(1 or 2) and (3 or 4 or 5)" is TRUE.

Because the aerobic soil metabolism half-life was not a strong predictor of leaching with his dataset, Johnson (1991) set a high value for this SNV to reduce its influence on the screening procedure. In contrast, the DTK model uses the aerobic soil metabolism half-life as a strong predictor. This may be a result of the different statistical procedures (no interactions in the Johnson models, versus interaction terms used here), or better half-life data made available since 1991. The DTK model results are compared with the SNV classification system by examining the classification accuracy using "Known Leachers" from California. Known leachers are commonly detected in ground water. For the comparison, any compound with a positive leaching percentage predicted by the DTK model at 1-m depth was classified as a potential leacher; a compound with a predicted leaching value of zero or less was classified as a non-leacher (negative values are artifacts of the regression model). This simple classification scheme produced an accuracy of 77% compared with 66% using the Boolean expression with the specific numerical values for solubility, KOC, and half-lives for hydrolosis, aerobic metabolism, and anaerobic metabolism. In the short term, the two coefficients of the DTK model can be derived for combinations of soils and climate for the five NAWQA agricultural study sites. In the longer term, a more robust transfer function can be developed by making the coefficients dimensionless and by including the variations of appearance and disappearance that may be expected from different production and disappearance half-lives.

Extrapolation of Field Results to the Watershed Using LEACHM and DTK Models
Variable transport and fate of metolachlor and its metabolites metolachlor ESA and metolachlor OXA as a result of distributed LEACHM modeling for the Morgan Creek watershed for the period 1995 to 2004 are shown in Fig. 9 . The spatial patterns are tied to interactions between controlling processes and variations in land use and unsaturated zone thickness. As expected, leaching of metolachlor and its metabolites was inversely related to unsaturated zone thickness. Residuals of metolachlor and metabolites in the soil column were predicted for areas in corn crop in 2004. For comparison, the DTK model was also applied to each hectare. Instead of assigning the value of the 1-m column to soils from 0- to 1-m depth, as was done to distribute the LEACHM results, the depths were allowed to vary from zero to the maximum thickness. Certain combinations of depth, half-life, and KOC can result in regression predictions of leaching less than zero; a value of zero was assigned in these cases as negative values are physically impossible but indicate extremely low potential for leaching. Overall, the masses of pesticides and metabolites predicted to leach by distributing the LEACHM point model results were similar to the masses predicted to leach by distributing the regression equation. Because the leaching predicted by the regression model was continuous for all depths, some leaching was predicted for the less mobile compounds applied where the unsaturated zone was less than 1 m thick, whereas no leaching was predicted by the 1 m deep LEACHM model (Fig. 10 ).


Figure 9
View larger version (62K):
[in this window]
[in a new window]

 
Fig. 9. Maps of the Morgan Creek watershed showing the total mass of metolachlor and metabolites applied or produced, transformed, leached, and residual in the unsaturated zone as of October 2004. A total of five corn-soybean rotations from 1995 through 2004 are simulated. All areas with corn–soybean rotation received a loading of 1.36 kg ha–1 of metolachlor preceding the corn crops. This resulted in a total load of 6.8 kg ha–1 for the 10 yr except for the intensive study field which received only four applications for a total of 5.44 kg ha–1. The pattern shows the location of corn and soybean fields, with greater leaching in areas with thinner unsaturated zones, and also greater residuals for areas that were in corn in 2004. Units are kg ha–1. White areas indicate values of zero.

 

Figure 10
View larger version (35K):
[in this window]
[in a new window]

 
Fig. 10. Pesticides and metabolites to leach to ground water as simulated by the distributed LEACHM model compared to that predicted by the two-term regression model. [Abbreviations: ATR, Atrazine; HYA, Hydroxyatrazine; DEA, Deethylatrzine; DIA, Deisopropylatrazine; DDA, Didealkylatrazine; MET, Metolachlor; MES, Metolachlor ESA; MOX, Metolachlor OXA; GLY, Glyphosate; AMP, AMPA; SIM, Simazine; and PAQ, Paraquat.]

 

    Conclusions, Implications, and Limitations
 TOP
 NOTES
 ABSTRACT
 INTRODUCTION
 Materials and Methods
 Results and Discussion
 Conclusions, Implications, and...
 REFERENCES
 
The fate and transport of pesticides in the unsaturated zone of Morgan Creek watershed was simulated with a combination of analytical and numerical models constructed to be consistent with observations. The LEACHM model was distributed throughout the watershed to account for variations in land use and unsaturated zone thickness as determined by MODFLOW. The calibration exercises indicated that water retention parameters described for soil textural classes by Meyer et al. (1997) were preferable to those estimated by Rosetta (Schaap et al., 2001) from measured soil texture and laboratory measurements of bulk density. Soil property statistics described by Meyer et al. (1997) are being used in forward simulations with pesticides of interest to predict mass loading for this and other NAWQA watersheds.

Little to no leaching of highly sorptive, short-lived compounds was predicted even in the model simulating a 1 m thick unsaturated zone. Some of these compounds, however, have persistent metabolites that were predicted to leach substantial masses to ground water. In particular, the metabolites metolachlor OXA and metolachlor ESA were simulated to leach 0.3 and 0.2% of the 15 metric tons of metolachlor applied in the watershed for the period 1995 to 2004 (Table 4 ). DEA, a metabolite of atrazine, was predicted to leach 0.2% of the 23 metric tons of atrazine applied. A regression model with interaction terms between depth, half-life, and organic carbon sorption coefficient provides a reduced form of the process-based model that can estimate pesticide leaching to ground water as a percent of parent compound loadings in areas of variable unsaturated zone thickness. The metamodel requires far fewer inputs but agrees well with LEACHM predictions and is potentially transferable to other pesticides applied in similar settings. Improvements to the models can make the coefficients dimensionless by relating depths, half-lives of production and decay, and KOC to some standardized residence time in the unsaturated zone. This application is a proof of concept to be applied to four other NAWQA-studied agricultural watersheds in California, Washington, Nebraska, and Indiana.


View this table:
[in this window]
[in a new window]

 
Table 4. Total masses of pesticides and metabolites applied, produced, degraded, retained, and leached for the 32-km2 Morgan Creek watershed for the period 1995–2004 as simulated with LEACHM and mapped to soybean-corn crop areas with variable unsaturated zone thickness. Pesticide abbreviations are defined in Table 1.

 

    ACKNOWLEDGMENTS
 
The authors would like to thank the farmers in the Morgan Creek Watershed for access and assistance in documenting the agricultural management practices and transport of agricultural chemicals through their landscape. We also are indebted to Dave Legates of the Univ. of Delaware and Gordon Heisler of the USDA Forest Service for providing additional climate data. John Hutson, the creator and developer of LEACHM, provided model guidance, source code and model documentation.


    NOTES
 TOP
 NOTES
 ABSTRACT
 INTRODUCTION
 Materials and Methods
 Results and Discussion
 Conclusions, Implications, and...
 REFERENCES
 
All rights reserved. No part of this periodical may be reproduced or transmitted in any form or by any means, electronic or mechanical, including photocopying, recording, or any information storage and retrieval system, without permission in writing from the publisher.


    REFERENCES
 TOP
 NOTES
 ABSTRACT
 INTRODUCTION
 Materials and Methods
 Results and Discussion
 Conclusions, Implications, and...
 REFERENCES
 




This article has been cited by other articles:


Home page
J. Environ. Qual.Home page
P. D. Capel, K. A. McCarthy, and J. E. Barbash
National, Holistic, Watershed-Scale Approach to Understand the Sources, Transport, and Fate of Agricultural Chemicals
J. Environ. Qual., May 1, 2008; 37(3): 983 - 993.
[Abstract] [Full Text] [PDF]


This Article
Right arrow Abstract Freely available
Right arrow Figures Only
Right arrow Full Text (PDF) Free
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Citing Articles
Right arrow Citing Articles via HighWire
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Webb, R. M. T.
Right arrow Articles by Linard, J.
Right arrow Search for Related Content
PubMed
Right arrow Articles by Webb, R. M. T.
Right arrow Articles by Linard, J.