Lake levels in a discontinuous permafrost landscape: Late Holocene variations inferred from sediment oxygen isotopes, Yukon Flats, Alaska

ABSTRACT

During recent decades, lake levels in the Yukon Flats region of interior Alaska have fluctuated dramatically. However, prior to recorded observations, no data are available to indicate if similar or more extreme variations occurred during past centuries and millennia. This study explores the history of Yukon Flats lake origins and lake levels for the past approximately 5,500 years from sediment analyses guided by previous work on permafrost extent, thermokarst, and modern isotope hydrology. Sediments dated by ²¹⁰Pb and AMS radiocarbon indicate stable chronologies following initial lake initiation. Subsequent lithology is autochthonous, and oxygen isotope ratios of endogenic carbonate reflect lake level change at multiple time scales. Sediment results indicate high lake levels between approximately 4000 and 1850 cal yr BP, which is interpreted to reflect wetter-than-modern conditions. Lower lake levels with short-lived high stands during the past approximately 800 years reflect generally arid conditions with brief wet intervals similar to the region’s moisture regime today. The millennial trend is one of increasing aridity and corresponds closely with fire reconstructions and regional paleoclimatic trends. We conclude that high-magnitude lake-level fluctuations and decadal scale trends occurred before the observational period and are persistent hydroclimatic features of the Yukon Flats region.

KEYWORDS:

• Alaska
• Holocene
• paleoclimate
• Yukon Flats
• lake sediment δ¹⁸O
• thermokarst

Introduction

Observations of lake-level declines between the late 1980s and 2014 in interior and northern Alaska subsequently generated efforts to identify spatial distributions, causes and mechanisms, and ecological impacts (Riordan, Verbyla, and McGuire 2006; Smith et al. 2005). Early studies focused on continuous permafrost regions and implicated warmer summer temperatures with an emphasis on the role of permafrost degradation and thermokarst draining (Hinzman et al. 2005). Subsequent assessments were broadened to include Subarctic boreal regions with discontinuous permafrost distributions, such as the Yukon Flats (YF) of interior Alaska (Figure 1). In addition to permafrost change, they highlighted the significance of surface–groundwater connections and direct climatic influences on lake hydrology and water budgets (Anderson et al. 2013; Jepsen et al. 2012; Lewis et al. 2014; Roach et al. 2011; Rover, Wylie, and Tieszen 2012).

Figure 1. Location of the Yukon Flats (inset) is shown on a digital elevation model (30 m), with surficial geology by Williams (1962) and the locations of Twelvemile and Track Lakes (large white circle symbols) and the locations of evaporite deposits (red-filled circle symbols) from Clautice and Mowatt (1981). Lakes from Anderson et al. (2013) are indicated as small open circle symbols for those with ¹⁸O-depleted water sources, such as snow and permafrost. Black-filled circle symbols indicate lakes sourced by precipitation, rivers, and groundwater (gw). The inset shows the locations of additional sites discussed in the text

The YF landscape and its Holocene history of thermokarst evolution appears to differ from other thermokarst-affected regions of Alaska (Edwards et al. 2016). Most notably, current lake-expansion rates are relatively low and thermokarst does not show broadly distributed evidence for extensive lake draining (e.g., Jones and Arp 2015) or multiple lake generations within one basin, such as observed in the northern Yukon Territory (Lauriol et al. 2009). Evidence suggests that following widespread early Holocene permafrost degradation and thermokarst initiation, the YF landscape has since been relatively stable, most likely because of the establishment of dense boreal forest cover and the insulating effect it has on permafrost.

Discontinuous permafrost distributions and widespread surface–groundwater connectivity are known to occur throughout the YF region (Hawkins 1995; Minsley et al. 2012; Williams 1970). To date, changes in YF permafrost have only been detected at shallow depths, and the subsequent influence by changes in groundwater-surface interactions on recent lake declines appears to be relatively minor (Jepsen et al. 2012, 2013; Riordan, Verbyla, and McGuire 2006).

Climate controls on YF lake levels are clearly indicated by substantial lake evaporation, as reflected by high dissolved ion concentrations and the wide occurrence of evaporite deposits (Figure 1; Clautice and Mowatt 1981; Heglund and Jones 2003). Changes in lake hydrologic balance in response to climate are also reflected by lake-water oxygen and hydrogen isotope ratios (δ¹⁸O and δ²H) as they evolve in response to changing hydrologic fluxes (Gibson, Birks, and Yi 2016), including lakes of thermokarst origin (MacDonald et al. 2017). Lake-water evaporation leads to enrichment in heavy isotopes (¹⁸O for oxygen and ²H for hydrogen) and is also promoted by restricted or closed-lake hydrology with lower lake levels.

Extensive YF lake-water isotope sampling (Figure 1) between 2007 and 2011 by Anderson et al. (2013) revealed that lake isotopic enrichment is closely associated with long-term lake-level declines, low relative humidity, increased temperatures, and reduced precipitation. Remote sensing inventories for 1985 through 2011 also indicate that hydrologically closed YF lakes are more vulnerable to declining levels. A higher proportion of closed-basin lakes have declined (15%, Chen et al. 2014) compared with the proportion of all lakes (open and closed) that have declined (3.4%, Rover, Wylie, and Tieszen 2012).

To develop a long-term context and better understanding of the significance of recent YF lake fluctuations, this study examines sediment records from two sites that have recently changed dramatically in response to climate: Twelvemile and Track Lakes (Figure 1). Both lakes underwent sustained water losses during the past two decades, with approximately 80 percent volume declines and subsequent enrichment in heavy oxygen and hydrogen isotopes in their lake water. Both lakes have also undergone full recovery to maximum lake levels within a thirteen-month period between 2014 and 2015 in response to higher than average precipitation. The lakes are located in similar settings in terms of surficial geology, vegetation, permafrost distributions, and surface–groundwater connectivity. Evidence presented here for stable radiocarbon chronologies and autochthonous lithology further indicates that following initial lake formation, minimal permafrost-induced changes in shorelines or basin morphometry has occurred.

Twelvemile and Track Lakes are alkaline and undergo endogenic carbonate (CaCO₃) precipitation and sedimentation. The sediment records therefore provide continuous calcium carbonate oxygen isotope ratios (δ¹⁸O_CaCO3). The δ¹⁸O_CaCO3 values reflect the values for past lake water (δ¹⁸O_H2O), which in turn reflects lake hydrology, lake levels, and climate (Anderson et al. 2016; Leng et al. 2005). Additional sediment measurements, lake-water δ¹⁸O_H2O observations, climate records, lake-level histories, and isotope mass balance model sensitivity experiments further inform the paleoclimatic interpretation of the Track and Twelvemile lakes δ¹⁸O_CaCO3 records. The combined data sets indicate that they accurately reflect lake-level fluctuations during the observational period and thereby provide a record of lake-level change in response to hydroclimate over long timescales beyond the period of monitoring. Our results show that recent YF lake-level dynamics are within the variation experienced during the past approximately 4,000 years and suggest that these lake systems are strongly resilient to climate-induced perturbations.

Study area

Track and Twelvemile lakes are located on alluvial lowlands within the Yukon Flats National Wildlife Refuge (YFNWR), which encompasses the confluence of the Yukon, Chandalar, Sheenjek, and Porcupine Rivers (Figure 1). Surrounding highlands include a dissected southern marginal upland 30–150 m higher than the lowlands. The YF basin is bounded by the Brooks Range to the north, Hodzana Highlands, White and Crazy Mountains, and Porcupine Plateau to the south, and Porcupine Plateau to the east. Alluvial lowland YF vegetation is characterized by a dynamic mosaic of open meadow, spruce forest (Picea glauca and P. mariana), birch forest (Betula neoalaskana and Populus tremuloides), and tall shrubs (salix and alnus; e.g., Jorgenson et al. 2010). Fire is the principal large-scale disturbance (Drury and Grissom 2008). Continuous and discontinuous permafrost distribution leads to abundant surface water such as bogs, fens, wetlands, and lakes. Lakes are principally of thermokarst, eolian, and fluvial origin because the area was unglaciated during the Pleistocene (Williams 1970).

The YFNWR is one of the most productive waterfowl breeding areas in North America, and lake levels are relevant to management. Preliminary YFNWR ground-based assessments by Rose and Guldager (2012) in the fall of 2010 explored relationships between lake levels, soils, groundwater, and permafrost and found a general lack of permafrost around two drying lakes. In contrast, permafrost was more frequently encountered around two expanding lakes and one lake with stable levels. In general, water-table depths during the growing season did not significantly differ among lake types.

Twelvemile Lake (official name, 66.44966°N, −145.54564°W, 115 m a.s.l.) is located approximately 12  miles (7 km) southeast of Fort Yukon on gently sloped Quaternary alluvium trending from the southeast to northwest (Figure 1). The lake is the terminus of a chain of lakes connected by surface and subsurface flow that are oriented parallel to the modern Yukon River floodplain. During high stands, streamflow is observed into and out of the basin. During low stands, the inflow and outflow stream channels are dry, and nearby spruce-covered permafrost plateaus with 5–20 m of relief define the watershed (Figure 2). Surface flows into the lake basin were observed and photographed by YFNWR pilots in April 2013 during an upstream ice jam flood event on the Yukon River. Distinct groundwater recharge localities were not observed, although permafrost surveys have identified that groundwater flow paths above and below permafrost are related to fluvial and lacustrine taliks, soils, and vegetation (Jepsen et al. 2013; Jorgenson et al. 2010; Minsley et al. 2012).

Figure 2. Twelvemile Lake is shown on a bare-earth digital terrain model derived from LIDAR (2 m) with locations for depth measurements (blue circle symbols) and modeled bathymetry (blue shading), core locations (white circle symbols), and shorelines from 2008 (light blue), 2013 (yellow), and 2016 (red)

Track Lake (unofficial name, 66.855°N, 145.170°W, 145 m a.s.l.) is located approximately 32 km north of Fort Yukon, between the Sheenjek and Christian rivers on gently sloped Quaternary alluvium trending from northeast to southwest (Figure 1). During high stands, stream flow is observed out of the basin but not in. During low stands, nearby spruce-covered permafrost plateaus with 5–20 m of relief define the watershed (Figure 3). Ice-jam floods of the Yukon River do not affect the Track Lake area. A prominent focal point of groundwater discharge is in the northeast basin. Groundwater paths are likely associated with fluvial taliks, soils, and vegetation similar to those documented at Twelvemile Lake.

Figure 3. Track Lake is shown on imagery for April 2014 and July 2016 (Digital Globe WV02), with coring location (white circle symbols), estimated bathymetry (dashed colored lines), and shorelines (solid lines) from 2009 (orange), 2014 (yellow), and 2016 (red)

The climate of the Yukon Flats is characterized as cold continental with large seasonal temperature contrasts, low humidity, and low precipitation. Mean annual temperature is approximately −6°C with recorded maxima of approximately  35°C and minima of −57°C. Subfreezing temperatures can occur in any month. The broad low-lying valley is strongly affected by cold-air drainage and winter temperature inversions. Mean annual total precipitation is approximately  250 mm with approximately  80 mm as snow (U.S. National Centers for Environmental Information, www.ncdc.noaa.gov). In general, precipitation is light and irregular with summer rainfall accounting for approximately 68 percent of annual totals (Figure 4; Shulski and Wendler 2007). Lakes are ice covered from seven to eight months of each year (September/October–April/May).

Figure 4. Fort Yukon SNOTEL monthly precipitation and rain/snow ratio (1997–2016) plotted against Yukon Flats (YF) relative humidity (NCEP/NCAR reanalysis percent anomaly from 1948 to 2013 mean for average of YF grid cells), Twelvemile and Track Lake surface areas (percent of maximum, diamond symbol), and sediment δ¹⁸O_CaCO3 (‰ VDPB, small circle symbol), living Chara δ¹⁸O (‰ VDPB, large filled circle symbol), and lake water δ¹⁸O_H2O (‰ VSMOW, open circle symbol) from approximately 1950 to 2016. MAP = mean annual precipitation

Anderson et al. (2013) obtained meteorological data from 1948 to 2011 from the U.S. National Centers for Environmental Prediction and National Center for Atmospheric Research (NCEP-NCAR) for an average of grid cells that encompass the YFNWR. The relative humidity record indicates that it was relatively dry between 1948 and the early 1970s, wetter between the 1980s and 1995, and dry again between 1995 and 2011 (Figure 4). For the period of record, wet periods with high relative humidity are associated with more snowfall and lower potential evaporation, whereas dry periods are associated with less snowfall and higher rates of potential evaporation (shown in Anderson et al. 2013).

Local precipitation measurements from 1997 to the present are available from the U.S. Department of Agriculture National Resource Conservation (USDA NRCS) snowpack telemetry network (SNOTEL) site at Fort Yukon. The twenty-year record shows anomalously high rainfall during July through September in 2014 (Figure 4). August 2015 had high rainfall, and September had two large snowfall events. Although mean annual 2016 precipitation was average, June and July were drier than normal and August was wetter than normal.

Methods

Bathymetry was determined in 2008 at Twelvemile Lake from 255 GPS-located sonar-measured water depths and 142 shoreline locations (Figure 2). The data were mapped on a 2 m grid cell bare earth digital terrain model derived from airborne light detection and ranging (lidar) data acquired between July 14 and September 3, 2009, by the U.S. Geological Survey and using a kriging method of interpolation with GIS. The bathymetry of Track Lake was estimated by assuming roughly evenly spaced contours between the GPS elevation measured of the lake’s surface in 2009 (145 m a.s.l.) and the location of the maximum measured depth (7.5 m, 138 m a.s.l) (Figure 3). Additional elevations were then estimated for the highest and lowest observed shorelines (144 m in 2011 and 147 m in 1985). Lake volumes and surface areas were calculated for both lakes at 0.5 m depth increments using a kriging interpolation method within ArcGIS 3D Analyst and are shown in Tables 1 and 2.

Table 1. Twelvemile Lake hypsometry estimates

Elevation (m)	Depth (m)	Area (m^2)	Volume (m^3)	Year	% Area	% Volume
131.0	14.5	6,253,126	17,652,919	1978	100.0	100.0
130.5	14	5,070,777	14,832,588	2016	81.1	84.0
130.0	13.5	4,156,839	12,534,595		66.5	71.0
129.5	13	3,361,714	10,665,409		53.8	60.4
129.0	12.5	2,691,604	9,159,641		43.0	51.9
128.5	12	2,228,292	7,940,001	2013	35.6	45.0
128.0	11.5	1,886,095	6,913,724		30.2	39.2
127.5	11	1,626,044	6,040,474	2008	26.0	34.2
127.0	10.5	1,407,851	5,287,076		22.5	30.0
126.5	10	1,204,527	4,637,423		19.3	26.3
126.0	9.5	1,072,631	4,070,142		17.2	23.1
125.5	9	975,879	3,559,049		15.6	20.2
125.0	8.5	894,762	3,091,817		14.3	17.5
124.5	8	820,414	2,663,136		13.1	15.1
124.0	7.5	751,622	2,270,263		12.0	12.9
123.5	7	679,124	1,912,361		10.9	10.8
123.0	6.5	613,448	1,589,193		9.8	9.0
122.5	6	553,160	1,297,531		8.8	7.4
122.0	5.5	495,306	1,035,464		7.9	5.9
121.5	5	430,496	803,804		6.9	4.6
121.0	4.5	352,258	606,856		5.6	3.4
120.5	4	276,727	450,003		4.4	2.5
120.0	3.5	232,381	323,092		3.7	1.8
119.5	3	190,316	217,275		3.0	1.2
119.0	2.5	142,500	133,562		2.3	0.8
118.5	2	100,087	74,309		1.6	0.4
118.0	1.5	70,355	31,767		1.1	0.2
117.5	1	30,546	6356		0.5	0.0
117.0	0.5	1,383	42		0.0	0.0

Table 2. Track Lake hypsometry estimates

Elevation (m)	Depth (m)	Area (m^2)	Volume (m^3)	Year	% Area	% Volume
147.0	9.5	375,351	842,226	2016	100.0	100.0
146.5	9	310,954	670,631	1985	82.8	79.6
146.0	8.5	241,351	533,243		64.3	63.3
145.5	8	201,205	422,929		53.6	50.2
145.0	7.5	164,378	331,353	2009	43.8	39.3
144.5	7	122,859	260,112	2010	32.7	30.9
144.0	6.5	99,390	205,244	2011, 2014	26.5	24.4
143.5	6	81,075	160,269		21.6	19.0
143.0	5.5	66,563	123,477		17.7	14.7
142.5	5	54,669	93,255		14.6	11.1
142.0	4.5	44,177	68,586		11.8	8.1
141.5	4	34,862	48,888		9.3	5.8
141.0	3.5	26,720	33,540		7.1	4.0
140.5	3	19,901	21,936		5.3	2.6
140.0	2.5	14,202	13,456		3.8	1.6
139.5	2	9,578	7,551		2.6	0.9
139.0	1.5	5,896	3,733		1.6	0.4
138.5	1	3,440	1,438		0.9	0.2
138.0	0.5	1,443	216		0.0	0.0

In July 2007, surface and water-column measurements were collected with a calibrated Eureka Amphibian™ and Manta sonde. Lake-water chemistry is presented by Halm and Guldager (2012). Near-surface water samples for isotope analyses were taken in 30 ml HDLP Nalgene bottles and sealed with no headspace. Water samples were prepared for isotope ratio analyses of oxygen and hydrogen by automated constant temperature equilibration with CO₂ and automated H/D preparation by chromium reduction coupled to an isotope ratio mass spectrometer. Water oxygen and hydrogen isotope results are reported in per mil (‰) as the relative difference of isotope ratios (δ) from the international measurement standard Vienna Standard Mean Ocean Water (VSMOW) defined by δ¹⁸O_H2O = ([¹⁸O/¹⁶O] _H2O/[¹⁸O/¹⁶O]_SMOW) − 1 and δ²H_H2O = ([²H/¹H] _H2O/[²H/¹H]_SMOW) − 1. Analytical uncertainties are within ± 0.5‰ and ± 0.05‰ for hydrogen and oxygen, respectively. Analyses of random replicate samples fell within the range of analytical error.

Sediment cores with an undisturbed sediment water interface were obtained from an anchored floatplane at Twelvemile Lake (core GL07, 36 cm length, 7.5 m water depth) and Track Lake (core A09, 38 cm length, 5 m water depth) using Glew (Glew, Smol, and Last 2001) and Aquatics Research™ gravity coring devices in July 2007 and 2009, respectively. Following vertical storage during the return flight to Fairbanks, cores were extruded at 0.5 cm increments. Additional modified Livingstone cores were obtained at Twelvemile Lake in September 2008 from a securely anchored floating platform in water depths of 10 m (core A08), 8 m (core B08), and 2.8 m (core C08). Black organic rich sediments in core A08 were disturbed after recovery by degassing and were not further analyzed.

At Track Lake, additional modified Livingstone “Bolivia” cores were obtained in July 2010 near the 2009 coring location from an anchored floatplane in 4–5 m water depth (core B10). Cores were stabilized with Zorbitrol Plus™ and were capped before transportation. Windy conditions in 2010 obstructed precise measurements of water and drive depths so that although every effort was made to begin coring just below the sediment–water interface, the dating of the core, discussed in detail further on, indicates that the uppermost depth was in fact between 40 cm and 50 cm below the sediment–water interface. We estimate that this left a gap of approximately 8–10 cm unattained sediment between the base of gravity core A09 to the top of Livingstone core B10.

Following whole-core magnetic susceptibility measurements at 1 cm increments with a Bartington susceptibility meter, the Livingstone cores were split, photographed, and visually logged for biogenic features, sedimentary structures, composition, and Munsell color. Dry bulk density was obtained from contiguous 0.5 cm or 1.0 cm increments. Percent organic and inorganic carbon was determined from dry, pulverized subsamples using total carbon (TC) and total inorganic carbon (TIC) measurements by a UIC Inc.™ carbon dioxide coulometer (Engleman, Jackson, and Norton 1985). Percent organic carbon was calculated as the difference between TC and TIC. The accuracy and precision for both TC and TIC is 0.1 percent.

Subsamples for organic carbon (OC) and nitrogen percentage and isotope analyses were acidified with 1 N HCl, rinsed to neutral, freeze dried, and pulverized for combustion to CO₂ and NO₂ by a Carlo Erba™ CN elemental analyzer coupled with a Finnigan Delta™ Advantage isotope ratio mass spectrometer. Organic carbon isotope ratios are reported in per mil as the relative isotope-ratio difference from the international standard Vienna Pee Dee Belemnite (VPDB), defined by δ¹³C_Org = ([¹³C/¹²C]_Org/[¹³C/¹²C]_PDB) − 1. Analytical uncertainties are ±0.2‰ and analyses of random replicate samples fell within the range of analytical error.

Subsamples for carbonate isotope analyses were obtained at 0.5 cm increments. Sediment carbonate particles include macroscopic Chara stem encrustations, ostracodes, gastropods, mollusks, and their shell fragments. Scanning electron microscopy and energy dispersive X-ray spectroscopy (SEM-EDS) indicated that the 32–63 µm sediment matrix contained mineral grains with endogenic carbonate growth structures. However, to remove as much of the larger biological carbonate shell fragments as possible, bulk sediments were wet sieved through nested screens (250, 125, 63, 32 µm), collected separately, and freeze dried. The finest fractions, either the less than 32 µm size or combined with the 32–64 µm sizes if the less than 32 µm sample was of insufficient mass, were visually examined for purity and powdered before subsampling for CO₂ extraction by a Kiel automated device coupled with a Finnigan Delta 252 isotope ratio mass spectrometer. Calcite oxygen and carbon isotope ratios are reported in per mil as the relative isotope-ratio difference from VPDB defined by δ¹⁸O_CaCO3 = ([¹⁸O/¹⁶O]_CaCO3/[¹⁸O/¹⁶O]_PDB) − 1 and δ¹³C_CaCO3 = ([¹³C/¹²C]_CaCO3/[¹³C/¹²C]_PDB) − 1. Analytical uncertainties are ±0.05‰ and 0.1‰ for oxygen and carbon, respectively, and analyses of random replicate samples fell within the range of analytical error.

The sediment core chronologies are based on ²¹⁰Pb analyses of bulk sediment and AMS radiocarbon measurements on identifiable terrestrial macrofossils, including charcoal. Terrestrial macrofossils were not present in sufficient quantities for radiocarbon measurements at most stratigraphic levels, and in some cases purified pollen samples were used. Both measured and calibrated radiocarbon ages (Calib 7.0.1; Stuiver and Reimer 1993) are reported, but only calibrated ages (cal yr BP) are used for discussion (Tables 3 and 4). ²¹⁰Pb was measured on 0.5 cm sediment slices at 1.5 cm intervals to 26 cm and 36 cm depth for Track and Twelvemile lakes, respectively (Figure 5). ²¹⁰Pb activities are low yet typical for high-latitude regions, with maxima of 9.2 pCi/g at Twelvemile Lake and 5.26 pCi/g at Track Lake. Down-core decreases are relatively smooth, and the break between sediments containing unsupported ²¹⁰Pb is well defined. A constant rate of supply model (CRS) was used to determine ²¹⁰Pb ages (Appleby 2001).

Table 3. Twelvemile Lake chronostratigraphic data

Core Depth (cm)	Material	^210Pb Age (AD)	^210Pb Age Error	AMS ^14C Lab #	^14C Age (^14C yr BP)	δ^13C (‰ VPDB)	Calibrated Age^1 (Cal yr BP)	1σrange^2 (Cal yr BP)
GL07
0	Surface	2007	1.3				−57
3.25	Lead-210	1996	1.3				−56.5
4.75	Lead-210	1990	1.3				−46.2
6.25	Lead-210	1984	1.4				−40.3
7.75	Lead-210	1976	1.6				−33.8
9.25	Lead-210	1968	1.9				−26.4
10.75	Lead-210	1960	1.3				−17.5
12.25	Lead-210	1952	1.4				−9.8
13.25	Lead-210	1947	1.4				−2.1
14.25	Lead-210	1942	1.6				2.7
15.25	Lead-210	1934	1.8				8.2
16.25	Lead-210	1923	2.3				16.0
17.25	Lead-210	1909	3.3				26.7
18.25	Lead-210	1893	4.7				40.6
18.75	Lead-210	1884	6.3				56.7
19.25	Lead-210	1877	7.7				66.2
19.75	Lead-210	1868	10.1				73.3
20.75	Lead-210	1846	19.5				82.2
25	Potamageton seed			WW-6805	304 ± 35		370	60
34	Vascular plant parts			WW-6806	544 ± 57		575	60
B08
41	Scirpus seeds			CURL-18248	2095 ± 15	−31.9	2060	50
48	Wood			CURL-10975	2235 ± 15	−25.1	2210	30
80	Wood			CURL-12032	2275 ± 20	−30.3	2280	65
134	Wood and plant fragments			CURL-12016	2565 ± 20	−26.1	2730	15
CO8
17.5	Wood			CURL-13036	1260 ± 15	−26.6	1215	40
38	Wood			WW-9128	2065 ± 30	−27.8	2050	55
49	Wood			CURL-13035	2445 ± 15	−26.8	2525	155

¹Calib 7.0.1, median of 1-sigma range.

²CRS age model.

Table 4. Track Lake chronostratigraphic data

Core Depth (cm)	Material	^210Pb Age (AD)	^210Pb Age Error	AMS ^14C Lab #	^14C Age (^14C BP)	δ^13C (‰ VPDB)	Calibrated Age^1 (Cal yr BP)	1σ range^2 (Cal yr BP)
ARG09
0	Surface	2009	1.3				−59
0.25	Lead-210	2007	1.3				−58
1.75	Lead-210	2001	1.3				−51
3.25	Lead-210	1994	1.3				−44
4.75	Lead-210	1986	1.3				−36
5.75	Lead-210	1980	1.4				−30
6.75	Lead-210	1974	1.5				−24
7.75	Lead-210	1967	1.6				−17
9.25	Lead-210	1956	1.7				−6
10.25	Lead-210	1947	2.1				3
11.75	Lead-210	1936	2.0				14
12.75	Lead-210	1913	3.3				37
13.75	Lead-210	1890	5.7				60
14.75	Lead-210	1868	10.4				82
15.75	Lead-210	1845	18.8				105
36.5	Pollen extract			CURL-18282	885 ± 15	−28.5	820	75
B10
8–9	Charcoal			CAMS-178263	1910 ± 90		1850	115
24.5	Pollen extract^3			CURL-18270	2915 ± 15	−29.4	3010	45
31.5	Seeds			CURL-15802	2465 ± 20	−26.6	2820	115
50	Wood			CURL-14989	3490 ± 20	−25.0	3770	55
71–74	Charcoal			CAMS-178263	4750 ± 70		5485	250

¹Calib 7.0.1, median of 1-sigma range.

²CRS age model.

³Not included in the age model.

Figure 5. Twelvemile core GL07 and Track core A09 ²¹⁰Pb total activity, CRS model ages, and sediment accumulation rates with error ranges

Age models were generated with ²¹⁰Pb and radiocarbon measurements using the Bayesian software Bacon (v2.2; Blaauw and Christen 2011 and 2013), which uses the most recent calibration curve, IntCal 13 (Reimer et al. 2013). Default prior assumptions were used, with the exception of thickness values that were set to 1. Bacon generated two equally plausible age models for Twelvemile Lake below the ²¹⁰Pb-dated depths that differ by approximately 400 years, which is an uncertainty greater than 50 percent of the length of the 600-year record. The cause is large 2-sigma calibrated age ranges (~ 200–300 years) for the seed and plant samples at 25 cm and 34 cm depth. Therefore, the Twelvemile Lake data are only shown by age for depths dated by ²¹⁰Pb, and interpretations of older data take the significant age uncertainty into consideration. The range of age errors for GL07 ²¹⁰Pb-dated depths (0.5–20.5 cm, 2008 to approximately AD 1870) is 1–10 yr (Figure 5, Table 3).

Individual age models for each Track Lake core, A09 and B10, avoid assumptions about the starting depth of core B10. Furthermore, core B10 data above the uppermost radiocarbon age at 9 cm depth of 1850 cal yr BP are not shown or interpreted. However, the age models indicate the strong likelihood of continuous sedimentation through the unattained sediments between the cores (Figure 6). Extrapolation of estimated sedimentation rates between the A09 ²¹⁰Pb-horizon and the lowermost radiocarbon age suggest that B10 has an uppermost depth of 40–50 cm, or between approximately 8 cm and 10 cm below the 36.5 cm depth of A09. The range of ²¹⁰Pb age uncertainties for A09-dated depths (0.5–16 cm, 2009 to approximately AD 1870) are 1–10 yr (Figure 5, Table 4). Below the ²¹⁰Pb-horizon, radiocarbon age-model uncertainty increases from 130 to 200 years. The range of B10 radiocarbon age-model uncertainty is from 100 to 500 years.

Figure 6. Bayesian age-depth models constructed with the age-depth modeling software Bacon v2.2 (Blaauw and Christen 2011 and 2013) for Track cores A09 and B10 on individual depth scales. Transparent blue symbols show calibrated radiocarbon dates and age-depth models are gray scale, with darker grays indicating more likely calendar ages. Stippled lines show 95 percent confidence intervals and the red dashed curve shows the best model fit based on the weighted mean age for each depth. Upper panels depict the Markov Chain Monte Carlo (MCMC) iterations (right panel) and prior (green) and posterior (gray histograms) distributions for the accumulation rate (middle) and memory (right panel). Extrapolated sedimentation rates suggest that unattained sediment gap was approximately 8–10 cm thick and spans from approximately 1850 to 800 cal yr BP

Results

Limnology

In July 2007, Twelvemile Lake had a maximum water depth of 10 m, and water-column measurements indicated a thermocline depth of 5–6 m with a 10°C temperature decline, no change in pH, a slight decline in dissolved oxygen (DO), and slight increase in specific conductivity (Spc). The high pH and Spc are greater than average for the YFNWR (Table 5). Alkalinity was in excess of dissolved magnesium (Mg), which exceeded calcium (Ca) by a factor of approximately 2. Dissolved sodium (Na) was below average. Inorganic carbon (DIC) concentrations were relatively high, and carbon isotope ratios (δ¹³C-DIC) were relatively low. The ice-free lake was oversaturated in dissolved carbon dioxide (CO₂) and methane (CH₄). Exposed shorelines examined in 2008 around the northern edges of the basin were silty sands delineated by strands of willow and aspen saplings. Southern shorelines were exposed organic-rich lake muds, and in some locations had thin crusts of white evaporite deposits. Submerged and emergent aquatic vegetation is abundant throughout the lake, such as pond weed (Potamogetonaceae sp.) and water lilies (Nymphaeaceae sp.), but Chara sp. algae was not observed.

Table 5. Aqueous geochemistry of Track and Twelvemile lakes¹

Parameter	Track^2	Twelvemile^3	YFNWR Mean^4
Field measurements
Temperature (°C)	19.5	10–17
pH	8.9	8.9	7.8
Specific conductivity (µS/cm)	600	500–700	442
Dissolved oxygen (mg/L	9.4	11.1
Dissolved ions (mg/L)
Ca	23.8	26–29	27.9
Mg	93.7	60–69	20.6
Na	40.4	18.4	45.9
K	19.7	23–27	10.0
Alkalinity	313	300–360	253
Cl	18.0	7.6–9.8	3.6
SO4	2.3	42–49	8.2
Dissolved carbon and gases
DIC (µmol/L)	5486	4,521–6,358
δ13C-DIC	4.6	−0.10–−0.16
CO2 (ppm)	391	180–1711
CH4 (ppm)	56.6	8.9–29.3

¹In Halm and Guldager (2012).

²Identified as Lake 159, sampled July 2009.

³Range of values between August 2007, March 2008, and May 2008.

⁴Summary statistics for shallow lakes (n = 129) from Heglund and Jones (2003).

In July 2009, Track Lake had a maximum water depth of 7.5 m. Water-column measurements were not obtained, but temperature stratification is considered likely based on a similar maximum depth, pH, and Spc as Twelvemile Lake. Alkalinity was in excess of Mg, which was higher than Twelvemile and exceeded Ca by a factor of from 3 to 4. Dissolved sodium was also higher. DIC was similar to Twelvemile, although δ¹³C values were approximately 4.5‰ more ¹³C-enriched. CO₂ was near equilibrium and CH₄ was oversaturated. The exposed shorelines examined in 2009 were organic and carbonate-rich lake muds. No evaporite crusts were observed. Submerged mats of Chara algae are abundant in addition to various pond weeds (Potamogetonaceae) and sedge (Scirpus sp.).

Isotope hydrology

Lake-water oxygen and hydrogen isotope ratios of YF lakes from Anderson et al. (2013) define a local evaporation line (LEL) with a slope of 4.89 (Figure 7). The intersection between the LEL with the global meteoric water line (GMWL) provides an estimate for the δ¹⁸O of mean annual precipitation (δ¹⁸O_p) between −23‰ and −24‰ that is similar to the δ¹⁸O values of the Yukon and Chandalar rivers and is also likely to reflect values for regional groundwater. A subset of YF lakes, which does not include Track and Twelvemile, plot below and parallel with the LEL. This group of lakes has a lower LEL-GMWL intersection value (less than −25‰) that reflects an isotopically depleted water source more dominated by snowmelt and/or permafrost thaw and active thermokarst.

Figure 7. Lake-water oxygen and hydrogen isotopes (a) from Twelvemile and Track Lakes (red diamond symbols) plotted with YF lake and river samples from Anderson et al. (2013) and the Global Meteoric Water Line (GMWL). Lakes sourced by precipitation, rivers, and groundwater are shown as either open or closed hydrology based on E/I ratios. Open hydrology lakes have E/I < 0.5 (square symbols) and closed hydrology lakes have E/I > 0.5 (cross symbols). Lakes sourced by snow and permafrost thaw are shown as a third group (open diamond symbols). A regression for lakes sourced by precipitation, rivers, and groundwater form a local evaporation line (LEL) with a GMWL intercept δ¹⁸O value of from −23‰ to −24‰. (b) Lake-water δ¹⁸O_H2O from Anderson et al. (2013) are plotted against dissolved sodium (Na⁺) concentrations

Lake-water enrichment in ¹⁸O and ²H isotopes along the LEL trajectory originating at the GMWL-LEL intersection reflects water loss by evaporation. Track and Twelvemile isotopic compositions on the LEL indicate that both lakes are dominantly sourced by precipitation delivered to the lakes via surface flow and groundwater. Using the isotope mass balance methods discussed in detail by Gibson, Birks, and Yi (2016), lake evaporation to inflow ratios (E/I) for Track and Twelvemile lakes during low stands are 1.29 and 1.12, respectively, indicating that the lakes had negative water budgets because of losses by evaporation (Anderson et al. 2013).

Thin crusts of white evaporite deposits surrounding lakes and within dry lake beds are observed throughout the YF, including the areas near Twelvemile and Track Lakes (Figure 1; Clautice and Mowatt 1981). The efflorescent salt films are less than 1 cm deep and are dominantly composed of calcite and dolomite, with some occurrences of trona (sodium bicarbonate). Brackish Mg-Na-HCO₃-Cl ponds in the YF form because evaporation exceeds inflow; notably, trona only occurs under extreme desiccation (Hawkins 1995; Heglund and Jones 2003). We note here a correspondence between elevated dissolved sodium with higher oxygen and hydrogen isotope compositions (Figure 7).

Recent lake-level fluctuations

Using Landsat imagery, Rover, Wylie, and Tieszen (2012) assessed YF lake-level trends between 1978 and 2012, including Twelvemile and Track Lakes, which we have supplemented here with additional aerial photography and Digital Globe WV02 imagery (Figures 2, 3, and 4; www.earthexplorer.usgs.gov). In 1953 and 1978, Twelvemile Lake’s levels were near overflow elevations. From 1985 to 2012, lake levels steadily declined by approximately 3 m, and the lake became hydrologically closed and occupied less than 25percent of its basin. Overbank flooding of the Yukon River caused by an ice jam in April 2013 renewed surface inflow, and by June 2013 lake levels had risen by approximately 1 m. The lake remained at this stage through April 2014. Between April 2014 and May 2016, lake levels steadily rose 2 m and reached overflow elevations.

Recent Track Lake trends were mostly similar to Twelvemile but with important differences. They were also near maximum between the 1950s and the early 1990s, but Track Lake’s levels subsequently fluctuated as Twelvemile’s steadily declined. Not until 2003 did Track steadily decline by approximately 3 m to occupy less than 25 percent of the basin by 2009. Both lakes remained near minima between 2009 and 2013, but Track’s levels did not rise in 2013. Between 2014 and 2016 both lakes rose, with Track’s levels rising by approximately 3 m to reach overflow elevations and representing a lake volume increase of approximately 75 percent.

As previously discussed by Anderson et al. (2013), lake-level declines observed throughout the YF during the early twenty-first century broadly correspond with climate trends. Higher temperatures, reduced precipitation, lower snow-water equivalent, and greater potential evaporation appear to have influenced lakes vulnerable to evaporative water loss, such as Twelvemile and Track. The rises at Track Lake when Twelvemile was in decline during the early 2000s are unexplained but illustrate the spatial heterogeneity that characterizes the YF region during drought. In contrast, during wet periods, lakes throughout the YF respond synchronously. Since 2014, high lake levels at both sites and the region as a whole correspond with increased summer and fall precipitation amounts.

Lake sediments

Calcite as a major phase of YF lake sediments was proposed by Hawkins (1995) because calcite precipitation drives water to reach the Mg-Na-CO₃ compositions required for trona formation. The Yukon Flat’s rivers are dominated by Ca and bicarbonate (HCO₃), and combined with bicarbonate formation in groundwater this leads to large alkalinity excess over calcite observed in YF’s lakes (Heglund and Jones 2003). In Track and Twelvemile, alkalinity is observed in excess over Ca and Mg, which suggests the possible occurrence of multiple authigenic carbonate phases, including calcite and aragonite. The small positive oxygen isotope fractionation between the polymorphs (+0.6‰) is relatively insignificant compared with the isotope variations and trends in the sediment records that are well over 2‰ (Shapley, Ito, and Donovan 2009 and references therein). Although dolomite was detected in evaporite crusts by Clautice and Mowatt (1981), it was not observed in lake sediments by Hawkins (1995) nor was it observed in Track and Twelvemile lakes.

Twelvemile Lake

Twelvemile basal sediments in core B08 have high magnetic susceptibility and bulk density, low percentages of OC and carbonate, and high rates of sedimentation with a maximum age of approximately 2700 cal yr BP (Figure 8). Sediments transition upcore to near equal parts OC (30–40%) and mineral matter (40–50%) with lower amounts of carbonate (5–25%). Basal sandy silts with similar ages were obtained from core C08 (Table 3). Magnetic susceptibility trends upcore to lower values as percentages of OC and carbonate trend to higher values. Nearshore core C08 has a similar but compressed sedimentary sequence as B10, with lower percentages of mineral material and higher percentages of OC and evidence for littoral sediment reworking (not shown). Oxygen isotope measurements from core B10 sediments produced prohibitively low CO₂ yields. Successful measurements were obtained from core GL07 to provide the 35 cm long Twelvemile δ¹⁸O_CaCO3 record that extends to approximately 600 cal yr BP, or approximately AD 1400 (Figure 8).

Figure 8. Core physical properties and geochemical properties including dry bulk density, magnetic susceptibility, weight percentages of calcium carbonate and organic carbon and mineral sediment, δ¹⁸O_CaCO3 and δ¹³C_CaCO3, C/N ratio, and organic carbon δ¹³C_Org of bulk sediment are plotted with the ²¹⁰Pb-horizon at AD 1850 and median calibrated radiocarbon ages from Twelvemile Lake

Isotopic equilibrium in Twelvemile Lake between the core-top carbonate and lake water was assessed using the calcite-water equilibrium equations of Epstein et al. (1953), Friedman and O’Niel (1977), and Tarutani, Clayton, and Mayeda (1969). Core-top δ¹⁸O_CaCO3 values of −11.2‰ and −12.5‰ are greater than 2‰ lower than the predicted equilibrium calcite values of −8.70‰ for the contemporaneously measured lake-water temperatures of 16.9°C. It is possible that endogenic carbonate production occurs in warmer littoral lake-water temperatures, which, if applied to the equations, slightly reduces the offset. It seems more likely that the offset is related to the fact that core-top δ¹⁸O_CaCO3 values integrate 2–7 years, which represent a wider range of lake-water temperatures and δ¹⁸O_H2O values. Equations by Leng et al. (2005) that use calcite-water isotope equilibria to predict a range of growing-season lake-water temperatures predict values between 15°C and 17°C, which are consistent with the measured lake-water temperature reported by Halm and Guldager (2012).

Twelvemile δ¹⁸O_CaCO3 values range from −10.5‰ to −13‰, and generally increase and show greater variability upcore. The δ¹³C_CaCO3 values range from −2.8‰ to 0.6‰, and generally decline before 600 cal yr BP until they rise within the uppermost 10 cm. The δ¹³C_Org values range from −27.5‰ to −30‰ and trend similarly upcore as δ¹³C_CaCO3. Higher δ¹³C_Org and δ¹⁵N_Org (not shown) may reflect some degree of active degradation of OC (Gälman et al. 2008). Relatively low C/N between 9.4 and 11.8 suggests dominantly aquatic sources. The δ¹⁸O_CaCO3 and δ¹³C_CaCO3 do not systematically covary prior to 1980 (Figure 9). After 1980, percentages of carbonate increased to 25 percent and percentages of OC declined to 40 percent as positive covariance between δ¹⁸O_CaCO3 and δ¹³C_CaCO3 strengthened (m = 0.87, R² = 0.73).

Figure 9. Twelvemile and Track δ¹⁸O_CaCO3 and δ¹³C_CaCO3 shown with linear regressions for the late twentieth and early twenty-first centuries

Comparison between historic Twelvemile Lake levels, sediment δ¹⁸O_CaCO3, and lake-water δ¹⁸O_H2O measurements indicate some consistent trends (Figure 4). Between 1978 and 2007, lake levels lowered and δ¹⁸O_CaCO3 values increased. High lake levels prior to 1978 are reflected by generally lower δ¹⁸O_CaCO3 values, although they are not the lowest values nor do they remain consistently low. Between 2007 and 2012, when lake levels were low, δ¹⁸O_H2O values remained high. After 2013, as lake levels rose, δ¹⁸O_H2O values declined by approximately 7.5‰.

Track Lake

Track Lake basal sediments are laminated dark gray brown silts of high magnetic susceptibility and density and low percentages of OC and carbonate (Figure 10). Sediments transition upcore to olive brown marl that is high in OC and carbonate (>70%), with lesser amounts of biogenic silica (10–20%) and mineral matter (<5%). Higher percentages of OC and carbonate occur as bulk density and magnetic susceptibility decline. Changes in the relative proportions between OC and carbonate are visually evident as centimeter-scale color banding. Darker shades of brown reflect higher OC, and lighter shades of olive reflect higher carbonate abundance. Aquatic macrofossils are abundant throughout the core, whereas terrestrial macrofossils were rarely found.

Figure 10. Core lithology, physical properties, and geochemical properties including dry bulk density, magnetic susceptibility, weight percentages of calcium carbonate and organic carbon, δ¹⁸O_CaCO3 and δ¹³C_CaCO3, C/N ratio, and organic carbon δ¹³C_Org of bulk sediment are plotted with the ²¹⁰Pb-horizon at AD 1850 and median calibrated radiocarbon ages from Track Lake

Living Chara is abundant in Track Lake, and isotopic equilibrium between living Chara calcite encrustation values δ¹⁸O (δ¹⁸O_Chara) and lake water at measured lake-water temperatures was assessed using the same equations that were described earlier and used for Twelvemile Lake. An estimated δ¹⁸O_Chara value of −7.92‰ at the lake’s measured water temperature of 19.47°C is similar to the measured δ¹⁸O_Chara value of −8.55‰. Track δ¹⁸O_CaCO3 values for core-top fine-grained sediment were lower than living Chara, between −11‰ and −13‰, but reflect a time-integrated signal and fall within the range of calcite-water isotopic equilibrium values predicted by the equations of Leng et al. (2005). Successful δ¹⁸O_CaCO3 measurements for cores A09/B10 provide the approximately 120 cm long record that extends to approximately 5500 cal yr BP.

Track δ¹⁸O_CaCO3 values vary between −10.5‰ and −15.5‰ and are shown with age-model uncertainty in Figure 11 as a proxy.ghost graph in Bacon. Following a decline of approximately 5‰ within basal sediments between approximately 5500 and 4000 cal yr BP, values gradually rose by approximately 2‰ between approximately 3500 and 1850 cal yr BP. Because of the gap in recovered sediment between the two cores, no data are available for the period between approximately 1850 and 800 cal yr BP. Values between approximately 800 cal yr BP and present were relatively high, between −13.5‰ and −15.5‰, with the lowest values of the period having occurred during the twentieth century. Comparison between lake-level fluctuations and the sediment δ¹⁸O_CaCO3 record indicate consistent patterns (Figure 4). From 1978 to the mid-1990s, δ¹⁸O_CaCO3 values were relatively low when lake levels were high. After approximately 1998, δ¹⁸O_CaCO3 values increased as lake levels declined.

Figure 11. Track δ¹⁸O_CaCO3 data plotted as a grayscale with Bacon v2.2 to display chronological uncertainties: the darker the gray, the greater the age certainty

Track δ¹³C_CaCO3 values range between 2‰ and −1.5‰. Values for δ¹³C_Org range between −32‰ and −27‰, and with C/N, increase upcore until 10 cm depth. The recent decrease in C/N in the upper 10 cm may reflect active OC degradation and possibly the Seuss effect, and may explain the recent decline in δ¹³C_Org. Covariance between δ¹⁸O_CaCO3 and δ¹³C_CaCO3 evolves upcore (Figure 9) from insignificance in basal sediments to a slightly negative relation between approximately 4000 and 1850 cal yr BP, and followed by strong positive covariance from AD 1945 to the core top (R² = 0.6 and m = 0.3).

Modeling climate-driven changes in lake water isotopes

Modeling lake-water isotopes provides helpful constraints on sediment oxygen isotope records even when knowledge of the hydrologic system is incomplete and true model calibration is not possible (Anderson, Finney, and Shapley 2011; Shapley, Ito, and Donovan 2008; Steinman and Abbott 2013). In this case we developed a model of Track Lake δ¹⁸O_H2O and conducted highly idealized experiments designed to explore plausible hydrologic controls that may explain prominent sediment δ¹⁸O_CaCO3 variations. Model runs were designed to test two specific questions: Within the limits of twentieth-century climate variability, is it possible to simulate (1) a decadal scale −3‰ excursion in δ¹⁸O_H2O and (2) the entire range of Track Lake δ¹⁸O_H2O values of 6‰?

The model is adapted from that of Benson and Paillet (2002) and is coded in Visual Basic in Microsoft Excel^TM, incorporating hypsometry, groundwater exchange, evaporation and precipitation fluxes, seasonal stratification, and precipitation seasonality, and is described in detail by Shapley, Ito, and Donovan (2008). Refinements for this study included seasonal snowpack and melt to simulate the effect of precipitation storage and seasonally pulsed runoff. Input climatic data are from the National Weather Service COOP station and SNOTEL at Fort Yukon. The nearest available monthly oxygen isotope compositions of precipitation (δ¹⁸O_p) are from Fairbanks, collected between 2007 and 2011 by B. Finney. Low δ¹⁸O_p values range from −25‰ to −30‰ in December and high values range from −11‰ to 18‰ from May to July. Monthly potential evaporation input is based on Thornthwaite estimates derived from NCEP/NCAR reanalysis climate parameters from Anderson, Finney, and Shapley (2011). Model performance was evaluated by comparison with modern lake volumes and measured δ¹⁸O_H2O values following adjustments of the most weakly constrained fluxes, such as monthly groundwater and evaporation rates (Figure 12). Modeled steady-state δ¹⁸O_H2O are within less than 1 percent of observed, and modeled volumes are within the range observed.

Figure 12. Time series of modeled monthly water fluxes (a) for Track Lake at steady state and (b) lake volume and isotope response to a step-change in climate from dry to wet conditions by +25 percent precipitation and −25 percent evaporation

To address question 1, step changes in a variety of water-balance parameters within the range of multiannual variability observed during the twentieth century were imposed, and time series of modeled lake volume and δ¹⁸O_H2O was evaluated. Results show that a 25 percent increase in precipitation and a 25 percent decrease in evaporation produce a sustained increase in lake volume over two decades (Figure 12). The δ¹⁸O_H2O response was twofold. First there was an approximately 3‰ decline within five years of the step shift in climate that was then followed by partial recovery over approximately thirty years to represent an overall 1.7‰ decline. This simulation suggests that the negative δ¹⁸O_CaCO3 excursion reflected by the sediment record during the early twentieth century could plausibly reflect a similar transient isotope response.

To address question 2, six model runs were carried out with combinations of twentieth-century climate extremes in precipitation amount and seasonal balances (Table 6). Simulated lake volumes ranged from approximately 1.4 to 6.1 × 10⁵ m³ and δ¹⁸O_H2O variations ranged from −3.7‰ to 10.7‰, which are similar to observed lake-volume variations of approximately 2.6 to 6.7 × 10⁵ m³ and sediment δ¹⁸O_Calcite values of −7‰ to −14.5‰. The results indicate that the entire 6‰ range in δ¹⁸O_H2O values recorded by the sediment data plausibly reflect steady-state responses to forcing by twentieth-century hydroclimate extremes.

Table 6. Model simulation experiments and lake volume and δ¹⁸O sensitivity

Simulation Scenario^1	Precipitation (cm)	Snowfall Proportion (%)	Evaporation (cm)	Lake Volume^2 (× 10^5 m^3)	Lake δ^18O (‰)	Lake δ^18O Difference from Base
Base	17	35	40	3.0	−6.5
Base precipitation, snow dominated	17	80	40	4.1	−8.7	−2.2
Wet precipitation, snow dominated	25	80	40	5.9	−10.7	−4.3
Dry precipitation, snow dominated	12	80	40	2.8	−8.0	−1.5
Base precipitation, rain dominated	17	20	40	2.8	−5.8	+0.6
Wet precipitation, rain dominated	25	20	40	4.9	−6.4	+0.1
Dry precipitation, rain dominated	12	20	40	2.0	−5.2	+1.3
Wet precipitation, snow dominated, low evaporation	25	80	30	6.1	−12.8	−6.4
Dry precipitation, rain dominated, high evaporation	12	20	55	1.4	−3.7	+2.8

¹High groundwater outflow model configuration with 0.15 surface runoff factor to maintain closed surface hydrology.

²Average of seasonally varying steady-state volume.

Discussion

The sediment cores from Track and Twelvemile lakes indicate that sedimentation following initial lake formation is low in lithic input and is dominated by endogenic carbonate and aquatic OC, which in turn reflects high levels of aquatic productivity and excess alkalinity. Bioinduced carbonate precipitation occurs within the water column as endogenic crystals and Chara encrustations. The data indicate that carbonate sediments are in isotopic equilibrium with lake water, which is consistent with a predominantly endogenic origin and insignificant contribution by detrital carbonate sources. Continuous sedimentation is indicated by the age models, which along with undisturbed sediment structures of aquatic origin, rules out the possibility of significant sediment reworking caused by active thermokarst. The simple lake bathymetries reflect stable taliks, which foster the groundwater–surface water connectivity indicated by the lake-water isotope compositions. The lakes evolve volumetrically and isotopically in response to evaporation and hydroclimate, which is unambiguously supported by ubiquitous lake-level lowering and rise in response to variations in aridity and precipitation. Although some changes in shallow-depth permafrost may have occurred at Twelvemile and Track lakes, the sediment core data indicate that this exerts a relatively minor influence on the hydrology compared to the rapid and readily detectable hydrologic response to changes in hydroclimate.

Oxygen isotope interpretation

The Track Lake isotope model illustrates that δ¹⁸O_H2O of lake water is influenced by watershed hydrology and hydroclimatic conditions, and is controlled by both (1) the combined δ¹⁸O of input waters from catchment runoff, precipitation falling directly on the lake, and groundwater inflow and (2) ¹⁸O enrichment caused by preferential evaporation of the ¹⁶O isotopes in lighter water molecules. Open lake systems with large catchment areas have low E/I, and the values and variations of δ¹⁸O_H2O reflect precipitation, δ¹⁸O_p. For closed-lake systems with small catchment areas and higher E/I, the values of δ¹⁸O_H2O are dominantly influenced by the isotopic effects caused by evaporation, which reflects effective moisture. Furthermore, for the closed systems presented here, the lake-water isotopic modification of approximately 10‰ in response to evaporation (Figure 7) far exceeds the ±0.5–1.5‰ variations that could potentially be attributed to changes in the relative proportions of calcite-aragonite mineralogy and/or lake-water temperatures.

When conditions are arid, high rates of evaporation lead to closed surface-water hydrology at Twelvemile and Track Lakes. High E/I and ¹⁸O-enriched lake water lead to high δ¹⁸O_CaCO3 values. Variations in δ¹⁸O_CaCO3 during arid periods when the lakes are closed to surface outflow are sensitive to changes in evaporation and effective moisture. When conditions are wet, or in the case of Twelvemile, the lake is affected by Yukon River floods, and the lakes are open to surface inflow they have lower E/I and more ¹⁸O-depleted lake water, which results in low δ¹⁸O_CaCO3 values. If the lakes reach overflow elevations and become open to surface inflow and outflow, the lakes become less sensitive to evaporation and δ¹⁸O_H2O values are more sensitive to changes in δ¹⁸O_p. Recent observations and isotope modeling suggest that Track and Twelvemile lakes are most often closed to surface outflow. The isotope model indicates that transition from a closed to open hydrology at both lakes is reflected by relatively large and rapid negative δ¹⁸O_H2O shifts of greater than 4–5‰.

Twelvemile Lake

Basal sediment ages indicate sand and silt deposition in the Twelvemile basin by fluvial processes at approximately 2500 cal yr BP, a result that is consistent with the northward lateral migration of the Yukon River across the Twelvemile basin to its present position, as suggested by three-dimensional geophysical models of permafrost (Minsley et al. 2012). Estimated downward freezing rates following northward migration suggest that the Twelvemile basin was isolated from the river by 900–1,800 years ago, which is broadly consistent with the onset of organic-rich lake sedimentation documented here by approximately 1500 cal yr BP.

Between approximately AD 1500 and 1700, the oxygen isotope record indicates a lake where δ¹⁸O_CaCO3 values varied independently of δ¹³C_CaCO3. Weak covariance between lake-water oxygen reflected by δ¹⁸O_CaCO3 and dissolved inorganic carbon isotopes reflected by δ¹³C_CaCO3 is typical when lake volumes remain stable (Li and Ku 1997). This is consistent with more frequent phases or events when Twelevemile hydrology is open to surface inflow in response to wetter climate and more frequent incursions of Yukon River floodwaters (Figure 8). The transition to higher δ¹⁸O_CaCO3 values by approximately AD 1800 would therefore reflect more frequent closure to surface inflow and outflow and indicate that conditions were drier. Since approximately AD 1800, higher and more variable δ¹⁸O_CaCO3 values reflect either a generally more arid climate with periodic wet periods or less frequent river incursions (Figure 13). High-frequency declines in δ¹⁸O_CaCO3 values during the 1950s and 1980s reflects observed wet periods. Therefore, we infer that high δ¹⁸O_CaCO3 values during the late 1800s and early 1900s reflect more arid conditions than during the observational period from the mid-twentieth century to the present.

Figure 13. Comparison between Track and Twevelmile δ¹⁸O_CaCO3 from AD 1850 to the present, with the summer temperature anomaly reconstruction at Firth River, northeast Alaska, from Picea glauca maximum latewood density, with uncertainty (Anchukaitis, D’Arrigo, and Andreu-Hayles 2013)

Track Lake

High δ¹⁸O_CaCO3 values that steadily declined between approximately 5500 and 4000 cal yr BP may reflect a period of significantly greater aridity relative to the present that trended toward wetter conditions (Figure 6). Basal sediments have high mineral silt, low OC, and relatively low OC δ¹³C values (δ¹³C_org) between −32‰ and −30‰ (Bade and Carpenter 2004). Our core may not have penetrated the full record, and the basal sediments, which are not structurally characteristic of thermokarst, may reflect a shallow lake with higher silt input from exposed shoreline erosion and eolian processes that are consistent with an arid climate. Low dissolved carbon δ¹³C values as reflected by organic carbon δ¹³C values (δ¹³C_org) can reflect respiration in groundwater sources and high rates of sediment organic matter decomposition (Finney et al. 2012). However, we note that declines in δ¹⁸O_CaCO3, if interpreted to reflect rising lake levels and an increase in moisture, appear to be in conflict with oxygen isotope records from non-thermokarst lakes in the Brooks Range, which indicates a region-wide dry period between 5 ka and 4 ka (Clegg and Hu 2010; Mann, Heiser, and Finney 2002).

Alternatively, the carbonate isotope evidence raises the possibility that basal Track sediments may reflect thermokarst processes in addition to hydroclimate. The relatively low carbonate δ¹³C values (δ¹³C_CaCO3) that were decoupled from high δ¹⁸O_CaCO3 values (Figure 9) suggest that a smaller and shallower lake developed with ¹³C-depleted dissolved carbon sources, which could have also been a product of methanogenic microorganisms (Bastiviken et al. 2003). In thermokarst settings, relatively low dissolved carbon δ¹³C values are associated with CO₂ contributions by higher rates of respiration and methanogenesis from degradation of permafrost carbon (Heslop et al. 2015; Sepulvida-Jauregui et al. 2015; Walter Anthony et al. 2014; Wooller et al. 2012). The highly elevated δ¹⁸O_CaCO3 values are most consistent with an early lake that was a small, highly evaporated water body, and if thermokarst was active the δ¹⁸O_CaCO3 declines between approximately 5.5 ka and 4 ka could reflect an increasing supply of ¹⁸O-depleted water that was potentially sourced from thawed permafrost during talik expansion (Figure 7). Unraveling the respective contributions by thermokarst and hydroclimate during the initiation of the Track Lake basin with more certainty will require additional investigation.

The development of endogenic marl sedimentation in Track Lake by approximately 4000 cal yr BP supports hydroclimate as the subsequent primary control on lake-water balance, which is further supported by the isotope model simulations (Figures 13 and 14). Low δ¹⁸O_CaCO3 values between approximately 4000 and 1850 cal yr BP indicate wetter conditions as relatively high lake levels and more frequently open hydrology. The subsequent δ¹⁸O_CaCO3 trend to gradually higher values suggests a multicentury trend to lower lake levels, which are inferred to reflect more frequent or more sustained arid periods. Persistently high δ¹⁸O_CaCO3 values since 800 cal yr BP indicate frequently closed hydrology in response to lower lake levels and even greater aridity. The prominently low δ¹⁸O_CaCO3 values during the 1950s through the 1980s reflect the observed multidecadal wet period, which the isotope modeling suggests may be an ongoing transient isotope response. The wet early twentieth century and wet conditions from 2014 to the present appear wetter than during the previous approximately 800 years. Below the gap in our record, similarly low δ¹⁸O_CaCO3 values were not previously recorded until approximately 2000 cal yr BP.

Figure 14. Comparison between the Track δ¹⁸O_CaCO3 with records of YF fire (Kelly et al. 2013), Alaska glacier activity (Barclay, Wiles, and Calkin 2009), Tangled Up Lake δ¹⁸O_Chara (Anderson, Abbott, and Finney 2001), Blue Lake varve thickness (Bird et al. 2009), and Firth River maximum late wood density summer temperature reconstruction (Anchukaitis, D’Arrigo, and Andreu-Hayles 2013). The normalized YF charcoal accumulation rate (CHAR) compilation shows trends on millennial (gray) and century (black) time scales. The red box around the Track Lake time series indicates the section of the sediment record that may have been influenced by thermokarst processes. The downward red arrow on the Tangled Up record indicates potential trends in regional δ¹⁸O_p discussed in the text

The role of precipitation oxygen isotopes

Although evaporative evolution is the major influence on Track Lake δ¹⁸O_H2O, the isotope model illustrates that isotopes in precipitation (δ¹⁸O_p), as stored snow and surface runoff, are also important (Figure 12). Variations in δ¹⁸O_p values in southern Alaska during the Holocene are known to be strongly influenced by changes in vapor source region and atmospheric circulation over the North Pacific in addition to the effects of temperature on condensation, topography, and continentality (Anderson et al. 2016). In contrast, variations in δ¹⁸O_p values in northern Alaska during the Holocene are currently unknown. Limited modern observations in Barrow, North Slope in situ vapor monitoring and surface water surveys, and measurements of glacier ice at McCall Glacier in the Brooks Range indicate the influence of Arctic Ocean vapor sources, which are in turn influenced by sea ice (Klein et al. 2015, 2016; Lachniet et al. 2016). Values of mean annual δ¹⁸O_p currently appear to range from approximately −20‰ to −23‰, with lower values reflective of colder temperatures at the vapor source area and at the time of condensation, longer transport, and more winter precipitation. Higher values reflect warmer temperatures at the vapor source area and at the time of condensation, shorter transport, and more summer precipitation.

Similar lake-water δ¹⁸O_H20 values were observed at Tangled Up Lake, located in the central Brooks Range, by Anderson, Abbott, and Finney (2001; Figure 1). Tangled Up is hydrologically open to surface inflow and outflow, and it was recognized that changes in the values of lake-water δ¹⁸O_H2O reflect changes in the values of δ¹⁸O_p. Although Anderson, Abbott, and Finney (2001) had suggested the temperature of condensation as a dominant control on δ¹⁸O_p values, we now recognize that vapor source and seasonality are likely as important as temperature and that the data are better understood and utilized as a proxy for δ¹⁸O_p.

The Tangled Up δ¹⁸O_p record is relatively coarse, with large age uncertainties (~250 yr), but data are continuous between 2000 cal yr BP and AD 1997 (Figure 14). Between 2000 and 1000 cal yr BP, values of δ¹⁸O_p varied considerably but declined overall from −20‰ to −23‰. Values subsequently fluctuated on the order of 3‰ from 1000 cal yr BP to the present. The implication for this study is that if a 3‰ decline in δ¹⁸O_p values between 2000 and 1000 cal yr BP also occurred in the Yukon Flats, then the trend of increasing δ¹⁸O_CaCO3 beginning in the upper portion of core B10, which appears to continue across the gap in the sediment record into the lower portion of core A09, reflects the effects of evaporative enrichment and the decline in δ¹⁸O_p and thus a relatively large Δδ¹⁸O decline. It is noteworthy that similar trends in Δδ¹⁸O during the late Holocene have also been observed for the interior Yukon Territory (Anderson et al. 2007). A better understanding of changes in precipitation isotope values across northern Alaska to further investigate this preliminary assessment will require the development of additional δ¹⁸O_p proxy records.

Late Holocene hydroclimate

The YF lake-level reconstructions based on Track Lake and Twelvemile Lake oxygen isotope data show hydroclimatic variations at time scales that range from millennia to decades for the past approximately 4,000 years. Furthermore, the lake levels reflect a balance between precipitation and evaporation, which in the YF region are both summer dominated. The millennial-scale hydroclimate trend at Track Lake is notably consistent with relatively low YF charcoal accumulation rates (CHAR) between approximately 4000 and 2000 cal yr BP and fire frequency (not shown), which gradually rise to a peak by 800 cal yr BP, based on a composite of fourteen lakes across southern marginal uplands (Kelly et al. 2013). Although millennial-scale increase in fire has been previously associated with the expansion of P. mariana (Higuera et al. 2009) and temperature trends (Clegg et al. 2011), Track Lake’s level reconstruction suggests a prominent role for effective moisture.

Higher effective moisture indicated by Track Lake between approximately 4 ka and 1.5 ka is broadly consistent with proxy records throughout the region (Figure 1) that indicate wetter conditions for sustained intervals during the same period, although with some subregional differences. Similarities include thicker varves at Blue Lake in the central Brooks Range interpreted to reflect a wetter regime from the end of the record at approximately 2 ka to 1.2 ka (Bird et al. 2009; Figure 14). Relatively low δ¹⁸O_CaCO3 values between approximately 4 ka and 2.5 ka at Takahula Lake in the south-central Brooks range are interpreted to reflect higher effective moisture (Clegg and Hu 2010). Farther west, relatively wet conditions between approximately 3.5 ka and 1.5 ka are indicated by reduced sand influx at Wolverine Lake in the Kobuk Sand Dunes (Mann, Heiser, and Finney 2002). Keche Lake, in the eastern Brooks Range, indicates more frequent periods of detrital carbonate between approximately 3.5 ka and 1.5 ka, which can be attributed to more surface runoff within the watershed due to wetter conditions (Chipman, Clegg, and Hu 2012). Marcella Lake, in the southwest Yukon, has relatively low δ¹⁸O_CaCO3 values between approximately 3 ka and 1.5 ka, which indicate higher effective moisture (Anderson et al. 2007).

The increase in YF CHAR between approximately 2000 and 800 cal yr BP suggests that greater biomass burning was coincident with higher aridity indicated by lower Track Lake levels. Although Track Lake data are not available for the period between approximately 1850 and 800 cal yr BP, a regional transition to drier conditions between approximately 1000 and 800 cal yr BP is reflected by a reduction in Blue Lake varve thickness across the period (Figure 14) and increased δ¹⁸O_CaCO3 values at Sevenmile and Marcella lakes (Anderson, Finney, and Shapley 2011). Anderson, Finney, and Shapley (2011) discussed aridity and “early Little Ice Age” glacial advances across Alaska’s major mountain ranges between approximately AD 1050 and 1300 (Figure 14; Barclay, Wiles, and Calkin 2009; Wiles et al. 2008). Track Lake levels were near modern with only minor variation between AD 1200 (800 cal yr BP) and 1850 and during the inception of the second major phase of Alaskan glacier advances, known as the “late Little Ice Age,” by approximately AD 1600.

Hydroclimatic trends during recent centuries indicate a cool-dry period between AD 1850 and 1900, as reflected by low Track Lake levels and relatively cool summer temperatures indicated by the maximum late wood density reconstruction from Firth River, located in the eastern Brooks Range (Figure 13; Anchukaitis, D’Arrigo, and Andreu-Hayles 2013). The first half of the twentieth century was increasingly warmer and wetter as reflected by a rise in lake levels at Track and Twelvemile and higher summer temperatures. The long-term lake-level data indicate that the observational period since the mid-twentieth century was generally wetter than earlier decades and the previous 800 years. However, the past 800 years were generally drier than the preceding 3,000 years, with the qualification that we cannot assess the period between approximately 1850 and 800 cal yr BP. The most recent time prior to the 20th century when conditions were as wet as the late twentieth century may have been approximately 2000 cal yr BP.

Conclusions

Paleoclimate histories that include surface hydrologic response provide a means to better understand ecosystem change and better inform land management. The new Track and Twelvemile sediment records suggest that the magnitude and frequency of lake-level fluctuations observed during the late twentieth and early twenty-first centuries are a persistent feature of the YF region. With these records, we can now place recent YF lake-level changes of the past few decades within a number of longer term contexts, which indicate that these lakes have previously experienced similar and more extreme hydroclimatic states during the late Holocene. Furthermore, the isotope model suggests that lake-level extremes, based on the entire isotopic range of sediment data, are plausibly driven by the range of twentieth-century hydroclimatic variations and extremes. Thus, the long-term perspectives developed here serve to highlight the resilient capacities of aquatic ecosystems in response to fluctuating lake levels (e.g., Lewis et al. 2014; Roach et al. 2011) and reinforce the notion proposed by Kelly et al. (2013) that the persistence of highly variable drought-induced fire regimes into the future is likely.

Future lake-level fluctuations will also continue to depend on permafrost stability, which reflects complex interactions among hydroclimate, topography, hydrology, soil, and vegetation (Jorgenson et al. 2010). However, while there is no evidence from this study to refute the proposition that future sustained lake-level declines are more likely to occur in response to permafrost loss, which boreal forest burning could promote (Jorgenson and Osterkamp 2005), consideration of a broad range of temporal perspectives is instructive. In this case, the longer perspectives on YF boreal forest, fire, and hydroclimate indicate that permafrost has remained relatively stable despite many centuries of high-severity fire and drought.

It is interesting to note the observation that shorter twentieth-century perspectives of YF fire frequency do not indicate clear relationships with climate (Drury and Grissom 2008), whereas millennial-scale reconstructions have a strong correspondence with the paleoclimatic data presented here. This conundrum is a consequence of persistent high-amplitude variations on multiannual to decadal time scales that will invariably plague the detection of long-term change based on short-term perspectives. As this study has shown, YF lake-level changes also show high-magnitude change and complex spatial responses on a range of timescales that should be assessed within a late twentieth- and early twenty-first-century context with caution. This study has shown that by framing recent lake-level fluctuations within a longer term context of multiple centuries to millennia, we can more confidently identify the occurrence of sustained trends, which will be essential to understand in terms of future impacts on wetland ecosystems, permafrost, and carbon.

Acknowledgments

This research is supported by the U.S. Geological Survey Land Change Science Research and Development program. Nikki Guldager and the Yukon Flats National Wildlife Refuge provided air transportation and field logistics. We thank David Dettman and Chris Eastoe at the University of Arizona, who provided carbonate isotope measurements, and Dan Engstrom of the Science Museum of Minnesota, who provided radiometric dating and age models. We thank Chad Wolack of the University of Colorado Institute for Arctic and Alpine Research, who provided macrofossil graphitization and AMS ¹⁴C ages, and Carl Jurkowski for his assistance. Paco VanSistine provided imagery and ARC-GIS; Colin Penn provided coulometry measurements. Jeremy Havens provided assistance with figures. We are grateful to Miriam Jones, Mary Edwards, and an anonymous reviewer for constructive criticisms and suggestions that served to significantly improve the manuscript. Any use of trade, firm, or product names is for descriptive purposes only and does not imply endorsement by the U.S. government.

Disclosure statement

No potential conflict of interest was reported by the authors.