Assessing the impact of regional rainfall variability on rapid pesticide leaching potential

Journal of Contaminant Hydrology 113 (2010) 56–65

Contents lists available at ScienceDirect

Journal of Contaminant Hydrology j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / j c o n h yd

Assessing the impact of regional rainfall variability on rapid pesticide leaching potential Gavan McGrath a,⁎, Christoph Hinz a, Murugesu Sivapalan b,c a

School of Earth and Environment, The University of Western Australia, M087, 35 Stirling Highway, Crawley 6009, Australia Departments of Geography and Civil and Environmental Engineering, University of Illinois at Urbana-Champaign, 220 Davenport Hall, 607 South Mathews Avenue, Urbana IL 61801, USA c Department of Water Management, Faculty of Civil Engineering and Geosciences, Delft University of Technology, Delft, Netherlands b

a r t i c l e

i n f o

Article history: Received 27 June 2009 Received in revised form 19 December 2009 Accepted 24 December 2009 Available online 6 January 2010 Keywords: Solute transport Preferential flow Runoff Risk assessment Pesticide Rainfall Seasonality Western Australia

a b s t r a c t The timing and magnitude of rainfall events are known to be dominant controls on pesticide migration into streams and groundwater, by triggering rapid flow processes, such as preferential flow and surface runoff. A better understanding of how regional differences in rainfall impact rapid leaching risk is required in order to match the scale at which water regulation occurs. We estimated the potential amount of rapid leaching, and the frequencies of these events in a case study of the southwest of Western Australia, for one soil type and a range of linearly sorbing, first order degrading chemicals. At the regional scale, those chemicals with moderate sorption and long half lives were the most susceptible to rapid transport within a year of application. Within the region, this susceptibility varied depending upon application time and seasonality in storm patterns. Those chemicals and areas with a high potential for rapid transport on average, also experience the greatest inter-annual variability in rapid leaching, as measured by the coefficient of variation. The timing and frequencies of rapid leaching events appeared to strongly relate to an area's relative susceptibility to rapid leaching. In the study region the results also suggested that frontal rainfall dominates rapid leaching along the western and southern coasts while convective thunderstorms play a greater role in the arid east. © 2009 Elsevier B.V. All rights reserved.

1. Introduction The regulation of pesticide use is often conducted at national and regional scales. In order to better manage the use of pesticides for the mitigation of environmental impact it is important to understand the spatial and temporal variabilities of the risk of off-site movement of these chemicals. Data in this regard is often lacking, particularly at regional scales, hence the need for models to assess leaching risk. It is often acknowledged that the timing and magnitude of rainfall events that have the potential to trigger rapid flow

⁎ Corresponding author. Tel.: +61 8 64883735; fax: +61 8 648810540. E-mail addresses: [email protected] (G. McGrath), [email protected] (C. Hinz), [email protected] (M. Sivapalan). 0169-7722/$ – see front matter © 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.jconhyd.2009.12.007

processes, such as surface runoff and preferential flow, are significant controls on pesticide transport (Flury, 1996; Kladivko et al., 2001). The antecedent soil moisture, rainfall intensity, storm duration and storm rainfall depth have also been shown to be controls on preferential flow (Beven and Germann, 1982; Edwards et al., 1993; Heppell et al., 2002) and on surface runoff (Sivapalan et al., 1987). The impact of their combined effect when considering multiple storms is less well understood. Recent research has begun to focus on the rainfall controls on rapid leaching risk, however, this work has largely assumed a statistically stationary climate (McGrath et al., 2008a, 2008b). At larger space and time scales, storm properties are likely to vary seasonally and also spatially. There is therefore the potential for significant differences in the temporal (seasonal) and spatial properties of rainfall that impact on the leaching potential of pesticides.

G. McGrath et al. / Journal of Contaminant Hydrology 113 (2010) 56–65

In the context of regional pesticide risk assessments (Corwin et al., 1997; Soutter and Musy, 1997; Stewart and Loague, 1999; Tiktak et al., 2006), the climate is often considered as just a model input as opposed to an intrinsic part of the system behaviour. At small spatial scales we know that rainfall is a significant control on pesticide leaching (Leonard, 1990; Flury, 1996; Kladivko et al., 2001). We also know that regional differences in rainfall are important. For example, Tiktak et al. (2004) found that pesticide leaching generally increased with increasing annual rainfall amount. However, areas with a high temporal variability of rainfall were also found to be associated with greater leaching. It is also known that pesticide loading to streams displays a strong seasonal component, due to both the seasonal nature of pesticide application, as well as the seasonality in rainfall and runoff (Larson et al., 1999). Given that rainfall patterns are expected to change with global climate change, there is clearly a need to better understand the role of our current climate in order to be able to infer the impact into the future. Current approaches to regional pesticide risk assessment use several layers of information that usually comprise pesticide–soil properties such as retardation, degradation, pesticide application, hydraulic and soil physical properties controlling water flow. Process based modelling then contributes to the predictions of concentrations or mass fractions leached (Loague et al., 1990; Soutter and Musy, 1997; Stewart and Loague, 1999; Tiktak et al., 2004; Tiktak et al., 2006). These models vary from the simple, yet versatile Attenuation Factor (Rao et al., 1985) which considers steady plug flow leaching, through to the dual domain MACRO model (Jarvis et al., 1994), which uses Richard's equation to simulate water flow in the soil matrix domain and a simple threshold to trigger flow in the macropore domain. Despite the existence of these more complex models, preferential flow is typically not considered in regional risk assessments, due in part to the difficulty of parameterisation and the computational burden (Rao et al., 1988). However, simulation of preferential flow is essential, as without it leaching can be significantly underpredicted (Herbst et al., 2005). There is also little consideration given to the variability of rainfall, or the temporal resolution required to capture the dominant transport processes (Soutter and Musy, 1997; Tiktak et al., 2004). At best it is treated as a historical time series for input to process modelling. Given that rapid flow events are often short lived processes, the temporal resolution of rainfall considered should be similarly fine scale (Blöschl and Sivapalan, 1995). We present an approach that specifically considers the statistical properties of fine temporal resolution rainfall and its impact on rapid transport. By doing this we add a first order control on the triggering of fast flow events that produce significant off-site pesticide movement (Leonard, 1990; Kladivko et al., 2001). We therefore add a prediction of the probabilities of off-site leaching and runoff events to occur in space and time. The probability is quantified here by the mean frequency of event triggering and this is compared to the rapid flow leaching potential. We denote this rapid flow leaching potential as the loading, which in essence is the amount of chemical released to fast flow pathways. Here, we assessed the impact of the seasonality in rainfall properties on the potential loading of pesticides to rapid flow

57

pathways, not the transport that results, as driven by rainfall variations within the Wheatbelt region of the southwest of Western Australia. We conducted numerical simulations of pesticide loading to a generic fast flow pathway. In the text we refer to runoff or preferential flow as examples, but we mean to refer only to an unspecified rapid transport route. While soil properties also vary significantly across the Wheatbelt, this analysis will focus on a single soil system in order to evaluate the magnitude of the impact of spatial differences in rainfall in isolation. We used a stochastic model of natural rainfall as a driver. This driver incorporates spatial and seasonal differences in the timing and magnitude of storm events. The objective of the simulations is to understand how the probabilistic nature of leaching risk is affected by chemical properties of pesticides in the context of the entire region, as well as spatio-temporal patterns of risk across the region due to differences in climate. The paper first presents the rainfall model used in the subsequent simulations. We then provide a simplification of a spatial regression derived previously (Hipsey et al., 2003) to parameterize a rainfall model. Because we focus on the impact of rainfall we present a simple water balance model for fast flow event triggering. We expect that other, perhaps more appropriate infiltration models could be easily substituted in the future, depending upon the researchers' needs. In the results we focus on how the pesticide retardation and degradation half-life impact upon the statistical properties of the regional loading as well as their spatial patterns. We also evaluated how the patterns of triggering vary spatially and seasonally. 2. Methodology 2.1. Rainfall modelling Rainfall is modelled as a series of independent storms with three basic properties: a storm duration tr [T]; an average intensity p [L T− 1]; and a time between each storm tb [T] (Fig. 1). The inter-storm duration however has a minimum time between storms of 7 h (Robinson and Sivapalan, 1997) and so we define an additional inter-storm time as t′b [T] such that tb = 7 + t′b. The storm properties tr, p and t′b are assumed to be independent and identically distributed random variables, following the exponential probability distribution

Fig. 1. Characterisation of rainfall events adapted from Struthers (2004). Storms have a duration tr, an average intensity p̄ and a time between events t′b =tb + 7. The average rainfall intensity is disaggregated to a variable intensity within each storm.

58

G. McGrath et al. / Journal of Contaminant Hydrology 113 (2010) 56–65

(pdf) in each case, such that their pdfs are given by (Robinson and Sivapalan, 1997; Hipsey et al., 2003): fTr =


 1 t Exp − r tr tr

ð1Þ


 1 p Exp − p p
 1 t′ = – ′ Exp − –b′

ð2Þ

fP = fTb′

t

b

t

ð3Þ

b

where tr̄ [T], p̄ [L T− 1] and t′̄ b [T] are the seasonally varying mean storm properties. The bounded random cascade approach (Menabde and Sivapalan, 2000) is used to distribute the rainfall randomly within the storm event at a temporal resolution of 6 min with parameters adapted to rainfall from the southwest of Western Australia (Hipsey et al., 2003). We assume that these cascade parameters are also constant in time and space though we acknowledge that currently this has not been assessed. Seasonality in storm properties is achieved by applying a sinusoidal function to generate particular mean storm properties for any day of the year via (Robinson and Sivapalan, 1997; Hipsey et al., 2003): ð4Þ

x = a cos½bðt−cÞ + d

where t [T] is the month of the year, x̄ is the mean storm variable (tr̄ , p̄ or tb̄ ) and a, b [T− 1], c [T] and d are empirical parameters derived from data by regression analysis. The units of x,̄ a and d are the same as the storm variable. Hipsey et al. (2003) conducted a spatial regression of these seasonality parameters for the southwest of Western Australia. The regression contained a number of predictor variables including latitude, distance from the coast, climatic aridity index (AI = Epann/Pann [−] where Pann [L T− 1] denotes the mean annual rainfall) and annual pan evaporation Epann [L T− 1]. In this paper we reduce the number of variables for this regression to just spatial coordinates by deriving spatial regressions for aridity index and annual pan evaporation. Monthly pan evaporation and rainfall data were obtained for 28 sites throughout the southwest of Western Australia from the Bureau of Meteorology (Bureau of Meterology, 2006). The distribution of these sites is shown in Fig. 2. For water balance modelling, seasonality in daily potential evaporation was also calculated at each of the above sites by fitting Eq. (4) to monthly variations in mean daily potential evaporation. Based upon the polynomial regression model used by Hipsey et al. (2003) we assumed a general spatial regression for Epann, AI and the seasonality parameters of daily potential evaporation (i.e. a, b, c and d) of the form: 2

2

2

z = c1 + c2 E + c3 N + c4 EN + c5 E + c6 N + c7 E N

2

ð5Þ

where E is the longitude in degrees east, N is the latitude in degrees north and z is the value of the variable under consideration at its particular spatial location. The one remaining variable, other than latitude and longitude, i.e. the distance from the coast, was determined numerically using latitude and longitude by calculating the minimum

Fig. 2. Modelled region of the southwest corner of Western Australia (a) and a close up (b) showing the simulated region (gray area), the location of climate stations used in the regression analysis, and the letters correspond to locations referred to in Fig. 4.

spherical distance between a point and a set of predefined curve segments which approximate the coastline of the southwest. 2.2. Pesticide transport simulation Preferential flow and solute loadings were modelled in a similar manner as in McGrath et al. (2008a). Rainfall enters a near surface mixing layer containing pesticides. From this the water and chemical are then partitioned into the soil matrix (the slow domain) and if triggered, a generic fast flow pathway, such as preferential flow and/or surface runoff. This partitioning is as far as we consider solute transport, as our interest is in the potential for rapid leaching. In the mixing layer chemicals are assumed to obey linear, equilibrium sorption with first order degradation. When the rainfall intensity exceeds the infiltration capacity of the soil matrix a proportion of rainfall and the same proportion of solute are released to the fast domain.

G. McGrath et al. / Journal of Contaminant Hydrology 113 (2010) 56–65

Mass balance equations for the slow domain and the mixing layer are given by the following: ð1−Wf Þθm zm

ð6Þ

dC Rθl zl l = −θl zl k−p½tCl dt

ð7Þ

where m denotes the slow domain, l the near surface mixing layer, Wf [L2L− 2] the areal fraction of soil occupied by the fast flow pathway, S [−] the fraction of the water filled porosity in storage, θl [L3L− 3] the saturated volumetric water content in the mixing zone, θm the saturated water content in the soil matrix, and zm [L] and zl [L] the depths of the matrix domain and the surface mixing layers, im [L T− 1] the infiltration rate, L [L T− 1] the rate of water losses from storage due to drainage and evaporation, p [L T− 1] the rainfall intensity, k [T− 1] the first order degradation rate, Cl [M L− 3] the solute concentration and finally R = 1 + Kdρb/θl [−] denotes the retardation factor, comprised in part by the equilibrium sorption coefficient (Kd [L3 M− 1]) and the soil bulk density (ρb [M L− 3]). After Struthers et al. (2007a,b), we adopted a simple modelling approach for describing the components of the matrix water balance. Losses from storage are described by: 8 s−sc > > ð1−Wf Þθm zm ep ½t + > > τ > < L½S = ep ½t > > > s−sr > > : ep ½t sm −sr

for Sc bS ≤1 for Sm bS ≤ Sc for S ≤ Sm


 1−Sc : im = min p; τ

ð9Þ

Infiltration into the fast domain is given by: ð10Þ

The mass of the chemical per unit area lost to fast flow during a rainfall event mp [M L− 2] was calculated via: tr

mp = ∫ C½tif dt: 0

Slow domain Field capacity storage Wilting point storage Residual storage Drainage response time Saturated water content Areal proportion Soil depth Near surface mixing cell Volumetric water content Mixing cell depth

Symbol

Units

Value

Sc Sm Sr τ θm Wf zm

– – – h mm3 mm− 3 m2/m2 mm

0.65 0.45 0.25 100 0.45 0.05 80

θl zl

mm3 mm− 3 mm

0.45 10

The parameters used in the simulation are adapted from Struthers et al. (2007a) to represent a moderately fast draining soil, and are summarized in Table 1. A thin soil depth, zm = 80 mm, was chosen in order to represent the highly variable water balance dynamics near the soil surface. Similarly small discretizations are commonly found in more complex modelling schemes. The choice of parameterization will of course impact on the results but we reiterate that the purpose here is to conduct a preliminary or first order analysis of the rainfall controls in isolation. 2.3. Limitations

ð8Þ

where Sc [−] describes the slow domain field capacity, τ [T] a representative drainage time scale, Sm [−] the wilting point, Sr [−] a residual storage and ep [L T− 1] the potential evaporation rate. During storms evaporation is assumed negligible. The infiltration rate into the slow domain is given by:

if = max½0; p−im 

Table 1 Model parameters used in the simulation. Variable

dS = im −L½S dt

59

ð11Þ

The above water balance model is intended to describe the partitioning of water between fast and slow flow pathways. One of the dominant soil profiles in the region is duplex soil, which is characterized by a fine to coarse soil layer overlying a compact clay layer (McArthur, 2004). The threshold trigger for rapid flow, inherent in the water balance model above, aims to represent the intermittent filling of this layer and the triggering of a generic fast flow process such as surface runoff or preferential flow.

There are a number of limitations with the above approach. The simple capacitance model cannot simulate upward water flow, which may be important for controlling the antecedent soil moisture, particularly during the warmer months in this semiarid climate. Additionally, the fixed infiltration threshold may be too simplistic to accurately capture the intermittent triggering of preferential flow for some soils. Finally, the model does not consider lateral redistribution of water through the landscape which may be relevant for the initiation of overland flow. We neglect a number of processes that may be important for modelling the actual concentration of pesticide leaching to streams and groundwater. Most significantly, we make no distinction between types of rapid flow processes, just that the chemicals are loaded to a generic fast flow pathway. We neglect details of the subsequent transport. So, for example, while it is known that in tile drain systems, surface runoff typically delivers chemicals at much higher concentrations than occurs via preferential flow (Kladivko et al., 2001), our model would describe equivalent loading. Non-equilibrium kinetic sorption and desorption processes may be important for the preferential transport of some chemicals, in particular soils. Additionally, we neglect soil temperature and soil moisture fluctuations in the mixing zone, which may in turn influence degradation rates, chemical solubility and sorption processes among others. As implied earlier, we also do not consider the details required to simulate runoff, and so slope, vegetation, soil texture or structure properties related to erosion are not factored into quantifying actual off-site chemical fluxes. However, for pesticide risk assessment the data required to identify and then parameterise these processes, not to mention their natural heterogeneity within the landscape, is often lacking. For that reason, the simplest models are often

60

G. McGrath et al. / Journal of Contaminant Hydrology 113 (2010) 56–65

deemed to have the greatest utility for regional risk assessment (Rao et al., 1988). 2.4. Chemical properties A total of 117 different combinations of retardation factor and degradation half-life, representing different pesticides were chosen for analysis to evaluate the variability of transport of a suite of chemicals within the region. Values for R were chosen from the range of 1 to 103. Values of k were chosen from the range of 10− 4 to 1 (day− 1). Pesticides were considered to be applied once a year although we acknowledge that in many circumstances these chemicals may actually be applied multiple times and in combination. This level of complexity is beyond the scope of this first order analysis. However we do consider the impact of different application times by modelling the date of application as occurring on either the first of March or September. The choice of application times is by no means a comprehensive list of the distribution of pesticide application times throughout the region. The March application date was chosen to represent a pre-sowing or knockdown type of application, while the September date represents a post-emergent spray, as is common practice in no-till cereal production.

Fig. 3. Comparison of modelled and observed aridity index.

Table 2 Results of spatial regressions. Variable

Symbol

r2

Mean square error

Aridity index Mean annual potential evaporation Seasonal variation in daily potential evaporation (Eq. (4))

AI Epann a b c d

0.97 0.99 0.97 0.84 0.88 0.99

0.23 2180 mm2 year− 2 0.042 mm2 day− 2 4 × 10−5 month− 2 6.8 × 10−3 month2 0.018 mm2 day− 2

2.5. Simulation and analyses A 1°× 1° grid, representing the southwest corner of Western Australia, was chosen for simulation. Fig. 2b shows the modelled grid. The latitude and longitude of the center of 32 of these grid cells were used as inputs to the spatial regressions (Eq. (5)) in order to generate the seasonality parameters (a, b, c, d) for the calculation of the monthly mean properties (Eq. (4)). Based on these means, random rain events were generated using Eqs. (1)–(3), and then disaggregated using the methodology of Menabde and Sivapalan (2000) to a resolution of 6 min. Each grid cell was assumed to be independent of the other and as such there was no accounting for spatio-temporal correlations of rainfall among cells as might be expected of actual rainfall. For each combination of location, type of chemical and time of application, 500 one year long simulations were conducted. Soil moisture was initialized by running the model for five years prior to the simulated pesticide application. Several variables from the modelled output were chosen for analysis. These include the frequency of preferential flow event triggering, statistics of the magnitude of annual loading to fast flow pathways, both at a point and regional totals, and the monthly distribution of loading. For simplicity we assumed equal areal contribution of each modelled grid to regional statistics.

Table 3 Regression equations for the seasonality parameters to estimate the mean storm statistics from Hipsey et al. (2003) and modifications since (Matthew Hipsey, 2006, pers. comm.) and additional spatial regressions to estimate seasonality in daily potential evaporation, annual potential evaporation and aridity index. Statistic Storm duration (tr̄ )

b

Storm intensity (p̄)

Inter-storm period (t̄b)

Daily potential evaporation

c d a b c d a b c d a b

3. Results and discussion

c

3.1. Spatial regressions There was generally excellent agreement between the spatial regression and the observed potential evaporation, the aridity index and the seasonality in daily potential evaporation. Fig. 3 for example shows comparisons of the regressed and observed aridity index. Two sided tests of all regressions were found to be significant at the p b 0.01 level. Table 2 gives details

Regressions a

d Aridity index

AI

Annual potential Epann evaporation

−2.791 × 10− 5C2 + 0.014C − 0.0281 for C N 55; and 0.0277C + 2.555 for C b 55 0.285 for C b 95; 0.325 for C N 160; and 0.3225 otherwise − 60 − 4.661 × 10− 5C2 + 0.0198C + 3.914 4.040 × 10− 3N2 + 0.239N + 3.775 0.507 for AI N 4.5; and 0.285 for AI b 4.5 − 60 0.0106AI2 + 0.01557AI + 0.402 1.342N2 − 69.1265N − 738.07 0.415 for N b 33.5; and 0.5(−6.472 × 10− 3N + 0.7381) for N N 33.5 0 0.0856Epann − 62.97 − 1156.7 + 23.19E − 0.1169E2 + 30.1N − 0.289NE + 0.1871N2 − 1.879 × 10− 5N2E2 49.15 − 0.6052E + 1.702 × 10− 3E2 + 1.478N − 1.203 × 10− 2EN + 1.137 × 10− 2N2 − 7.537 × 10− 7E2N2 625.3 − 7.308E + 1.763 × 10− 2E2 + 21.29N − 1.763EN + 1.527N2 − 1.045 × 10− 5E2N2 − 2429 + 35.84E − 1.272E2 − 38.57N + 0.3296EN − 0.3083N2 + 2.137 × 10− 5E2N2 12,424 − 122.4E + 0.1494E2 + 680.3N − 5.754EN + 5.491N2 − 3.947 × 10− 4E2N2 − 629,727 + 9423.8E − 33.41E2 − 9584N + 86.14EN − 74.39N2 + 5.612 × 10− 3E2N2

N = latitude (degrees north); E = longitude (degrees east); C = distance from coast (km); Epann = annual potential evaporation (mm/year); AI = aridity (dryness) index = Mean Epann / Mean annual rainfall.

G. McGrath et al. / Journal of Contaminant Hydrology 113 (2010) 56–65

of the goodness of fit and mean square errors of the adopted model. Table 3 summarises the regressions used to model the spatio-temporal variation in rainfall and potential evaporation. While these regressions have been conducted to make the (Hipsey et al., 2003) storm model more operational, there are a number of limitations of this approach that should be noted. The first is that we have not evaluated the impact of our regressions on the resulting storm statistics. Given the excellent agreement between data and model we can be reasonably confident that our simplification does not impact the original model's predictions significantly. However, as can be seen (Fig. 3) a negative aridity index is predicted for one site, Albany, which is the most southerly climate station in the data set. Clearly this is not physically realistic and therefore, care should be taken in applying the regression model at the extremities of the region. Problems with non-physical parameters were not found for any of the locations used in the water balance modelling. Fig. 4 shows simulated mean storm properties and daily potential evaporation for selected locations in the modelled domain. These properties are expected of the region in general (Bureau of Meterology, 2006). It can be seen that in the north and east the mean storm intensity is highest in the summer months, peaking in February, while storm duration, inter-storm duration and evaporation tend to peak in early January. Along the mid-western to south-central coast the mean storm intensity is highest in June. Throughout the region the mean storm duration is of a similar magnitude and changes very little seasonally, while there are significant increases in potential

61

evaporation to the north and east. The mean storm intensity and the mean inter-storm duration display a similar seasonal pattern and also change spatially in a fashion similar to potential evaporation. 3.2. Regional loading In order to assess how chemical properties impact pesticide loading at a regional scale, for each chemical the mean annual load at the 32 sites was used to determine the regional average and the spatial coefficient of variation of annual loads. For the suite of pesticides used in the simulations Fig. 5 shows contours of these statistics as a function of the retardation factor and degradation rate for a March application. Similar types of relationships with chemical properties were found for other application times. It can be seen that pesticides with a retardation factor R ∼ 101.5 and long degradation half-lives (low k) pose the highest risk of rapid leaching on a regional average basis. In contrast, rapidly degrading (high k) and weakly sorbing (low R) chemicals pose the lowest risk of rapid transport. The peak in rapid leaching potential as a function of sorption is consistent with the results of Larsson and Jarvis (2000). They found a similar peak in simulated leaching with preferential flow relative to simulations conducted without it. Predictions of leaching potential made by the Attenuation Factor (AF) (Rao et al., 1985) and many other matrix-only flow models (Herbst et al., 2005) would suggest instead that weak sorption and long half-lives pose the greater risk of leaching to groundwater and that this risk decreases with

Fig. 4. Spatial variability of the seasonality in simulated storm properties. Shown are the seasonal variations in mean daily potential evaporation and mean storm properties including, storm duration, storm intensity, and inter-storm duration for 8 sites (a to h) corresponding to the locations shown on the map in Fig. 2b.

62

G. McGrath et al. / Journal of Contaminant Hydrology 113 (2010) 56–65

The spatial coefficient of variation shows a similar relationship with chemicals properties as did the regional average. The greatest spatial variation in annual loads occurs for a range of moderate to strongly sorbing solutes and low degradation rates. The least spatial variation occurs for weakly sorbing, rapidly decaying solutes. This result is interesting because it suggests that those chemicals which are more susceptible to leach on average are also the ones which are likely to display the greatest variation in loading throughout the region. There is ample evidence to support the hypothesis that the source zone for rapid pesticide transport is often near the soil surface (Leonard, 1990; Allaire-Leung et al., 2000). Therefore, not only do soil and chemical characteristics control rapid leaching, but rainfall, which directly impacts the surface layer, will play a significant role too. That role is best expressed for moderately adsorbed chemicals which tend to show greater variability in rapid leaching, due again to their relative flexibility. The randomness in the timing and magnitude of rainfall triggering rapid flow events coupled with the rainfall infiltrating the slow domain causes the higher flexibility chemicals to also display greater variability. 3.3. Spatial patterns of loading

Fig. 5. Effect of pesticide properties on regional annual loads. Shown are contours of the regional average (a) and the regional coefficient of variation (b) of the mean annual pesticide loading Log10(loading as a fraction of that applied) as a function of chemical properties, the retardation factor R and the degradation rate k (Log10 day− 1).

increasing sorption strength. It should be remembered however, that we are only considering the rapid leaching potential and not slow flow through the soil matrix. The peak with respect to sorption, found by Larsson and Jarvis (2000), was expressed by them as the ratio of total leaching from both slow and fast domains to total leaching when no preferential flow was simulated. We have previously argued (McGrath et al., 2008a, 2009) that rather than just being the result of the occurrence of preferential flow the peak in loading is caused by a balance between the retention of a chemical at the place where it is available for rapid transport and the release of the chemical to rapid flow processes, when they occur. Weakly adsorbed chemicals are primarily lost to slow domain transport prior to the triggering of runoff events. Very strongly adsorbed chemicals are instead held so tightly that they are largely not available for transport via either pathway (particle-facilitated transport neglected). This leaves moderately sorbing chemicals with greater flexibility to rapidly leach on average.

The patterns of average annual loading and the inter-annual coefficient of variation of loading as a function of chemical properties are reflected in Fig. 6. This figure shows the spatial distribution of rapid leaching statistics for March and September application times and for nine example pesticides (i.e., nine combinations of R and k). Pesticides with intermediate sorption strengths and low degradation rates show greater average annual rapid leaching potential (see Fig. 6a, c) than those which more strongly or more weakly sorb. Additionally, those chemicals which are more prone to rapid leaching also have a greater coefficient of variation in annual loads, both at a point and throughout the region (see Fig. 6b, d). The impact of different application times further emphasises the interaction of chemicals' properties and climate. A March application occurs at the onset of a cooler, wetter period with relatively frequent frontal rainfall which extends throughout the entire region. A September application on the other hand, occurs at a time when the climate warms and dries, and intense thunderstorms become more prevalent in the interior. For a March application, loading is generally elevated along the coast as well as throughout the interior. In comparison, a September application has lower loading along the central and northern coast as well as the southeast coast. There is also a greater loading potential of weakly to moderately sorbing chemicals along the southwest coast for a September application. A band of lower loading potential is also evident in a March application, which extends through the center of the region, from the northwest through to the southeast. This band appears to broaden and shift southwest for a September application. These patterns can be explained by the balance between retention of chemical in the mixing layer and its release to rapid flow as moderated by climate. The release is in part determined by the frequency of rapid flow events which we show in Fig. 7. The reduction in the loading potential along the central and northern coasts can be explained by the relatively few such events in the spring and summer months. Along the south coast there are still a significant number of

G. McGrath et al. / Journal of Contaminant Hydrology 113 (2010) 56–65

63

Fig. 6. Spatial patterns of the mean ((a) and (c)) and coefficient of variation ((b) and (d)) of the annual pesticide load log10(fraction of applied) to fast flow pathways for a March and September application as a function of pesticide properties the retardation factor R and the degradation rate k (day− 1).

events in spring (September through to November). Therefore, rapid flow events occur often enough for weakly sorbing chemicals to be available for rapid transport at this time. For the March application, however, the small amount of rainfall, which occurs before the soil is primed to initiate preferential flow, is sufficient to leach the chemical into the soil matrix and away from the location where it can be rapidly mobilised. The band of low loading is consistent with a band of low event frequencies that remain low throughout the year. To the west of this band, events occur frequently during the cooler

months. To the east, events occur more frequently than the band in winter and there is also a greater probability of summer events to occur. In the eastern Wheatbelt rainfall is more conducive to rapidly leaching weakly sorbing solutes as it occurs infrequently but when it occurs it is often higher in intensity (see Fig. 4). The importance of the seasonality in rainfall on the distribution of pesticide loading is further illustrated for the most southerly site in the region (Fig. 8). Shown is the monthly loading for a suite of chemicals applied in March, differentiated

64

G. McGrath et al. / Journal of Contaminant Hydrology 113 (2010) 56–65

Fig. 7. Spatial distribution of the monthly frequency of fast flow events, the x axis of each plot denotes the month of the year beginning with January, and the y axis the mean number of rapid flow events per month. The vertical scale ranges from 0 to 10 and is cropped for clearer presentation.

by their strength of sorption. As expected, weakly sorbing solutes are rapidly dissipated and become unavailable for rapid transport within a few months after application. As the retardation factor increases there is a tendency for the peak in loading to be shifted from March toward July. While this is true too for the most strongly sorbing solutes, they tend to leach more similar amounts throughout the year. The timing of the peak in loading suggests that streams and groundwater are susceptible to relatively high concentrations of a number of different chemicals at the same time. This is consistent with observations by Larson et al. (1999) who noted the occurrence of seasonal pulses of high concentrations of several pesticides in the months following typical pesticide application. These results suggest that rapid leaching of weakly sorbing and very strongly sorbing chemicals is not necessarily strongly impacted by seasonal variations in climate. Rather, weakly sorbing chemicals will be strongly influenced by the structure of rainfall shortly after application, while strongly sorbing chemicals will be influenced by the annual average rate of leaching events. Intermediately sorbing chemicals, on the other hand, appear to be strongly controlled by seasonality in rainfall and evaporation. It is these seasonal differences in the timing and magnitude of rainfall events which may be responsible for the greater coefficient of variation of loading for these chemicals.

Fig. 8. Mean monthly pesticide load as a function of the retardation factor for the most southerly location simulated. The gray scale denotes Log10 of the average proportion of the applied mass loaded to fast flow pathways for a March application with k = 10−2 day− 1.

G. McGrath et al. / Journal of Contaminant Hydrology 113 (2010) 56–65

4. Conclusions Our results suggest that the timing and frequency of rapid flow events may provide a useful measure of the rapid transport potential of surface applied chemicals. Given the added uncertainty of predicting the actual solute flux, this may provide a reasonable measure of how “hot” a particular location is in terms of its susceptibility to leaching for example. Additionally, the timing and frequency can also give an indication of when “hot” leaching moments are likely to occur. Further research is required to assess this measure as an alternative to flux prediction of leaching risk assessment. The slow flow domain is expected to be the dominant pathway for weakly adsorbing, persistent pesticides. Our approach could be easily extended with existing models to additionally show the spatial and temporal variabilities of transport from the slow domain. In order to evaluate the overall environmental effect however, more specific assumptions about the mechanics of preferential flow and surface runoff would need to be made than we have considered here. Finally, despite the development of a number of process based models for regional leaching risk assessment, there remains a significant absence of validating predictions to actual data collected at these scales. This is essential in order to progress in this area. There is therefore the need to develop cost-effective sampling and analysis methodologies as well as infrastructure in order to fulfill this need. The lack of monitoring programs in many parts of the world, including the southwest of Western Australia, is one clear indication of this need. Acknowledgements The authors would like to acknowledge financial support by the Australian Research Council's Linkage Projects funding scheme (project number LP0211883), Water Corporation of Western Australia and Centre for Groundwater Studies. The third author is grateful to Delft University of Technology for providing facilities and financial support during a stay as visiting professor, during which time this work was completed. References Allaire-Leung, S.E., Gupta, S.C., Moncrief, J.F., 2000. Water and solute movement in soil as influenced by macropore characteristics 1. Macropore continuity. J. Contam. Hydrol. 41 (3–4), 283–301. doi:10.1016/S0169-7722(99)00079-0. Beven, K., Germann, P., 1982. Macropores and water flow in soils. Water Resour. Res. 18, 1311–1325. doi:10.1029/WR018i005p01311. Blöschl, G., Sivapalan, M., 1995. Scale issues in hydrological modelling—a review. Hydrol. Proc. 9, 251–290. doi:10.1002/hyp. 3360090305. Bureau of Meterology, 2006. URL http://www.bom.gov.au/, last accessed 1 June 2009. Corwin, D., Vaughan, P., Loague, K., 1997. Modeling non-point source pollutants in the vadose zone. Environ. Sci. Technol. 31 (8), 2157–2175. doi:10.1021/ es960796v. Edwards, W., Shipitalo, M., Owens, L., Dick, W., 1993. Factors affecting preferential flow of water and atrazine through earthworm burrows under continuous no-till corn. J. Environ. Qual. 22, 453–457. Flury, M., 1996. Experimental evidence of transport of pesticides through field soils — a review. J. Environ. Qual. 25 (1), 25–45. Heppell, C.M., Worrall, F., Burt, T.P., Williams, R.J., 2002. A classification of drainage and macropore flow in an agricultural catchment. Hydrol. Process. 16 (1), 27–46. doi:10.1002/hyp. 282. Herbst, M., Fialkiewicz, W., Chen, T., Pütz, T., Thiéry, D., Mouvet, C., Vachaud, G., Vereecken, H., 2005. Intercomparison of flow and transport models applied

65

to vertical drainage in cropped lysimeters. Vadose Zone J. 4, 240–254. doi:10.2136/vzj2004.0070. Hipsey, M.R., Sivapalan, M., Menabde, M., 2003. A risk-based approach to the design of rural water supply catchments across Western Australia. Hydrol. Sci. J. 48 (5), 709–727. doi:10.1623/hysj.48.5.709.51456. Jarvis, N.J., Stahli, M., Bergstrom, L., Johnsson, H., 1994. Simulation of dichlorprop and bentazon leaching in soils of contrasting texture using the MACRO model. J. Environ. Sci. Health A 29 (6), 1255–1277. doi:10.1080/ 10934529409376106. Kladivko, E.J., Brown, L.C., Baker, J.L., 2001. Pesticide transport to subsurface tile drains in humid regions of north America. Crit. Rev. Environ. Sci. Tech. 31 (1), 1–62. doi:10.1080/20016491089163. Larson, S., Gilliom, R., Capel, P., 1999. Pesticides in streams of the United States: initial results from the national water-quality assessment program. Tech. rep., United States Geological Survey, Water Resources Investigations Report, pp. 98–4222. Larsson, M.H., Jarvis, N.J., 2000. Quantifying interactions between compound properties and macropore flow effects on pesticide leaching. Pest Manag. Sci. 56, 133–141. doi:10.1002/1526-4998(200002)56:2 b133::AIDPS103N3.0.CO;2-N. Leonard, R., 1990. Movement of pesticides into surface waters. In: Cheng, H.H. (Ed.), Pesticides in the Soil Environment. Soil Science Society of America, Madison, USA, pp. 303–349. Loague, K., Green, R., Giambelluca, T., Liang, T., 1990. Impact of uncertainty in soil climatic and chemical information in a pesticide leaching assessment. J. Contam. Hydrol. 5, 171–194. doi:10.1016/0169-7722(90)90004-Z. McArthur, W., 2004. Reference soils of south-western Australia. Department of Agriculture, Western Australia, Perth, Australia. McGrath, G.S., Hinz, C., Sivapalan, M., 2008a. Modelling the impact of withinstorm variability of rainfall on the loading of solutes to preferential flow pathways. Eur. J. Soil Sci. 59, 24–33. McGrath, G.S., Hinz, C., Sivapalan, M., 2008b. Modeling the effect of rainfall intermittency on the variability of solute persistence at the soil surface. Water Resour. Res. 44, W09432. doi:10.1029/2007WR006652. McGrath, G.S., Hinz, C., Sivapalan, M., 2009. A preferential flow leaching index. Water Resour. Res. 45, W11405. doi:10.1029/2008WR007265. Menabde, M., Sivapalan, M., 2000. Modeling of rainfall time series using bounded random cascades and Levy-stable distributions. Water Resour. Res. 36 (11), 3293–3300. doi:10.1029/2000WR900197. Rao, P.S.C., Hornsby, A.G., Jessup, R.E., 1985. Indices for ranking the potential for pesticide contamination of groundwater. Soil Crop Sci. Soc. Florida Proc. 44, 1–8. Rao, P.S.C., Jessup, R.E., Davidson, J.M., 1988. Mass flow and dispersion. In: Grover, R. (Ed.), Environmental Chemistry of Herbicides. CRC Press, Boca Raton, Florida, pp. 21–43. Robinson, J., Sivapalan, M., 1997. Temporal scales and hydrological regimes: implications for flood frequency scaling. Water Resour. Res. 33 (12), 2981–2999. doi:10.1029/97WR01964. Sivapalan, M., Beven, K., Wood, E.F., 1987. On hydrologic similarity 2. A scaled model of storm runoff production. Water Resour. Res. 23 (12), 2266–2278. doi:10.1029/WR023i012p02266. Soutter, M., Musy, A., 1997. Coupling 1D Monte-Carlo simulations and geostatistics to assess groundwater vulnerability to pesticide contamination on a regional scale. J. Contam. Hydrol. 32, 25–39. doi:10.1016/S0169-7722 (97)00075-2. Stewart, I., Loague, K., 1999. A type transfer function approach for regional scale pesticide leaching assessment. J. Environ. Qual. 28, 378–387. Struthers, I., 2004. Climate and soil controls upon the water balance of mine tailings landscapes. Ph.D. thesis, University of Western Australia, 35 Stirling Highway, Crawley, 6009, Western Australia. Struthers, I., Sivapalan, M., Hinz, C., 2007a. Conceptual examination of climate– soil controls upon rainfall partitioning in an open-fractured soil: I. Single storm response. Adv. Water Resour. 30 (3), 505–517. doi:10.1016/j. advwatres.2006.04.006. Struthers, I., Sivapalan, M., Hinz, C., 2007b. Conceptual examination of climate– soil controls upon rainfall partitioning in an open-fractured soil II: response to a population of storms. Adv. Water Resour. 30 (3), 518–527. doi:10.1016/ j.advwatres.2006.04.005. Tiktak, A., de Nie, D.S., Piñeros Garcet, J., Jones, A., Vanclooster, M., 2004. Assessment of the pesticide leaching risk at the pan-European level: the EuroPEARL approach. J. Hydrol. 289, 222–238. doi:10.1016/j.jhydrol.2003.11.030. Tiktak, A., Boesten, J.J.T.I., van der Linden, A.M.A., Vanclooster, M., 2006. Mapping ground water vulnerability to pesticide leaching with a processbased metamodel of EuroPEARL. J. Environ. Qual. 35 (4), 1213–1226. doi:10.2134/jeq2005.0377.