Sensitivity analysis of an updated bidirectional air – surface exchange model for elemental mercury vapor

A box model for estimating bidirectional air– surface exchange of gaseous elemental mercury (Hg 0) has been updated based on the latest understanding of the resistance scheme of atmosphere–biosphere interface transfer. Simulations were performed for two seasonal months to evaluate diurnal and seasonal variation. The base-case results show that water and soil surfaces are net sources, while vegetation is a net sink of Hg 0. The estimated net exchange in a domain covering the contiguous US and part of Canada and Mexico is 38.4 and 56.0 Mg as evasion in the summer and winter month, respectively. The smaller evasion in summer is due to the stronger Hg 0 uptake by vegetation. Modeling experiments using a two-level factorial design were conducted to examine the sensitivity of flux response to the changes in physical and environmental parameters in the model. It is shown that atmospheric shear flows (surface wind over water and friction velocity over terrestrial surfaces), dissolved gaseous mercury (DGM) concentration, soil organic and Hg content, and air temperature are the most influential factors. The positive effect of friction velocity and soil Hg content on the evasion flux from soil and canopy can be effectively offset by the negative effect of soil organic content. Significant synergistic effects are identified between surface wind and DGM level for water surface, and between soil Hg content and friction velocity for soil surface, leading to ∼ 50% enhanced flux compared to the sum of their individual effects. The air–foliage exchange is mainly controlled by surface resistance terms influenced by solar irradiation and air temperature. Research in providing geospatial distribution of Hg in water and soil will greatly improve the flux estimate. Elucidation on the kinetics and mechanism of Hg(II) reduction in soil/water and quantification of the surface resistances specific to Hg species will also help reduce the model uncertainty.


Introduction
Mercury (Hg) is a persistent, bioaccumulative pollutant released into the atmosphere from a variety of anthropogenic and natural sources.The anthropogenic release (2000-2400 Mg yr −1 ) primarily comes from fossil fuel combustion, waste incineration, metal smelting and cement production (Pacyna et al., 2003(Pacyna et al., , 2006;;Streets et al., 2005Streets et al., , 2009;;Pirrone et al., 2010).The natural sources include biomass burning, volcanic activities, geological weathering from Hg enriched substrates and so-called re-emission (i.e., the emission of previously deposited Hg) at the atmosphere-biosphere interface (Gustin et al., 2008;Mason and Sheu, 2002).While the anthropogenic emissions have been estimated and continuously updated with reasonable consistency since the 1990s, the estimates for natural emissions have been highly uncertain (1500-5207 Mg yr −1 ), primarily due to a lack of understanding in the air-surface exchange of Hg 0 .Since the natural release can account for up to two-thirds of global mercury input to the atmosphere (Friedli et al., 2009;Pirrone et al., 2010), better quantification of the mass input is critical Published by Copernicus Publications on behalf of the European Geosciences Union.
X. Wang et al.: Sensitivity analysis of an updated bidirectional air-surface exchange model in assessing the global biogeochemical cycling of mercury (Lindberg et al., 2007).
Air-surface exchange is an important component in atmospheric mercury modeling for estimating Hg 0 evasion and deposition over soil, water and vegetation.For terrestrial surfaces, the soil Hg evasion has been calculated using the statistical relationships obtained from the measured Hg 0 flux and observed environmental factors such as temperature, solar irradiance, leaf area index and Hg content (Bash et al., 2004;Gbor et al., 2006;Lin et al., 2005;Shetty et al., 2008;Xu et al., 1999;Selin and Jacob, 2008;Smith-Downey et al., 2010).Such an approach oversimplifies the role of environmental factors in the exchange process because Hg 0 flux was measured in a limited number of locations where the environmental parameters (e.g., soil properties and meteorology) are specific to those sites.In addition, most of these models treat vegetation as a net evasion source of Hg 0 , which is inconsistent with later assessments that suggest vegetation to be a net sink (Gustin et al., 2008;Hartman et al., 2009).Recent isotopic tracer studies have shown that plant roots serve as a barrier that prevents translocation of inorganic Hg in soil to other parts of plants (Cui et al., 2014).It has also been suggested that Hg absorbed on foliage can be transported to stems and roots (Yin et al., 2013).
To better represent the surface process, algorithms parameterizing the transport resistances at soil and foliage interfaces have been developed to calculate the multilayered, bidirectional flux by the gradient between an ambient Hg level and a "compensation point" inferred from the surface characteristics (Bash, 2007(Bash, , 2010;;Scholtz et al., 2003;Zhang et al., 2009a;Sutton et al., 2007).Such an approach is more scientifically sound and mathematically robust.The model results were also more consistent with those from isotopic tracer studies (Bash, 2010).However, the complicated parameterization makes it difficult to understand the relative importance of each parameter on the simulated flux.There are also multiple model assumptions on the model variables that lack field data for verification.Although the model results can be constrained by air concentration and wet deposition, the assumptions could increase the uncertainty of model estimates and limit the improvement of model algorithms.
The objectives of this study are to present an updated Hg 0 air-surface exchange model and to quantitatively examine the relative importance of the physical and environmental variables implemented in the model.Coupled with the latest understanding in the partitioning and mass transfer at different atmosphere-biosphere interfaces, we integrated the bidirectional air-surface exchange model (Bash, 2007(Bash, , 2010) ) and the surface resistance schemes of Hg dry deposition and photochemical reaction (Zhang et al., 2003(Zhang et al., , 2009a;;Lin et al., 2006) for quantifying the air-surface exchange of Hg 0 .Two monthly (summer and winter) simulations were performed to investigate the seasonal and diurnal variability of the modelestimated flux.A systematic set of sensitivity simulations using a multi-step factorial design of experiments was per-formed to investigate the effect of significant model parameters and their interconnections.Based on the sensitivity results, processes that control Hg 0 air-surface exchange over different natural surfaces are discussed and research needs for future model improvement are proposed.

Model description
The total air-surface exchange is the sum of Hg 0 fluxes from water, soil (including bare lands and soil under the canopy) and foliage surfaces.The direction (evasion or deposition) of the flux is driven by the gradient between atmospheric Hg 0 concentration and a compensation point that represents the Hg 0 concentration at the interface between the atmosphere and the surface of interest.The magnitude of the flux is determined by the ratio of concentration gradient to surface resistance (for terrestrial surfaces) or by the product of overall mass transfer coefficient and concentration gradient (for water surfaces).The nomenclature and dimension of the entire set of model variables are detailed in Table 1.The parameterization of each model component is briefly described below.

Air-water exchange
The flux over fresh water and oceanic surfaces, F w , is calculated using a two-film mass transfer model with the transfer rate limited by the diffusion in the water boundary layer (Poissant et al., 2000): where K w is the overall mass transfer coefficient estimated by the wind speed at 10 m above water surface and the mass transfer ratio of CO 2 / Hg across the air-water interface (Shetty et al., 2008), C w is the dissolved gaseous mercury (DGM) concentration in surface water, H w is the dimensionless Henry's law constant.K w and H w are calculated using formulation described earlier (Poissant et al., 2000;Lin and Tao, 2003).

Air-terrestrial exchange
The terrestrial system is divided into two categories: the canopy biomes (leaf area index, LAI > 0) and the bare lands (LAI = 0, referring to sparsely vegetated land, bare tundra and snow/ice surfaces).The total flux from the canopy is made up of the air-soil exchange flux and the air-foliage exchange flux (air-stomata and air-cuticle).Over the canopy system, a multi-layer canopy resistance scheme modified after Bash (2010) and Zhang et al. (2003) was applied (Fig. 1).
Compared to the earlier mechanistic schemes, this model also (1) includes the foliage storage effect, (2) considers photochemical reduction on foliage, and (3) updates the Fraction of stomatal blocking under wet condition dimensionless a Value for base-case simulation, Xu et al. (1999); b for 0-20 cm topsoil, the bulk density is 1.1-1.3g cm −3 and organic carbon content is 3.3 kg m −2 in the US (Calhoun et al., 2001;Guo et al., 2006), so assuming f oc is 2 % in the 0-5 cm topsoil; c value for base-case simulation, Bash (2010); d based on the negligible solubility (Henry's constant = 0.139 M atm −1 ) and chemical inertness (Zhang et al., 2009b(Zhang et al., , 2012 resistance terms.The flux over canopy biomes, F cnp , is estimated as where t is the time duration, R a is the aerodynamic resistance, R b is the quasi-laminar layer resistance and C atm is the atmospheric Hg concentration.R a and R b are calculated according to Marsik et al. (2007).χ cnp is the overall compensation point parameterized as a weighted average at the air-cuticle, air-stomata and air-soil interfaces as illustrated in Fig. 1 (Bash, 2010;Zhang et al., 2009a): where χ c is the cuticular compensation point, χ s is the stomatal compensation point, χ g is the soil compensation point, R c is the cuticular resistance, R s is the stomatal resistance, R g is the soil diffusion resistance, and R ac is the in-canopy aerodynamic resistance.The individual compensation points are described by Eqs. ( 6), ( 9) and (15).

Air-soil exchange
In the absence of vegetation (when LAI = 0), the flux from bare lands (F bls ) can be estimated as In the presence of vegetation (when LAI > 0), the flux from soil under canopy (F g ) is calculated as where R ac accounts for the resistance of gas diffusion from ground to the lower canopy and is assumed to be common for all gaseous species (Zhang et al., 2002b).The compensation point at air-soil interface (χ g ) can be expressed as (Bash, 2010) where Hg 0 sl is the concentration of Hg 0 bound to soil, calculated as a reduction product of Hg(II) using soil Hg content and a pseudo-first-order rate constant related to solar irradiance (Gustin et al., 2002).H is Henry's constant parameterized following Andersson et al. (2008).f oc is the fraction of organic carbon in surface soil (0-5 cm).K oc is the partition coefficient of Hg 0 between soil organic carbon and water.
R g is the Hg 0 diffusion resistance over a ground surface (soil, ice/snow) (Zhang et al., 2002b): where R g(SO 2 ) and R g(O 3 ) are the diffusion resistances of SO 2 and O 3 , α Hg 0 is the Hg 0 scaling factor based on SO 2 , and β Hg 0 is Hg 0 scaling factor based on O 3 .The formulation of R g(SO 2 ) and R g(O 3 ) has been described previously (Zhang et al., 2003).

Air-cuticle exchange
Air-cuticle exchange flux is calculated as (Bash, 2010) where LAP denotes the leaf-air partitioning coefficient for Hg 0 (Rutter et al., 2011) and [Hg 0 c ] is the concentration of Hg 0 bound to foliar cuticular surface, calculated as the photoreduction product of a fraction of newly deposited Hg(II) on foliar interfaces (Graydon et al., 2009): where [Hg II+ c,DD ] is the concentration loading of total dry deposited Hg(II) on cuticle, Hg II+ c is the concentration of the deposited Hg(II) residing on cuticular surfaces, [Hg II+ w ] is the concentration of Hg(II) that can be washed off from leaves, f rxn is the fraction of Hg(II) that can be photoreduced, f fixed is the fraction of Hg(II) fixed into tissue and not available for re-emission or wash-off, and T l is the leaf thickness.f rxn and f fixed are parameterized following Smith-Downey et al. (2010).R c is the cuticular resistance calculated as (Zhang et al., 2002b) 1

Air-stomata exchange
The air-stomata exchange flux is estimated as (Bash, 2010) It is assumed that the uptake of Hg species through stomata is predominantly Hg 0 due to its abundance in the atmosphere (Capiomont et al., 2000;Millhollen et al., 2006;Stamenkovic and Gustin, 2009).As such, the dissolved Hg 0 in the stomatal compartment, [Hg 0 s ], can be formulated as where [Hg 0 s,DD ] is the concentration of newly deposited Hg 0 stored in the stomatal compartment.The overall stomatal resistance is calculated as (Zhang et al., 2002b) where R st is the resistance associated with stomata, R me is resistance associated with mesophyll reservoir, and W st is the fraction of stomatal blocking under wet condition.The detailed formulation of R st and R me and W st can be found elsewhere (Zhang et al., 2002b(Zhang et al., , 2003(Zhang et al., , 2012)).

Modeling experiments for sensitivity analysis
To explore the sensitivity of model response to the changes in model parameters, the model results obtained by varying the input values of model parameters were compared systematically.This was accomplished by a technique called two-level factorial design of experiments, which estimates the change of model response (i.e., the simulated air-surface exchange flux of Hg 0 ) caused by the changes in model parameters at two levels (i.e., a high and a low typical value).Additional discussion on the application of factorial design of experiments is provided in the Supplement.In short, the values of model parameters were varied individually (i.e., changing one input value of a model parameter at a time) and in combinations (i.e., simultaneously changing the input values of two or more model parameters) at the two selected levels.Then the Hg 0 exchange fluxes at the changed input values were calculated and compared to estimate the sensitivity.This method is statistically robust, and therefore the synergistic and antagonistic interactions among model parameters can be estimated with indications of statistical significance.The studied parameters include both physical and environmental factors.Their respective experimental levels are shown in Tables 2-4.The principle of factor sparsity (Myers et al., 2009) states that the main factor effects and lowerorder interactions dominate most system responses (i.e., it is unlikely for three or more predominant factors influencing the system response simultaneously).Therefore, the effect of interaction terms higher than second order was not considered.
For the water surface, there are four parameters (factors) driving the model simulation (Table 2).Therefore, there are 16 (2 4 ) possible combinations of input parameter values at two levels.In this case, the model was run 16 times (so-called "full factorial design" because all possible experimental runs are analyzed), and the data were compared for identifying the significant factors.For the bare lands, there are 2048 possible combinations at two levels (2 11 for 11 parameters).In this case, we first selected 32 runs (so-called 2 11−6 "fractional design" because only a fraction of all possible experimental runs is analyzed) that allowed for statistical analysis of all single-factor effects and two-factor interactions in order to screen out the significant factors.After this initial screening, a two-level full factorial design was applied to the five most significant factors based on a 95 % confidence level to assess the sensitivity (the results of the 2 11−6 design are shown in the Supplement).For the canopy ecosystem, there are 15 model factors (Table 4).In this case, a 2 15−9 fractional design (64 runs) was first utilized to identify 11 predominant factors.Then a successive 2 11−6 design was applied to obtain the five most significant factors, whose sensitivity was assessed using a 2 5 full factorial design (the results of the 2 15−9 and 2 11−6 designs are shown in the Supplement).The sensitivity results were obtained from the final full factorial design for watersheds, bare lands and canopy ecosystems.The data analysis of the factorial experiments was conducted using Minitab ® 16 professional version.

Model configuration and data
The modeling domain is in Lambert conformal projection covering mainly the contiguous United States (CONUS), with 156 × 118 grid cells at 36 km spatial resolution.Hourly meteorological data were prepared using the Weather Re A base-case simulation was performed in a summer and a winter month (August and December 2009) to evaluate the seasonal and diurnal variability of the air-surface exchange.The base case refers to the modeling utilizing the values listed in Table 1 with the meteorological parameters extracted from WRF output.In the simulation, the atmospheric Hg 0 concentration retrieved from the output of the Hg extension of Community Multi-scale Air Quality modeling system (CMAQ-Hg) version 4.6 for the same modeling period was applied to represent the air concentration of Hg 0 .The simulation does not directly incorporate the feedback of the air-surface exchange to the air concentration.However, for a regional model domain (CONUS), natural evasion and deposition of Hg 0 does not significantly modify the ambient concentration (Lin et al., 2005;Gbor et al., 2006) since the time required for air turnover is relatively short (typically 3-4 days) and the air concentration of Hg is mainly controlled  by the boundary conditions (Pongprueksa et al., 2008).In the model experiments, the concentration of Hg 0 was tested as a sensitivity parameter.

Results of base-case simulations
The model estimates a net emission of 38.4 Mg in the summer month (16.6 Mg from water, 45.0 Mg from soil and −23.2 Mg from foliage) and 56.0 Mg in the winter month (33.9 Mg from water, 29.5 Mg from soil and −7.4 Mg from foliage) for the entire domain.The evasion from water bodies accounts for ∼ 50 % of the total natural emission (the cumulative net release of Hg 0 caused by the air-surface exchange process) because of the large water areal coverage in the domain (59 %).Vegetation represents a net sink, which is different from earlier estimates using the evapotranspiration approach (Bash et al., 2004;Shetty et al., 2008) but consistent with recent observational studies (Gustin et al., 2008;Stamenkovic and Gustin, 2009).For the terrestrial system, the total emission is 43.9 Mg in 2 months.Assuming the annual emission is 5-6 times of the 2-monthly sum and excluding the emission from Canada, Mexico and Caribbean lands, the model-estimated annual emission in the contiguous US is 118-141 Mg yr −1 , comparable to the recent estimates (95-150 Mg yr −1 ) using flux scaling methods (Ericksen et al., 2006;Hartman et al., 2009;Zehner and Gustin, 2002).

Air-water exchange
Over water surface, the mean simulated flux is 1.6 and 3.1 ng m −2 h −1 in the summer and winter month (Figs.2a and 3a), respectively.Water bodies in the domain are net sources, producing fluxes typically in the range of 1-4 ng m −2 h −1 , similar to earlier measurements (Mason et al., 2001a;Andersson et al., 2011).The spatial distribution is primarily driven by the surface wind speed.Temperature, air Hg 0 and DGM concentration play a much less significant role because a constant DGM was assumed (40 ng m −3 ) and the Hg 0 level over water was in a narrow range (1.4-1.8 ng m −3 ).The Pearsons's correlation coefficient (r) be-tween flux and wind speed is much stronger than the value between flux and temperature (0.56 vs. 0.18).The flux in the winter month is greater because of stronger winds in the northeastern corner of the domain.The emission flux does not show clear diurnal variation in both months because wind speed is the most dominant factor (Fig. 4a).

Air-soil exchange
Soil surfaces have been suggested to be a net source of Hg (Gustin et al., 2008;Hartman et al., 2009), which is also shown in the base-case model results (Figs. 2 and 3).The mean flux from bare lands (0.7 and 0.6 ng m −2 h −1 in the summer and winter month) is lower than the value from soil under the canopy (4.3 and 2.7 ng m −2 h −1 ) because of the land use classification.The bare lands in the domain include sparsely vegetated land, bare ground tundra and snow/ice lands.The flux contribution from such land use types is largely from the southern portion of the domain.The simulated flux from soil under canopy is comparable to those reported at background sites, −0.1-7 ng m −2 h −1 (Ericksen et al., 2006;Kuiken et al., 2008a, b;Carpi and Lindberg, 1998).
The simulated Hg 0 flux from soil under canopy is controlled by the degree of vegetation coverage (LAI), air temperature, friction velocity, air Hg concentration and solar irradiation.In the summer month, the flux in eastern US is lower due to heavy vegetation coverage that increases the incanopy aerodynamic resistance (R ac ) (Zhang et al., 2002a).Higher flux occurs in the central and western US because of the smaller LAI and higher air temperature (Figs.2c and S8 in the Supplement).In the winter month, the higher air temperature and longer sunlit hours cause the higher flux in the south (Fig. 3c, S8 in the Supplement).Among the environmental parameters, LAI has the greatest influence on the estimated flux (r = 0.45).The spatially average soil flux for the entire domain shows a typical diurnal variation caused by air temperature and solar irradiance (Gabriel et al., 2006).The detailed impact of the model variables is discussed in the sensitivity analysis.

Air-foliage exchange
Vegetation represents a net sink of Hg 0 in the base-case simulations.The mean simulated air-foliage exchange is   irradiance (Rutter et al., 2011).The overall diurnal variation in the model domain exhibits the feature of air-foliage exchange (Fig. 4d).
The simulated flux from soil under canopy and foliar surfaces is highly dependent on the resistance terms.Presently, the values of cuticular (R c ), stomatal (R g ) and soil (R s ) resistances of Hg are not well understood (Holmes et al., 2011) and have been estimated by relating to the measured resistance of O 3 , SO 2 and H 2 O (Bash, 2010;Scholtz et al., 2003;Zhang et al., 2003).There have been experimental efforts to determine R c and R s based on Fick's law by introducing isotopic Hg tracer to plants grown in an environmentally controlled chamber (Rutter et al., 2011).The resistances were found to depend on temperature, solar irradiance and Hg species, with reported R c and R s ranging from 150 to 50 000 m s −1 at 0-35 • C and 0-170 W m −2 (Millhollen et al., 2006;Rutter et al., 2011).The simulated flux in the base case applied similar resistance values in the model.However, the lack of deterministic relationships between the resistance terms and environmental parameters still represents an uncertainty and there is a need to better quantify the resistance for Hg 0 .

Sensitivity of exchanges over water bodies
Figure 5 shows the change of air-water flux due to the change of model variables from the low to the high experimental level (Table 2).Individually, wind speed is the most significant parameter (p = 0.003), followed by DGM (p = 0.004) and surface temperature (p = 0.059).On average, increasing wind speed from 0.001 to 20 m s −1 enhanced the flux by 7.6 ng m −2 h −1 (p = 0.003); increasing the DGM from 15 to 240 ng m −3 increases the flux by 7.0 ng m −2 h −1 (p = 0.004).A higher air Hg 0 concentration slightly decreases the evasion flux.There is a significant synergistic effect caused by wind speed and DGM concentration (p = 0.004).Increasing both variables simultaneously from the low to high level (Table 2) causes an additional 48 % increase in the evasion flux.The wind speed and surface temperature also have a synergistic effect, although not as significant (p = 0.059), followed by the effect enhanced by DGM concentration and surface temperature (p = 0.076).The effects of higher DGM concentration and air Hg 0 concentration offset each other, leading to a nearly zero effect on flux (p = 1.000).
In the base case, a uniform DGM concentration was assumed.The spatially constant DGM level represents a  2 (water bodies).T denotes air temperature at water surface, GEM denotes air Hg concentration, DGM denotes dissolved gaseous Hg concentration in surface water, W denotes wind speed, and " * " denotes the interaction effects.significant uncertainty since other environmental parameters such as temperature, wind speed can be estimated reliably through meteorological simulations at a high spatial resolution.The mechanism leading to the net DGM formation in surface water is complex and not fully understood (Qureshi et al., 2010).It has been suggested that dissolved organic matter (Amyot et al., 1994(Amyot et al., , 1997)), hydroxyl radicals (Zhang and Lindberg, 2001) and oxyhalide radicals (e.g., OCl − , OBr − ) (Lalonde et al., 2001) can participate in the sunlight-induced processes that produce DGM.Data on measured net DGM concentration over vast water bodies are not readily available because of a limited number of cruise campaigns (Mason et al., 1998(Mason et al., , 2001b;;Andersson et al., 2011).Strode et al. (2007) and Soerensen et al. (2010) estimated the global distribution of DGM in sea water and showed that accurate representation of DGM concentration is key for calculating air-water exchange.More knowledge on the temporal and spatial distribution of net DGM concentration in surface water can greatly reduce the model uncertainty.Experimental investigation to better understand the chemical pathways leading to net DGM formation will also help constrain the model estimate.

Sensitivity of exchange over bare lands
Figure 6 illustrates the model response to the model variables over bare lands (Table 3).Soil Hg content, friction velocity, air temperature and the scaling factor β Hg 0 (Eq.7) have a positive effect on the simulated Hg flux, while the soil organic content has a negative effect.On average, increasing soil Hg content from 50 to 1000 ng g −1 soil enhances the flux by 55.3 ng m −2 h −1 (p = 0.013); increasing friction velocity from 0.0001 to 1 m s −1 increases the flux by  3 to isolate the significant factors.T denotes surface air temperature, foc denotes fraction of organic carbon in soil, UST denotes friction velocity, SM denotes soil mercury content, b denotes scaling factor for Hg reactivity (β Hg 0 ), and " * " denotes interaction effects.54.8 ng m −2 h −1 (p = 0.014).On the other hand, increasing the soil organic content from 0.6 to 10 % reduces the flux by 54.2 ng m −2 h −1 (p = 0.015).There are several notable interactions among the model variables.First, the positive effects of soil Hg content and friction velocity can be completely offset by soil organic content (Fig. 6).An increase in soil organic content substantially decreases the soil Hg compensation point (Eq.6), suggesting the significant role of soil organic matter in preventing Hg from evading (p = 0.025).There is a strong synergistic effect between friction velocity and soil Hg content (p = 0.022), leading to an additional 46 % increase compared to the sum of the two individual effects (Fig. 6).Quasi-laminar layer resistance (R b ) and aerodynamic resistance (R a ) both decrease with increasing friction velocity.Coupled with the increased soil Hg compensation point at higher soil Hg content (Eq.6), the flux is greatly enhanced (Fig. 6).Overall, friction velocity, soil Hg and organic content are the most influential parameters for Hg exchanges over bare lands.Other parameters including temperature, Hg scaling factor (β Hg in Eq. 7) and other interaction terms have less significant effects.4).For comparison, the sensitivity results for air-soil exchange under canopy are also shown.The forcing of air-canopy exchange is dominated by the air-soil exchange under canopy at the two experimental levels.This resembles the Hg 0 emission characteristics  (Zhang et al., 2003); j presence of precipitation events has effects on ground and cuticular resistance terms (Zhang et al., 2003); k soil moisture ≤ 20 % represents low level, and > 20 % represents high level (Zotarelli et al., 2010).observed in a gas exchange system, which suggested that the evasion from soils is much greater than the emission from the plants grown in the chamber (Frescholtz and Gustin, 2004;Frescholtz et al., 2003).After the factor prescreening step (Figs.S2-S7 in the Supplement), the simulated flux is particularly sensitive to five parameters.Friction velocity (positive effect, p = 0.020), soil Hg content (positive effect, p = 0.028) and soil organic content (negative effect, p = 0.030) are the most significant model parameters (Fig. 7).These effects are similar to the sensitivity results of air-soil exchange over bare lands (Figs. 6 and 7) but slightly weaker based on the p values because of the "shielding" of vegetation coverage that modifies the values of the resistance terms (R b and R ac ) (Zhang et al., 2002a).Highly moist soil (soil moisture content > 20 %, Table 4) has a negative effect because it effectively increases soil diffusion resistance (R g ) (Zhang et al., 2003), although the effect is less significant (p = 0.289).Air temperature also has a positive effect, as anticipated (p = 0.180).The synergistic effect caused by friction velocity and soil Hg content is significant for the aircanopy exchange (p = 0.028, Fig. 7), enhancing the evasion flux by 47 % (77.8 ng m −2 h −1 ).Both soil organic content  4 to isolate the significant factors.The overall flux at canopy level is dominated by the soil flux under the canopy (the sensitivity of foliar exchange is shown in Fig. 8).T denotes surface air temperature, foc denotes fraction of organic carbon in soil, UST denotes friction velocity, SM denotes soil Hg content, MC denotes soil moisture, and " * " denotes interaction effects.and highly moist soil condition can offset the positive effects caused by higher friction velocity, soil Hg content and air temperature at different degrees (Fig. 7), with the soil organic content being more influential.Higher soil organic content at high soil moisture (> 20 %) yields a weak positive effect (p = 0.340), which is interpreted as the combined negative effect of the two parameters is smaller than the sum of the two individual effects.Overall, these characteristics resemble the air-soil exchange because the air-canopy exchange is dominated by the air-soil exchange under canopy.

Sensitivity of exchange over canopy
Atmospheric mercury can deposit on the surface of cuticle or be accumulated in leaves through stomatal uptake (Fig. 1).For cuticular exchange, air temperature has a significant positive effect (Fig. 8).Since air-cuticle exchange is mainly deposition (negative flux), this means that a higher air temperature leads to smaller deposition or greater evasion (p < 0.001).Friction velocity has a strong negative effect (i.e., higher deposition at higher friction velocity, p < 0.001) on the simulated flux.Higher soil organic content (p = 0.009) and highly moist (> 20 %) soil (p = 0.194) increase the simulated flux (i.e., weaken the deposition) by decreasing the canopy compensation point (χ c in Eq. 8).Under this circumstance, Hg deposits preferentially to soil and therefore there is a reduced deposition on cuticle.Higher soil Hg content decreases the flux (p = 0.008) by increasing the overall compensation point (χ cnp in Eq. 8), suggesting greater deposition on cuticle at higher soil Hg content.For stomatal exchange, the trend of the single-factor effect is the same as that of cuticular exchange.4 to isolate the significant factors.T denotes surface air temperature, foc denotes fraction of organic carbon in soil, UST denotes friction velocity, SM denotes soil Hg content, MC denotes soil moisture, and " * " denotes interaction effects.
Several notable interaction effects are observed for foliar exchanges.For cuticular exchange, the deposition is reversed from deposition to evasion at the high air temperature level, leading to the overall positive interaction effect for air temperature and friction velocity (Fig. 8, p < 0.001).The positive effect of soil organic content significantly offsets the negative effect of friction velocity (p = 0.010) and soil Hg content (p = 0.016).For stomatal exchange, the only significant interaction effect is between soil organic and Hg content, which is more strongly dominated by soil organic content.Overall, the foliar exchange is primarily controlled by air temperature and friction velocity because the resistance terms can be affected by the two variables.This is in contrast to the evapotranspiration approach where soil Hg content plays a predominant role in simulated Hg 0 evasion flux (Bash et al., 2004;Gbor et al., 2006).
In this analysis, the effect of solar irradiance is not as significant as the selected parameters under the resistance model scheme and has been ruled out during the prescreening for the model variables (Sect.2.4 and Figs.S2-S6).In the model, solar irradiation can influence the flux in three ways: (1) through modifying the rate constant of Hg(II) reduction in soils and foliage (Eqs.6, 10 and 16), (2) through forcing the change of aerodynamic resistance (R a and R ac ) and (3) through forcing the change of cuticular and stomatal resistance terms (R c and R st ).For air-soil exchange, the effect of solar irradiance on the reduction rate constant is the most sensitive process (Eqs.6 and 10).The photoreduction of Hg(II) in soils has been suggested to be responsible for the increased soil flux observed under sunlit condition (Gustin et al., 2002).There have been kinetic studies showing that increasing UV-A intensity by 75 % approximately doubles the photoreduction rate in the aqueous phase (Qureshi et al., 2010).However, the effect of lights on the kinetics of Hg(II) reduction in soils is poorly understood.In this modeling, the photoreduction rate constant was set to a mean value (Eq.6).This limits a full examination of the true impact of solar irradiation on the simulated Hg flux.Results from experimental studies on Hg(II) photoreduction rates will help reduce this model uncertainty.For foliar exchange, solar irradiation has a weak positive effect on the flux (i.e., slightly weakens deposition, Supplement Fig. S4) but has a significant positive effect on the stomatal exchange (p = 0.004, Supplement Fig. S5).

Conclusions
An updated model for estimating the bidirectional airsurface exchange of Hg is presented based on the current understanding of surface resistance schemes.From the basecase results, water and soil surfaces are net sources and vegetation is a net sink of Hg 0 .Each natural surface exhibits a different diurnal and seasonal variation.Sensitivity analysis of model variables using a two-level factorial design of experiments shows that atmospheric shear flows (surface wind over water and friction velocity of terrestrial surfaces), dissolved gaseous mercury (DGM) concentration, soil organic and Hg content, and air temperature are the most influential factors controlling the magnitude of the atmosphere-biosphere exchange of Hg 0 .However, the positive effect of friction velocity and soil Hg content on the evasion flux from soil and canopy can be greatly offset by the negative effect of soil organic content.Significant synergistic effects are identified between surface wind and DGM level for water surface, and between soil Hg content and friction velocity for soil surface, leading to ∼ 50 % enhanced flux in the combined effect compared to the sum of their individual effects.The airfoliage exchange is mainly controlled by surface resistance terms controlled by environmental parameters such as solar irradiation and air temperature.
The uncertainty in this modeling assessment is primarily from the lack of knowledge in (1) the spatial distribution of organic and Hg content in soil and DGM concentration in water, (2) the reduction mechanism and kinetics of Hg(II) in soil and water and (3) the values of resistance terms over different natural surfaces.More research in providing geospatial distribution of Hg in water and soil will greatly improve the model estimate.Further elucidation on the interaction of Hg and organic carbon in top soil and surface water, as well as quantification of the surface resistance terms specific to Hg species, will also help improve the model scheme.Recent field and experimental investigations have suggested that organic carbon in soil potentially shapes the distribution of Hg in forest at continental scales (Obrist et al., 2011) and that the long-term Hg evasion from soil is highly related to the Hg and organic carbon interactions (Smith-Downey et al., 2010).
Given the predominance of soil organic content in reducing soil Hg evasion flux using the mechanistic approach in this study, soil organic content is likely the controlling factor determining the intensity of air-soil Hg 0 exchange.
The Supplement related to this article is available online at doi:10.5194/acp-14-6273-2014-supplement.

Figure 1 .
Figure 1.Resistance scheme implemented in the air-surface exchange model following Sutton et al. (2007) and Zhang et al. (2009a).
search and Forecasting (WRF) model version 3.4 with the Noah land surface model.The model algorithms were coded in FORTRAN 90 and Network Common Data Form (NetCDF) version 4.1.The gridded model results were visualized by the Visualization Environmental for Rich Data Interpretation (VERDI) version 1.4.

Figure 3 .
Figure 3. Monthly mean of the simulated Hg 0 flux (ng m −2 h −1 ) in the winter month: (a) flux from water bodies, (b) flux from bare lands, (c) flux from soil under the canopy and (d) flux from foliage.

Figure 4 .
Figure 4. Diurnal variation of mean simulated Hg 0 for the entire model domain (UTC-7): (a) flux from 727 water body, (b) total flux from soils (soil under the canopy and bare lands), (c) flux from foliage, and (d) 728 flux for the total domain.729

Figure 4 .
Figure 4. Diurnal variation of mean simulated Hg 0 for the entire model domain (07:00 local time): (a) flux from water bodies, (b) total flux from soils (soil under the canopy and bare lands), (c) flux from foliage and (d) flux for the total domain.

Figure 5 .Figure 5 .
Figure 5. Sensitivity analysis based on the 2 4 factorial design shown in Table 2 (water body).T denotes 731 air temperature at water surface; GEM denotes air Hg concentration; DGM denotes dissolved gaseous Hg 732 concentration in surface water; W denotes wind speed."*" denotes the interaction effects.733 Figure 5. Sensitivity analysis based on the 2 4 factorial design shown in Table2(water bodies).T denotes air temperature at water surface, GEM denotes air Hg concentration, DGM denotes dissolved gaseous Hg concentration in surface water, W denotes wind speed, and " * " denotes the interaction effects.

Figure 6 .Figure 6 .
Figure 6.Sensitivity analysis based on the 2 5 factorial design for bare lands after pre-screening model 735 variables shown in Table 3 to isolate the significant factors.T denotes surface air temperature; foc 736 denotes fraction of organic carbon in soil; UST denotes friction velocity; SM denotes soil mercury 737 content; b denotes scaling factor for Hg reactivity (βHg0)."*" denotes interaction effects.738 Figure 6.Sensitivity analysis based on the 2 5 factorial design for bare lands after prescreening model variables shown in Table3to isolate the significant factors.T denotes surface air temperature, foc denotes fraction of organic carbon in soil, UST denotes friction velocity, SM denotes soil mercury content, b denotes scaling factor for Hg reactivity (β Hg 0 ), and " * " denotes interaction effects.

Figure 7
Figure 7 illustrates the sensitivity of simulated Hg flux over canopy to the model variables (Table4).For comparison, the sensitivity results for air-soil exchange under canopy are also shown.The forcing of air-canopy exchange is dominated by the air-soil exchange under canopy at the two experimental levels.This resembles the Hg 0 emission characteristics

Figure 7 .Figure 7 .
Figure 7. Sensitivity analysis based on the 2 5 factorial design at canopy level after pre-screening model 740 variables shown in Table 4 to isolate the significant factors.The overall flux at canopy level is dominated 741 by the soil flux under the canopy (the sensitivity of foliar exchange is shown in Figure 8).T denotes 742 surface air temperature; foc denotes fraction of organic carbon in soil; UST denotes friction velocity; SM 743 denotes soil Hg content; MC denotes soil moisture."*" denotes interaction effects.744 Figure 7. Sensitivity analysis based on the 2 5 factorial design at canopy level after prescreening model variables shown in Table4to isolate the significant factors.The overall flux at canopy level is dominated by the soil flux under the canopy (the sensitivity of foliar exchange is shown in Fig.8).T denotes surface air temperature, foc denotes fraction of organic carbon in soil, UST denotes friction velocity, SM denotes soil Hg content, MC denotes soil moisture, and " * " denotes interaction effects.

Figure 8 .Figure 8 .
Figure 8. Sensitivity analysis based on the 2 5 factorial design for foliar exchange after pre-screening 746 model variables shown in Table 4 to isolate the significant factors.T denotes surface air temperature; foc 747 denotes fraction of organic carbon in soil; UST denotes friction velocity; SM denotes soil Hg content; 748 MC denotes soil moisture."*" denotes interaction effects.749 Figure 8. Sensitivity analysis based on the 2 5 factorial design for foliar exchange after prescreening model variables shown in Table4to isolate the significant factors.T denotes surface air temperature, foc denotes fraction of organic carbon in soil, UST denotes friction velocity, SM denotes soil Hg content, MC denotes soil moisture, and " * " denotes interaction effects.

Table 1 .
Model variables and units in the base-case simulation.

Table 2 .
Examined model variables and the experimental levels of factorial design for air-water exchange.

Table 3 .
(Jones et al., 2004ables and the experimental levels of factorial design for air-soil exchange over sparsely vegetated land, bare ground tundra and snow/ice surfaces.Kwun and You (2009); b suggested default value for modeling of volatilized contaminant to air by USEPA (2004); c upper limit of the forest soils(Jones et al., 2004); d Akkarappuram and Raman (1988) e Carpi and Lindberg (1998); f has effect on ground and cuticular resistance, Zhang et al. (2003); g Zhang et al. (2012) ; h Zhang et al. (2009a) ; i air temperature below dew point represents low level and higher than dew point represents high level a