Environ. Sci. Technol. 2008, 42, 2964–2969
Evaluation of Two Current Approaches for the Measurement of Carbon Dioxide Diffusive Fluxes from Lentic Ecosystems NICOLAS SOUMIS,* RENÉ CANUEL, AND MARC LUCOTTE GEOTOP/Institut des sciences de l’environnement, Université du Québec à Montréal, C.P. 8 888, Succ. Centre-Ville, Montréal, Québec, Canada, H3C 3P8
Received September 19, 2007. Revised manuscript received January 17, 2008. Accepted January 25, 2008.
BLE and SC depend on methodological and environmental factors or a combination of both. This implies that, under favorable combinations, BLE and SC can provide accurate estimates of fCO2. In this paper, we describe an experiment conceived by Dr. Richard Carignan (University of Montréal, Québec, Canada) that enabled us to evaluate the accuracy of fCO2 measurements performed with the BLE from Cole and Caraco (9) and with a particular model of SC. Both approaches have been used routinely by our research team over the last years (4, 11, 15). Furthermore, this opportunity allowed us to identify factors responsible for the biases related to these approaches. This experiment (DISE) consisted in adding dry ice to a lake and monitoring the subsequent evasion of CO2. Some information reported here was collected by other researchers and the sources of this information are hereafter identified with their name and affiliation.
Experimental Section The dry ice sowing experiment (DISE) consisted in adding dry ice to a lake and monitoring the subsequent evasion of carbon dioxide (CO2). DISE allowed us to evaluate two approaches commonly used for measuring aquatic CO2 diffusive fluxes: the boundary layer equation (BLE) from Cole and Caraco (1998) and a particular model of static chamber (SC). CO2 evasion measurements with both approaches were compared to CO2 mass budgets as a relative reference to define their recovery coefficients (F). F for the BLE and the SC over the whole measurement period were 101 ( 14% and 115 ( 56%, respectively. Results from discrete sampling intervals revealed that the BLE generally provided estimations in good agreement (80–130%) with the mass budgets during both daytime and nighttime. Variations in F for the BLE were related to wind speed and, consequently, piston velocity (k600). The SC overestimated CO2 evasion during daytime (149 ( 39%), and underestimated it during nighttime (57 ( 18%). Variations in F for the SC were related to k600, stemming mainly from the alteration of the air/ water temperature gradient.
Introduction Most unproductive lakes (1, 2) and reservoirs (3, 4) are supersaturated with carbon dioxide (CO2), thus acting as sources of carbon for the atmosphere. Because this gas contributes to climate change (5, 6), accuracy in the assessment of CO2 evasion from lentic ecosystems is crucial. Boundary layer equations (BLE) (e.g., see refs 7-9) and static chambers (SC) (e.g., see refs 4, 10, and 11) are commonly used for measuring CO2 diffusive fluxes (fCO2). Whereas their simplicity and low cost contributed to their widespread use, the accuracy of these approaches is criticized because they alter or disregard certain processes involved in gas exchange (12–14). Previous studies have compared BLE and SC to eddy covariance (EC) and surface renewal models (SRMs) (13, 14). EC and SRMs are generally preferred to BLE and SC because they do not affect gas exchange and better integrate the biogeochemical processes involved, therefore providing more insight into the latter (14). Interestingly, results from BLE and SC are either consistent (14) or inconsistent (13) with those from EC and SRMs, suggesting that biases related to * Corresponding author phone: (514) 987-3000, ext. 6791; fax: (514) 987-3635; e-mail:
[email protected]. 2964
9
ENVIRONMENTAL SCIENCE & TECHNOLOGY / VOL. 42, NO. 8, 2008
Description of the Measurement Approaches. BLE theoretically consider the air/water interface as a thin two-layer film between the gaseous and aqueous phases. The main body of each fluid is assumed to be well mixed, and the principal resistance to CO2 flux between interfacial layers is molecular diffusion (8, 16). Consequently, Fick’s law applies; f CO2 ) kCO2 × ([CO2]surf - (p CO2atm × KH))
(1)
where fCO2 depends on the difference between the actual CO2 concentration at the water surface ([CO2]surf) and the CO2 concentration that water would have if it were in equilibrium with the atmosphere (pCO2atm × KH, where pCO2atm is the atmospheric CO2 partial pressure and KH is Henry’s constant corrected for water temperature). The difference between the actual CO2 concentration in the water and its concentration in equilibrium with the atmosphere was determined empirically by measuring [CO2]surf and pCO2atm (four replicates each). Replicates with over 25% of variation were rejected from calculations. The CO2 piston velocity constant in eq 1 (kCO2) was calculated according to eq 2. Also known as gas transfer velocity or gas exchange coefficient, piston velocity refers to the rate at which a gas crosses the air/water interface. KCO2 ) k600 ×
( ) 600 ScCO2
0.67
(2)
In this experiment, the piston velocity normalized to a Schmidt number (Sc) of 600 (k600) was calculated according to Cole and Caraco (9) (eq 3), where U10 is the wind speed at 10 m above the water surface. k600 ) 2.07 + (0.215 × U101.7)
(3)
The Schmidt number for CO2 corrected for temperature in eq 2 (ScCO2) is calculated according to eq 4 (7), where TW° is the water surface temperature: ScCO2 ) 1911.1 - (118.11TW°) + (3.4527TW°2) -
(0.04132TW°3)
(4)
We also evaluated a modified version of the SC developed by Duchemin et al. (4). The present model (Figure 1) weighs about 5 kg and consists in an acrylic rectangular prism opened at the lower end (29 cm L × 30 cm W × 65 cm H; total volume: 56.6 L; surface area sampled: 870 cm2). Walls of the SC are covered with aluminized polyester sheeting in order to reduce 10.1021/es702361s CCC: $40.75
2008 American Chemical Society
Published on Web 03/07/2008
1 m above the water surface. Because eq 3 is valid for wind measurements performed at 10 m above the water surface, wind speed data (U1) were extrapolated from 1 (Z1) to 10 m (Z10) according to eq 5 (17). U10 ) U1 ×
FIGURE 1. Cross section of the static chamber (SC) model evaluated during DISE (not to scale). any warming effect inside the enclosure. A polystyrene float is fitted around the SC to ensure an air headspace of 34.8 L above the water surface and to maintain its lower end at 25 cm below the surface, this immersed part acting as a skirt. During positioning of the SC, a vent located on the top is opened to prevent pressurization inside the enclosure. The vent is then shut as soon as the SC is positioned. Gas samples (30 mL) were collected from the headspace every 15 min over a 1 h series with 60 mL syringes using a 3 m long tube. Assuming that the variation of CO2 concentrations inside the headspace is a first-order function over the sampling period, fCO2 values were calculated using a linear regression. Acceptance of results is based upon three criteria: (1) the initial CO2 concentration inside the headspace must be ( 10% of that measured in the atmosphere; (2) the correlation coefficient (r2) must be >0.80; and (3) the regression slope must be significantly different from zero (p < 0.05). DISE and Evaluation of Measurement Approaches. DISE was held in August of 2001 in Lake Croche (Laurentian Biological Station, University of Montréal), a small temperate oligotrophic lake consisting of two basins connected by a shallow stream, which limits water exchange between them. DISE was circumscribed to one basin of Lake Croche (area: 63 340 m2; volume: 36 619 m3; mean depth: 4.8 m). A total of 525 ( 25 kg of dry ice was evenly spread with shovels along the 3.2 m bathymetric curve. The amount of dry ice was calculated to mimic the CO2 concentration in Lake Croche during the spring overturn. The 3.2 m curve was chosen because it corresponded to the depth of the thermocline at the beginning of DISE (R. Carignan, University of Montréal). By doing so, we attempted to enrich the base of the epilimnion, while limiting the amount of CO2 reaching the hypolimnion. Dry ice sowing was performed under low wind conditions to maximize the amount of CO2 initially kept by the basin. Three littoral and three pelagic sites were sampled during the two days preceding the addition of dry ice in order to characterize background CO2 emissions produced by respiration and photomineralization. Postsowing sampling was restricted to a littoral and a pelagic site. On each site, water samples were collected at every meter with a 3 L Kemmerer bottle fitted with a gastight spigot in order to compute epilimnetic CO2 mass budgets. Three SC (as replicates) were also deployed. Concomitantly, air and water samples were collected with 60 mL syringes for ƒCO2 measurements using the BLE. Wind speed was measured with an anemometer
( ) Z10 Z1
m
(5)
The exponent m is the vertical extrapolation factor for wind speed, which depends on the surrounding surface roughness, site features, flow characteristics, stability, and height (18). It can be defined empirically through long-term averages of measurement series performed at two different heights (17, 18). However, it is generally assumed that m varies between 0.143 (1/7) for terrain with low relief and 0.200 (1/5) for forested or hilly terrain. So as to depict the moderate morphometry of Lake Croche and its watershed, we set an m value of 0.167 (1/6) as a compromise between both extreme types of terrain (17). To evaluate the BLE and the SC, an assessment of CO2 evasion from the basin was performed with both approaches during 13 (BLE) and 11 (SC) sampling intervals of various duration (26–396 min) over a 31 h period following the onset of an uniform horizontal distribution of the CO2 sublimated from the dry ice. We hereafter refer to this period as the postsowing measurement period (PSMP), its beginning corresponding to time zero. The CO2 evasion rate for a given sampling interval (E0-i) according to both approaches was determined using eq 6, E0-i )
f CO20 + f CO2i 2
× t0-i
(6)
where fCO20 and fCO2i are CO2 diffusive fluxes measured by either the BLE or the SC at the beginning and the end of the sampling interval, respectively; and t0-i is the duration of the sampling interval. E0-i from each approach were converted into total CO2 losses by multiplying the results of eq 6 by the area of the basin. Values for E0-i were then compared to epilimnetic CO2 mass budgets (EMB) computed for the corresponding sampling intervals. EMB were used as a relative reference for determining the accuracy of time-integrated measurements from the BLE and the SC, and were computed from 13 CO2 concentration profiles alternately performed at the pelagic and the littoral sampling sites. Each profile was assumed to be representative of the whole basin at the moment of sampling and was converted to a total amount of CO2 according to hypsometric data (R. Carignan, University of Montréal). Then, cumulated CO2 loss during each sampling interval was determined by the absolute difference between the total amount of CO2 at the beginning and at the end of each sampling interval. Finally, a regression model of the cumulated CO2 loss vs time allowed us to determine a specific amount of CO2 evaded during each sampling interval. Measurements from the BLE and the SC represent the summation of CO2 evasion from both dry ice sublimation and background CO2 emissions, whereas EMB solely account for the difference between an initial and a final CO2 bulk concentration. Gross total evasion of CO2 (GTE) from the BLE or the SC was therefore corrected by subtracting background emissions (BE) in order to reflect the net loss of CO2 attributable to the sublimation of dry ice. BE values were estimated from the average of measurements obtained with both approaches. Recovery coefficients for each approach (F) were finally determined according to eq 7 for each sampling interval. F)
GTE - BE × 100 EMB
(7)
Analyses CO2 analyses were conducted by gas chromatography (GC). A Varian Star-3400 system was equipped with a 1 mL sampling VOL. 42, NO. 8, 2008 / ENVIRONMENTAL SCIENCE & TECHNOLOGY
9
2965
FIGURE 2. CO2 cumulated evasion for the basin of Lake Croche during DISE, according to CO2 concentration profiles and hypsometric data. The amount of CO2 evaded during a given sampling interval (from t0 to ti) is obtained by computing the difference in the cumulated CO2 loss (ti - t0). The area on the graph circumscribed by the two vertical dashed lines corresponds to the PSMP, during which the sampling intervals were performed (time zero corresponds to the first postsowing measurement). The error bars correspond to the uncertainties on the epilimnetic CO2 mass budgets. Circles: profiles from the pelagic sampling site; squares: profiles from the littoral sampling site; open symbols: daytime measurements; solid symbols: nighttime measurements. loop, a steel packed Hayesep-Q column (200 mm long, 3 mm in diameter, 80/100 mesh, isothermal oven at 50 °C), and a thermal conductivity detector. Carrier gas was helium (30 mL min-1). Samples were analyzed within 3 h of their collection. Gaseous samples (SC and atmospheric samples) were directly injected into the GC. Aqueous samples (gas profiles and BLE samples) were first allowed to reach ambient temperature prior to analysis. While in syringes, they were then equilibrated with 30 mL of ultrapure nitrogen by vigorous shaking, following the headspace technique described by McAuliffe (19). Finally, water was discarded and the gas was injected into the GC. CO2 concentrations were calculated according to the solubility coefficient corrected for temperature (20).
Results At time zero of the experiment, more than 60% of the dry ice was remaining in the epilimnion of the basin (R. Carignan, University of Montréal). CO2 profiles indicated that CO2 enrichment was mainly restricted to the first 4 m of the water column. Also, depth-time evolution of isotherms showed that dry ice addition caused a slight erosion of the thermocline, which passed to 3.5 m during the entire PSMP (R. Carignan, University of Montréal). This value was therefore used as the lower depth limit in order to compute the EMB. The best regression fit describing the cumulated CO2 loss within the basin during the PSMP is a second-order function (Figure 2). We calculated average BE values of 0.37 ( 0.23 and 0.20 ( 0.07 mmol CO2 m-2 h-1 during the daytime and the nighttime, respectively. Assuming that the basin metabolism was unaffected by the addition of dry ice, which was confirmed by measurements of photosynthesis and respiration rates (D. Planas, UQAM), BE for the whole basin during the PSMP represented an efflux of 1.03 ( 0.65 kg CO2 h-1 during the daytime and 0.56 ( 0.19 kg CO2 h-1 during the nighttime, which is not taken into account by the EMB. Table 1 provides the recovery coefficients (F) of both approaches for each sampling interval. When considering 2966
9
ENVIRONMENTAL SCIENCE & TECHNOLOGY / VOL. 42, NO. 8, 2008
the sampling intervals separately, F values of the BLE varied between 80 and 130%, whereas F values of the SC varied between 41 and 210%. Over the whole PSMP, average F values of the BLE and the SC were 101 ( 14% and 115 ( 56%, respectively, with no statistical difference (Tukey-Kramer HSD test). In addition, the average F values of the BLE (twosided t test: DF ) 12; p ) 0.8109) and the SC (two-sided t test: DF ) 10; p ) 0.3846) did not significantly differ from a theoretical F value of 100%. During the whole PSMP, the average measurement ratio (SC/BLE) was 1.1. This ratio varied from 0.46 to 1.86 with increasing flux magnitude, indicating that the SC underestimated CO2 evasion at low flux values and overestimated it at high flux values as compared to the BLE. Both approaches provided fairly similar estimates of k600 for U10 ranging between 1.5 and 4.49 m s-1 (Figure 3). However, below 1.5 m s-1, k600 estimates from both approaches increasingly diverged, resulting in distinct y-intercept values. Equation 3 and the regression fit for the SC data indicate respective k600 values of 2.07 and 1.31 cm h-1 in the absence of wind. Another important observation from Figure 3 is that the k600 computed from fCO2 measurements with the SC were highly variable for a given wind speed.
Discussion Our results show that, under the set of environmental conditions encountered during DISE (U10: 0.41–4.49 m s-1; scarce and mild rain events; penetrative convection occurring on most sampling occasions), both approaches allowed for accurate measurements of CO2 evasion when a series of diurnal and nocturnal data are averaged. However, one must not neglect the fact that F varied substantially among the sampling intervals, especially in the case of the SC (Table 1). To explain this situation, we hereafter assess the potential influence of some factors on F.
Accuracy of the BLE Values of F for the BLE showed a slight gradual decline during the PSMP (Figure 4A), as did U10 and k600 (not shown). Figure 4B reveals that F was significantly related to k600, being optimal when k600 was 2.46 cm h-1. Figure 4C shows that wind speed also significantly affected F, this coefficient being optimal when U10 was 1.30 m s-1. The significant relations in Figures 4B and 4C suggest that the BLE evaluated here may not accurately assess the effect of wind on k600, and consequently fCO2. Cole and Caraco (9) obtained a noticeable scatter in their relation of k600 vs U10 (r2 ) 0.61), which strongly suggests that phenomena other than wind speed also influence k600. These include rain-driven turbulence (21, 22) and penetrative convection (14, 23, 24) Although determining k600 on the sole basis of wind speed can lead to important biases, because other processes affecting gas exchange are disregarded (14), it appears that such was not the case during DISE. It is unlikely that rain affected gas exchange during the PSMP since it briefly (ca. 15 min.) occurred twice at a very low rate (4 mm h-1; R. Carignan, University of Montréal). Penetrative convection, which occurred during most sampling, was likely and adequately accounted for by eq 3, because Figure 4D shows no relation between the air/water temperature gradient (∆T ° ) Tair° - Twater°) and F. ∆T ° acts as an indicator for penetrative convection, which is known to promote gas transfer through thermally induced turbulence at the air/ water interface (24, 25). Penetrative convection is assumed to occur when the air is colder than the water surface (∆T ° < 0). When the thin upper layer of the water column loses heat to the colder air above, it sinks and is replaced at the surface by a water layer whose gas concentration is different from that in the atmosphere (23). Despite the
TABLE 1. Recovery Coefficients of the Measurement Approaches for Each Sampling Interval sampling interval (experiment time in d) A (0.091 @ 0.160) B(0.160 @ 0.239) C (0.239 @ 0.305) D (0.305 @ 0.402) E (0.402 @ 0.478)e F (0.478 @ 0.561)e G (0.561 @ 0.641)e H(0.641 @ 0.659)e I (0.659 @ 0.844)e J (0.844 @ 0.909)e K (0.909 @ 1.054) L (1.054 @ 1.144) M (1.144 @ 1.268) A (0.016 @ 0.179) B (0.179 @ 0.245) C (0.245 @ 0.315) D (0.315 @ 0.433) E(0.433 @ 0.502)e F (0.502 @ 0.586)e G (0.586 @ 0.861)e H (0.861 @ 0.923)e I (0.923 @ 1.069) J (1.069 @ 1.169) K (1.169 @ 1.294)
BE (kg)b
EMB (kg)c
G (%)d
Boundary Layer Equation (BLE) 1.66 7.52 ( 0.31 1.90 7.26 ( 0.43 1.58 6.59 ( 1.55 2.33 8.74 ( 2.04 1.82 5.67 ( 0.25 1.99 5.97 ( 0.14 1.92 5.16 ( 0.25 0.43 1.17 ( 0.09 4.44 12.41 ( 0.99 1.56 4.25 ( 0.23 3.48 10.16 ( 0.59 2.16 5.94 ( 1.22 2.98 7.57 ( 1.60
1.70 ( 1.07 1.94 ( 1.22 1.62 ( 1.02 2.39 ( 1.50 1.01 ( 0.34 1.11 ( 0.37 1.07 ( 0.36 0.24 ( 0.08 2.47 ( 0.83 0.87 ( 0.29 3.57 ( 2.24 2.22 ( 1.39 3.05 ( 1.92
4.48 5.02 4.11 5.90 4.51 4.81 4.52 1.00 9.95 3.35 7.21 4.29 5.67
130 ( 25 106 ( 26 121 ( 45 108 ( 43 103 ( 9 101 ( 8 91 ( 10 93 ( 12 100 ( 13 101 ( 11 91 ( 32 87 ( 43 80 ( 44
Static Chambers (SC) 26.35 ( 4.71 8.00 ( 1.83 9.36 ( 1.67 11.95 ( 2.04 4.27 ( 0.72 3.85 ( 1.00 9.77 ( 2.96 2.31 ( 0.36 9.92 ( 1.53 9.88 ( 1.56 10.69 ( 1.18
4.01 ( 2.52 1.62 ( 1.02 1.72 ( 1.08 2.90 ( 1.83 0.92 ( 0.31 1.12 ( 0.38 3.67 ( 1.23 0.83 ( 0.28 3.59 ( 2.26 2.46 ( 1.55 3.08 ( 1.94
10.66 4.18 4.35 7.14 4.06 4.83 14.94 3.18 7.22 4.73 5.66
210 ( 50 152 ( 50 175 ( 46 127 ( 38 82 ( 19 56 ( 22 41 ( 21 47 ( 14 88 ( 38 157 ( 47 134 ( 40
duration (h)
3.91 1.58 1.68 2.83 1.66 2.02 6.60 1.49 3.50 2.40 3.00
GTE (kg)a
a Gross total CO2 evasion (as measured by the SC or the BLE). b Background CO2 emissions (as determined during the presowing period). c Epilimnetic CO2 mass budget. No associated error was calculated for this parameter. d Recovery rate. The absolute error associated to the recovery rate considers the error on the net loss of CO2 evaluated by each measurement approach and the error on the background CO2 emissions. e Indicates nighttime measurements.
FIGURE 3. Relations between piston velocity (k600) and wind speed at 10 m (U10) for the approaches evaluated during DISE. The relation from Cole and Caraco (9) is represented by the solid curve. The dashed line is the linear regression fit obtained with k600 data computed from the flux measurements by the SC (circles: presowing measurements; squares: postsowing measurements; open symbols: daytime measurements; solid symbols: nighttime measurements). absence of a relation between ∆T° and F, Figure 4D suggests that eq 3 was particularly prone to inaccuracy when ∆T° g -2 °C (i.e., when penetrative convection was low or absent). The fact that Cole and Caraco (9) conducted their experiment under conditions where penetrative convection occurred most of the time may constitute an explanation as to why eq 3 failed to correctly estimate k600 when penetrative convection was low/absent during DISE. However, we cannot explain why both under- and overestimation occurred in this situation. Although underestimation and overestimation of fCO2 occurred with the BLE, it is remarkable that this approach
FIGURE 4. Factors influencing the recovery rate (G) of the BLE from Cole and Caraco (9) during the PSMP. (A) Evolution of G during the PSMP. (B) Influence of k600 on G. (C) Influence of wind speed (U10) on G. (D) Influence of air/water temperature gradient (∆T ° ) Tair° - Twater°) on G. Open circles: daytime measurements; solid circles: nighttime measurements. The intersection of dashed lines indicates the x value at which G is optimal (i.e., 100%) according to regression fits. generally provided fairly good estimates of CO2 evasion during DISE. However, the fact that DISE was performed under calm conditions and that Cole and Caraco (9) based their BLE on data from low-wind sites brings up an important question regarding either large lakes with great fetches or hydroelectric reservoirs in which vertical advection is frequent due to both VOL. 42, NO. 8, 2008 / ENVIRONMENTAL SCIENCE & TECHNOLOGY
9
2967
FIGURE 5. Factors influencing the recovery rate (G) of the SC during the PSMP. (A) Evolution of G during the PSMP. (B) Influence of k600 on G. (C) Influence of wind speed (U10) on G. (D) Influence of air/water temperature gradient (∆T ° ) Tair° Twater°) on G. Open circles: daytime measurements; solid circles: nighttime measurements. The intersection of dashed lines indicates the x value at which G is optimal (i.e., 100%) according to regression fits. large fetches and huge water outflows. In these environments, it is possible that the BLE evaluated here underestimates the actual k600 and the ensuing fCO2.
Accuracy of the SC The evolution of F for the SC followed a diel cycle during the PSMP (Figure 5A), and the significant difference between diurnal (149 ( 39%) and nocturnal (57 ( 18%) values of F is substantiated by a t-test (DF ) 9; p ) 0.0016). Furthermore, one-sided t-tests corroborate that the SC overestimated fCO2 during the day (high side; DF ) 6; p ) 0.0076) and underestimated them during the night (low side; DF ) 3; p ) 0.0086). Figure 5B indicates a positive relation between k600 and F, and that the best estimates of fCO2 occurred when k600 was 2.44 cm h-1. Interestingly, this value is very similar to the optimal k600 for the BLE (2.46 cm h-1; Figure 4B). No relation was found between U10 and F (Figure 5C). The absence of a wind effect on SC accuracy during DISE contradicts the findings from Matthews et al. (16), who identified wind as a possible source of error while estimating diffusive fluxes with their SC. Contradictions between both experiments probably stem from the distinct weights of the SC tested. It is likely that our 58.5 L SC, made of thick acrylic walls, is heavier than that used by Matthews et al. (16), which is made from a 16.4 L polycarbonate water bottle. The heavier weight of our SC possibly reduces wind-induced movements of the SC on water, at least under the low wind conditions encountered during DISE. For instance, Matthews et al. (16) noted that their SC were constantly dragging the water surface at wind speeds
