Introduction

Levels of carbon dioxide (CO2) not only can serve as an indicator of autotrophic or heterotrophic nature of an aquatic environment1, but also can be an important parameter to elucidate calcification and potential pH changes in a water body2,3,4. Due to increasing human consumption of fossil fuels and other anthropogenic activities, atmospheric CO2 level has increased from 280 μatm in the 1800 s to levels above 400 μatm in 2015 5. As the most important greenhouse gas on earth, inventory and fluxes of CO2 across different reservoirs and its role in the global carbon cycle and thus climate and environmental changes have been the focus of many recent research programs6,7,8.

Unlike the open ocean which contains abundant dissolved inorganic carbon (DIC, usually >2000 μmol/kg) and typically serves as a sink of CO27,9,10,11, lakes and inland waters could serve as a source of CO2 to the atmosphere12,13,14 and may play an essential role to the local, regional and even global carbon cycles and climate. Recently, carbon cycling and magnitude of CO2 fluxes across the air-water interface in lakes and inland waters have received increasing attention15,16,17,18, although studies on carbon dynamics in the Great Lakes remain scarce17,19.

The Laurentian Great Lakes are the largest freshwater system on Earth and receive vast amounts of organic and inorganic carbon from surrounding terrestrial ecosystems. Over the past decades, the Great Lakes have experienced significant ecological and environmental changes due to the introduction of invasive species, notoriously nonindigenous quagga mussels (Dreissena rostiformis bugensis), leading to a decrease in primary production and increasing water clarity in the Great Lakes20,21,22,23,24,25. In addition to impacts on foodweb structure and ecological functions, changes in biogeochemical cycles of nutrients have also been documented26,27,28. However, specific changes in carbon dynamics, as well as the impacts and biogeochemical consequences after the colonization of invasive quagga mussels in the Laurentian Great Lakes remain poorly understood. The direction and magnitude of air-lake CO2 fluxes after dreissena introduction in the Great Lakes are largely unquantified, especially in Lakes Michigan, Huron and Ontario17,29,30,31.

To examine the response of carbon dynamics to the introduction of invasive quagga mussels and linkages between invasive species and changes in biogeochemical cycling, open lake water samples were collected from all of the Laurentian Great Lakes, including Lake Superior, Lake Michigan, Lake Huron, Lake Erie and Lake Ontario, during August 2013 (Fig. 1). Furthermore, seasonal water samples were also collected from Lake Michigan between 2013 and 2015. In addition to water isotopic composition (δ2H and δ18O), total alkalinity (TA), dissolved inorganic carbon (DIC), and pH were measured prior to the evaluation of the partial pressure of CO2 (pCO2) and air-water CO2 fluxes.

Figure 1
figure 1

Sampling locations in the Laurentian Great Lakes during August 2013 (The map is created from Ocean Data View (Schlitzer, R., Ocean Data View, odv.awi.de, 2015) and edited on Microsoft PowerPoint 14.5.8 (https://products.office.com/en-us/home).

Results

Variations in pH, TA and DIC in the Great Lakes

The Laurentian Great Lakes can be characterized as a high-pH and high-carbonate ecosystem although Lake Superior had a relatively lower pH and carbonate abundance compared to the other Great Lakes (Fig. 2 and Table 1). For example, pH values in Great Lakes waters were typically higher than 8 except for two sampling sites in Lake Superior (Table 1). Lakes Michigan, Erie and Ontario had higher pH values, averaging 8.22 ± 0.04, 8.19 ± 0.08 and 8.20 ± 0.04, respectively, followed by Lake Huron (8.07 ± 0.06) and Lake Superior (7.98 ± 0.05).

Table 1 Surface water pH, total alkalinity (TA) and dissolved inorganic carbon (DIC), as well as partial pressure of CO2 (pCO2) and air-lake CO2 fluxes (FCO2) from the Laurentian Great Lakes based on TA-pH and DIC-pH data, respectively.
Figure 2
figure 2

Spatial distributions of surface water pH, total alkalinity (TA), dissolved inorganic carbon (DIC), dissolved organic carbon (DOC) and water isotope composition (δ2H and δ18O) in the Laurentian Great Lakes during the summer of 2013.

The figure is created from MATLAB R2013 (http://www.mathworks.com/products/matlab/index.html?s_tid=gn_loc_drop) and modified by Microsoft PowerPoint 14.5.8 (https://products.office.com/en-us/home).

Concentrations of DIC and TA typically exceeded the threshold of 1000 μmol/kg except for Lake Superior (Table 1). Values of TA varied from 830 to 2200 μmol/kg for all the Great Lakes, with an average of 843 ± 10 μmol/kg in Lake Superior, 2168 ± 22 μmol/kg in Lake Michigan, 1600 ± 59 μmol/kg in Lake Huron, 1872 ± 22 μmol/kg in Lake Erie and 1815 ± 19 μmol/kg in Lake Ontario (Table 1). Similarly, different DIC concentrations were observed among the five Great Lakes, showing the highest value in Lake Michigan (average of 2065 ± 24 μmol/kg) and lowest values in Lake Superior (average of 793 ± 26 μmol/kg). Generally, these DIC and TA values are higher than those in most freshwater systems (e.g., refs 32, 33 and 34) and comparable to those of global ocean basins35 (especially data of Lake Michigan, Table 1S). Our data here are also similar to the long-term data from previous studies in the Great Lakes36. High DIC and TA in Great Lakes waters are related to the spatial distribution of limestone and carbonate weathering in the Great Lakes basin. It is the high carbonate abundance that allows the thriving and rapid colonization of the mussel community in the Great Lakes37, which in turn has “re-engineered” the lake ecosystem during the past decades22,23,24,25.

The surface distribution (Fig. 2) showed that TA, pH and DIC concentration all increased consistently from the upper Great Lakes (i.e., Lake Superior) to Lake Huron and then to the lower Great Lakes (i.e., Lakes Erie and Ontario), similar to the increasing δ2H and δ18O values in surface waters along the Great Lakes. This is consistent with the general water transport pathway from the upper to the lower Great Lakes, suggesting an accumulative effect on water chemistry in the lower Great Lakes from surrounding riverine inputs38. Among the five Great Lakes, Lake Michigan had notably high values for each of the measured parameters including pH, TA, DIC, and water isotopes (Fig. 2 and Table 1), which is related to its high abundance of limestone and the fact that it is a semi-closed basin which allows the accumulative effect of evaporation in Lake Michigan.

Variations in pCO2 along the Great Lakes

As shown in Table 1 and Fig. 3, our results consistently showed an evident CO2 supersaturation (>398 μatm) for all the Great Lakes, with the highest pCO2 in Lake Michigan (762 ± 88 μatm) and Lake Huron (774 ± 92 μatm), followed by Lakes Erie (725 ± 126 μatm) and Ontario (647 ± 47 μatm). Even in Lake Superior, with the lowest TA and DIC abundance (both <1000 μmol/kg) and lowest pCO2 values (461 ± 77 μatm), CO2 supersaturation was still evident (Fig. 3). This indicated that surface waters of all the Laurentian Great Lakes can be a net CO2 source to the atmosphere, consistent with those observed for other global lakes16,18.

Figure 3
figure 3

The averaged partial pressure of CO2 (pCO2) in the Laurentian Great Lakes during summer months from 1983–2006 and 2013.

LS, LM, LH, LE and LO denote Lakes Superior, Michigan, Huron, Erie and Ontario, respectively. The horizontal line represents the averaged atmospheric pCO2 value in 2013 (398 μatm, from NOAA at http://www.esrl.noaa.gov/gmd/). Data from 1983–2006 are calculated from surface water pH and alkalinity obtained from the EPA’s biannual Great Lakes surveillance program, which begin from 1983 except for Lake Superior beginning from 1992 (see Urban and Desai (2009)30). The pCO2 values for Lakes Superior, Huron and Erie were derived from the averaged values from previous studies17,30,39, excluding outlier data points.

Additionally, our seasonal pCO2 data from Lake Michigan between 2013 and 2015 also showed a significant CO2 supersaturation during all sampling seasons. For example, average pCO2 was 762 ± 88 μatm in August 2013, 699 ± 17 μatm in May 2014, 769 ± 124 μatm in October 2014, 815 ± 48 μatm in October 2015, and 944 ± 83 μatm in December 2015 (Fig. 4). Furthermore, CO2 supersaturation existed throughout the whole water column of open Lake Michigan although pCO2 values decreased slightly in the middle water column due to biological uptake (see examples of vertical profiles in October 2015, Fig. 4).

Figure 4
figure 4

Seasonal variations in partial pressure of CO2 (pCO2) in surface waters of Lake Michigan (upper panel) between 2013 and 2015 and the vertical profiles of pCO2 and other hydrological data (pH, dissolved inorganic carbon, DIC) in open Lake Michigan during October 2015 (lower panels).

The horizontal line represents the averaged atmospheric pCO2 value (398 μatm) during 2013–2015.

Discussion

The invasive species, especially filter feeders including zebra and quagga mussels, have been shown to strongly affect the foodweb structure and biogeochemical cycles of nutrients after their introduction to the Great Lakes21,22,23,24,25. Thus, the equilibrium and interaction of CO2 between lake and atmosphere may have been altered during the past decades, particularly after quagga mussels became predominant. As shown in Fig. 3, the pCO2 values in the Great Lakes observed in 2013 are significantly higher than the median values of pCO2 observed during the same sampling season between 1983–200617,30,39 (Table 1), especially for Lakes Michigan (p < 0.001), Huron (p < 0.001), Erie (p = 0.004) and Ontario (p = 0.035). Compared to the increase in atmospheric CO2 levels from 1983 to 2013 (~60 μatm), the increase in pCO2 in the water column of the Great Lakes was considerably larger. Taking Lakes Michigan and Huron as examples, during 1983–2006, the average ∆pCO2 value, the difference in pCO2 between atmosphere (averaging ~362 μtm during that time periods) and lake waters, was about 118 and 22 μatm for Lakes Michigan and Huron30,39, respectively. In comparison, ∆pCO2 values in these two lakes during summer 2013 reached a new high of ~400 μatm, averaging 365 ± 89 and 377 ± 92 μatm for Lake Michigan and Lake Huron, respectively (Table 1), showing an increase of 4 times in ∆pCO2 values in August between 1983–2006 and 2013 (Fig. 3). For Lakes Erie and Ontario, the increase in pCO2 values between 1983–2006 and 2013 was also significant, although large errors existed in historical CO2 data30,39 for the shallower lakes with high primary production and anthropogenic influences including seasonal hypoxia and algal bloom, especially in Lake Erie40,41.

Moreover, our seasonal CO2 data and its vertical distribution in Lake Michigan (Fig. 4) further demonstrated the increased CO2 supersaturation in the water column. For example, pCO2 in Lake Michigan during October 2015 ranged from ~600 μatm in subsurface waters to ~770 μatm in deeper waters (Fig. 4). Elevated pCO2 in deeper waters resulted from the decomposition of organic matter and decreased biological uptake also implied a potential higher CO2 emission flux from lake waters when the CO2-enriched deeper water is upwelled to surface during winter mixing seasons42. Furthermore, our seasonal data clearly show a consistent increase in pCO2 from summer to fall and then to winter, reaching as high as 944 ± 83 μatm during December 2015; an increase by ~30% compared to other seasons (Fig. 4). Thus, the pCO2 data obtained during summer 2013 (Fig. 3) might represent the lower limit of CO2 supersaturation in the Great Lakes, although other factors could have an influence, such as partial ice-cover during winter.

In contrast to Lake Michigan and other lower Great Lakes, the increase in pCO2 in Lake Superior during summer between 1983–200617,30 and 2013 does not seem to be significant (p = 0.665, Fig. 3, Table 1). Nevertheless, increased CO2 supersaturation in Lake Michigan and other lower Great Lakes coincides with the distribution and density of invasive quagga mussels among all the Great Lakes. For example, high quagga mussel densities (e.g., >19,000-mussels/m2) have been reported for Lake Michigan22,25,27,43, but a very low density for Lake Superior. Although zebra mussels arrived into the Great Lakes before quagga messuls, they have much weaker filter ability, and changes in foodweb structure and ecological function were more evident only after quagga mussels became predominant25,44,45,46. In addition, quagga mussels did not start blanketing offshore regions of the Great Lakes, such as Lake Michigan until 200525,44. Therefore, CO2 data collected during 1983–200617,30,39 can be considered representative of the time period prior to or during the early stage of the colonization of invasive quagga mussels in the Great Lakes. Linking the spatiotemporal variations of quagga mussel population with the changes in the extent of CO2 supersaturation during summer sampling time, we hypothesized that increasing CO2 supersaturation in the Great Lakes during the past decade has been induced by the colonization of invasive quagga mussels and the subsequent biogeochemical response to that colonization.

Possible pathways and mechanisms to support our hypothesis include (1) decreasing primary production in the Great Lakes after the introduction of invasive quagga mussels20,21,22,23,47,48; (2) increase in water clarity24 which enables the light penetration into deeper waters and consequently enhances the photo-degradation of natural organic matter in the water column; and (3) metabolic processes of vast quantities of invasive quagga mussels blanketing lake floor27,49, although other processes, such as changes in sediment-water processes, coastal erosion and water chemistry after the colonization of quagga mussels, are also important. All these changes are inter-correlated and are the result of the colonization of invasive quagga mussels and all these processes would directly or indirectly favor the accumulation and release of CO2 in the water column, resulting in increased CO2 abundance in the Great Lakes, especially those heavily infested with quagga mussels (e.g., Lake Michigan25). The close coupling between invasive species and carbon dynamics elucidated here in the Great Lakes clearly shows how small invasive quagga mussels could have caused basin-scale changes in carbon dynamics and biogeochemical cycling.

The increase in pCO2 levels after the colonization of invasive quagga mussels in the Great Lakes evidently exceeded the increase in average atmospheric CO2 level from 1983 to 2013 (~60 μatm). Thus, enhanced CO2 emission fluxes from lake waters to the atmosphere can be expected. Based on these summer pCO2 data, the daily CO2 emission fluxes ranged from 8.37 to 14.48 mmol-C/m2/d with an average of 10.38 ± 2.52 mmol-C/m2/d for Lake Michigan and ranged from 7.69 to 13.60 mmol-C/m2/d, averaging 10.71 ± 2.61 mmol-C/m2/d for Lake Huron (Fig. 5 and Table 1). These CO2 emission fluxes were 3 to 4 times higher than those in Lakes Michigan and Huron at the same sampling time during 1983–2006, based on their differences in ∆pCO2 values. Lakes Erie and Ontario had CO2 emission fluxes of 9.35 ± 3.59 mmol-C/m2/d and 7.13 ± 1.33 mmol-C/m2/d respectively, while an emission flux of 1.84 ± 2.17 mmol-C/m2/d was estimated for Lake Superior. Overall, the emission fluxes of CO2 in all five Great Lakes follow the order of Huron ≥ Michigan > Erie > Ontario > Superior during summer 2013 (Fig. 5 and Table 1), distinct from the order of Superior > Michigan > Huron > Ontario > Erie observed prior to or during the early stage of quagga mussel colonization30,50 and further attests to the strong impact of invasive species on carbon cycles in the Laurentian Great Lakes.

Figure 5
figure 5

The averaged daily CO2 emission fluxes (FCO2), and lake-wide integrated CO2 emission fluxes in all Laurentian Great Lakes during the summer of 2013, derived from pH-dissolved inorganic carbon (pH-DIC) data pair.

LS, LM, LH, LE and LO denote Lakes Superior, Michigan, Huron, Erie and Ontario.

Although CO2 emission fluxes from the Great Lakes were comparable to or lower than those from small lakes (e.g., refs 51, 52 and 53), the large surface areas still made the Great Lakes a strong CO2 source to the atmosphere. We estimated the lake-wide integrated CO2 fluxes for each Great Lake and those estimated flux values are 1.5 × 108 mol-C/d for Lake Superior, 6.0 × 108 mol-C/d for Lake Michigan, 6.4 × 108 mol-C/d for Lake Huron, 2.4 × 108 mol-C/d for Lake Erie, and 1.4 × 108 mol-C/d for Lake Ontario (Fig. 5). Collectively, up to (7.7 ± 1.0) × 1012 g-C/yr or 7.7 ± 1.0 g-C/yr in the form of CO2 can be released to the atmosphere from the Laurentian Great Lakes, which is comparable to the annual DIC flux to the ocean from the Mississippi River and many other world rivers54,55. This annual CO2 emission flux from the Laurentian Great Lakes is also higher than the annual organic carbon burial rate in all inland seas16,56, and makes up 1% of the annual CO2 emission from global deforestation57. It should be noted that our estimated CO2 emission flux could be the lower limit since it was calculated from summer sampling when pCO2 levels were lower relative to winter seasons42 (e.g., Fig. 4). Consequently, the Laurentian Great Lakes certainly serve as a significant CO2 source to the atmosphere especially after the colonization of invasive species, and play an essential role in regional and global CO2 budgets and cycling. This is a paradox in that the Great Lakes with similar carbonate abundance and pH as those in global oceans, should be expected to absorb CO2 from the atmosphere and become acidified in the face of the accelerating rise of atmospheric CO258. However, the Great Lakes have responded actively to the changing ecosystem and are making a positive feedback to climate change, with decreased whiting events (or CaCO3 saturation/precipitation) during summer24, increased CO2 supersaturation and emission fluxes after the introduction of invasive quagga mussels and resultant decrease in primary production, enhanced degradation of natural organic matter due to increased water clarity, and direct CO2 release from quagga mussel respiration27,49.

Materials and Methods

Sampling

During August 2013, surface water samples were collected from open lake stations on all of the Great Lakes: including Lake Superior, Lake Michigan, Lake Huron, Lake Erie and Lake Ontario (LS, LM, LH, LE and LO). Specific sample locations are shown in Fig. 1 and Table 1. Additionally, seasonal transect sampling, including a vertical profile at an open lake station (43°11.5362′N; 87°39.9566′W, 104 m depth), was also conducted in open Lake Michigan during 2013 to 2015. Detailed sampling locations can be found in ref. 28. Surface water samples were taken from a hull-mounted, all-Teflon and stainless steel, high-speed pumping system (2 m depth). Hydrographic parameters were recorded with a Hydrolab Datasonde 5A and coordinated with ship’s positioning systems, such as temperature (°C) and chlorophyll-a fluorescence (V). Vertical profiles of hydrographic characteristics (SeaBird SBE 25Plus) and Niskin bottle-collected water samples were taken in open Lake Michigan at the offshore station. After collection, samples were saved for pH and total alkalinity (TA) and water isotopes, and filtered through pre-rinsed syringe filters for the measurements of dissolved inorganic carbon (DIC) and dissolved organic carbon (DOC).

Analysis

Concentrations of DOC were measured on a Shimadzu TOC-L analyzer using the high temperature combustion method59. Total dissolved carbon (TDC) was measured on the same TOC analyzer without acidification and sparging. The concentrations of DIC were calculated by the difference between TDC and DOC concentrations60. Community consensus seawater reference from the University of Miami for DOC and the certified reference seawater (CRS) from the Scripps Institution of Oceanography for DIC53, as well as working standards were measured as samples to ensure data quality. After calibration with three standard solutions (pH = 4.00, 7.00 and 10.00), the pH electrode (Sartorius PB-11) was used to measure the pH. The precision and accuracy of pH was ±0.01. Based on the Gran titration procedure61, the sample (~30 mL) for TA was titrated by the CRS-calibrated HCl solution (~0.02 M) to an endpoint pH of 4.5. The results of DIC, DOC and TA are all reported in μmol/kg. Stable hydrogen and oxygen isotopic composition (δ2H and δ18O) of lake waters was measured on a Picarro cavity ring down spectrometer (L2130-i). Standard mean ocean water (Kona water) was used as a standard for both δ2H and δ18O with a precision of ±0.02‰.

Calculations of pCO2 and flux of CO2

pCO2 was calculated using the CO2 program62 and measured pH and DIC concentrations or TA under the freshwater option. The CO2 fluxes (F) were calculated based on the one-dimension flux model:

where k (cm/h) is the gas transfer velocity of CO2, and the Ko (mol/m3/atm) is the solubility coefficient of CO2 at in situ temperature and salinity. Values of k and Ko were calculated based on the method of ref. 63 and ref. 64, respectively, and the monthly average wind speed in the study area was used for the calculation. The pCO2-water and pCO2-air are the partial pressures of CO2 (in μatm) in surface waters and air, respectively. The global annual averaged surface air pCO2 for 2013 from NOAA (http://www.esrl.noaa.gov/gmd/) was used as the pCO2-air value (398 μatm). When the difference in partial pressure of CO2 between water and air (∆pCO2) is higher than zero, the emission of CO2 from surface waters to the atmosphere occurs. In contrast, a negative value of ∆pCO2 indicates uptake of atmospheric CO2 in surface waters.

Uncertainties of pCO2 and statistical analysis

Both DIC-pH and TA-pH pairs were used to derive the pCO2 values in lake waters. Uncertainties of pCO2 derived from pH ( ± 0.01) were estimated to be  ± 16 μatm. For the discussion on CO2 fluxes, pCO2 data derived from pH-DIC were used, while data derived from pH-TA were used to evaluate possible contribution from non-carbonate alkalinity65. As shown in Table 1, pCO2 values calculated from pH-TA pair were, on average, 48.8 ± 14.4 μatm higher than those derived from pH-DIC pair. In other words, overestimation of pCO2 from non-carbonate alkalinity could be up to 6.5 ± 1.4% in Great Lake waters even under the generally low DOC concentrations (Fig. 2). A one-way ANOVA analysis was performed to determine the significance of differences in pCO2 between different sample groups.

Additional Information

How to cite this article: Lin, P. and Guo, L. Do invasive quagga mussels alter CO2 dynamics in the Laurentian Great Lakes? Sci. Rep. 6, 39078; doi: 10.1038/srep39078 (2016).

Publisher's note: Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.