Environmental Pollution 187 (2014) 182e192

Contents lists available at ScienceDirect

Environmental Pollution
journal homepage: www.elsevier.com/locate/envpol

Mercury and methylmercury stream concentrations in a Coastal Plain

watershed: A multi-scale simulation analysis
C.D. Knightes a, *, H.E. Golden b, C.A. Journey c, G.M. Davis a, P.A. Conrads c,
M. Marvin-DiPasquale d, M.E. Brigham e, P.M. Bradley c
a
US Environmental Protection Agency, Office of Research and Development, Ecosystems Research Division, Athens, GA 30324, USA
b
US Environmental Protection Agency, Office of Research and Development, Ecological Exposure Research Division, Cincinnati, OH, USA
c
US Geological Survey, South Carolina Water Science Center, Columbia, SC, USA
d
US Geological Survey, Menlo Park, CA, USA
e
US Geological Survey, Mounds View, MN, USA

a r t i c l e i n f o a b s t r a c t

Article history: Mercury is a ubiquitous global environmental toxicant responsible for most US fish advisories. Processes
Received 28 May 2013 governing mercury concentrations in rivers and streams are not well understood, particularly at multiple
Received in revised form spatial scales. We investigate how insights gained from reach-scale mercury data and model simulations
29 December 2013
can be applied at broader watershed scales using a spatially and temporally explicit watershed hydrology
Accepted 31 December 2013
and biogeochemical cycling model, VELMA. We simulate fate and transport using reach-scale (0.1 km2)
study data and evaluate applications to multiple watershed scales. Reach-scale VELMA parameterization
Keywords:
was applied to two nested sub-watersheds (28 km2 and 25 km2) and the encompassing watershed
Mercury
Exposure
(79 km2). Results demonstrate that simulated flow and total mercury concentrations compare reasonably
Watershed to observations at different scales, but simulated methylmercury concentrations are out-of-phase with
Modeling observations. These findings suggest that intricacies of methylmercury biogeochemical cycling and
Scale transport are under-represented in VELMA and underscore the complexity of simulating mercury fate
and transport.
Published by Elsevier Ltd.

1. Introduction biomagnifies in food webs (Bloom, 1992). Exposure to MeHg may

cause severe human health effects, including immune system
Mercury (Hg) is a ubiquitous global environmental pollutant. In suppression, neurodevelopmental delays in children, compromised
the United States (US), Hg is responsible for most toxicant fish cardiovascular health in adults (Mergler et al., 2007), and wildlife
advisories. In 2010, 66,400 km2 of lakes and 1,770,000 km of rivers effects, including behavioral, neurochemical, hormonal, and
were under advisory due to Hg (USEPA, 2010). Over 4700 water reproductive changes (Scheuhammer et al., 2012, 2007).
bodies in the US are 303(d) listed as impaired due to Hg (http:// Mercury fate and transport studies have advanced under-
water.epa.gov/lawsregs/lawsguidance/cwa/tmdl/mercury/index. standing of processes regulating Hg dynamics in aquatic ecosys-
cfm). Atmospheric deposition is a primary source of Hg to surface tems, including atmospheric deposition to vegetation and soils,
waters, via direct deposition or watershed transport. A large transformations of Hg in the landscape, watershed hydrological
portion of deposition is often derived from sources beyond the transport, and Hg transformations in the water and sediments
jurisdiction of individual states or regions (Fitzgerald et al., 1998; (Branfireun et al., 1996; Chen et al., 2008; Gao et al., 2006; Golden
MPCA, 2007; NEIWPCC, 2007). Biological Hg hotspots occur in re- et al., 2012; Harris et al., 2007; Henry et al., 1995.; Hudson et al.,
gions with no direct industrial impact (Evers et al., 2007). After 1994; Kelly et al., 1995; Knightes et al., 2009; Kotnik et al., 2002;
deposition, Hg may be converted to methylmercury (MeHg) (Benoit Larssen et al., 2008; Macleod et al., 2005; Shanley et al., 2008;
et al., 2002; Eckley and Hintelmann, 2006; Gilmour and Henry, Vette et al., 2002; Wang et al., 2009). Mercury stream concentra-
1991; Hintelmann et al., 2000), which bioaccumulates and tions, however, often relate to out-of-channel processes rather than
processes within the aquatic habitat (Bradley et al., 2011; Brigham
et al., 2009). Recent studies have moved towards simulating Hg fate
* Corresponding author.
and transport at the watershed scale (Futter et al., 2012; Golden
E-mail address: [email protected] (C.D. Knightes). et al., 2012; Golden and Knightes, 2011). Reach-scale studies have

0269-7491/$ e see front matter Published by Elsevier Ltd.

http://dx.doi.org/10.1016/j.envpol.2013.12.026
 C.D. Knightes et al. / Environmental Pollution 187 (2014) 182e192 183

demonstrated that filtered total Hg and filtered MeHg fluxes and 1999; USEPA, 1996, 2001); total suspended sediment (TSS) (Ullrich et al., 2001);
concentrations each have positive correlations with flow and with dissolved organic nitrogen (DON) (Fishman, 1993; Patton and Kryskalla, 2003);
dissolved inorganic nitrogen (DIN) (Fishman, 1993; Patton and Kryskalla, 2003); and
dissolved organic carbon (DOC) concentrations (Bradley et al., DOC (Aiken et al., 1992) were measured, as cited. Soil concentrations for THg, MeHg,
2013; Brigham et al., 2009). We must understand and adequately and total organic carbon (TOC) at the focused reach location (near New Holland site)
represent the watershed-scale governing processes in mechanistic were sampled across the O, A, and B horizons in hilltop, hillslope, riparian zone, and
models to link atmospheric Hg deposition to in-stream Hg expo- wetland soils (Woodruff et al., 2010).
sure concentrations effectively.
Pollutant management strategies are often implemented at
2.3. Watershed model description
broad spatial scales (i.e., state, regional, and national levels). While
The Visualizing Ecosystems for Land Management Assessment (VELMA) model
research is often performed at the field scale, regional and regula-
is a spatially-distributed, eco-hydrological differential mass balance model that
tory actions are taken at larger scales that are too extensive to fully simulates soil water infiltration and redistribution; evapotranspiration, surface and
characterize with limited time and resources. Scaling up is partic- subsurface runoff; C and N cycling in plants and soils; and the transport of DOC, DIN,
ularly important for Hg since limited information is available at and DON from the terrestrial landscape to streams (Abdelnour et al., 2013, 2011).
broader watershed and regional scales (Ali et al., 2013; Bradley VELMA is designed to simulate and evaluate the impacts of changes in climate, land
use, land cover, and forest harvest, regrowth, and fire on hydrological, ecological,
et al., 2013; Clark et al., 2009; Sivapalan, 2003). Mechanistic
and biogeochemical watershed processes (Abdelnour et al., 2013, 2011). Mercury
models can be used to address scaling, because models simulate cycling dynamics were added to simulate divalent inorganic Hg (Hg(II)) and MeHg
future scenarios and improve process understanding. Because (Golden et al., 2012). Previous field studies have demonstrated a correlation of water
model complexity and sensitivity are connected, model develop- column filtered THg and filtered MeHg concentrations with water column DOC
concentrations and streamflow (Bradley et al., 2013; Brigham et al., 2009; Kolka
ment requires a balance between system complexity and under-
et al., 1999; Schuster et al., 2008; Shanley et al., 2008). The Hg module simulates
standing, parameter parsimony, and model usability watershed Hg dynamics, including methylation, demethylation, reduction, sorption
(Lindenschmidt, 2006; Snowling and Kramer, 2001). The combi- on soils and organic matter, and DOC-Hg complexation. VELMA does not simulate
nation of mechanistic modeling with field data improves our un- erosion of soils, so simulated VELMA output is for the flux and concentrations of
derstanding of Hg fate and transport. dissolved Hg(II) and MeHg. VELMA simultaneously simulates hydrology, N, C, and
Hg, thus enhancing VELMA’s utility for simulating Hg exposure concentrations and
We investigate the applicability of reach-scale Hg data and
associated climate and land use/cover change impacts.
model simulations to broader watershed scales. We apply the VELMA’s modeling domain is a three-dimensional matrix covering the topo-
Visualizing Ecosystems for Land Management Assessment (VELMA) graphical surface and four soil layers. It uses a distributed soil column framework to
model, which couples watershed hydrology and biogeochemical simulate lateral and vertical movement of water, heat, N, C, and Hg. The soil column
cycling (nitrogen (N), carbon, (C), and Hg) (Abdelnour et al., 2013, model has three linked sub-models: (1) a hydrological model simulates vertical and
lateral movement of water within the soil, including a variable source area flow
2011) and simulates daily flow, total Hg (THg) and MeHg fluxes, component similar to TOPMODEL (Feaster et al., 2010; Wolock, 1993) and losses of
and stream and soil THg and MeHg concentrations. Our goals were water from the soil and vegetation to the atmosphere; (2) a soil temperature model
to (1) parameterize, calibrate, and simulate hydrology and Hg fate simulates daily ground soil layer temperatures from surface air temperature; and (3)
and transport processes at the reach scale (0.1 km2) and (2) assess a biogeochemistry model simulates C, N, and Hg dynamics. The soil column model is
placed within a catchment framework to create a spatially-distributed model
the utility and ease of application of the reach-scale model to larger
applicable to watersheds and landscapes. Adjacent soil columns communicate
watersheds. We define “parameterize” as the process of assigning through down-gradient water transport. Surface and subsurface lateral flow are
values to model parameters and “calibrate” as the process of routed by a multiple flow direction method. A 30 m digital elevation model (DEM) is
adjusting parameter values to improve model performance. used to determine flow direction and compute flow contribution area.
Without further calibration, the reach-scale model was applied to The governing equations for Hg soil concentration, CHg(II) [g Hg(II) per m3 soil] or
CMeHg [g MeHg per m3 soil] for each soil layer in each grid cell, at a given time within
three watersheds: two headwater sub-watersheds (28 km2 and the catchment, are given by equations (1) and (2), which are solved using forward
25 km2) and the encompassing McTier Creek watershed (79 km2). Euler finite difference approximation with a daily time-step.

2. Methods dCHgðIIÞ
¼ DHgðIIÞ þ dHgðIIÞ þ kd;e CMeHg  km;e CHgðIIÞ  kr CHgðIIÞ þ Jin;HgðIIÞ  Jout;HgðIIÞ
dt
2.1. Site description
(1)
McTier Creek is an Edisto River headwater tributary in the Sand Hills Region of
the Upper Coastal Plain physiographic province in South Carolina, US dCMeHg
¼ DMeHg þ dMeHg  kd;e CMeHg þ km;e CHgðIIÞ þ Jin;MeHg  Jout;MeHg (2)
(Fig. 1) (Bradley et al., 2011, 2010, 2012; Feaster et al., 2010; Golden et al., 2012, 2013). dt
Mercury burdens in largemouth bass (Micropterus salmoides) in the Edisto River
D is the throughfall deposition of Hg(II) or MeHg [g m2 d1]; d is the release of
basin are among the highest in top predator fish in the US (Scudder et al., 2009).
Hg(II) or MeHg via litter decay [g m3 d1]; kd is the effective demethylation rate
McTier Creek watershed has a drainage area of 79 km2 and consists primarily of
constant [d1]; km is the effective methylation rate constant [d1]; kr is the rate
forests with extensive riparian wetlands (Bradley et al., 2011, 2010; Feaster et al.,
constant for reduction and immediate volatilization [d1]; Jin is the net Hg(II) or
2010; Golden et al., 2012; NLCD, 2001; Scudder Eikenberry et al., 2012). The hy-
MeHg influx [g m3 d1]; and Jout is the respective Hg(II) or MeHg outflux
drologic setting of McTier Creek Watershed includes low-gradient surface water
[g m3 d1]. Release of Hg via decay is simulated given the ratio of Hg(II) to C and
channels, perennially inundated wetlands characterized by standing or slow flowing
MeHg to C.
water, and seasonally inundated wetlands; all replenished by groundwater seepage
The demethylation (transformation of MeHg to Hg(II)) and methylation (trans-
(Bradley et al., 2011, 2012). The shallow subsurface soils consist primarily of coarse
formation of Hg(II) to MeHg) processes are biological processes, which effective rate
clean sand; rain infiltrates quickly with little overland runoff. The depth to bedrock
constants are functions of temperature and soil moisture, calculated as:
varies from 0 m to 15 m, characterized by bedrock outcrops in the headwaters and
an increasing depth to bedrock gradient to the watershed outlet.
kd;e ¼ kd $QdTTd $q (3)

2.2. Data compilation

TTm
km;e ¼ km $Qm $q (4)
The New Holland, South Carolina US Geological Survey (USGS) gage 02172305
1
(NWIS, 2010) is located at the outlet of the McTier Creek watershed. Upstream of the where kd is the base demethylation rate constant [d ]; km is the base methylation
New Holland gage is the Monetta, South Carolina USGS gage 02172300 (NWIS, 2010). rate constant [d1]; Qd is the temperature coefficient for demethylation (non-
The two collectively measured streamflow data in our study from 13 June 2007 to 30 dimensional); Qm is the temperature coefficient for methylation (non-dimensional);
September 2009. Surface water samples were collected over a range of flow con- Tm is the base methylation temperature [ C]; Td is the base demethylation temper-
ditions at the New Holland, Monetta, and Gully Creek sampling locations (Fig. 1). ature [ C]; and q is the soil moisture content [mwater
3
m3
soil]. Reducing Hg(II) to Hg(0)
Ultra-trace-level clean-sampling and analysis procedures were used (Bradley et al., with subsequent volatilization is simulated as a first-order process; Hg(0) is not
2010; Brigham et al., 2009; Lewis and Brigham, 2004). Concentrations of THg explicitly modeled in this system. Fluxes of Hg(II) and MeHg are governed by
(Brigham et al., 2009; Olson and DeWild, 1999; USEPA, 1996; USGS, variously dated); transport of dissolved species (freely dissolved and DOC-complexed) and linked to
MeHg (Bloom, 1989; DeWild et al., 2002; Horvat et al., 1993; Olson and DeWild, lateral and vertical flow.
184 C.D. Knightes et al. / Environmental Pollution 187 (2014) 182e192

Fig. 1. Study site: McTier Creek Watershed, South Carolina, US. Streamflow is measured for the entire McTier Creek Watershed ends at the New Holland Gage (02172305) and for the
upper west sub-watershed at the Monetta Gage. Stream sample locations include the Monetta Gage, Gully Creek at Old Shoals Road, and the New Holland Gage. The Focused Reach
study area is delineated by the lower right figure with soil sample locations.

Partitioning and complexation are modeled assuming instantaneous equilib- 2.4. Watershed model setup
rium, with Hg(II) and MeHg partitioned between different phases using linear dis-
The McTier Creek VELMA model uses a three-dimensional grid of 30 m  30 m
tribution coefficients for soil, soil organic carbon (SOC), DOC, and the aqueous phase.
cells as the model domain. A watershed boundary was delineated to the New
The distribution fraction, fia , of Hg species, i, to a given phase, a, is given by:
Holland Gage using VELMA’s pre-processor. The areal extents of the two McTier
headwater sub-watersheds were determined by assigning sub-watershed outlets.
Kda
fia ¼ i
½a (5) The hydrologic flow path was determined by a multiple flow direction method,
KdDOC
i
½DOC þ KdSOC
i
½SOCKdSoili ð1  fÞrb þ 1 where flow is fractionally allocated among grid cells via terrain slope.

where a is DOC, SOC, water, or soil; i is Hg(II) or MeHg; Kda i [L kg1] is the partition
coefficient between i and a; f [m3 m3] is porosity, and rb [kgsoil m3] is soil bulk 2.4.1. Model parameterization and calibration
density. Further details on VELMA’s hydrology and biogeochemistry modules can be We assessed how well VELMA represents the McTier Creek system at its current
found elsewhere (Abdelnour et al., 2013, 2011) level of complexity. We performed calibration in a step-wise manner: (1) flow, (2)
 C.D. Knightes et al. / Environmental Pollution 187 (2014) 182e192 185

soil Hg concentrations, and (3) stream Hg concentrations. The simulation period was Table 1
1 January 2001 to 31 December 2009. We used site-specific observations in com- Model parameters and values for the zones and layers of the McTier Creek VELMA
bination with literature values as our parameter starting values. Calibration was Model.
through visual analysis to capture trends and patterns in observed data, while
adjusting parameters within ranges of observed and literature values, emphasizing Parameter Description Zone Layer Value
mechanistic inference and system understanding. Parameters and their values are km Methylation rate Hilltop 1,2,3 0.001 d1
presented in Table 1. Air temperature and precipitation data are from Feaster et al. constanta Hillslope 1,2,3 0.007 d1
(2010). Thicknesses and hydraulic conductivities of the four soil layers were based Riparian 1,2,3 0.01 d1
on the field observed soil layer types and associated literature values and were All 4 0 d1
adjusted to best match observed flows (Table 1) (Bradley et al., 2012). Further hy- kd Demethylation Hilltop 1,2 0.003 d1
drology simulation and evaluation details are presented in Golden et al. (2012). rate constanta 3 0.006 d1
The Hg module was calibrated for the Focused Reach using the New Holland site, Hillslope 1,2 0.015 d1
based on Bradley et al. (2012). We used observed values for initial soil concentrations 3 0.03 d1
of Hg(II) and MeHg (Woodruff et al., 2010) and adjusted the partition and Riparian 1,2 0.015 d1
complexation coefficients based on the range of values presented in the literature 3 0.03 d1
(Allison and Allison, 2005) to maintain constant soil concentrations and to best All 4 0 d1
match observed THg stream concentrations. Methylation and demethylation base kr Reduction rate All 1 0.0001 d1
rate constants and associated parameters were based on field experiments per- constantb All 2,3,4 0
formed on McTier Creek soil and sediment samples and were adjusted to best match Q10,m Methylation All All 1.14
observed stream concentrations at the New Holland site. temperature coefficientc
Mercury loading to watershed soils is simulated as throughfall and litter decay. Q10,d Demethylation All All 1.04
We assume that annual throughfall is equal to annual observed wet deposition temperature coefficientc
(10 mg Hg/m2/yr), based on nearby Congaree Swamp MDN station observations. Tm Base temperature All All 15
Some studies have suggested that throughfall may be higher than open area MDN for methylationc
observations, which could change the values of calibrated parameters, increase loss Td Base temperature All All 20
to the subsurface, or increase simulated Hg flux and stream concentrations. for demethylationc
We divided deposition into 98.5% Hg(II) and 1.5% MeHg (Graydon et al., 2008). Ks,v Saturated vertical All 1 4.37 m d1
Litterfall decay releases C, N, Hg(II), and MeHg to the top soil layer, based on mass- hydraulic conductivityd 2 1.978 m d1
per-mass relationships (Table 1) (Abdelnour et al., 2013; Obrist et al., 2011). Once the 3 3.243 m d1
Focused Reach calibration was complete, we simulated two McTier Creek sub- 4 0.163 m d1
watersheds (Monetta and Gully Creek) and the encompassing McTier Creek Ks,h Saturated horizontal All 1 4.1 m d1
watershed (New Holland) without further adjustment or calibration. Parameter hydraulic conductivityd 2 0.1 m d1
values, references, and descriptions are presented in Table 1. 3 4.1 m d1
4 0.074 m d1
2.4.2. LOADEST: load estimator model z Soil layer thicknesse All 1 0.61 m
We use S-LOADEST (Brigham et al., 2009; Journey et al., 2012; Runkel et al., 2 0.61 m
2004) to compute upper and lower confidence limits (95% UCL/LCL) for MeHg and 3 1.22 m
THg concentrations at the New Holland site to additionally evaluate VELMA simu- 4 Variable, based on
lations. S-LOADEST employs statistical regression analyses based on the rating-curve depth to bedrock
method (Cohn et al., 1992, 1989). The regression model computes daily loads using n Porosityf All All 0.390
relations between constituent load and any of seven explanatory variables that are qfc Field Capacityf All All 0.180
functions of streamflow and time; time functions are represented as increasing or qwp Wilting Pointf All All 0.055
decreasing trends or as seasonal changes. Streamflow and constituent load are log- FHg(II) Hg(II) annual depositiong All 1 9.85 ug m2 yr1
transformed to improve fit. To account for bias produced in back-transformation, the All 2,3,4 0
computed daily loads are adjusted with the minimum variance unbiased estimator FMeHg MeHg annual depositiong All 1 0.15 ug m2 yr1
for cases of no censored data (Cohn et al., 1989) or the adjusted maximum likelihood All 2,3,4 0
method (AMLE) (Cohn, 2005). The load equation used here is: C:N Carbon to Nitrogen Ratioh All All 50
Hg(II):C Hg(II) to Carbon Ratioi All All 38.961 ng g1
L ¼ ßo þ ß1 ln Q þ ß2 ln Q 2 þ ß3 sinð2pTÞ þ ß4 cosð2pTÞ (6) MeHg:C MeHg to Carbon Ratioi All All 0.039 ng g1
Kd,DOC,Hg(II) DOC complexation All All 106.2 L kg1
1
where L is the natural logarithm (ln) of the estimated load [mg d ], Q is the daily coefficient for Hg(II)j
mean streamflow [ft3 s1], T is centered time [decimal years], and bn is the estimated Kd,SOC,Hg(II) SOC partition All All 104.9 L kg1
coefficients for each variable. coefficient for Hg(II)j
Kd,soil,Hg(II) Soil partition All All 104 L kg1
coefficient for Hg(II)j
3. Results Kd,DOC,MeHg DOC complexation All All 106 L kg1
coefficient for MeHgj
The VELMA simulated results and observed data are presented Kd,SOC,MeHg SOC partition All All 103.9 L kg1
coefficient for MeHgj
in Figs. 2e5 for the Focused Reach (0.1 km2), New Holland Kd,soil,MeHg Soil partition All All 102.5 L kg1
(79.4 km2), Monetta (28 km2), and Gully Creek (25 km2) sites, coefficient for MeHgj
respectively. Each figure presents four plots for the period from 1 a
Methylation and demethylation rate constants were based on radiotracer Hg203
June 2007 to 1 October 2009: (a) precipitation, observed stream- and stable isotope Hg200 experiments using McTier Creek, SC soils and sediments.
flow, and simulated streamflow; (b) observed and simulated THg Rates ranged from 0.0003 to 0.0860 d1, mean 0.029 d1 for radiotracer experiments,
stream concentrations; (c) observed and simulated MeHg stream and 0.0003 to 0.011 d1, mean 0.003 d1 for stable isotope experiments. These pa-
concentrations; and (d) observed and simulated DOC, DON, and rameters were calibrated to observed soil concentrations and literature (USEPA, 1997).
b
Reduction rate constant adjusted and based on literature (USEPA, 1997).
DIN stream concentrations. The New Holland gage sampling site c
Temperature coefficients and base temperatures based on site specific studies in
02172305 (NWIS, 2010) was used for the Focused Reach and New SC, GA, and FL.
d
Holland observed Hg concentrations. LOADEST THg and MeHg Hydraulic conductivities based on observed soil types and structures, literature
predictions are presented for New Holland and the Focused Reach values (Freeze and Cherry, 1979).
e
Soil layer thicknesses based on site-specific measurements.
to estimate confidence bounds; there were not enough observa- f
Porosity, field capacity, and wilting point based on default value for soil type
tions for Monetta and Gully Creek to use LOADEST (requires at least (Abdelnour et al., 2013, 2011).
20 observations). Table 2 presents summary statistics for observed g
Annual flux based on observed MDN data at the Congaree Swamp station, MeHg
and simulated data, including root mean square error (RMSE) for flux based on 1.5% of total flux (Graydon et al., 2008).
h
streamflow, THg, and MeHg, and R-squared (R2) for streamflow. The i
Carbon to nitrogen ratio based on VELMA parameterization (Abdelnour et al., 2013).
Hg(II) and MeHg to C ratio (Graydon et al., 2008; Obrist et al., 2011).
simulated data statistics are presented for the entire period of the j
Complexation and coefficients based on literature (Allison and Allison, 2005).
186 C.D. Knightes et al. / Environmental Pollution 187 (2014) 182e192

Fig. 2. Simulated and observed results for the Focused Reach for June 2007 to October 2009: a) simulated streamflow and precipitation; b) THg stream concentration for simulated,
observed and LOADEST Mean, UCL, and LCL; c) MeHg stream concentrations for simulated, observed, and LOADEST Mean, UCL, and LCL; d) simulated DOC, DON, and DIN.

study, as the sampling strategy for the observed data was pur- (9.06 ng L1), maximum (9.88 ng L1) and minimum (0.82 ng L1).
posefully designed to represent the range of environmental The simulated maximum occurs during a single event in the winter
conditions. of 2009, when the maximum is higher than the UCL of LOADEST.
For simulated MeHg (Fig. 2c), VELMA does not capture the
observed temporal pattern, and there are periods when it exceeds
3.1. Focused reach simulations
the LOADEST UCL. The simulated mean (0.136  0.0526 ng L1) is
significantly higher (p ¼ 0.039, df ¼ 47) than the observed mean
For the Focused Reach simulations, VELMA generally predicted
(0.121  0.0581 ng L1), and the range (0.261 ng L1) is slightly
THg concentrations (Fig. 2b) within the LOADEST confidence limits.
higher than the observed range (0.200 ng L1). The simulated
During August 2008 to October 2008, simulated THg results often
minimum (0.060 ng L1) and maximum MeHg (0.321 ng L1) are
exceed the LOADEST UCL and from February to June of 2008 and
higher than the corresponding observed (0.040 ng L1 and
2009, the simulated results are generally lower than the LOADEST
0.240 ng L1).
LCL. The simulated mean THg (1.78  1.26 ng/L) is significantly
lower (p ¼ 0.001, df ¼ 43) than the observed mean THg
(2.78  1.76 ng L1). VELMA simulates several spikes in observed 3.2. New Holland, Monetta, and Gully Creek simulations
THg concentration that arise during precipitation events. High
concentrations of THg are observed in JanuaryeFebruary 2009, We applied the Focused Reach VELMA model parameterization
when neither LOADEST nor VELMA results match observations. The and calibration to VELMA simulations of the larger watersheds,
simulated range (7.10 ng L1), maximum (7.65 ng L1) and mini- New Holland (79.4 km2), Monetta (28 km2), and Gully Creek
mum (0.551 ng L1) THg are less than the observed range (25 km2), without any adjustments, to assess how well processes at
 C.D. Knightes et al. / Environmental Pollution 187 (2014) 182e192 187

Fig. 3. Simulated and observed results for New Holland for June 2007 to October 2009: a) simulated streamflow and precipitation; b) THg stream concentration for simulated,
observed and LOADEST Mean, UCL, and LCL; c) MeHg stream concentrations for simulated, observed, and LOADEST Mean, UCL, and LCL; d) simulated and obeDOC, DON, and DIN.

the reach-scale apply to multiple watershed scales. The simulated February 2009. For both Monetta and Gully Creek, the simulations
daily hydrology compared adequately to the observed, but with captured the spring and summer 2008 and June to August 2007
better results for the smaller watershed (all p < 0.05), Monetta observations, but over-predicted December 2007 and October
(R2 ¼ 0.62), than for New Holland (R2 ¼ 0.46), with similar RMSEs 2008. The simulated mean THg stream concentrations for Gully
for each (Table 2). This compares well to hydrology results for the Creek (2.74  1.75 ng L1) and New Holland (2.68  1.726 ng L1)
New Holland gage for two other watershed models: GBMM are not statistically different (p ¼ 0.10, df ¼ 5; p ¼ 0.37, df ¼ 45) than
(R2 ¼ 0.54) and TOPMODEL (R2 ¼ 0.64) (Golden et al., 2012). Each the observed means (2.28  1.76 ng L1, 2.78  1.028 ng L1). The
case had a delayed rise and fall in the simulated hydrograph mean THg stream concentration for Monetta (2.74  1.74 ng L1)
compared to the observed. For Monetta and New Holland, both was statistically higher (p ¼ 0.001, df ¼ 7) than the observed
simulated mean streamflows (0.729 mm d1 and 0.744 mm d1, (1.58  0.507 ng L1). For simulated MeHg, model simulation was
respectively) for the study period were higher than the observed comparable for the larger watersheds and the Focused Reach. For
mean streamflow (0.63 mm d1, 0.66 mm d1). The simulated New Holland, the simulated stream MeHg concentrations are over-
Monetta maximum streamflow (6.04 mm d1) is lower than the predicted compared to the observed. Given the smaller number of
observed (7.45 mm d1). observations, results are less clear for Gully Creek and Monetta;
For simulated THg, model simulations for the larger watersheds October 2008 (Gully Creek) and December 2007 were over-
compared reasonably well with that for the Focused Reach. For New predicted. The mean simulated MeHg stream concentrations for
Holland, the simulated THg stayed mostly in the bounds of the New Holland (0.206  0.0957 ng L1) and Monetta
LOADEST UCL and LCL. As with the Focused Reach model, the (0.217  0.106 ng L1) are statistically higher (p < 0.001, df ¼ 45;
simulations over-predicted from August 2007 to December 2007 p ¼ 0.024, df ¼ 6), than the observed (0.121  0.0581 ng L1 and
and did not capture the high concentrations from January 2009 to 0.127  0.0987 ng L1), while Gully Creek had no significant
188 C.D. Knightes et al. / Environmental Pollution 187 (2014) 182e192

Fig. 4. Simulated and observed results for Monetta for June 2007 to October 2009: a) simulated streamflow and precipitation; b) THg stream concentration for simulated and
observed; c) MeHg stream concentrations for simulated and observed; d) simulated DOC, DON, and DIN.

difference (p ¼ 0.08, df ¼ 5) between simulated cover similar ranges, with the observed having a higher maximum
(0.210  0.0996 ng L1) and observed (0.142  0.0898 ng L1) MeHg MeHg/THg. Interestingly, the observed MeHg/THg covers a con-
stream concentrations. tinuum of data (0.008e0.14 ng L1), while the simulated fall into
two discrete ranges (0.017e0.045 ng L1and 0.09e0.11 ng L1).
3.3. Simulated versus observed THg and MeHg
4. Discussion
Fig. 6 shows simulated versus observed for THg, MeHg and the
ratio MeHg/THg for New Holland, Gully Creek, and Monetta. The We investigate how well VELMA model structure simulates THg
model’s tendency to over-predict is evident from the large number and MeHg concentrations in associated streams and assess the
of data points above the y ¼ x line, however there are under- potential for model scale-up. Our goal is to improve understanding
predicted results. The high observed concentrations, with associ- of processes that govern Hg cycling dynamics and exposure con-
ated low simulated values for THg in February 2009, are evident in centrations at multiple spatial scales and evaluate their simulation
the lower right corner of Fig. 6a. The simulated THg covers a wider in VELMA. The success of streamflow simulations at two water-
range of values than the observed values - except for the pair of sheds of different size suggests the potential of a single parame-
high observations. For MeHg, many observed and simulated data terization for hydrological questions at larger scales in this region.
pairs are grouped around the y ¼ x, between approximately 0.1 and VELMA simulations of THg concentrations showed comparable
0.25 ng L1, with a cluster of high simulated MeHg concentrations results across different scales. For all scales, simulated THg con-
(0.3e0.45 ng L1) with associated low observed MeHg concentra- centrations did not follow the time varying sinusoidal pattern that
tions (0.04e0.15 ng L1). The observed and simulated MeHg/THg LOADEST predicts, and neither VELMA nor LOADEST completely
 C.D. Knightes et al. / Environmental Pollution 187 (2014) 182e192 189

Fig. 5. Simulated and observed results for Gully Creek for June 2007 to October 2009: a) simulated streamflow and precipitation; b) THg stream concentration for simulated, and
observed; c) MeHg stream concentrations for simulated and observed; d) simulated and observed DOC, DON, and DIN.

capture the irregular pattern of THg concentrations. The observed wetlands and hydrological flowpaths not captured in VELMA. The
THg concentrations are within the bounds of the LOADEST UCL and simulated mean MeHg concentrations were statistically higher
LCL. VELMA-simulated mean THg concentrations are not statisti- than the observed for all scales except Gully Creek. The model’s
cally different than the observed mean THg concentrations for all methylation rate constants came from isotope experiments, which
except Monetta, were the simulated is higher. Neither VELMA nor may represent ideal conditions rather than environmental condi-
LOADEST predicted the high THg concentrations of January to tions, and may thus be higher than field values, which could result
February 2009, which suggests an atypical event that neither in model over-prediction of MeHg. For all scales, VELMA-simulated
model was able to capture. These high concentrations are not mean concentrations were no more than twice the observed values
preceded by a precipitation event, a high flow event, or high and captured the shape and periodicity of MeHg concentrations,
instream DOC concentrations. However, a nearby beaver dam but not the timing.
breach substantially changed the area of the adjacent wetland The phase-shifted results of the MeHg simulations suggest that
during this period. Studies have shown that rapid changes in underlying mechanisms governing MeHg transport to and con-
wetted area may elevate THg concentrations (Bodaly et al., 1997; centrations in the stream are not adequately represented in the
Snodgrass et al., 2000; St. Louis et al., 2001). current model. This suggests that temperature is not the limiting
The VELMA simulations of MeHg concentrations showed an factor in the methylation rate as assumed herein. It is feasible that
increase in RMSE for Gully Creek, Monetta and New Holland sim- there is a substrate that is being depleted that is resulting in a
ulations compared to the Focused Reach, suggesting that the decrease in methylation rate in the summer. Methylation is affected
parameterization and calibration of the Focused Reach for MeHg by many factors: the sulfur cycle (sulfate stimulates methylation,
was not as transferrable to broader spatial scales. This may suggest while sulfide reduces Hg(II) bioavailability for methylation), pres-
that different processes control MeHg production at diverse scales, ence of wetlands (wetlands are known to be a zone of increased
such as signals from zones of MeHg sources and sinks, including methylation), pH (acidity has been correlated to increased MeHg
190 C.D. Knightes et al. / Environmental Pollution 187 (2014) 182e192

Table 2
Summary statistics for streamflow, THg, and MeHg.

Focused reach Gully Creek Monetta New Holland

Observed Simulated Observed Simulated Observed Simulated Observed Simulated

1
Streamflow [mm d ]
Mean 0.678 0.741 0.63 0.729 0.66 0.744
Maximum 3.87 6.60 7.45 6.04 5.02 6.45
Minimum 0.172 0.124 0.11 0.119 0.08 0.129
Range 3.70 6.47 7.34 5.92 4.94 6.32
Standard deviation 0.452 0.732 0.73 0.689 0.60 0.7230
RMSE 0.046 0.051 0.050
R2 0.51 0.62 0.46
THg [ng L1]
Mean 1.78 2.28 2.75 1.58 2.74 2.78 2.68
Maximum 7.65 3.44 8.67 2.38 8.66 9.88 8.66
Minimum 0.551 1.04 0.779 0.96 0.765 0.82 0.764
Range 7.10 2.40 7.89 1.42 7.90 9.06 7.90
Standard deviation 1.26 1.028 1.75 0.507 1.74 1.76 1.726
RMSE 10.6 9.6 12.3 12.3
MeHg [ng L1]
Mean 0.136 0.142 0.210 0.127 0.217 0.121 0.206
Maximum 0.321 0.250 0.485 0.280 0.514 0.240 0.463
Minimum 0.060 0.070 0.083 0.040 0.084 0.040 0.0824
Range 0.261 0.180 0.403 0.240 0.431 0.200 0.381
Standard deviation 0.0526 0.0898 0.0996 0.0987 0.106 0.0581 0.0957
RMSE 0.0093 0.050 0.052 0.23

bioaccumulation), organic matter (affects complexation and perform when scaled up to sub-watershed and watershed scales.
bioavailability of Hg and MeHg), iron (may impact net methylation Two perspectives helped to evaluate these goals: (1) how well the
rate), anoxia (increased methylation), and flooding (linked to model performed compared to observed data, and (2) how well it
increased methylation) (Munthe et al., 2007). These processes are performed when applied at larger scales without further
not currently incorporated in the VELMA algorithms to simulate calibration.
MeHg. The VELMA model was designed to correlate MeHg con- From these applications, we found that the simulated stream-
centrations with DOC. At McTier, however, correlation of MeHg flow model performed comparably at both the smaller (Monetta,
with DOC concentrations is not as strong as at other sites (Bradley 28 km2) and the larger (New Holland, 79.4 km2) watershed scales.
et al., 2013; Brigham et al., 2009), including one farther down- We conclude that within this Coastal Plain setting, VELMA has
stream in the Edisto basin (Bradley et al., 2013). This suggests that potential for scaled-up hydrologic simulation. It performed
the MeHg correlation with hydrology and DOC are not enough to comparably across the scales examined here for THg, but perfor-
fully represent dynamics at our study’s scale. mance decreased with watershed size. VELMA simulates THg and
MeHg concentrations in the range of observations, but does not
5. Conclusions capture the daily temporal signal. Specifically for MeHg, the timing
of the observed sinusoidal time series is not adequately simulated.
We presented a recently developed spatially-explicit, mecha- These results suggest that the existing VELMA THg model could be
nistic, differential mass balance model, VELMA, for simulating hy- adequate for some applications. VELMA may suffice if annual or
drology, N, C, and Hg cycling in McTier Creek, a Coastal Plain seasonal concentrations are required rather than daily concentra-
watershed in South Carolina, US. Our goals were to (1) parame- tions or if one needs to simulate within a factor of 2.
terize and simulate hydrology and Hg fate and transport processes The observed MeHg concentrations exhibit a strong sinusoidal
at the Focused Reach scale using field study data and (2) assess how pattern, captured well by LOADEST, but the sinusoidal pattern does
well the reach-scale parameterization and process formulations not follow temperature. Observed MeHg concentrations increase in

Fig. 6. Predicted versus Observed for all sites: a) THg stream concentrations; b) MeHg stream concentrations; c) MeHg/THg.
 C.D. Knightes et al. / Environmental Pollution 187 (2014) 182e192 191

mid-winter and peak in late spring, while temperature peaks in Bloom, N., 1989. Determination of picogram levels of methylmercury by aqueous
phase ethylation, followed by cryogenic gas chromatography with cold vapour
August. There may be a different limiting factor, such as sulfate,
atomic flourescence detection. Can. J. Fish. Aquat. Sci. 46, 1131e1140.
organic carbon, or other electron acceptors, all of which could su- Bloom, N.S., 1992. On the chemical form of mercury in edible fish and marine
persede the effect of temperature. Additional laboratory and invertebrate tissue. Can. J. Fish. Aquat. Sci. 49, 1010e1017.
fieldwork linked with process modeling would assist in improving Bodaly, R.A., St. Louis, V.L., Paterson, M.J., Fudge, R.J.P., Hall, B.D., Rosenberg, D.M.,
Rudd, J.W.M., 1997. Bioaccumulation of mercury in the aquatic food chain in
the current formulations. Watershed transport of THg and MeHg newly flooded areas. Mercury Eff. Environ. Biol. 34, 259e287.
from the watershed to receiving streams may occur via erosion Bradley, P., Burns, D., Murray, K., Brigham, M., Button, D., Chasar, L., Marvin-
(MeHg sorbed to soil particles), direct runoff, or via association with DiPasquale, M., Lowery, M., Journey, C.A., 2011. Spatial and seasonal variability
of dissolved methylmercury in two stream basins in the eastern United States.
DOC the flows into the stream. VELMA models MeHg concentra- Environ. Sci. Technol. 45, 2048e2055.
tions based on an assumed strong correlation between DOC and Bradley, P., Journey, C., Brigham, M., Burns, D., Button, D., Riva-Murray, K., 2013.
MeHg, which is not evident at McTier Creek but is present at the Intra- and inter-basin mercury comparisons: importance of basin scale and
time-weighted methylmercury estimates. Environ. Pollut. 172, 42e52.
basin scale (Bradley et al., 2013). Our modeling results indicate Bradley, P., Journey, C., Chapelle, F., Lowery, M., Conrads, P., 2010. Flood hydrology
additional factors control MeHg stream concentrations, which are and methylmercury availability in coastal Plain Rivers. Environ. Sci. Technol. 44,
not captured by correlations with streamflow or DOC alone. 9285e9290.
Bradley, P., Journey, C.A., Lowery, M., Brigham, M., Burns, D., Button, D., Chapelle, F.,
All VELMA grid cells are simulated as soil layers; in-stream Lutz, M., Marvin-DiPasquale, M., Riva-Murray, K., 2012. Shallow groundwater
processes, such as settling and resuspension, photo-oxidation, mercury supply in a Coastal Plain Stream. Environ. Sci. Technol. 46, 7503e7511.
and methylation in the sediment active layer, are not explicitly Branfireun, B.A., Heyes, A., Roulet, N.T., 1996. The hydrology and methylmercury
dynamics of a Precambrian Shield headwater peatland. Water Resour. Res. 32,
represented. This could be overcome by linking VELMA to a water
1785e1794.
quality model or by incorporating a stream network. VELMA Brigham, M.E., Wentz, D.A., Aiken, G.R., Krabbenhoft, D.P., 2009. Mercury cycling in
currently does not simulate erosion or sediment transport and as- stream ecosystems. 1. Water column chemistry and transport. Environ. Sci.
sumes uniform land-use/land-cover across the watershed. Technol. 43, 2720e2725.
Chen, C.W., Herr, J.W., Goldstein, R.A., 2008. Model calculations of total maximum
A challenge in modeling Hg dynamics at these and management daily loads of mercury for drainage lakes. J. Am. Water Resour. Assoc. 44, 1295e
scales is to adequately capture system processes while simulta- 1307.
neously balancing model complexity, transparency, and applica- Clark, M., Rupp, D., Woods, R., Tromp-van Meerveld, H., Peters, N., Freer, J., 2009.
Consistency between hydrological models and field observations: linking pro-
bility, in addition to minimizing the number of parameters. Many cesses at the hillslope scale to hydrological responses at the watershed scale.
factors control MeHg concentrations in streams. Therefore, further Hydrol. Process. 23, 311e319.
advances concerning key processes to integrate into watershed Cohn, T., 2005. Estimating contaminant loads in riversdan application of adjusted
maximum likelihood to type 1 censored data. Water Resour. Res. 41, W07003.
models at multiple scales e considering both computational and Cohn, T.A., Caulder, D.L., Gilroy, E.J., Zynjuk, L.D., Summers, R.M., 1992. The validity
scientific factors e are needed. of a simple statistical model for estimating fluvial constituent loads e an
VELMA provides a framework for incorporating improvements empirical study involving nutrient loads entering Chesapeake Bay. Water
Resour. Res. 28, 2353e2363.
in Hg watershed process science. Simulating the cycles of N, C, and Cohn, T.A., DeLong, L.L., Gilroy, E.J., Hirsch, R.M., Wells, D.K., 1989. Estimating con-
Hg and hydrology permits investigation of how management stituent loads. Water Resour. Res. 25, 937e942.
strategies for each will impact the others. Our work used results of a DeWild, J.F., Olson, M.L., Olund, S.D., 2002. Determination of Methyl Mercury by
Aqueous Phase Ethylation, Followed by Gas Chromatographic Separation with
focused reach study to parameterize VELMA and assessed its utility
Cold Vapor Atomic Fluorescence Detection. U. S. Geological Survey, Middleton,
for headwater watershed scales. Based on our results, there is po- WI, p. 14.
tential to scale-up to large basin and regional scales where regu- Eckley, C.S., Hintelmann, H., 2006. Determination of mercury methylation potentials
latory and management decisions are made. in the water column of lakes across Canada. Sci. Total Environ. 368, 111e125.
Evers, D., Han, Y.-J., Driscoll, C.T., Kamman, N., Goodale, M.W., Lambert, K.F.,
Holsen, T.M., Chen, C.Y., Clair, T.A., Butler, T., 2007. Biological mercury hotspots in
Acknowledgments the northeastern United States and southeastern Canada. BioScience 57, 29e43.
Feaster, T., Golden, H., Odom, K., Lowery, M., Conrads, P., Bradley, P., 2010. Simula-
tion of streamflow in the McTier Creek watershed, South Carolina. In:
We appreciate helpful suggestions from Katie Price, Roy Sidle, Survey, U.S.G. (Ed.), Scientific Investigations Report. U.S. Geological Survey,
and the insightful comments of two anonymous reviewers. p. 61.
Fishman, M.J., 1993. Methods of Analysis by the U.S. Geological Survey National
Mention of trade names, products, or services does not convey, and Water Quality LaboratorydDetermination of Inorganic and Organic Constitu-
should not be interpreted as conveying, official EPA approval, ents in Water and Fluvial Sediments. U.S. Geological Survey, Denver, Colo, p. 217.
endorsement or recommendation. The views expressed in this Fitzgerald, W.F., Engstrom, D.R., Mason, R.P., Nater, E.A., 1998. The case for atmo-
spheric mercury contamination in remote areas. Environ. Sci. Technol. 32, 1e7.
paper are those of the authors and do not necessarily represent the Freeze, R.A., Cherry, J.A., 1979. Groundwater. Prentice-Hall, Inc., Englewood Cliffs, NJ.
views or policies of the USEPA or USGS. Futter, M.N., Poste, A.E., Butterfield, D., Dillon, P.J., Whitehead, P.G., Dastoor, A.P.,
Lean, D.R.S., 2012. Using the INCA-Hg model of mercury cycling to simulate
total and methyl mercury concentrations in forest streams and catchments. Sci.
References Total Environ. 424, 219e231.
Gao, N., Armatas, N.G., Shanley, J.B., Kamman, N.C., Miller, E.K., Keeler, G.J.,
Abdelnour, A., McKane, R., Stieglitz, M., Pan, F., Cheng., Y., 2013. Effects of harvest on Scherbatskoy, T., Holsen, T.M., Young, T., McIlroy, L., Drake, S., Olsen, B., Cady, C.,
carbon and nitrogen dynamics in a Pacific Northwest forest catchment. Water 2006. Mass balance assessment for mercury in Lake Champlain. Environ. Sci.
Resour. Res. 49. http://dx.doi.org/10.1029/2012WR012994. Technol. 40, 82e89.
Abdelnour, A., Stieglitz, M., Pan, F., McKane, R., 2011. Catchment hydrological re- Gilmour, C.C., Henry, E.A., 1991. Mercury methylation in aquatic systems affected by
sponses to forest harvest amount and spatial pattern. Water Resour. Res. 47, acid deposition. Environ. Pollut. 71, 31e169.
W09521. http://dx.doi.org/10.1029/2010WR010165. Golden, H., Knightes, C., Conrads, P., Davis, G., Feaster, T., Journey, C.A., Benedict, S.,
Aiken, G.R., McKnight, D.M., Thorn, K.A., Thurman, E.M., 1992. Isolation of hydro- Brigham, M., Bradley, P., 2012. Characterizing mercury concentrations and
philic organic acids from water using nonionic macroporous resins. Org. Geo- fluxes in a Coastal Plain watershed: insights from dynamic modeling and data.
chem. 18, 567e573. J. Geophys. Res. Biogeosci. 117, 17.
Ali, G., Oswald, C.J., Spence, C., Cammeraat, E.L.H., McGuire, K.J., Meixner, T., Golden, H., Knightes, C., Conrads, P., Feaster, T., Davis, G., Benedict, S., Bradley, P.,
Reaney, S.M., 2013. Towards a unified threshold-based hydrological theory: 2013. Climate change and watershed mercury export: variable responses to
necessary components and recurring challenges. Hydrol. Process. 27, 313e318. diverse projections and models. Environ. Toxicol. Chem. 32, 2165e2174.
Allison, J., Allison, T., 2005. Partition Coefficients for Metals in Surface Water, Soil, Golden, H.E., Knightes, C.D., 2011. Simulated watershed mercury and nitrate flux
and Waste. U.S. Environmental Protection Agency, Athens, GA. responses to multiple land cover conversion scenarios. Environ. Toxicol. Chem.
Benoit, J.M., Gilmour, C.C., Heyes, A., Mason, R.P., Miller, C.L., 2002. Geochemical and 30, 773e786.
biological controls over methylmercury production and degradation in aquatic Graydon, J.A., St Louis, V.L., Hintelmann, H., Lindberg, S.E., Sandilands, K.A.,
ecosystems. In: ACS Symposium Series. Biogeochemistry of Environmentally Rudd, J.W.M., Kelly, C.A., Hall, B.D., Mowat, L.D., 2008. Long-term wet and dry
Important Trace ElementsAmerican Chemical Society, Washington, DC, deposition of total and methyl mercury in the remote boreal ecoregion of
pp. 262e297. Canada. Environ. Sci. Technol. 42, 8345e8351.
192 C.D. Knightes et al. / Environmental Pollution 187 (2014) 182e192

Harris, R.C., Rudd, J.W.M., Amyot, M., Babiarz, C.L., Beaty, K.G., Blanchfield, P.J., Patton, C.J., Kryskalla, J.R., 2003. Methods of Analysis by the U.S. Geological Survey
Bodaly, R.A., Branfireun, B.A., Gilmour, C.C., Graydon, J.A., Heyes, A., National Water Quality Laboratory e Evaluation of Alkaline Persulfate Digestion
Hintelmann, H., Hurley, J.P., Kelly, C.A., Krabbenhoft, D.P., Lindberg, S.E., as an Alternative to Kjeldahl Digestion for Determination of Total and Dissolved
Mason, R.P., Paterson, M.J., Podemski, C.L., Robinson, A., Sandilands, K.A., Nitrogen and Phosphorus in Water. U.S. Geological Survey, Denver, Colorado,
Southworth, G.R., St. Louis, V.L., Tate, M.T., 2007. Whole-ecosystem study shows p. 33.
rapid fish-mercury response to changes in mercury deposition. Proc. Natl. Acad. Runkel, R.L., Crawford, C.G., Cohn, T.A., 2004. Load Estimator (LOADEST)da
Sci. 104, 16586e16591. FORTRAN Program for Estimating Constituent Loads in Streams and Rivers. U.S.
Henry, E.A., Dodge-Murphy, L.J., Bigham, G.N., Klein, S.M., Gilmour, C.C., 1995. Total Geological Survey, p. 69.
mercury and methylmercury mass balance in an alkaline hypereutrophic urban Scheuhammer, A.M., Basu, N., Evers, D.C., Heinz, G.H., Sandheinrich, M.B.,
lake (Onondaga Lake, NY). Water Air Soil Pollut. 80, 509e518. Bank, M.S., 2012. Ecotoxicology of mercury in fish and wildlife: recent advances.
Hintelmann, H., Keppel-Jones, K., Evans, R.D., 2000. Constants of mercury methyl- In: Bank, M.S. (Ed.), Mercury in the Environment: Pattern and Process. Uni-
ation and demethylation rates in sediments and comparison of tracer and versity of California Press, Berkley, CA.
ambient mercury availability. Environ. Toxicol. Chem. 19, 2204e2211. Scheuhammer, A.M., Meyer, M.W., Sandheinrich, M.B., Murray, M.W., 2007. Effects
Horvat, M., Bloom, N.S., Liang, L., 1993. Comparison of distillation with other current of environmental methylmercury on the health of wild birds, mammals, and
isolation methods for the determination of methyl mercury compounds in low fish. AMBIO J. Hum. Environ. 36, 12e19.
level environmental samplesdpart 1. sediments. Anal. Chim. Acta 281, 135e Schuster, P.F., Shanley, J.B., Marvin-Dipasquale, M., Reddy, M.M., Aiken, G.R.,
152. Roth, D.A., Taylor, H.E., Krabbenhoft, D.P., DeWild, J.F., 2008. Mercury and
Hudson, R.J.M., Gherini, S.A., Watras, C.J., Porcella, D.B., 1994. Modeling the organic carbon dynamics during runoff episodes from a northeastern USA
biogeochemical cycle of mercury in lakes: the Mercury Cycling Model (MCM) watershed. Water Air Soil Pollut. 187, 89e108.
and its application to the MTL study lakes. In: Watras, C.J., Huckabee, J.W. Scudder, B.C., Chasar, L.C., Wentz, D.A., Bauch, N.J., Brigham, M.E., Moran, P.W.,
(Eds.), Mercury Pollution: Integration and Synthesis. Lewis Publishers, Boca Krabbenhoft, D.P., 2009. Mercury in Fish, Bed Sediment, and Water from
Raton. Streams Across the United States, 1998e2005. U.S. Geological Survey, p. 74.
Journey, C.A., Burns, D.A., Riva-Murray, K., Brigham, M.E., Button, D.T., Feaster, T.D., Scudder Eikenberry, B., Riva-Murray, K., Smith, M., Bradley, P., Button, D., Clark, J.,
Petkewich, M.D., Bradley, P.M., 2012. Fluvial Transport of Mercury, Organic Burns, D., Journey, C.A., 2012. Environmental Settings of Selected Streams
Carbon, Suspended Sediment, and Selected Major Ions in Contrasting Stream Sampled for Mercury in New York and South Carolina, 2005e09. U.S. Geological
Basins in South Carolina and New York, October 2004 to September 2009. U.S. Survey, p. 36.
Geological Survey Scientific Investigations Report, p. 125. Shanley, J.B., Mast, M.A., Campbell, D.H., Aiken, G.R., Krabbenhoft, D.P., Hunt, R.J.,
Kelly, C.A., Rudd, J.W.M., St. Louis, V.L., Heyes, A., 1995. Is total mercury concen- Walker, J.F., Schuster, P.F., Chalmers, A., Aulenbach, B.T., Peters, N.E., Marvin-
tration a good predictor of methyl mercury concentration in aquatic ecosys- Dipasquale, M., Clow, D.W., Shafer, M.M., 2008. Comparison of total mercury
tems. Water Air Soil Pollut. 80. and methylmercury cycling at five sites using the small watershed approach.
Knightes, C.D., Sunderland, E.M., Barber, M.C., Johnston, J.M., Ambrose, R.B., 2009. Environ. Pollut. 154, 143e154.
Application of ecosystem-scale fate and bioaccumulation models to predict fish Sivapalan, M., 2003. Process complexity at hillslope scale, process simplicity at the
mercury response times to changes in atmospheric deposition. Environ. Toxicol. watershed scale: is there a connection? Hydrol. Process. 17, 1037e1041.
Chem. 28, 881e893. Snodgrass, J.W., Jagoe, C.H., Bryan Jr., A.L., Brant, H.A., Burger, J., 2000. Effects of
Kolka, R.K., Grigal, D.F., Verry, E.S., Nater, E.A., 1999. Mercury and organic carbon trophic status and wetland morphology, hydroperiod, and water chemistry on
relationships in streams draining forested upland/peatland watersheds. mercury concentrations in fish. Can. J. Fish. Aquat. Sci. 57, 171e180.
J. Environ. Qual. 28, 766e775. Snowling, S.D., Kramer, J.R., 2001. Evaluating modelling uncertainty for model se-
Kotnik, J., Horvat, M., Jereb, V., 2002. Modelling of mercury geochemical cycle in lection. Ecol. Model. 138, 17e30.
Lake Velenje, Slovenia. Environ. Model. Softw. 17, 593e611. St. Louis, V.L., Rudd, J.W.M., Kelly, C.A., Hall, B.D., Rolfhus, K.R., Scott, K.J.,
Larssen, T., de Wit, H.A., Wiker, M., Halse, K., 2008. Mercury budget of a small Lindberg, S.E., Dong, W., 2001. Importance of the Forest Canopy to fluxes of
forested boreal catchment in southeast Norway. Sci. Total Environ. 404, 290e methyl mercury and total mercury to boreal ecosystems. Environ. Sci. Technol.
296. 35, 3089e3098.
Lewis, M.E., Brigham, M.E., 2004. Chapter A5, section 5.6.4.B, low-level mercury. In: Ullrich, M., Tanton, T.W., Abdrashitova, S.A., 2001. Mercury in the aquatic envi-
Wilde, F.D., Radtke, D.B., Gibs, J., Iwatsubo, R.T. (Eds.), National Field Manual for ronment: a review of factors affecting methylation. Crit. Rev. Environ. Sci.
the Collection of Water-quality Data. Technol. 31, 241e293.
Lindenschmidt, K.-E., 2006. The effect of complexity on parameter sensitivity and USEPA, 1996. Method 1669-Sampling Ambient Water for Trace Metals at EPA Water
model uncertainty in river water quality modeling. Ecol. Model. 190, 72e86. Quality Criteria Levels. U.S. Environmental Protection Agency, Washington, D.C.,
Macleod, M., McKone, T.E., Mackay, D., 2005. Mass balance for mercury in the San p. 35
Franciso Bay Area. Environ. Sci. Technol. 39, 6721e6729. USEPA, 1997. Mercury Study Report to Congress. In: Fate and Transport of Mercury
Mergler, D., Anderson, H.A., Chan, L.H.M., Mahaffey, K.R., Murray, M., Sakamoto, M., in the Environment, vol. III. U.S. Environmental Protection Agency, Washington,
Stern, A.H., 2007. Methylmercury exposure and health effects in humans e a DC.
worldwide concern. Ambio 36, 3e11. USEPA, 2001. Method 1630dmethyl Mercury in Water by Distillation, Aqueous
MPCA, 2007. Minnesota Statewide Mercury Total Maximum Daily Load. Minnesota Ethylation, Purge and Trap, and Cold Vapor Atomic Fluorescence Spectrometry.
Pollution Control Agency. U.S. Environmental Protection Agency, Washington, D.C., p. 41
Munthe, J., Bodaly, R.A., Branfireun, B.A., Driscoll, C.T., Gilmour, C.C., Harris, R., USEPA, 2010. Biennial National Listing of Fish Advisories. U.S. Environmental Pro-
Horvat, M., Lucotte, M., Malm, O., 2007. Recovery of mercury-contaminated tection Agency, Office of Water, Washington, DC. USGS, variously dated. Na-
fisheries. Ambio 36, 33e44. tional field manual for the collection of water-quality data. In: Survey, U.S.G.
NEIWPCC, 2007. Northeast Regional Mercury Total Maximum Daily Load. New (Ed.), U.S. Geological Survey Techniques of Water-resources Investigations. U.S.
England Interstate Water Pollution Control Commission. Geological Survey, Reston, VA.
NLCD, 2001. National Land Coverage Data-2001. Multi-Resolution Land Character- Vette, A.F., Landis, M.S., Keeler, G.J., 2002. Deposition and emission of gaseous
istics Consortium (MRLC). mercury to and from Lake Michigan during the Lake Michigan mass balance
NWIS, 2010. National Water Information System: USGS Water Data for the Nation. study (July, 1994eOctober, 1995). Environ. Sci. Technol. 26, 4525e4532.
Obrist, D., Johnson, D.W., Lindberg, S.E., Luo, Y., Hararuk, O., Bracho, R., Battles, J.J., Wang, Z., Zhang, X., Xiao, J., Zhijia, C., Yu, P., 2009. Mercury fluxes and pools in three
Dail, D.B., Edmonds, R.L., Monson, R.K., Ollinger, S.V., Pallardy, S.G., subtropical forested catchments, southwest China. Environ. Pollut. 157, 801e
Pregitzer, K.S., Todd, D.E., 2011. Mercury distribution across 14 U.S. Forests. Part 808.
I: spatial patterns of concentrations in biomass, litter, and soils. Environ. Sci. Wolock, D., 1993. Simulating the variable-source-area concept of streamflow gen-
Technol. 45, 3974e3981. eration with the watershed model TOPMODEL. In: Survey, U.S.G. (Ed.), Water
Olson, M.L., DeWild, J.F., 1999. Techniques for the collection and species-specific Resources Investigations Report. U.S. Geological Survey.
analysis of low levels of mercury in water, sediment, and biota. In: Woodruff, L.G., Cannon, W.F., Knightes, C., Chapelle, F., Bradley, P., Burns, D.,
Morganwalp, D.W., Buxton, H.T. (Eds.), U.S. Geological Survey Toxic Substances Brigham, M., Lowery, M., 2010. Total Mercury, Methylmercury, and Selected
Hydrology Program, Proceedings of the Technical Meeting, Charleston, South Elements in Soils of the Fishing Brook Watershed, Hamilton County, New York,
Carolina, March 8e12, 1999dVolume 2 of 3dContamination of Hydrologic and the McTier Creek watershed, Aiken County, South Carolina, 2008. In: U.S.
Systems and Related Ecosystems. U.S. Geological Survey, pp. 191e200. Geological Survey Data Series 516. U. S. Geological Survey, p. 10.