Abstract
Arctic neoglaciation following the Holocene Thermal Maximum is an important feature of late-Holocene climate. We investigated this phenomenon using a transient 6000-year simulation with the CESM-CAM5 climate model driven by orbital forcing, greenhouse gas concentrations, and a land use reconstruction. During the first three millennia analyzed here (6–3 ka), mean Arctic snow depth increases, despite enhanced greenhouse forcing. Superimposed on this secular trend is a very abrupt increase in snow depth between 5 and 4.9 ka on Ellesmere Island and the Greenland coasts, in rough agreement with the timing of observed neoglaciation in the region. This transition is especially extreme on Ellesmere Island, where end-of-summer snow coverage jumps from nearly 0 to virtually 100% in 1 year, and snow depth increases to the model’s imposed maximum within 15 years. This climatic shift involves more than the Milankovitch-based expectation of cooler summers causing less snow melt. Coincident with the onset of the cold regime are two consecutive summers with heavy snowfall on Ellesmere Island that help to short-circuit the normal seasonal melt cycle. These heavy snow seasons are caused by synoptic-scale, cyclonic circulation anomalies over the Arctic Ocean and Canadian Archipelago, including an extremely positive phase of the Arctic Oscillation. Our study reveals that a climate model can produce sudden climatic transitions in this region prone to glacial inception and exceptional variability, due to a dynamic mechanism (more summer snowfall induced by an extreme circulation anomaly) that augments the traditional Milankovitch thermodynamic explanation of orbitally induced glacier development.
Introduction
Due to greater boreal summer insolation, Arctic climate during the early-mid-Holocene was typically warmer with contracted glacial extent compared with the late-Holocene (e.g. Marcott et al., 2013; Wanner et al., 2008). A recent study by McKay et al. (2018) summarized these changes and concluded that glaciers were ‘uniformly retreated or absent throughout the Arctic’ from 8 ka until 6.5 ka, then expanded to near-modern conditions by 4.5–2 ka (all dates in this paper refer to years before AD 1950). This period of formation and growth of glaciers is known as the Holocene ‘Neoglaciation’ (Porter, 2007). Other studies have confirmed the approximate timing of this transition from relatively warm to cold conditions in the Canadian Arctic. The Ellesmere Ice Shelf began forming shortly after 5.5 ka (England et al., 2008), and the Agassiz Ice Cap on Ellesmere Island indicates a gradual long-term cooling trend since its early Holocene thermal peak (Lecavalier et al., 2017). Evidence of neoglaciation on Baffin Island around 4.5–5 ka has been identified in several studies (Miller et al., 2005, 2013; Pendleton et al., 2019), and glacial responses to regional climate controls were ‘near-instantaneous’ and ‘synchronous throughout eastern Canadian Arctic and possibly eastern Greenland’ (Margreth et al., 2014). However, McKay et al. (2018) do not find support for a synchronous Arctic-wide neoglaciation, which suggests that the expansion of high-latitude glaciers following the Holocene Thermal Maximum resulted from multiple regional cooling episodes and varying regional ice-growth thresholds.
The Holocene Neoglaciation is a recent example of a large Arctic cooling anomaly sufficient to produce widespread glacier expansion, but the Canadian Arctic is known to be a region susceptible to glacial inception generally and was the presumed source of the Laurentide Ice Sheet (Andrews and Mahaffy, 1976; Dong and Valdes, 1995; Williams, 1978). Furthermore, this region and Greenland have experienced very rapid climatic fluctuations on timescales of decades or even shorter (Alley et al., 1993; Bourgeois et al., 2000; D’Andrea et al., 2011). A variety of physical processes can promote glacial inception, which is ultimately driven by waning summertime insolation (Milankovitch, 1941). Amplifying feedbacks include reduced ocean heat transport (Calov et al., 2005; Khodri et al., 2001), expanded sea ice cover (Kageyama et al., 2004; Vavrus, 1995), and vegetation feedbacks (Gallimore and Kutzbach, 1996; Yoshimori et al., 2002). The basis of most explanations for glacial inception is that reduced summertime snow melt, rather than enhanced snowfall, is the key factor (Jochum et al., 2012; Otieno and Bromwich, 2009).
Our study revisits this traditional perspective by describing the pivotal role played by increased (summer) snowfall and the associated short-term regional weather patterns that promoted an inferred and abrupt glaciation in the Canadian Arctic in a transient global climate model simulation of the late-Holocene (past 6000 years). Climate models have struggled to replicate instances of glacial inception in the paleo record (e.g. Dong and Valdes, 1995; Otieno and Bromwich, 2009; Rind et al., 1989), due to their numerous simplifications such as muted topography (Birch et al., 2017; Vavrus et al., 2010) and inadequate representation of the positive physical feedbacks identified above. Here, we investigate the simulated rapid buildup of snow cover on Ellesmere Island and attribute its timing and abruptness to two consecutive summers around 4.9 ka featuring extremely favorable atmospheric circulation patterns for generating cold weather and exceptional snowfall. We suggest that this mechanism may play an important but unrecognized role in glacial inceptions.
Methods
We ran version 1.2.2.1 of the Community Earth System Model (CESM), a fully coupled global climate model whose atmospheric component is the Community Atmospheric Model version 5.3 (CAM5; Hurrell et al., 2013; Neale et al., 2012). The horizontal resolution of the atmosphere and land components used here is 1.875° latitude × 2.5° longitude, which is coarser than the standard CESM-CAM5, and the atmosphere consists of 30 levels. The sea ice component and 60-level ocean module use a variable grid spacing that averages 0.47° meridionally and is a uniform 1.125° zonally. CESM-CAM5 has been widely used in modeling studies of past, present, and future climates (Hezel et al., 2014; Kay et al., 2015; Otto-Bliesner et al., 2015; Xiao et al., 2014). It has one of the most realistic simulations of present-day global climate and Arctic climate/sea ice cover (Jahn et al., 2016; Knutti et al., 2013). In addition, this model’s predecessor versions (CCSM3, CCSM4) have been used in a number of studies of glacial inception (Colleoni et al., 2012; Jochum et al., 2012; Vavrus et al., 2010, 2018).
Our simulation of the late-Holocene climate starts at year 6000 BP and ends at year 1850 CE, the commonly used termination point of the pre-industrial era (Gent et al., 2011; Otto-Bliesner et al., 2015). In this paper we describe the first half of the simulation (6–3 ka), within which is the abrupt buildup of snow cover on Ellesmere Island. The boundary conditions for our transient experiment are the well-established orbital parameters (eccentricity, obliquity, and precession) from Berger (1978); greenhouse gas (GHG) concentrations of CO2, CH4, and N2O measured in Greenland ice cores (Loulergue et al., 2008; MacFarling Meure et al., 2006; Monnin et al., 2001) and the land cover reconstruction of Kaplan et al. (2010). Short-term variations in solar irradiance and volcanic forcing were not included. The model was spun up from the control simulation of the Last Millennium Ensemble (Otto-Bliesner et al., 2015) and was integrated for 500 years using constant 6 ka boundary conditions, which resulted in a quasi-steady state. We then integrated the model for 3000 years to 3 ka with annual updates to orbital parameters, GHG concentrations, and land cover composition.
A major limitation of CESM-CAM5 is that, like most state-of-the-art climate models, it does not explicitly represent glacial dynamics or even the formation of glaciers. In addition, the model imposes a maximum snow depth of 1 m snow water equivalent (SWE) to avoid an unrealistic accumulation of ocean salinity in the absence of explicit glacial calving. The nucleation of glaciers and ice caps must therefore be inferred, which we do here by assuming glaciation when a grid point’s end-of-summer (taken to be August) snow cover reaches a concentration of 1.0 and/or a depth of 1 m SWE. Presumably, simulated snow cover that survives an entire summer of melt would form the nucleus of a glacier and would eventually spread laterally under sufficiently cold regional conditions. Our identification of glaciated locations therefore serves as a very conservative estimate.
Results
Simulated mean-Arctic (60°–90°N) snow depth exhibits an upward trend during the late-Holocene, rising 23% from an annual mean of 9.5 cm at 6 ka BP to 11.7 cm by 3 ka BP (Figure 1a). This increase is presumably driven by the Milankovitch-based expectation of cooler summers from waning summer insolation, because greenhouse forcing is positive during this period (e.g. CO2 concentration rises from 263 to 273 ppm). The most conspicuous feature of this time series is the abrupt rise in Arctic-mean snow depth between 5 and 4.9 ka, centered around year 4931 BP. There was no prior indication of such a rapid change during the first millennium of the simulation; in fact, the regionally averaged snow depth was stable during the four centuries beforehand. Overall, the rate of snow depth increase between years 6000 and 4931 BP was 0.0535 cm century−1. Following the abrupt jump was a relaxation period of about two centuries, when the snow depth dropped from 10.6 to 10.2 cm before resuming its rise for the remaining 1700 years. This second positive trend of 0.074 cm century−1 is 38% larger than the initial trend through 4931 BP.

(a) Time series of Arctic (60°–90°N) mean snow depth (cm) between 6 and 3 ka and (b) linear trend in August SWE snow depth (m century−1) during the century of rapid snow cover buildup (5–4.9 ka). The Ellesmere Island study point (78°N, 80°W) is marked by an arrow.
A closer examination of the century with maximum Arctic snowpack buildup (5–4.9 ka) reveals that the rapid increase was confined to Ellesmere Island and coastal Greenland (Figure 1b), although the associated cooling encompassed the entire Arctic and led to a large expansion of sea ice cover (Supplemental Figure S1, available online). The largest snow accumulations occurred in far northern and northwestern Greenland and the highlands of Ellesmere Island. The snow buildup on Ellesmere was not only large but also extremely sudden, as shown for an island grid point identified in the figure (78°N, 80°W). The time series of August snow cover fraction and SWE during that century (Figure 2a) shows that end-of-summer snow cover was nearly non-existent during the decade preceding the transition but jumped to virtually complete coverage within a single year (4937 BP). This abrupt and sustained lateral expansion was followed immediately by a very large increase of snow depth, whose gridbox-averaged value was virtually zero in the years preceding the transition but rose to the model’s maximum of 1 m SWE within two decades. This location appears to have been on the brink of a similar abrupt change in the earlier part of the century, when August snow fraction also reached 1 for five consecutive years and snow depth reached almost 0.4 m SWE before less favorable summers resumed in the ensuing decades. Conversely, during the first millennium of the simulation (6–5 ka), August snow depth never reached the 1 m maximum, and August snow concentration remained below 1 in over 95% of the years (not shown).

August mean (a) snow cover fraction in solid circles and snow depth (SWE, m) in open circles and (b) surface temperature (°C) on Ellesmere Island (78°N, 80°W) during the century of rapid snow buildup from 5 to 4.9 ka. The melting point is indicated in (b).
These snow cover features are consistent with the simulated local mean-August surface temperatures on Ellesmere Island, which hovered on either side of the melting point until the transition (Figure 2b). The aborted snowpack buildup during the second decade was fostered by predominantly sub-freezing temperatures for several summers, but in most subsequent decades the surface remained above the melting point during August. In fact, this was consistently the case during the entire decade prior to the abrupt snowpack increase. Beginning in 4937 BP, however, the August temperature fell modestly below freezing (−2.5°C) and then dropped even more in the following few summers to a minimum of 9.4°C below freezing in 4932 BP and remained below freezing for the remainder of the century while snow concentration was 100%. This regime shift was promoted by a positive snow-albedo feedback, as the August surface albedo rose from 0.58 during the years immediately preceding the transition to around 0.80 in subsequent years with complete snow coverage.
Zeroing in on the decade when the abrupt transition occurred (4941–4932 BP) helps to explain the cause of the abrupt snow cover buildup. The annual cycle of monthly snow depth illustrates a very regular seasonal buildup and meltoff of the snowpack on Ellesmere Island during the first several years of the decade (Figure 3a). Accumulation occurred during every winter and spring, culminating in an annual maximum in June that was followed by rapid and complete meltoff during July and August. During 4937 BP, however, this cycle was disrupted by severely muted summer melting, such that the end-of-summer snow depth remained over 100 mm. During the following summer the expected melting was nearly non-existent and resulted in the snowpack remaining at 300 mm SWE. The normal annual cycle then resumed during the following years but at a much higher level, with August snow depths of around 400 mm or greater.

Monthly mean (a) snow depth (SWE, mm) and (b) water-equivalent snowfall (mm day−1) on Ellesmere Island (78°N, 80°W) during the decade of 4941–4932 BP. Each July is marked by a solid circle for reference.
The reason for the pair of consecutive years with a short-circuited melt phase is not only the cold weather during these two summers that reduced melting (Figure 2b) but also the abundance of summertime snowfall that boosted accumulation (Figure 3b). The amount of water-equivalent snowfall in the summers of 4937 and 4936 BP was double the long-term average during the decades prior to the abrupt transition: 1.05 and 1.12 mm day−1, respectively, compared with a mean of 0.54 mm day−1 in the preceding summers that century (years 5000–4938 BP). Thus, the total accumulated snowfall for these two summers was about 100 mm (water equivalent), or around 50 mm above the summer average before the transition, which represents a substantial portion of the normal snow pack depletion during previous summers (Figure 3a).
This excessive snowfall on Ellesmere Island was fostered by cold conditions and a very favorable weather pattern for precipitation (Figure 4). Anomalously cyclonic circulation prevailed locally and over most of the central Arctic during these two summers. The sea level pressure (SLP) during the first summer was 2 hPa below the century-long June–August average and was centered directly above Ellesmere Island, followed by a much more intense and widespread cyclonic anomaly of up to 9 hPa centered over the Arctic Ocean during the following summer (Figure 4a and b). The deviation during the second summer (4936 BP) represents an extremely positive phase of the Arctic Oscillation pattern (Thompson and Wallace, 1998), as reflected by the very high correlation (0.9) between the observed summertime Arctic Oscillation and North Pole SLP. In fact, the SLP at the North Pole that summer was the lowest (by 2.5 hPa) of any other summer that century. It was also a deeper pressure anomaly than in any summer during the entire first millennium of the simulation (from 6 to 5 ka) and nearly the very deepest over the entire 3000-year integration (Supplemental Figure S2, available online). The extremely cyclonic circulation during summer 4936 BP represents more than a three standard deviation departure from the long-term average, and it runs counter to the increasing secular trend of summertime North Pole SLP from 6000 to 3000 BP (1.5 hPa increase). The associated snowfall maps illustrate the impact of these very anomalous circulation patterns during summers 4937 and 4936 BP, consisting of above-normal regional snowfall over Ellesmere Island and western Greenland, as well as much of the Arctic Ocean (Figure 4c and d).

Summer (June–August) anomalies in years 4937 and 4936 BP of (a and b) sea level pressure (SLP, hPa) and (c and d) water-equivalent snowfall (mm day−1).
Following the transitional century from 5 to 4.9 ka, the simulation ran another 1900 years, during which time Ellesmere Island remained in a generally cold regime. August surface temperatures initially exceeded the melting point only occasionally for a few centuries and then stayed consistently below freezing during nearly the final 1500 years of the simulation (Supplemental Figure S3, available online). After the abrupt transition in the 4930s BP decade, the August snow coverage on Ellesmere Island remained at 100% for 88 consecutive years and then near the maximum for all but nine out of the next 136 years (the snow depth was at or close to the 1 m SWE depth limit for about the first 60 of these years; Supplemental Figure S4, available online, label A). After that, the model produced about 2.5 centuries of oscillatory behavior (Supplemental Figure S4, available online, label B) with a wide range of end-of-summer snow concentrations, but the majority of years had full coverage. Around 4470 BP, a less abrupt but more permanent transition into a fully glaciated state ensued, with consistent August snow concentrations of 100% and snow depths of 0.8 m or more in virtually every summer for the remainder of the simulation (Supplemental Figure S4, available online, label C).
Discussion and conclusion
The findings of this study lead to three primary conclusions. First, they show that a climate model can simulate dramatically rapid climatic transitions in a region where glacial inception is known to occur and where Greenland ice cores have documented similarly rapid transitions (Alley et al., 1993; Andrews and Mahaffy, 1976). Second, they demonstrate a complementary mechanism to the traditional Milankovitch thermodynamic perspective that cooler summers cause less snow melt and therefore a buildup of snow cover. Our simulation indicates that cooler summers can also produce more snowfall and thereby promote a rapid thickening of the snow pack. Third, the simulated atmospheric circulation patterns during the identified pair of transitional summers provide a logical explanation for the unusual amount of summer snowfall that occurred on Ellesmere, and both of these summers – especially the second – feature a well-defined mode of circulation variability (a strongly positive phase of the AO).
Our conclusions need to be interpreted in the context of several factors and unresolved questions, perhaps the most important of which is the constraint imposed by CESM’s treatment of snow cover. This limitation includes both the absence of an explicit glaciation process and the model’s snow-capping parameterization, which restricts SWE to 1 m depth regardless of the local snow budget. These simplifications short-circuit the expected evolution of surviving summertime snow cover into glaciers and thus suggest that our inferences of neoglaciation are conservative estimates. In particular, this model constraint calls into question the realism of the multi-centennial transitional period on Ellesmere Island following the rapid snow cover buildup between 5 and 4.9 ka that preceded the eventual permanently glaciated state (phase B in Supplemental Figure S4, available online). Another uncertainty is that the direction of causality is not clear: did the regional summertime cooling pattern around Ellesmere Island during the 4930s BP decade cause more snowfall or vice versa? Therefore, we cannot easily isolate the relative importance of temperature versus precipitation in contributing to the rapid development and maintenance of the year-round snow pack. An additional question is the cause of the extreme, short-term atmospheric circulation anomalies conducive to the summertime snowfall and snow cover buildup during years 4937 and 4936 BP. The slightly positive trend in central Arctic summer SLP during the first 1100 years of the simulation suggests that this pair of summers were not a response to the overarching orbital driver that operates on much longer timescales. Because our model integration contains no short-term external forcings such as volcanic eruptions, the results therefore suggest that an extreme internal climate variation triggered the abrupt transition event. These findings resemble those of Drijfhout et al. (2013) and Kleppin et al. (2015), who simulated rapid regional cooling around Greenland initiated by stochastic SLP variations. However, those abrupt changes required a positive feedback loop involving coupled atmosphere-ocean-sea ice interactions on decadal time scales, whereas our rapid change mechanism depends more simply on a short-term atmospheric circulation anomaly.
Without additional model simulations with slightly perturbed initial conditions, we cannot know whether the region would have experienced neoglaciation around the same time without the contribution from the pair of exceptional summers. Nevertheless, our findings contribute to a better understanding of Holocene Neoglaciation in particular and more broadly to insights into possible dynamical mechanisms responsible for rapid Arctic climate change and glacial inception evident in the paleoclimate record.
Supplemental Material
Figure_S1 – Supplemental material for Rapid neoglaciation on Ellesmere Island promoted by enhanced summer snowfall in a transient climate model simulation of the middle-late-Holocene
Supplemental material, Figure_S1 for Rapid neoglaciation on Ellesmere Island promoted by enhanced summer snowfall in a transient climate model simulation of the middle-late-Holocene by Stephen J Vavrus, Feng He, John E Kutzbach and William F Ruddiman in The Holocene
Supplemental Material
Figure_S2 – Supplemental material for Rapid neoglaciation on Ellesmere Island promoted by enhanced summer snowfall in a transient climate model simulation of the middle-late-Holocene
Supplemental material, Figure_S2 for Rapid neoglaciation on Ellesmere Island promoted by enhanced summer snowfall in a transient climate model simulation of the middle-late-Holocene by Stephen J Vavrus, Feng He, John E Kutzbach and William F Ruddiman in The Holocene
Supplemental Material
Figure_S3 – Supplemental material for Rapid neoglaciation on Ellesmere Island promoted by enhanced summer snowfall in a transient climate model simulation of the middle-late-Holocene
Supplemental material, Figure_S3 for Rapid neoglaciation on Ellesmere Island promoted by enhanced summer snowfall in a transient climate model simulation of the middle-late-Holocene by Stephen J Vavrus, Feng He, John E Kutzbach and William F Ruddiman in The Holocene
Supplemental Material
Figure_S4 – Supplemental material for Rapid neoglaciation on Ellesmere Island promoted by enhanced summer snowfall in a transient climate model simulation of the middle-late-Holocene
Supplemental material, Figure_S4 for Rapid neoglaciation on Ellesmere Island promoted by enhanced summer snowfall in a transient climate model simulation of the middle-late-Holocene by Stephen J Vavrus, Feng He, John E Kutzbach and William F Ruddiman in The Holocene
Footnotes
Acknowledgements
We would like to acknowledge high-performance computing support from Cheyenne (doi:10.5065/D6RX99HX) and computing resources (doi:10.5065/D6RX99HX) provided by NCAR’s Computational and Information Systems Laboratory, sponsored by the National Science Foundation and other agencies.
Funding
The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work was funded by National Science Foundation grants AGS-1602771 and AGS-1602967.
Supplemental material
Supplemental material for this article is available online.
References
Supplementary Material
Please find the following supplemental material available below.
For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.
For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.
