Feedback Between Deglaciation, Volcanism, and Atmospheric CO2

An evaluation of the historical record of volcanic eruptions shows that subaerial volcanism increases globally by two to six times above background levels between 12 ka and 7 ka, during the last deglaciation. Increased volcanism occurs in deglaciating regions. Causal mechanisms could include an increase in magma production owing to the mantle decompression caused by ablation of glaciers and ice caps or a more general pacing of when eruptions occur by the glacial variability. A corollary is that ocean ridge volcanic production should decrease with the rising sea level during deglaciation, with the greatest effect at slow spreading ridges. CO 2 output from the increased subaerial volcanism appears large enough to in ﬂ uence glacial/interglacial CO 2 variations. We estimate subaerial emissions during deglaciation to be between 1000 and 5000Gt of CO 2 above the long term average background ﬂ ux, assuming that emissions are proportional to the frequency of eruptions. After accounting for equilibration with the ocean, this additional CO 2 ﬂ ux is consistent in timing and magnitude with ice core observations of a 40 ppm increase in atmospheric CO 2 concentration during the second half of the last deglaciation. Estimated decreases in CO 2 output from ocean ridge volcanoes compensate for only 20% of the increased subaerial ﬂ ux. If such a large volcanic output of CO 2 occurs, then volcanism forges a positive feedback between glacial variability and atmospheric CO 2 concentrations: deglaciation increases volcanic eruptions, raises atmospheric CO 2 , and causes more deglaciation. Such a positive feedback may contribute to the rapid passage from glacial to interglacial periods. Conversely, waning volcanic activity during an interglacial could lead to a reduction in CO 2 and the onset of an ice age. Whereas glacial/interglacial variations in CO 2 are generally attributed to oceanic mechanisms, it is suggested that the vast carbon reservoirs associated with the solid Earth may also play an important role.


Introduction
Volcanism fluxes carbon from Earth's interior to its exterior fluid envelope and links among volcanism, the carbon cycle, and climate are well known to operate at long time scales (≥10 6 yr) (Walker et al., 1981). Short-term (10 0 yr) changes in climate and weather also result from volcanic eruptions; for example, the 1991 Mount Pinatubo eruption injected enough aerosols into the atmosphere to decrease Earth's surface temperature by half a degree Celsius for about a year (Hansen et al., 1992). Evidence has also accumulated that short-term changes in weather and environment influence volcanism: the subtle variations associated with Earth's tides (Johntson and Mauk, 1972;Hamilton, 1973;Sparks, 1981), daily variations in atmospheric pressure and temperature (Neuberg, 2000), seasonal changes in water storage (Saar and Manga, 2003;Mason et al., 2004), and other short-term changes in the environment (Kennett and Thunell, 1975;Rampino et al., 1979;Dzurisin, 1980) have all been implicated as influencing the timing of volcanic eruptions.
In this paper we explore the relationships among glacial loading, volcanic eruptions, and climate on the intermediate time scales pertinent to glacial/interglacial variations (10 3 -10 5 yr). If subtle changes in weather and climate are sufficient to influence the timing of volcanism, even larger effects could be expected from the massive changes associated with deglaciation (Hall, 1982). Indeed, deglaciation is observed to coincide with increased volcanism in Iceland (Sigvaldason et al., 1992;Gudmundsson, 2000;Sinton et al., 2005;Licciardi et al., 2007), France and Germany (Nowell et al., 2006), the Western United States (Jellinik et al., 2004;Bacon and Lanphere, 2006), and Chile (Best, 1992;Gardeweg et al., 1998). Of these, the clearest demonstration comes from Iceland, where dates from lava flows (Sigvaldason et al., 1992), sulphate concentrations in Greenland ice cores (Zielinksi et al., 1997), and table mountains (Licciardi et al., 2007) are all consistent with increased volcanism during or following deglaciation. The effect on Iceland has also been modeled to result from decompression melting of the mantle caused by removal of an approximately 2 km thick ice sheet during deglaciation (Jull and Mckenzie, 1996;Maclennan et al., 2002). Furthermore, volcanic Earth and Planetary Science Letters 286 (2009) [479][480][481][482][483][484][485][486][487][488][489][490][491] eruptions are more likely to occur when the confining pressure associated with magma reservoirs and fluid transport in confining rock is decreased by ice removal.
Observations and some theory thus suggest that deglaciation increases volcanism. In this paper we assess the global extent of changes in volcanism during the last deglaciation and consider the consequences for climate, particularly with regard to changes in the concentration of atmospheric CO 2 .

Volcanic activity through the last deglaciation
To explore the global extent and magnitude of increased volcanism during the last deglaciation, we combine together two datasets (Siebert and Simkin, 2002;Bryson et al., 2006) comprising the date (Fig. 1) and location (Fig. 2) of eruptions over the last 40 ky. Although the dataset from Siebert and Simkin (2002) covers only the Holocene, it is more complete over this interval and gives an indication of the size of the eruption using the Volcanic Explosivity Index (VEI), which will be discussed in more detail later. Redundant events between the datasets are removed, as are events without age estimates or which are not bracketed between certain dates. We also exclude small events (VEI ≤ 2) because their consistent identification is less likely (Siebert and Simkin, 2002), leaving a total of 5352 volcanic events during the last 40,000 yr. The data are included in the electronic supplement.
The dates of most volcanic events are uncertain (Siebert and Simkin, 2002), and we use a probability distribution to describe when each occurred. The reported calendar age uncertainties are used for each event or, in the case of radiocarbon, the dates and their uncertainties are adjusted to calendar ages using the CALIB 5.0 program (Stuiver et al., 2005). Ages without a reported uncertainty are assumed to have a normal probability distribution with a standard deviation of 10% of the age, which is large relative to most dating uncertainties. In this manner, the dataset of 5352 events and their uncertain ages is transformed into an equal number of probability distributions spanning the interval from 40,000 yr ago to the present (Fig. 1). The data show a marked observational bias with 80% of the dated eruptions occurring in the last 1000 yr. This temporal bias presents a major challenge to assessment of the amount and distribution of volcanism.

Mapping of volcanic activity
In its most basic form, the hypothesis of deglacial triggering of volcanism implies that volcanism should increase when and where ice decreased. To test this prediction, we map the average frequency of eruptions during three distinct epochs: the glacial (g, 40-20 ka), the deglacial (d, 17-7 ka), and the late Holocene (h, 5-0 ka) as a function of latitude, ϕ, and longitude, θ. Mapping is accomplished using a timeand space-weighted average, f ðϕ; θÞ = ∑ N i = 1 P i;j λ 2 = ðλ + r i ðϕ; θÞÞ 2 , where P i,j is the probability that event i occurred during interval j. P i,j is multiplied by a spatial-weighting term that depends on the distance, r, between each point on the globe (ϕ, θ) and the location of each volcanic event, i, and has a smoothing length scale, λ, of 500 km.
The frequency maps are, to first order, consistent with the deglacial triggering hypothesis. The deglacial frequency of eruptions has maxima at high latitudes, whereas the glacial and late Holocene maps show less distinct features (Fig. 2). Owing to the temporal bias, however, it is more informative to consider relative changes in frequency between epochs (Fig. 3). We will refer to the ratio of deglacial to glacial eruption frequencies, taken on a point-by-point basis, as d/g (Fig. 3a) and the ratio of deglacial to Holocene eruption frequencies as d/h (Fig. 3c).
Values of d/g are greatest in those regions that experienced deglaciation (Fig. 3a,b). In particular, the Southern Andes (McCulloch et al., 2000), Alaska and the Aleutian Arc (Denton and Hughes, 1981;Yu et al., 2008), the Cascades and Cordillera regions (Denton and Hughes, 1981;Dyke, 2004), Iceland (Licciardi et al., 2007), Western Europe (Boulton et al., 2001), and Eastern Russia (Grosswald, 1998;Bigg et al., 2008) all experienced significant deglaciation and exhibit higher ratios. Conversely, volcanic regions that experienced little ice unloading exhibit low d/g ratios (Fig. 3c,d), to include Southeast Asia, Africa, and the tropical Americas. Furthermore, the d/g and d/h patterns resemble one another, having a cross-correlation of 0.7 when sampled at the sites of each of the 5352 eruptions in our database. One major difference between the d/g and d/h volcanic frequency ratios is that only d/h registers increased values in Northern New Zealand and Hawaii, both of whose volcanoes experienced deglaciation (Alloway et al., 2007;Blard et al., 2007).
A concern that arises in mapping eruption frequency is that differential preservation of volcanic evidence in glaciated and unglaciated regions will lead to biases. Processes such as ash deposited atop ice or tephra emplaced and then scoured by ice would tend to destroy evidence of eruptions during the glacial period. However, glacial destruction of evidence would tend to increase d/g ratios and lower d/h ratios. The consistent d/g and d/h patterns suggest that the signal is not a preservational artifact. Furthermore, the ratio of glacial to late Holocene volcanism (g/h, not shown), which ought to be most susceptible to preservational bias, shows no systematically high or low values in glaciated regions.

Deglaciation versus volcanic activity
To better quantify the relationship between volcanism and deglaciation, we compare the distribution of volcanic frequency ratios to an independent estimate of the magnitude of deglaciation. The dimensions of the great ice sheets are relatively well documented for the Last Glacial Maximum, but no comprehensive map of the extent of  (Siebert and Simkin, 2002;Bryson et al., 2006) for (a) the period prior to 2 ky bp and (b) the period after 2 ky bp. The continuous line indicates the most probable year of an eruption and the horizontal bars indicate the 95% confidence interval. Events are ordered from the oldest to youngest most probable year. Note the observational bias toward recent years-77% of the observed events occur within the last 1000 years. mountain glaciers and ice caps exists. To obtain a global indication of changes in glaciation, we use estimates of modern accumulation relative to ablation, assuming that regions with a less negative ice volume balance today are likely to have been more glaciated in the past. Accumulation is estimated from the NCEP/NCAR reanalysis of precipitation (Kalnay et al., 1996). Annual ablation is estimated from the reanalysis of temperatures at two meter elevation using a positive degree day approach, where daily average temperatures exceeding freezing are multiplied by 0.15 mm/(day°C) and summed over the year (e.g. Braithwaite, 1984). All ice melt is assumed to run off. Positive ice volume balances are partially suppressed in this estimate because mountainous regions-with their high elevation, substantial shading, and concentrated precipitation-are poorly resolved and because the modern climate is warm relative to the temperatures over the last 40,000 yr. Reassuringly, however, regions predicted to have the greatest modern ice balance tend to be currently glaciated (Cogely, 2003;Fig. 4) and to have been more glaciated in the past (Denton and Hughes, 1981;Grosswald, 1998;Dyke, 2004).  5 shows the independently estimated indices of ice volume balance and deglacial eruption frequency plotted against one another for each volcanic site. The d/g eruption frequency is approximately one in regions with insignificant deglaciation but stretches out to values of six for glaciated regions (Fig. 5a). For example, the relative frequency of volcanic events along the Western Pacific Rim increases by more than a factor of five between Indonesia and Kamchatka. The relationship between the d/h eruption frequencies and ice volume balance (Fig. 5b) is similar to the d/g relationship, albeit somewhat more diffuse and perhaps having two branches. These strong function-like relationships between independently estimated indices support the hypothesis that deglaciation triggers widespread volcanic activity.

A quantitative estimate of global rates of volcanism
The foregoing discussion documents increased volcanism in regions experiencing deglaciation, but quantifying the magnitude of this increase requires accounting for the strong temporal reporting bias. Fig. 5 indicates a uniform rate of volcanism in unglaciated regions, and this suggests a normalization method. We assume that regions not associated with deglaciation erupted with a constant frequency, c o . The value of c o is estimated from the average rate of eruptions in unglaciated regions over the last 2 ky. The observed frequency of eruptions, c′, is never greater than c o and the two can be related using an observational bias term, c′ t = c o b t . The bias term, b t , appears to follow a power-law (Siebert and Simkin, 2002), though our analysis does not depend on its exact form.
For glaciated regions we have x′ t = x t b t , where the true frequency, x t , is time dependent. Assuming that the sampling bias is consistent between both groups permits us to use the frequency of unglaciated eruptions as a control by which to estimate the glaciated frequency, x t = c o x′ t / c′ t . Global volcanic frequency is then f t = c o (1 + x′ t / c′ t ), which we present as fractional deviations from the most recent 2 ky, f t /f o = (1 + x′ t / c′ t )/(1 + x o / c o ). The key term in this equation is the ratio between the observed eruption rates, x′ t /c′ t , with the other terms acting to scale the result.
It is uncertain whether the observational bias is consistent between the glaciated and unglaciated groups. Whereas evidence of volcanic eruptions in glaciated regions may be disturbed or destroyed by the action of the ice itself, exposure in glaciated regions is often excellent, making identification easier. In addition, Western European  Fig. 2. Contours indicate (a) deglacial/glacial eruption frequency ratios and (b) deglacial/Holocene eruption frequency ratios. (b,d) are the averages across longitude. Temporal reporting bias causes the ratios in (a,b) to average higher than one, while those in (c,d) are lower than one. Nonetheless, both (a,b) and (c,d) show consistent increases in frequency at high latitudes, indicating that the spatial pattern is not an artifact of observational bias. Fig. 4. Ice volume balance for the modern climate system. Ice volume balance (contours) is calculated as the number of meters of precipitation per year minus the meters of ice ablation, both estimated from the NCEP/NCAR reanalysis of meteorological data (Kalnay et al., 1996) (see text for details). Also indicated are the locations of volcanoes (dots) and regions containing modern day glaciers (squares) (Cogely, 2003). The bold contours at −6 and − 9 m/yr indicate the two divisions used to distinguish glaciated and unglaciated regions in subsequent calculations. Note that the orography associated with volcanoes tends to lead to higher volume balances. and North American volcanic systems have also been more thoroughly studied, whereas the historical record in important areas such as Indonesia is poorly known. There is much room for further work in quantifying the global volcanic record.
We divide volcanoes into groups that did and did not experience deglaciation using the modern ice volume balance (see Fig. 5). Volcanoes with a strongly negative ice volume balance today are unlikely to have been glaciated during the late Pleistocene. Plausible cutoffs are −9 m/yr (giving an approximately equal frequency of events between the glaciated and unglaciated groups during the last 2 ky, x o /c o~1 ) and −6 m/yr (based on where eruption frequencies most markedly shift toward higher values). Regions included for each value can be seen in the contours of Fig. 4. The −9 m/yr cutoff includes the high peaks of the Mexican volcanic belt, the central and southern Andes, and northern New Zealand, all regions known to have been glaciated, while the -9 m/yr cutoff does not.
Either cutoff indicates that the ratio of glaciated to unglaciated eruption frequencies, x′ t /c′ t , is lowest during the last glacial, increases sharply near 12 ka, peaks near 7 ka, and subsides toward glacial levels in the recent past (Fig. 6). The −6 m/yr cutoff indicates a doubling of the global frequency of volcanic events between 12 and 7 ka relative to the last 2 ky, while − 9 m/yr gives a factor of six increase-we expect that the true value is between these bounds.
To test the statistical significance of the increase in global volcanism during the deglaciation (Fig. 6), we adopt the null-hypothesis that volcanic systems are independent of deglaciation. If the null-hypothesis holds, random reassignment of volcanoes to the glaciated (g′) or unglaciated (u′) classes should not influence the results. Of the 10 5 shuffling trials conducted, 99% result in global rates of volcanism limited to 0.5 and 1.5 times modern levels for all time periods, inconsistent with the doubling or more observed during the deglaciation. Further analysis in which the spatial covariance of the volume balance field is preserved during shuffling increases the width of the confidence interval, but the doubling is then still significant at the 95% level.
While volcanism begins to increase with the deglacial rise in sea level near 18 ka, it is intriguing that the large uptick in volcanism appears later, at 12 ka (Fig. 6b). We favor the explanation that deglaciation began with the eastern Laurentide and Antarctic ice sheets and, initially, did not much influence volcanic regions. Indeed, the highly active volcanic regions covered by the western Laurentide appear to have become more glaciated near 18 ka and not to have undergone significant deglaciation until 12 ka (Dyke, 2004). Likewise, the volcanic regions in Alaska (Yu et al., 2008) and Iceland (Maclennan et al., 2002;Licciardi et al., 2007) experienced the most pronounced deglaciation near 12 ka. Mauna Kea in Hawaii experienced deglaciation after 15 ka (Blard et al. 2007). South American deglaciation is less well constrained but appears to have proceeded in a series of steps between 17 and 11 ka (McCulloch et al., 2000). Thus, the discrepancy in timing between the initial rise in sea level and increased volcanism may reflect the differing regional histories of deglaciation.
Bryson and co-authors have used earlier versions of their dataset (Bryson et al., 2006) to estimate global rates of volcanism. They assume that the observational bias follows a power-law form (Bryson, 1989) or some other specified distribution (Bryson and Bryson, 1998) and use anomalies from these assumed distributions as indicators of changes in global volcanism. Apart from methodology and use of the updated database (Bryson et al., 2006), the present results also differ from Bryson (1989) and Bryson and Bryson (1998) in that they are supported by the geography of the increased volcanism, include many more observations over the last 12 ky (Siebert and Simkin, 2002), show a more distinct increase in volcanism during the deglaciation and subsequent decline during the Holocene, and show little variability during the last glacial. However, all reconstructions agree in showing anomalously large volcanism during the last deglaciation.
Our estimate of the time history of global volcanism is also broadly consistent with an independent, albeit more regional, estimate from the Greenland Ice Sheet Project II core (Zielinksi, 2000). The Greenland record of volcanism, which depends on excess SO 4 concentrations in the ice to identify volcanic events, indicates that the greatest frequency of volcanic events during the last 40ky occurred between 15 and 8 ka and that the largest eruptions occurred between 13 and 7 ka, consistent with our reconstruction (Fig. 6c). Greenland concentrations also reflect proximity of the volcanic events, atmospheric transport and deposition of sulphate, and changes in ice accumulation rates. Thus, the excess sulphate is not directly comparable to the eruption ratio that we calculate, though the correspondence between these independent reconstructions does support the reliability of both as indicators of volcanic activity.
Reconstructions of volcanic activity from Antarctic ice cores show conflicting results. A reconstruction of volcanism over the last 40 ky using volcanic SO 4 from Epica Dome C, Antarctica, indicates little change in volcanism (Castellano et al., 2004), whereas an estimate from Siple Dome, based on the optical properties of dust, indicates heightened volcanism near 10 ka (Bay et al., 2004). The Antarctic signal, or lack thereof, may reflect the fact that there are fewer subaerial volcanoes at high southern than at high northern latitudes. Siebert and Simkin (2002) identify some 435 volcanoes north of 40°N, and only 70 south of 40°S.

Physical mechanisms relating deglaciation and volcanism
The results presented above rely upon inferences drawn from existing data. Of course, they can be questioned with regard to the completeness of the volcanic record. In this section we discuss how the foregoing inferences are consistent with the physical consequences expected from deglaciation.

Depressurization and magma production
The fraction of melt in the Earth's mantle is sensitive to pressure changes. In particular, experimental data and analyses show that the amount of melt in the mantle, in regions above the solidus, increases by about 1% for each 100 MPa of pressure decrease (McKenzie, 1984;Langmuir et al., 1992). For example, unloading 1 km of ice leads to a 10 MPa depressurization and a 0.1% increase in the melt percentage, which equates to 60 m of melt production in a melt region 60 km thick. For convergent margins there is uncertainty in the thickness of the melting region and the distribution of melt between the mantle wedge and arc. The modeling results of Cagnioncle et al. (2007) suggest ã 60 km thick mantle melting region beneath arcs, a value that we adopt here. There are also questions regarding how the pressure exerted by surface loading will be distributed at depth. Whereas detailed modeling of the pressure field exerted by ice atop arcs would be useful, the large spatial extent associated with glaciers and ice cap systems during the Last Glacial Maximum leads to the possibility that ice loading influences much of the thickness of the underlying magmatic column. For example, during the Last Glacial Maximum, southern Alaska appears to have been covered beneath contiguous ice extending from the Alexander Archipelago in the East to the Aleutian Islands in the West, and extending northward from the Pacific hundreds to thousands of km inland (Manley and Kaufman, 2002); the Cascades of northwest North America lay beneath the Western lobe of the Laurentide ice sheet (Dyke, 2004); and southern South America was covered beneath extensive ice fields (Denton and Hughes, 1981).
The volume of mountain glaciers and small ice caps are estimated to have decreased from 1.9 million km 3 during the Last Glacial Maximum to 0.1 million km 3 today (Denton and Hughes, 1981). If 15% of this ice volume loss influences magmatic production (consistent with unloading ice of 100 m thickness from a 100 km swath along 25,000 km of convergent margin), we anticipate 16,000 km 3 of magma productionroughly equivalent to doubling global subaerial volcanism for 5000 yr. Or, if the equivalent of half the glacier ablation is involved, perhaps then also including the ablation of the Cordilleran ice sheet in western North America, magma production equates to 35,000 km 3 , or more than a five-fold increase in subaerial volcanism over 5000 yr. These crude estimates are consistent with the lower and upper bounds on deglacial volcanic activity inferred from the observed record of volcanic eruption.
The prediction of increased magmatic emplacement permits a test of our hypothesis. The deglacial magma production on Iceland has been estimated at 3000 km 3 (Jull and Mckenzie, 1996;Maclennan et al., 2002), or~10% of the total increase which we estimate, but we are more concerned with continental magmatism. Crisp (1984) estimates that the ratio of subsurface to surface emplacement is~10:1 for continental magmatism, and this leads us to expect the equivalent of 2 m of tephra being emplaced across a 100 km swath and along 25,000 km of convergent margin, where we have accounted for tephra being about half the density of the original melt. How broadly material is distributed will depend on the explosivity of individual eruptions and the amount of transport subsequent to emplacement.
Mount Mazama, in Northwest North America, is one of the few volcanoes sufficiently well mapped and dated to permit estimation of the volume of erupted material during deglaciation. Bacon and Lanphere (2006) identify 60 km 3 of material erupted at Mount Mazama during the last deglaciation, almost half of the total mapped at that site over the last 400 ky, corresponding to a roughly 60 m thickness over a 1000 km 2 region. In addition, Bacon and Lanphere (2006) suggest that several times this amount was emplaced beneath Mount Mazama. Although an isolated study, this work indicates that the magma products from one volcano contributed about 1% (4× 60 km 3 /30,000 km 3 ) of the total melt  Fig. 4 for a map indicating which volcanoes are in the glaciated and unglaciated groups. Note that the y-axis is logarithmic. (b) Estimates of the global frequency of volcanic eruptions using the − 6 m/yr (thin line) and − 9 m/yr (bold line) volcano groupings, where global frequency is normalized relative to the last 2 ky. The 99% confidence interval for the null-hypothesis of no systematic difference between glaciated and unglaciated events (dotted lines) indicates that the deglacial increase in eruptions is significant. (c) An eruption index based on volcanic SO 4 from a Greenland ice core (Zielinksi, 2000). production that we estimate. However, Hildreth et al. (2003) suggest that Mount Mazama had the largest eruptive volume in the Cascades. Further, Mount Mazama had its most recent caldera forming eruption at 5.9 ka, accounting for 28 km 3 of emplaced material, which is rather late in the deglaciation.
A lack of detailed mapping makes it unclear whether the volumes emplaced at Mount Mazama are isolated examples or globally indicative. The equivalent of one hundred volcanoes such as Mount Mazama would have had to be active during the last deglaciation to be consistent with our estimate of global magma output. More work on the magmatic emplacement of individual volcanoes spanning the last 40 ka is needed to test whether global volcanism substantially increased during the last deglaciation.

Other possible contributors to increased volcanism
The previous section considered only magma production in the mantle as a contribution to increased volcanism, but there are other effects associated with glaciation and deglaciation that may encourage eruption. Eruptability can be viewed as a balance between the forces generated by magmatic melt and gas production within the volcanic edifice on the one hand and the confining pressure and integrity of the enclosing rock on the other (e.g. Jellinik et al., 2004). Removing ice reduces the confining pressure and could induce an eruption. Ongoing glacial erosion throughout a glacial period may further reduce confining pressure as well as the structural integrity. Finally, far field effects, such as the unloading of the continents and rising sea level, may encourage volcanism by opening passageways or altering the pressure in magma chambers (Nakada and Yokose, 1992;McGuire et al., 1997).
The interesting possibility also arises that glacial/interglacial cycles come to pace the timing of certain volcanic eruptions (for a general discussion of volcanic pacing see Jupp et al., 2004). For pacing to occur, the volcanic systems should have an intrinsic time scale similar to that of the~100 ky glacial/interglacial variations (e.g. Strogatz, 1994). However, it is difficult to ascertain the number and size of volcanic systems exhibiting periods of more than about 1 ky because of the diminishing completeness of the volcanic record with age. Estimates do exist, however, for the scaling of global eruption frequency with the size of the eruptions (Siebert and Simkin, 2002;Mason et al., 2004). Size is generally categorized according to the Volcanic Explosive Index (VEI). A VEI 1 eruption involves 10 4 -10 6 m 3 of tephra and eruptions having a VEI of x N 1 involve 10 x + 4 -10 x + 5 m 3 of tephra. Mason et al. (2004) estimate that the global frequency of eruptions decreases with VEI almost as rapidly as 1/VEI for VEIs between 2 and 7, but that the total tephra mass flux increases more rapidly than VEI. For example, they estimate that VEI 2 eruptions occur every few years while VEI 7 eruptions occur once every thousand years, yet the net volume erupted by VEI 7 eruptions actually contributes more than three times that of all VEI 2 eruptions. This scaling suggests that a small number of large, low-frequency volcanic systems could substantially contribute to the total volcanic flux.
Large eruptions appear to have occurred with a quasi-100 ky period at Mount Mazama (Bacon and Lanphere, 2006), Western Europe (Nowell et al., 2006), and the South Eastern United Sates (Jellinik et al., 2004), indicating that at least some low-frequency volcanic systems are candidates for being paced by the glacial cycles. Increased eruptions during deglaciation, either directly forced by deglacial processes or paced to coincide with deglaciation, would lead to greater CO 2 emissions if they tap deep, pressurized magmatic reservoirs that still retain their CO 2 or if the eruption leads to depressurization of magmatic chambers and encourages the generation of new melt.

Submarine volcanism
Melting ice on land has the downstream consequence of raising sea level, and the glacial mechanisms explored for increasing subaerial volcanic eruption have complementary implications for decreasing submarine volcanism through pressure increases caused by rising sea level. Glaciation and deglaciation are expected to lead to anticorrelated variations in subaerial and submarine volcanism, and it is necessary to assess the relative magnitudes of these effects.
Virtually no historical data exists for submarine eruptions, so we are restricted to theoretical estimates based on expected melt production beneath ridges. Given the factor of three difference in density between water and mantle peridotite, the~130 m deglacial rise in sea level is equivalent to suppressing δz = 45 m of mantle upwelling. This effect is independent of spreading rate, but the fractional change in upwelling depends on the geometry of the melting region and the spreading rate. We approximate the geometry of the melt region as an isosceles triangle with its vertex at the spreading center and with its base oriented perpendicular to the spreading axis and parallel to the sea floor. If the area of the triangular melt region is conserved, the mean upwelling rate at the base, u, can be related to the full spreading rate of the ridge, s, as u = (s/2) tan θ, where θ is the angle between the sea floor and the legs of the triangular melt region.
For this triangular geometry, the fractional reduction in melt production at ridges, δz/z, is a function of θ and s, where z is the thickness of mantle upwelled at a ridge over 10 ky of deglaciation. Assuming θ of 45°, the rapidly spreading southern East Pacific Rise, which has an s of about 200 m/ky, upwells 1 km of material during the deglaciation, yielding a small fractional reduction in melt production, δz/z = 45 m/1000 m~5%. The slower spreading mid-Atlantic ridge (s~30 m/ky ) yields a 30% reduction. For very slow spreading ridges, such as the Gakkel Ridge and Southwest Indian Ridge (s b 15 m/ky), complete cessation of magma production could occur. Note that rapid rises in sea level could lead to short, global cessations of ridge melt production. A global average spreading rate of 50 m/ky is obtained by dividing the~3400 km 2 /ky of global mean production of oceanic crust (Rowley, 2002) by 65,000 km of ridge length. This spreading rate implies a~20% reduction in submarine magma production over thẽ 10 ky course of deglaciation. Lund and Asimow (2008) recently presented a similar scenario using a more complete melt model, and considered the implications of changes in melt production on the geochemistry of seawater in the East Pacific.
Pacing of eruption by glacial changes in sea level appears unlikely to operate at ridges because eruptions tend to occur at years to decades at fast-spreading ridges (Perfit et al., 1994) and are not expected to reach the~100 ky time scales associated with glacial/interglacial changes in sea level, even at slow spreading regions. Ocean ridge volcanism is characterized by a large number of small flows (Perfit et al., 1994), rather than the small numbers of large eruptions that dominate subaerial volumes, a difference consistent with the thin lithosphere characteristic of the ridge environment.
If the glacial/interglacial variations in sea level do lead to modulations in the production of melt at mid-ocean ridges, a signature of these variations might be expected in the bathymetry surrounding spreading centers. Quasi-periodic variability in the sea floor elevation on length scales of kilometers have been documented, and given spreading rates on the order of 1-10 cm/yr, the variability is plausibly consistent with the time scale of the glacial cycles. However, our first order analysis of high-resolution bathymetry at the fast-spreading East Pacific Rise (Carbotte et al., 2004) did not reveal significant correlation between bathymetry and estimates of past changes in sea level. We nonetheless suspect that analysis of ridge bathymetry will eventually reveal some imprint of glacial-interglacial variations in sea level, once sufficient data are available. Diminished eruptive activity during deglaciation may also be detected by dating of young lavas from the seafloor, should this become more viable.

Implications for the carbon cycle
The data and some theory indicate that deglaciation drives a widespread increase in subaerial volcanism and decline in submarine volcanism. While both of these are important effects on global volcanism, the subaerial effects for CO 2 are substantially larger, as we discuss below. We postulate that elevated volcanism during deglaciation contributes to the rise in atmospheric CO 2 during deglaciation, with the ensuing warming constituting a positive feedback upon the deglaciation. Conversely, waning volcanic activity during the Holocene would contribute to cooling and reglaciation, thus tending to suppress volcanic activity and promote the onset of an ice age. This hypothesis depends on the amount of CO 2 emitted from volcanoes, and how it is partitioned between the atmosphere and other carbon reservoirs. Thus, in principle, such a calculation depends upon nearly all aspects of the climate system, including parts of the solid earth. Here we seek first order estimates.

CO 2 emissions from subaerial volcanoes
First, we discuss the modern subaerial volcanic CO 2 flux. The most direct approach is to rely on data from currently active volcanoes (Williams et al., 1992) or estimates derived from CO 2 / 3 He (Sano and Williams, 1996;Marty and Tolstikhin, 1998;Hilton et al., 2002). These estimates cluster near 0.1 Gt CO 2 /yr. However, a recent simulation of arc volcanism combined with observational studies (Gorman et al., 2006) suggests that while the range of emissions across these studies is plausible, the upper end of the range,~0.14 Gt CO 2 /yr, is most likely. These estimates have the potential short-coming that the very short historical snapshot of CO 2 emissions may not be representative of longer periods.
An alternative approach is to multiply arc magma production rates by their average primary CO 2 concentration. Reymer and Schubert (1984) estimate long term crustal production to be 20-40 km 3 per km of arc length per Ma, but this estimate has been criticized as too low by a factor of two (Dimilanta et al., 2002), and any such estimates, which are based on the amount of crust that persists through time, are minima with respect to magma additions because they are the net of production after losses due to erosion. A value of 80 km 3 /km/Ma and 35,000 km total arc length gives a magma production rate of approximately 3 km 3 /yr, in accord with other estimates (Dimilanta et al., 2002).
Primary CO 2 concentrations cannot be determined directly because CO 2 is almost entirely degassed prior to eruption. Instead, we adopt an indirect approach using element ratios. Carbon isotope data and CO 2 / 3 He ratios both indicate that the mantle contributes only 10-20% of the total CO 2 at arc volcanoes (Fischer et al., 1998;Marty and Tolstikhin, 1998). The concentration of carbon in the mantle can be estimated by multiplying an average CO 2 /Nb ratio of~500 (Saal et al., 2002;Cartigny et al., 2008) by an average Nb content of~3 ppm in arc basalts, yielding a 0.15% CO 2 mantle contribution in primary magmas. Since the mantle contribution is only 10-20% of the total CO 2 , primary CO 2 concentrations would be 0.75-1.5%. Multiplying this concentration by the magma production rate given above, assuming a density of 3 Gt/km 3 , yields an average global emission rate of 0.1 Gt/yr for convergent margins. Note that multiplication of the mantle CO 2 /Nb ratio by the Nb content in arc basalts depends on the assumption that Nb behaves consistently during mantle melting at arcs, ridges, and islands. This is likely, given that Nb is not readily transferred from the down-going slab. Nonetheless, were Nb held in a residual phase in the mantle wedge during the melting that forms arc magmas, our estimates of mantle Nb concentration and mantle-derived CO 2 would be biased low.
It is harder to parse emissions from non-convergent margin subaerial volcanoes, but they likely add another 0.05 Gt/yr (Marty and Tolstikhin, 1998;Hilton et al., 2002), leading to total subaerial emissions of 0.15 Gt/yr. All of these estimates are thus consistent with subaerial volcanic emissions between 0.1 and 0.15 Gt CO 2 /year. Submarine volcanic emissions are also estimated to be about 0.1 Gt CO 2 /year (Marty and Tolstikhin, 1998;Saal et al., 2002;Cartigny et al., 2008), consistent with CO 2 concentrations being ten times less and eruption volumes being eight times greater than at arcs. To estimate the time history of CO 2 fluxes, we multiply current subaerial volcanic emissions by the ratio between past and present eruption frequencies, assuming proportional changes between the frequency of volcanic events and CO 2 emissions. To take into account uncertainties in the eruption response to deglaciation, the time history of volcanic CO 2 emissions is estimated using random draws from a uniform distribution of ice volume balance cutoffs (bounded between −6 m/yr and − 9 m/yr) and modern subaerial fluxes of CO 2 (bounded between 0.1 and 0.15 Gt of CO 2 per year). Repeating this procedure many times provides an ensemble of plausible time histories of volcanic CO 2 emissions. The ensemble average indicates that volcanoes emitted 3000 Gt of CO 2 during the last deglaciation above a baseline scenario of current emissions, and 90% of all time histories fall between 1000 and 5000 Gt of CO 2 emissions. These are large numbers. By way of comparison, 3000 Gt of volcanic CO 2 emissions corresponds to roughly a century of anthropogenic emissions, at current rates.
It is also necessary to consider the reduction in CO 2 emissions from ocean ridge volcanism. As discussed above, we estimate that global submarine volcanism decreased by about 20% for 10 ky. That leads to a reduction of 200 Gt CO 2 , or an order of magnitude lower than the increased emissions expected from subaerial volcanism.
Before proceeding to the implications of these emissions, it is important to note that the relationship between increased volcanism and increased CO 2 is likely to be complex and is, in any case, poorly constrained. Depressurization associated with deglaciation is expected to increase the melt fraction. However, if all the CO 2 in the source is already in the silicate melt, then the amount of CO 2 brought to the surface may not depend on changes in the amount of melt at all, since CO 2 would simply be diluted by the additional volume of magma. On the other hand, if deglaciation increases melt production at depth or increases the volume of the source region, we expect greater CO 2 emissions.
It is also unclear when and how CO 2 in crustal magma reservoirs reaches the atmosphere. CO 2 solubility is low enough that much of it is released at depth, and this gas might escape at times and locations not directly associated with eruption. There is evidence that intrusive emplacement of magma accounts for high CO 2 fluxes, for example, at Mount Etna (Allard et al., 1994), Mount Vesuvius (Iacono-Marziano et al., 2009), and Yellowstone (Werner and Brantley, 2003), where the CO 2 is derived from the cooling magma as well as interaction between the magma and sedimentary sources of CO 2 . Although not accounted for here, it also seems an interesting further possibility that the Laurentide ice sheet formed an impermeable layer atop passive emission regions, perhaps in Yellowstone or the Cordillera region, that inhibited CO 2 fluxes into the atmosphere during the glacial.

Ocean carbonate compensation
Having an estimate of the volcanic contribution of CO 2 during the deglaciation, we now turn to its partitioning between atmosphere and ocean. A 1000 to 5000 Gt CO 2 release from subaerial volcanoes might be expected to increase ocean acidity and lead to a shoaling of the oceanic carbonate saturation horizon. For example, injection of~3000 Gt of CO 2 into the ocean, accompanied by a 4°C ocean warming (Schrag et al., 1996) and 100 ppm increase in atmospheric CO 2 concentration coming out of the last glacial, would cause the carbonate saturation horizon to shoal by about 1 km (Fig. 7). (Note that we compute the change in the saturation horizon accounting for the influence of pressure upon dissolution but not the more uncertain and less important vertical gradients in temperature or salinity.) Such a shoaling is consistent with observations of carbonate dissolution in the Pacific (Farrell and Prell, 1989;Thunell et al., 1992) but not the Atlantic (Thunell et al., 1992).
The above approach, however, ignores the increased terrestrial carbon storage expected to accompany warming climates and the retreat of ice on land. Indeed, organic land storage likely accounts for the~0.3‰ increase in ocean δ 13 C observed between the glacial and Holocene (Curry et al., 1988). Normally the change in δ 13 C is postulated to result from a biospheric uptake of~1500 Gt of CO 2 in isotopically light organic matter, but an uptake of~2000 Gt would then also compensate for a 3000 Gt of volcanic CO 2 emissions having an isotopic ratio of −2.8 ± 1.2‰ (Sano and Marty, 1995;de Leeuw et al., 2007). When biospheric uptake is taken together with mean ocean warming, the expected change in the carbonate saturation horizon is indistinguishable from zero (Fig. 7), particularly given further uncertainties associated with the carbonate system such as coral reef building (Vecsei and Berger, 2004). A small change in saturation is consistent with estimates of carbonate ion concentration during the last glacial (Broecker and Clark, 2001;Anderson and Archer, 2002), which suggest changes in the distribution of water masses but little change in overall concentration. We also note that the additional flux of CO 2 from volcanoes is consistent with inferences that no one oceanic mechanism is capable of explaining the glacial/interglacial changes in atmospheric CO 2 (e.g. Kohfeld et al., 2005;Marchitto et al., 2005).
The above considerations lead us to adopt a first order approach that assumes changes in carbonate compensation play no significant role in determination of the atmospheric CO 2 concentration coming out of the last glacial.

A simple time variable model of atmospheric CO 2
A simple atmosphere/ocean box model, similar to that of Kheshgi (2004), is adopted to represent the time variable volcanic influence upon atmospheric CO 2 , Here C a and C o are the amounts of inorganic carbon in the atmosphere and ocean, measured in Gt of CO 2 . The atmosphere-ocean flux is F = (C′ a − C′ o (1 − q) / q) / τ, where the primes indicate anomalies away from equilibrium. The fraction of volcanic carbon remaining in the atmosphere once in equilibrium with the ocean, q, is taken to be between 10% and 15% (Montenegro et al., 2007). Estimates of the equilibration time scale, τ, range from~300 yr (Archer, 2005) tõ 1800 yr (Montenegro et al., 2007), or longer (Wunsch and Heimbach, 2008), and we assign wide bounds on τ of 300-2000 yr.
V is the volcanic flux of CO 2 into the atmosphere, as estimated in the foregoing section. This volcanic flux of carbon is assumed to be in balance with uptake by a constant silicate weathering, W o , over the course of a 100 ky glacial cycle. The average CO 2 emissions between 100 and 40 ka, where we have no eruption data, are assumed to equal the rates between 40 and 20 ka, where we do have data.
This model is simplistic because, among other reasons, it does not account for the uptake of carbon by the biosphere, nor does it account for increases in atmospheric CO 2 resulting from ocean warming or other mechanisms (e.g. Toggweiler and Russell, 2008). However, similar results are obtained using our volcanic emission estimates in a more sophisticated carbon box model (F. Joos, personal communication) that includes a representation of ocean carbonate compensation and the biosphere (Joos et al., 1996(Joos et al., , 2004, indicating that our model gives plausible results. To explore the range of atmospheric CO 2 scenarios, parameters are drawn from uniform distributions between the previously discussed bounds: −9 to −6 m/yr for the ice volume cutoff used to distinguish glacial and non-glacial volcanoes, 0.1-0.15 Gt CO 2 /yr for the modern volcanic flux, 300-2000 yr for the ocean equilibration time, and 10-15% for the airborne CO 2 fraction. Model integrations are initialized at 40 ka with the atmosphere and ocean in equilibrium. We take the mean of an ensemble of 10 4 runs as the best estimate and use the range that covers 90% of all model runs as an indicator of the uncertainty (Fig. 8).

Atmospheric CO 2 concentrations
The model estimate of the time history of atmospheric CO 2 can be compared with atmospheric CO 2 observations obtained from the Dome C (Monnin et al., 2001) and Taylor Dome (Indermühle et al., 2000) Antarctic ice cores. We consider four distinct intervals (Fig. 8), 1. Glacial, 40-18 ka: The model results indicate atmospheric CO 2 decreases by 10 ppm (5 to 20 ppm, 90% confidence interval), marginally consistent with the observed 20 ppm decrease, suggesting that the trend toward lower atmospheric CO 2 levels during glaciation is partly attributable to excess weathering relative to volcanic emissions. 2. First half of the deglaciation, 18-13 ka: CO 2 increases modestly bỹ 10 ppm (5 to 40 ppm, 90% c.i.), whereas the observations show a 50 ppm rise, highlighting that factors independent of volcanism exert influence on glacial-interglacial variations in CO 2 (e.g. Broecker and Peng, 1982). Indeed, changes in ocean circulation, continental biomass, methane release in deglaciating regions, weathering, and marine biology including coral reefs all likely influence CO 2 changes, but to uncertain degrees. 3. Second half of the deglaciation, 13-7 ka: The reconstruction indicates a 40 ppm (15 to 70 ppm, 90% c.i.) increase in volcanic CO 2 , consistent with observation, particularly with respect to the sharp uptick near 12 ka. 4. Late Holocene, 7-0 ka: Volcanic CO 2 contributions wane owing to lower volcanic activity and on-going equilibration with the oceans and weathering, while observations instead indicate rising CO 2 levels during this interval. It appears that this divergence between the modeled volcanic CO 2 and observations is peculiar to the Holocene, as opposed to prior interglacials, a point we will take up elsewhere.
Excess volcanic emission of CO 2 coming out of the last glacial is estimated to account for roughly half the deglacial increase in atmospheric CO 2 . This contribution can be roughly tallied along with other sources and sinks of CO 2 . An ocean warming of 4°C would give another~40 ppm increase in atmospheric CO 2 concentration during the deglaciation. Marine δ 13 C values in combination with the expected volcanic influence on δ 13 C suggest that land carbon storage increased by 2000 Gt CO 2 (Section 4.2), giving a~30 ppm drawdown in atmospheric CO 2 . Tallying these sources and sinks, and assuming temperature and land influences are evenly distributed over the deglaciation, suggests an Fig. 7. Modeled change in the carbonate saturation horizon between the last glacial and present interglacial for various volcanic inputs of CO 2 into the ocean-atmosphere system. One scenario considers only changes in CO 2 and temperature (dashed line), and the other also includes organic storage of carbon on the continents (solid line). Our estimate of 3000Gt of volcanic CO 2 emissions requires little change in the carbonate saturation horizon.
atmospheric CO 2 budget of 10 + 20 − 15 = 15 ppm for the first half of the deglaciation and 40 + 20 − 15 = 45 ppm for the second half from volcanoes, ocean warming, and land carbon uptake, respectively. Increased ventilation of the deep ocean, or other contributions, can then be invoked to contribute another~35 ppm during the first half of the deglaciation in order to account for the full deglacial increase. While there are significant remaining uncertainties associated with volcanic emissions of CO 2 , as well as with other sources, it appears that coupling between the deep earth and atmosphere/ocean system plays a potentially important role in the determination of glacial/interglacial changes in atmospheric CO 2 concentrations.

Carbon isotopes
Volcanic emissions of CO 2 are radioactively inert, and thus have implications for atmospheric radiocarbon activity. The average rate of decline in atmospheric Δ 14 C over the last 40 ka is − 15‰/ka, as indicated by the estimates of Hughen et al. (2004). This rate of decrease appears more than twice as fast (− 33‰/ka) during the last deglaciation, 18-7 ka, and Broecker and Barker (2007) have called particular attention to the rapid decline between 17.5 and 14.5 ka (−70‰/ka), the source of which remains unknown. Our estimates place the bulk of the volcanic emissions later in the deglaciation, overlapping with a yet more rapid drop in atmospheric Δ 14 C between 12.5 and 10 ka, − 90‰/ka, which is intriguing.
However, in the idealized circumstance of 3000 Gt of volcanic CO 2 and instantaneous mixing between the atmosphere and ocean, the atmospheric radiocarbon activity is expected to decrease by onlỹ 20‰. While the transient atmosphere and surface ocean radiocarbon anomalies may be larger, it appears that these levels of volcanic emissions are a minor influence upon the time history of atmospheric radiocarbon values. Similarly, the expected change in atmospheric δ 13 C from volcanic emissions appears small relative to the observed variability (Smith et al., 1999). Note that the fractional increase in atmospheric CO 2 concentration is about ten times that of the fractional decrease in radiocarbon concentration primarily because, unlike for carbon isotopes, ocean-atmosphere equilibration of the partial pressure of CO 2 is strongly influenced by changes in ocean pH (as described by the Revelle factor, see e.g. Broecker and Peng, 1982).

Further discussion
The foregoing discussion has dealt with volcanic eruptions and CO 2 emissions since the last glacial. Here we consider the glacio-volcano-CO 2 hypothesis in a somewhat larger context by briefly discussing three topics: the aerosol cooling effect associated with explosive eruptions, the CO 2 flux we estimate relative to Antarctic and Southern Ocean climate indicators, and the interactions between volcanism and glaciation over the Plio-Pleistocene and earlier periods of the Earth's history.

Warming, cooling, or both?
In addition to increasing atmospheric CO 2 , volcanism also increases atmospheric aerosol loading and thereby cools the climate (e.g. Rampino et al., 1988). Mount Pinatubo injected about 17 Mt of SO 2 into the atmosphere in 1991 and had a peak radiative cooling effect of 4 W/m 2 at the surface, causing surface temperatures to cool by about 0.5°C (Hansen et al., 1992), which diminished with an e-folding time scale of approximately one year. By way of comparison, we estimate volcanism contributes~40 ppm to the deglacial rise in atmospheric CO 2 , leading to a~1 W/m 2 increase in radiative forcing. In this rough view, volcanic CO 2 forcing would offset a Pinatubo-like eruption every four years. Large eruptions such as Pinatubo occur about once a century during the late Holocene, so that even if global volcanism increased by a factor of five during the deglaciation, we still only expect a Pinatubo size eruption every twenty years, suggesting that the aerosol cooling is smaller than the CO 2 warming. Furthermore, Pinatubo caused a strong cooling, in part, because it injected aerosols into the stratosphere near the equator-which then reflect more sunlight and remain aloft longer-whereas deglaciation is expected to drive primarily high latitude volcanism.
The short-term cooling associated with volcanism may, however, also be important. In particular, it is tempting to speculate that the spike in volcanism near 12 ka may have contributed to the resumption of glacial-like conditions in and around the North Atlantic during the Younger Dryas. The competing influences of volcanic CO 2 and aerosol emissions may be likened to the case of the tortoise and the hare, with CO 2 steadily warming the climate and aerosols driving relatively brief intervals of cool conditions.  (Fairbanks and Peltier, 2006). (b) Number of volcanic events per ky for glaciated (solid lines) and unglaciated (dashed lines) volcanoes, where the y-axis is logarithmic. (c) Estimated global frequency of volcanic events (solid line), normalized so that the frequency during the last 2 ky equals one. Also shown is the 99% coverage interval for the null-hypothesis of no systematic difference between glaciated and unglaciated events (dashed lines). (d) The contribution to atmospheric CO 2 from volcanic activity. The solid lines in (b), (c), and (d) represent the averages across 10,000 different estimates made using random draws from the range of parameter values described in the text, while the gray shading in (c) and (d) indicate the interval within which 90% of those estimates fall. (e) CO 2 concentrations from the Dome C (Monnin et al., 2001) and Taylor Dome (Indermühle et al., 2000) Antarctic ice cores, placed on a consistent timescale (Monnin et al., 2004) (dots), and a smoothed version using a 2 ky window (solid line). Also shown is the residual atmospheric CO 2 after subtracting the volcanic CO 2 contribution (dash-dot line). The vertical shaded bar is between 12 and 7 ka, when volcanic frequency appears greatest.

Relationship with changes in polar climates
The correspondence between CO 2 and Antarctic temperature through the last deglaciation has been interpreted as evidence for Southern Ocean control over atmospheric CO 2 concentrations (Monnin et al., 2001) and might be construed to rule out a volcanic contribution. However, it may be that Antarctic air temperature-as well as southern sea surface temperature (Barrows et al., 2007)respond to atmospheric CO 2 concentrations without regard to the source of that CO 2 . Such a view point is supported by a model simulation that reproduces Antarctic warming through the last deglaciation by prescribing atmospheric CO 2 values, along with smaller contributions from changes in ice sheet elevation and Earth's orbital configuration (Timmermann et al., 2009), without specification of the CO 2 source. Perhaps Southern Hemisphere climate remains nearly in equilibrium with atmospheric CO 2 concentrations through the last deglaciation, unlike northern climate, because there were no large changes in the distribution of southern continental ice volume (e.g. Huybers and Denton, 2008).
If our estimates are correct that volcanism accounts for much of the rise in atmospheric CO 2 during the latter half of the deglaciation, there still remains the outstanding question as to what drove the initial rise in CO 2 . Anderson et al. (2009) recently showed that the concentration of opal increased at three separate sites in the Southern Ocean during the last deglaciation, which they interpret as indicating increased upwelling of silicate. Upwelling is then also expected to increase out-gassing of CO 2 into the atmosphere. While one of the records presented by Anderson et al. (2009) indicates maxima corresponding with the two major rises in CO 2 during the deglaciation, all three indicate that the strongest upwelling occurred during the first half of the deglaciation, raising the possibility that an initial contribution of CO 2 from the Southern Ocean helped drive deglaciation and, thereby, volcanic CO 2 contributions during the second half of the deglaciation.
Another issue is why the Bölling/Alleröd warming documented at Greenland and many other northern sites is not associated with an increase in volcanism. The answer appears to be that the Bölling/ Alleröd warming was primarily a winter warming associated with decreased winter sea ice (e.g. Denton et al., 2005). Warmer winters would still be too cold to melt ice and may actually serve to increase moisture availability and cause ice caps and glaciers to grow.

Volcanism prior to the last glacial
The feedback between glaciation and volcanism suggests that deglaciation would not only pace eruptions but possibly also be paced by eruptions. That is, the conditions exist for volcanic eruptions and glacial cycles to mutually influence one another's timing so as to become synchronized (e.g. Strogatz, 1994). It may be that the progression of Pleistocene climate towards larger and more asymmetric glacial cycles (e.g. Huybers, 2007) can, in part, be understood as the synchronization and attendant amplification of the feedback between volcanic systems and glacial variability.
It is also useful to consider the general conditions that would give rise to the glacio-volcano-CO 2 feedback outlined here. One necessary condition is for ice and volcanoes to be in proximity (e.g. Fig. 4). At a basic level, a volcano's orography tends to promote precipitation and its elevation helps to retain that precipitation as ice. Furthermore, the current plate configuration may be peculiarly conducive to generating glaciated volcanoes, as it places many volcanoes at high latitudes and along the western margins of continents, which are then well situated to capture precipitation from moisture laden westerlies and to retain it as ice. A second condition for evoking the glacial-volcano-CO 2 feedback is sensitivity of glacier mass to changes in atmospheric CO 2 . Thus, for example, limited ice volume during the Paleocene and Eocene would limit the feedback.
Turning to a contrasting environment, the massive ice unloading postulated to occur at the termination of a snowball Earth episode would presumably lead to a dramatic increase in volcanism (Hoffman et al., 1998). However, the high atmospheric CO 2 conditions thought to accompany deglaciation of a snowball would serve to minimize the influence of additional volcanic CO 2 because radiative forcing scales nearly logarithmically with CO 2 concentration. (Volcanism may, however, play a role in the termination of a snowball through deposition of tephra leading to a decrease in ice albedo.) It thus seems that conditions during the Pleistocene-wherein the Earth has been precariously poised between glaciated and unglaciated states, atmospheric CO 2 concentrations have been modest, and the plate configuration places volcanoes in cold and wet climates-make this epoch particularly well-suited to evoke the glacial-volcano-CO 2 feedback.

Conclusions
A global reconstruction indicates that volcanism increased two to six times above background levels during the last deglaciation, with the signal driven by those regions that underwent deglaciation. This result is consistent with regional studies (e.g. Gardeweg et al., 1998;Jellinik et al., 2004;Bacon and Lanphere, 2006;Nowell et al., 2006;Licciardi et al., 2007), earlier inferences drawn from the Bryson and Bryson (1998) database of volcanic eruptions, and the volcanic index reconstructed from the Greenland ice core (Zielinksi et al., 1997). Depressurization associated with ablation of glaciers and ice caps could cause decompressional melting of the mantle sufficient to sustain a factor of five increase in volcanic output for 5 ky. In addition, glacial variability may pace the timing of low-frequency eruptions so as to coincide with deglaciation. Similar logic suggests ocean ridge production should decrease with rising sea level during deglaciation, although direct evidence from the ridge environment is lacking.
The most firm conclusion from this work is that subaerial volcanoes respond to deglaciation, but there may also be a feedback. We suggest that deglacially induced anomalies in volcanic activity cause imbalances in the atmospheric carbon budget which accumulate through deglaciations and persist into interglacials. Factor of two to six increases in the rate of volcanic emissions, persisting for thousands of years, are estimated to increase atmospheric CO 2 concentrations by 20-80 ppm, with the majority of the emissions occurring during the latter half of the deglaciation. While multiple other mechanisms can cause glacial/interglacial CO 2 variability, this increase in volcanism is expected from the effects of deglacial unloading and coincides with the observed secondary deglacial rise of atmospheric CO 2 . Thus, the deglacial rise in atmospheric CO 2 can, in part, be understood as a feedback induced by the deglaciation itself and mediated by volcanic activity. By similar logic, the glacial drawdown in CO 2 may partly owe to a deficit in volcanic emissions relative to CO 2 drawdown by weathering and other processes.
This analysis suggests that the Earth system is deeply coupled. As long as the climate and continental configuration engender colocation of volcanoes and ice, we expect interactions between the Earth's interior, surface, and atmosphere to amplify and modify the cycling between glacial and interglacial climates.