A volumetric water budget of Devils Lake (USA): non-stationary precipitation–runoff relationships in an amplifier terminal lake

ABSTRACT

Devils Lake, a terminal lake in eastern North Dakota, rose more than 9 m between 1992 and 2013, producing a 286% increase in lake area, and causing more than US$1 billion in direct damages. An annual volumetric lake water budget is developed from monthly hydroclimatological variables for the period 1951–2010 to investigate the rapid lake expansion. The lake is an amplifier terminal lake in which long-term climatic changes are amplified by positive feedback mechanisms, causing the lake to transition from a precipitation-dominated to a runoff-dominated water budget. Factors specific to the Devils Lake Basin further amplify this positive feedback relationship. These include principles of fill–spill hydrology that operate between individual sub-basins within the closed basin, and between the innumerable wetland complexes within each sub-basin. These factors create a pronounced non-stationary precipitation–runoff relationship in the basin during both long-term wetting and drying phases.

KEYWORDS:

• Devils Lake
• terminal lake
• volumetric lake water budget
• amplifier lake
• non-stationarity
• precipitation–runoff relationship

Editor A. Castellarin Associate editor H. Lins

1 Introduction

A large portion of the surface of Earth drains to inland lowland depressions that form terminal lakes. These closed basin lakes are normally saline and are often referred to as inland saline (salt) lakes. The total volume of water contained within saline lakes is comparable to that of freshwater lakes, and they form an important global water resource (Williams 1996, 2002).

Terminal lakes are dynamic hydrological systems, and their record of historical lake level fluctuations is considered a sensitive indicator of climate change (Mason et al. 1994). Such lakes can experience large and rapid changes in lake level, lake area, and lake volume due to the integrated effects of multiple climatic forcings (Szesztay 1974, Mason et al. 1994). Lake level fluctuations are driven by short-term and long-term variations in climate, and also by lake bathymetry such that the lake area–lake volume relationship of a saline lake can have preferred states in lake level (Street 1980, Mohammed and Tarboton 2011, Haghighi et al. 2016). Terminal lake level fluctuation is also extremely sensitive to human modification of basin hydrology, including human abstraction or augmentation of surface or groundwater inflow, diversion of tributaries, and alteration of basin land use/land cover (Williams 1996, 2002, Coe and Foley 2001, Micklin 2007).

Devils Lake is a terminal saline lake in the Devils Lake closed basin in eastern North Dakota (Fig. 1). The closed basin has an area of 9868 km², of which 8599 km² drains to Devils Lake, and 1269 km² drains to Stump Lake, an adjacent saline lake to the east. The rise in lake level, and increase in lake area and lake volume of Devils Lake has been large and sudden since 1993 (Fig. 2). From 30 September 1992 to 30 September 2013 the lake level (L_L, m) rose from 433.62 to 442.75 m (ΔL_L = 9.13 m), lake area (L_A, km²) expanded from 179.9 to 694.6 km² (ΔL_A = 286%), and lake volume (L_V, km³) increased from 0.704 to 4.003 km³ (ΔL_V = 469%) (Fig. 2). Devils Lake began spilling into Stump Lake in 1999, and the two lakes now form a single lake system. At a lake level of 444.40 m Devils Lake spills naturally through the Tolna Coulee into the Sheyenne River, and becomes part of the Red River of the North drainage basin.

Figure 1. Location of Devils Lake Basin, North Dakota (USA). Lake extent shown is approx. 1999. Large, lake open water surfaces are shown in black; small, wetland open water surfaces are shown in grey. Source: Sethre et al. (2005).

Figure 2. Lake volume (L_V, km³) time series for Devils Lake; lake volumes are shown to 31 January 2014. Source: USGS (2013, 2014).

The lake has experienced catastrophic flooding over the past two decades (Larson 2012). The lake flooding that resulted required the raising and extension of levees to protect the City of Devils Lake, the raising of roads and bridges to preserve critical infrastructure, the destruction or relocation of more than 500 structures, and the inundation of more than 60 000 hectares of prime agricultural farmland (Zheng et al. 2014). Direct federal and state payments for flood damages have exceeded US$1 billion, although the total direct/indirect and tangible/intangible flood losses are certainly much higher (Smith and Katz 2013).

The abrupt and dramatic change in the hydrological regime at Devils Lake has not been fully explained. Basic hydrological investigations by Ryan and Wiche (1988), Wiche (1992), Wiche and Pusc (1994), and Vecchia (2008) provide an overview of the lake hydrology, while Hoerling et al. (2010) present a preliminary investigation of the climatological basis of the recent lake rise. Todhunter and Knish (2014) examined the relationship between the frequency of weather types and lake volume changes in the period 1965–2010, and found statistically significant changes in the frequency of selected weather types over time that indicated a trend toward increased advection of more humid weather types consistent with the historical rise in lake level. Their results and those of Hoerling et al. (2010) support a climatic explanation for the historical lake rise at Devils Lake, but the results are not as strong nor as profound as might have been anticipated given the unprecedented rate of lake rise.

The purpose of this study is to investigate the hydrological basis of the sudden rise in Devils Lake. The first objective is to construct an annual volumetric water budget for the lake for water years 1951–2010 that considers both hydroclimatic variability and the effects of lake bathymetry. The second objective is to quantify the relationship between precipitation input to the basin and runoff inflow to the lake for the same study period.

2 Background

2.1 Study area

Devils Lake Basin lies within the Drift Prairie physiographic region which was traversed by a continental ice sheet during the last glaciation (Bluemle and Clayton 1984). Devils Lake includes a series of separate bays arranged in a linear sequence, separated from one another by a swale of higher elevation, and underlain by the Spiritwood Aquifer. The bays were formed by glacial thrusting of subglacial materials caused by elevated pore water pressures beneath the advancing ice sheet (Bluemle and Clayton 1984).

The closed basin includes a series of sub-basins, each drained by a major coulee that flows south and empties into a small lake (Ryan and Wiche 1988). These sub-basin lakes form a natural chain of lakes; each lake is located within a glacial depression, separated from each of the adjacent lakes in the chain by a swale with a unique spill elevation. Runoff from one sub-basin and lake does not directly contribute to the next lower sub-basin until the swale elevation separating the two sub-basins is reached and the lake can spill into the lower lake in the chain. Runoff along the chain of lakes proceeds in an east to west direction, and enters into northwest Devils Lake through Big Coulee. Prior to 1979, all runoff produced within the chain of lakes moved east to west through Big Coulee and into Devils Lake. In 1979 a levee was constructed along the natural outlet from Dry Lake, and Channel A was constructed connecting Dry Lake to Sixmile Bay on Devils Lake. Detailed maps of all channels, sub-basins, and the chain of lakes are available in Wiche and Pusc (1994) and from the North Dakota State Water Commission.¹

Surface runoff now enters Devils Lake through three sources: the eastern sub-basins that enter Devils Lake through Channel A, the western sub-basins that enter Devils Lake through Big Coulee, and a limited amount of local runoff from areas adjacent to Devils Lake. The upper basin lands consist of clay-rich surficial deposits underlain by Cretaceous shale; both possess very low permeability, and groundwater flow rates are extremely small (Wiche and Pusc 1994).

A great number of topographically closed depressions of variable size and depth were left by the continental ice sheet in the prairie landscape (Zhang et al. 2009, Shaw et al. 2013). Water storage in these prairie potholes is normally highest in the spring following the melt of the winter snow cover (Hayashi et al. 2016). A large portion of the upper basin drains internally to these pothole depressions, which frequently do not fill to a high enough level to spill to adjacent lower level pothole depressions. Surface water drainage between adjacent depressions is thus normally disconnected. The region lacks a well-defined surface drainage system, and the fraction of the area contributing to the surface stream network varies on a seasonal and inter-annual basis. The volume of water stored within prairie potholes varies dynamically across time, which results in great temporal variation in the number, average depth, and permanence of prairie potholes in the landscape. Pothole abundance, depth, and permanence are greatest in the early spring and during deluge periods, and decrease through summer and during drought periods (Zhang et al. 2009). Regional runoff totals are primarily driven by snowmelt and rain-on-snow processes (Berghuijs et al. 2016).

This complex natural hydrological setting has been extensively modified by human action over time. Extensive historical land cover change from native prairie grasses to cultivated crops has occurred, which would have resulted in decreased surface infiltration and increased surface runoff (Van Der Kamp et al. 2003). The region has experienced moderate land cover change between 1973 and 2000 (Auch et al. 2011, Drummond et al. 2012), including extensive wetland drainage (Johnson 2013).

2.2 Late Holocene and historical lake level fluctuations

Geological field investigations show that within the past 4000 years Devils Lake has gone completely dry four times, and spilled over to the Red River of the North three times (Bluemle et al. 1999). The natural condition for the lake is to be rising toward overflow to the neighbouring Stump Lake or Sheyenne River, or falling toward desiccation as a dry lake bed (Bluemle et al. 1999). Widely fluctuating lake levels are the normal condition for the lake, while long periods of stable lake levels are uncommon. Multiple paleolimnological studies of Holocene climate variability from northern Great Plains saline lakes show a similar pattern of lake level instability (Laird et al. 2003, Shapley et al. 2005, Grimm et al. 2011). These studies reveal a coherent regional pattern during the late Holocene of large oscillations between longer periods of high salinity–dry conditions, and shorter periods of low salinity–wet conditions. The transitions between these two modes appear to be relatively rapid, suggesting large and abrupt changes in the long-term hydroclimatic drivers controlling surface moisture availability, and rapid response of the lake system to such changing inputs (Laird et al. 2003).

2.3 Hydroclimate modes

Studies of the climatological basis of the historical rise in lake level at Devils Lake identified 1980 as one such transition point between drier and wetter hydroclimatic modes (Vecchia 2008, Hoerling et al. 2010, Todhunter and Fietzek-DeVries 2016). These low-frequency variations in regional hydrology must be related to persistent, large-scale, low-frequency circulation anomalies (Hirschboeck 1987), but the precise mix and phase of climate indices driving the two modes has yet to be identified.

The two hydroclimatic modes are associated with distinctive hydrological responses, producing what has been referred to as a drought to deluge cycle (Winter and Rosenberry 1998, Todhunter and Rundquist 2004, Zhang et al. 2009). Todhunter (2016) developed a hydroclimatic profile of the two climate modes that included: climate variables (precipitation, mm; air temperature, °C), lake hydrology variables (lake evaporation, mm; net lake evaporation, mm), closed lake variables (inflow factor, %; climate parameter, fraction), selected climatic water balance variables (aridity index, fraction; precipitation minus potential evapotranspiration, mm; potential evapotranspiration minus actual evapotranspiration, mm; mean annual soil moisture, mm; dryness index, fraction; moisture index, fraction), and hydrology variables (surface inflow, km³; runoff ratio, %). Mode 1 was a cooler and drier mode, while Mode 2 was a warmer and wetter mode. Mean annual precipitation onto the lake increased from 466 to 579 mm from Mode 1 to Mode 2 (113 mm or 24% increase). This modest but sustained increase in the primary hydrological driver was amplified into much larger changes in regional hydrology between the two modes. Mean annual surface inflow to the lake increase from 0.041 to 0.198 km³ (383% increase) and mean annual runoff ratio increased from 1.08 to 4.12% (282% increase) from Mode 1 to Mode 2.

2.4 Regional hydrology

Devils Lake has experienced a precipitous increase in lake volume since 1993 (Fig. 2). This increase would appear to require additional causation than the transition from a drier and cooler hydroclimate mode to a wetter and warmer hydroclimate mode, and the basin’s location along a sharp hydrotone gradient (Todhunter 2016). In the absence of human drivers of landscape change, the relationship between precipitation input and streamflow output in most environments is assumed to be linear and stationary (Gupta et al. 2015).

There are three reasons to question this assumption for the Devils Lake Basin. First, even in basins where the annual precipitation–runoff relationship varies linearly across the normal range of precipitation variations, a nonlinear relationship may occur with extreme high and low precipitation variations due to the effect of basin hydrological memory, and the manner by which basin antecedent wetness conditions control the partitioning of climatic input to hydrological output (Risbey and Entekhabi 1996). Fundamental changes in internal basin dynamics during multiyear dry periods can impact the annual precipitation–runoff relationship (Saft et al. 2016a). The climate mode transition that occurred in 1980 makes such a non-stationary relationship probable (Todhunter and Fietzek-DeVries 2016).

Second, principles of fill–spill hydrology operate within the multiple wetland complexes within each sub-basin (Spence and Woo 2003). These principles have been studied for wetlands in the Canadian Prairie Pothole Region, and must certainly hold for the North Dakota Prairie Pothole Region (Spence 2007, Shook and Pomeroy 2011, Ehsanzadeh et al. 2012, Shaw et al. 2012, 2013, Dumanski et al. 2015, Shook et al. 2015). Snowmelt is the dominant source of runoff, and the transformations of snowmelt to runoff, and runoff to streamflow are nonlinear (Shook et al. 2015). Local depressions within a sub-basin do not contribute runoff to that sub-basin’s lake until its spill elevation is reached, and the wetland complex is hydrologically connected to its lake. Sub-basin contributing area is spatially and temporally dynamic, with the contributing area of each sub-basin expanding and contracting based upon the wetland storage levels and wetland connectivity experienced across the multitude of sub-basin wetland complexes.

Third, the sub-basins in the upper Devils Lake basin follow similar principles of fill–spill hydrology (Ryan and Wiche 1988). Each sub-basin drains to the next lowest sub-basin only when its lake level exceeds the elevation of the topographic swale separating the two adjacent lakes and the two adjacent sub-basins become hydrologically connected. The entire upper basin contributes runoff into Devils Lake only when every lake in the chain is above its spill elevation and flowing over its swale to the next lowest lake in the chain. If any lake along the chain is below spill elevation its sub-basin area, and that of all sub-basins above it in the chain, is not contributing runoff to Devils Lake (Ryan and Wiche 1988). The entire upper basin is hydrologically connected to Devils Lake only when multiple thresholds of lake level are exceeded. Thus, the area contributing surface inflow to the lake (R_i) is dynamic through time, expanding during wet phases, and contracting during dry phases, and producing a non-stationary contributing area to Devils Lake. It is hypothesised that all three effects contribute to the non-stationary runoff regime apparent in the hydrology of Devils Lake (Saft et al. 2015, Shook et al. 2015).

3 Data and methods

3.1 Monthly lake water budget data

Monthly water budget data for Devils Lake were provided by the US Geological Survey (USGS) North Dakota Water Science Center (A. Vecchia, personal communication, 14 October 2011). These included monthly totals of lake precipitation (P_L, mm), lake evaporation (E_L, mm) and inflow (R_i, km³) for January 1950 to December 2010. Net lake evaporation (E_L − P_L, mm) was determined as a secondary variable from the primary data. These data are considered the best available estimate of the lake hydrology, and served as the basis for a lake level frequency analysis for Devils Lake (Wiche and Vecchia 1996), and two versions of a stochastic simulation model used to evaluate the impact of pumping from proposed outlets upon future lake levels (Vecchia 2002, 2008). The total inflow term (R_i) includes the sum of inflow from the USGS streamgauges at Big Coulee and Channel A, and local ungauged inflow from the peripheral area of Devils Lake. Full details concerning the development of the monthly P_L, E_L, and R_i time series are given in Wiche and Vecchia (1996) and Vecchia (2002, 2008). A standard water year running from 1 October to 30 September is examined to capture winter snow accumulation, spring snowmelt production, and summer thunderstorm precipitation within a single year.

3.2 Lake water budget

For a large lake with high lake volume inertia it is appropriate to analyse the lake water budget over an annual time step. The change in lake volume (ΔL_V) for an annual time step is given by the water balance equation for a terminal lake:

ΔL_V = R_i + P_LL_A + G_i – E_LL_A – G_o (1)

where R_i is surface runoff into the lake (km³), P_L is direct precipitation onto the lake (mm), L_A is lake area (km²), E_L is evaporation from the lake (mm), and G_i and G_o are groundwater inflow and outflow from the lake (Mason et al. 1994). Shallow groundwater inflow from glacial till deposits and deep groundwater contributions to the lake are estimated to be only 0.0037 km³ year⁻¹, and do not have an important effect upon annual lake fluctuations (Wiche and Pusc 1994). In the summer of 2005 construction began on a series of outlets to pump and remove water from the Devils Lake system to mitigate flood damages (North Dakota State Water Commission 2014). When operative, the two outlets have a pumping capacity of 17.0 m³ s⁻¹. The last year of minimal pumping volume was 2010, which marks the ending date of the analysis (W. Schuh, personal communication, North Dakota State Water Commission, 22 June 2015).

By assuming that net groundwater flow (G_i − G_o) is negligible and can be dropped from the analysis (Wiche and Pusc 1994), and ending the analysis in 2010 when human abstraction by pumping became an important water balance term, Equation (1) can be simplified:

ΔL_V = R_i + P – E (2)

where ΔL_V is the change in lake volume (km³), and surface inflow to the lake (R_i) is the water year total of the monthly values provided by the USGS. Lake precipitation (P, km³) and lake evaporation (E, km³) are discretized annual totals obtained by multiplying the water year total precipitation on the lake (P_L, mm) and evaporation from the lake (E_L, mm) by the lake area (L_A, km²) on 30 September (Mohammed and Tarboton 2012).

3.3 Historical lake level and lake volume time series

Historical observations of lake level for Devils Lake were obtained from USGS Water Data for the Nation (USGS (United States Geological Survey) 2013). This dataset is a composite of infrequent 19th century survey elevations, annual observations to 1930, discontinuous daily observations to 1965, and daily observations to the present.

Advances in satellite remote sensing have led to recent improvements in our ability to measure lake level and lake area, and thus to estimate lake volume variations (Gao 2015). The complete lake hypsometry defining the relationships between lake level, lake area, and lake volume is normally not known for closed lakes, but is available for Devils Lake. Tables of L_L–L_A and L_L–L_V relationships for lake levels in the ranges 426.7–445.6 m (1400.0–1462.0 feet) were obtained from the USGS (2014), which allowed for the creation of a L_V time series from the historical L_L series.

4 Results

4.1 Annual volumetric lake water budget totals

Direct precipitation onto the lake (P_L), evaporation from the lake (E_L), and lake net evaporation (E_L – P_L), in units of depth (mm), are the standard measurements available for examining the hydrology of terminal lakes (Gao 2015). Surface runoff into the lake (R_i), is usually measured as a flow rate (m³ s⁻¹) but is easily converted to an annual surface runoff volume (km³). A plot of the annual totals for P_L (mm), E_L (mm), and E_L − P_L (mm) for water years 1951–2010 is available in Todhunter (2016). Table 1 gives the mean (X), standard deviation (SD) and coefficient of variation (CV) for each variable for the 60-year study period.

Table 1. Summary statistics for annual water budget variables for standard lake measurements (mm) at Devils Lake. X: mean; SD: standard deviation; CV: coefficient of variation. Source: USGS.

Period	Variable	PL (mm)	EL (mm)	EL – PL (mm)
POR 1951–2010	X	523	763	241
	SD	110	64	142
	CV (%)	21	8	59
Phase 1 1951–1968	X	471	736	266
	SD	101	37	126
	CV (%)	21	5	47
Phase 2 1969–1992	X	515	794	279
	SD	117	73	144
	CV (%)	23	9	52
Phase 3 1993–12 010	X	585	750	165
	SD	80	55	130
	CV (%)	14	7	79

Annual precipitation onto the lake surface (P_L) averages 523 mm, and has a coefficient of variation (CV) of 21% (Table 1). Annual lake evaporation (E_L) averages 763 mm, and is much less variable, with a CV of only 8% (Table 1). The lake freezes over during the cold season, so lake evaporation normally only occurs between April and November. The maximum lake depth is currently about 15.0 m, with an average lake depth of about 5.8 m. The shallow lake depth, large wind fetch, and active regional wind regime produce a well-mixed lake with a single thermal layer that experiences minimal lake heat storage, and a maximum average monthly lake evaporation of 153 mm occurring in August.

Lake net evaporation (E_L − P_L) is the only natural source of water loss unless the lake level exceeds the ordinary high water mark of 444.4 m. Average lake net evaporation is 241 mm (Table 1), with annual totals ranging between −73.4 and 586.7 mm; negative values occurred in only three of the 60 years. The CV for lake net evaporation is 59%, the second highest CV for the four lake water budget variables (Table 1). A weak but statistically significant inverse linear relationship occurs between P_L and E_L (Fig. 3, P = 0.031, SAS for Windows v.9.4), which leads to high inter-annual variability in lake net evaporation (Table 1). Evaporation from shallow freshwater lakes is controlled by the lake energy budget (Schertzer 1997, Lenters et al. 2005). Net radiation is the dominant driver of inter-annual variation of lake evaporation, and net radiation and lake evaporation have a strong positive relationship in this region (Lenters et al. 2005). Moderate inverse relationships are found between lake evaporation and the lake–air temperature difference and Bowen ratio (Lenters et al. 2005). Years with increased annual precipitation are associated with higher levels of cloud cover, decreased solar radiation transmission, reduced surface net radiation, and reduced lake evaporation (Lenters et al. 2005).

Figure 3. Plot of annual lake evaporation (E_L, mm) against annual lake precipitation (P_L, mm) for Devils Lake for water years 1951–2010.

Annual totals for the lake water budget variables determined from the volumetric analysis for water years 1951–2010 are shown in Figure 4 and summarized in Table 2. These are the annual totals for P_L, E_L, and E_L − P_L (mm) shown in Todhunter (2016), but here converted to volumetric totals for P, E, and E − P (km³). Direct precipitation onto the lake (P, km³) exceeded the surface inflow to the lake (R_i, km³) in 44 of the 60 years; the volume of surface inflow to the lake was greatest for the remaining 16 years. Total water input to Devils Lake over the 60-year period was 15.73 km³, of which 8.33 km³ (53.0%) was by direct precipitation on the lake, 7.18 km³ (45.7%) was by surface inflow, and 0.22 km³ (1.3%) by groundwater inflow. Total water output from Devils Lake was 11.56 km³ over the 60-year period, of which lake evaporation contributed 11.47 km³ (99.2%). Pumping from the lake outlets totalled only 0.09 km³ (0.8%) and was inconsequential through 2010, but has become an important output variable since 2010. Pumping removed 0.063 km³ of water in 2010, which is comparable to the average annual net evaporation (E − P) of 0.054 km³ (Table 2). Pumping can now potentially double the average lake net evaporation, and will have an important role in the future lake water budget.

Table 2. Summary statistics for annual volumetric water budget variables (km³) at Devils Lake. Source: USGS.

Period	Variable	P (km^3)	E (km^3)	E – P (km^3)	Ri (km^3)
POR 1951–2010	X	0.139	0.193	0.054	0.120
	SD	0.106	0.122	0.045	0.169
	CV (%)	76	63	84	141
Phase 1 1951–1968	X	0.054	0.084	0.030	0.012
	SD	0.014	0.012	0.014	0.021
	CV (%)	27	14	47	177
Phase 2 1969–1992	X	0.102	0.155	0.054	0.062
	SD	0.031	0.027	0.030	0.070
	CV (%)	31	18	56	113
Phase 3 1993–2010	X	0.274	0.351	0.078	0.305
	SD	0.092	0.097	0.066	0.197
	CV (%)	34	28	85	65

Figure 4. Annual volumetric water budget totals (km³) for Devils Lake for water years 1951–2010.

Surface inflow to the lake (R_i) is the most variable lake water budget component, and is the lake input most responsible for the dramatic rise in lake level at Devils Lake. Mean annual R_i is 0.120 km³ (Table 2), but annual totals range between 0.0 and 0.720 km³, producing a CV of 141%. Surface inflow to the lake is dominated by spring snowmelt production, with peak monthly R_i totals occurring in April. Winter precipitation (Dec–Feb) is only 9% of annual precipitation, but spring surface inflow (Mar–May) is 69% of annual surface inflow to the lake. Summer runoff (Jun–Aug) to the lake is only 24% of annual R_i although summer precipitation is 46% of annual P_L. A preliminary examination of annual R_i in Figure 4 reveals three distinct clusters. First, a period of minimal R_i from the beginning of the period of analysis through the late 1960s; second, a period of moderate R_i beginning in the late 1960s and lasting through the early 1990s; finally, a period of much greater R_i beginning in the early 1990s and continuing through the end of the period of analysis.

The net change in lake volume (ΔL_V) obtained from the discretized water budget shown in Figure 4 is 3.42 km³, while the measured change in lake volume between 30 September 1950 and 30 September 2010 is 4.17 km³, for a closure error of −0.75 km³ (−18%). This small closure error is most likely attributable to the assumptions and uncertainties associated with development of the lake water budget database (Vecchia 2008), and the use of an annual time step in discretizing the volumetric lake water budget.

Several additional variables used to characterize terminal lakes are given in Table 3 (Szesztay 1974, Street 1980, Mason et al. 1994). These include the inflow factor (IF, %):

(3)

Table 3. Summary statistics for selected annual volumetric water budget totals for Devils Lake.

Period	Runoff ratio (%)	Inflow factor (%)	Precipitation factor (%)	LA/A (%)
POR 1951–2010	2.6	30.2	69.8	2.8
Phase 1 1951–1968	0.3	12.3	87.7	1.3
Phase 2 1969–1992	1.5	30.2	69.8	2.2
Phase 3 1993–2010	6.4	48.1	51.9	5.2

A precipitation factor (PF, %) can also be determined as a complementary variable:

(4)

Neglecting the minimal groundwater input, surface inflow contributes 30.2% of lake volume input for the study period, with precipitation on the lake providing the remaining 69.8% (Table 3).

4.2 Ratio of lake area to closed basin area

The area of a closed basin (A, km²) is the sum of the lake area (L_A, km²) and the basin area (B_A, km²). Figure 5 shows lake area as a percentage of closed basin area for 1901–2014. The percentage averages 2.8% over the study period, but varies over 1.1 < L_A/A < 9.0%. This variable has an important control over the relative contribution of P_L and R_i to lake input, and has been found to range between less than 1.0% and nearly 50% for closed basins throughout the world (ILEC 2005). Even after the spectacular rise in lake level the percentage is still less than 9% for the Devils Lake basin, placing it toward the lower end of the distribution among global closed basins.

Figure 5. Lake area (L_A, km²) of Devils Lake as a percentage of the Devils Lake closed basin area (A, km²) for 15 June 1901 to 31 January 2014.

4.3 Precipitation–runoff relationship

The precipitation–runoff relationship is a statistical model that shows drainage basin runoff response across a range of years of variable wetness (Saft et al. 2015). It aggregates drainage basin runoff response to precipitation input at an annual time step, and is a simple yet effective characterization of drainage basin precipitation–runoff processes (Saft et al. 2015). If human influences upon runoff are minor and precipitation is the main driver of runoff variations, the plot of annual precipitation against annual runoff will assume a linear pattern in most hydrological environments (Gupta et al. 2015). For cold regions, however, Dumanski et al. (2015) and Shook and Pomeroy (2012) have shown that changes in precipitation phase, and the duration of precipitation events can also influence hydrological response even with no change in annual precipitation amount. In a water-limiting environment a nonlinear relationship may arise due to human influences upon drainage basin processes, or due to climatic drivers other than annual precipitation totals (Saft et al. 2015). The scatter plot of annual lake precipitation (P_L, mm) versus annual surface inflow (R_i, km³) for the study period reveals no coherent aggregate linear or nonlinear relationship (Fig. 6). For example, for an annual P_L total of 600 mm, annual runoff (R_i) ranges from 0.0 km³ to more than 0.7 km³.

Figure 6. The annual precipitation–runoff relationship for Devils Lake for water years 1951–2010.

The runoff ratio (RR, %) is defined as the volume of runoff to the lake (R_i) divided by the volume of precipitation falling on the sub-basin (P_B, mm). It is determined by:.

(5)

where P_B is assumed equal to P_L. The time series of RR for the 60-year study period is shown in Figure 7. The average RR value is 2.6% (Table 3), but varies for individual years between 0.0 and 14.3%. Only one year had no runoff (RR = 0.0), although many years experienced minimal runoff with RR values close to zero. The inter-annual variation in RR is extreme, with a CV of 136%, and RR ranging by more than three orders of magnitude. A nine-term binomial filter has been applied to the annual RR time series in Figure 7 to help identify long-term patterns of variation (Richard Heim, Jr., personal communication, Meteorologist, NOAA National Centers for Environmental Information, 5 March 2018).

Figure 7. Runoff ratio (%) for Devils Lake basin for water years 1951–2010. Individual water years are shown with open triangles; a nine-term binomial filter (solid curve) has been applied to the time series to clarify temporal patterns. Horizontal (red) lines show the average values for Phases 1–3.

The use of double-mass curves to check the consistency of the relationship between hydrological variables is a standard method in hydrology (Searcy and Hardison 1960). In Figure 8 cumulative surface runoff to the lake (R_i, km³) is plotted against cumulative precipitation on the lake (P_L, mm) for water years 1951–2010. If the relationship between the two variables is a fixed ratio the two variables should plot as a linear trend, with scatter about the line due to short-term basin memory effects. A break in the plotted double-mass curve would indicate a change in the relationship between precipitation input and runoff generation over time (Searcy and Hardison 1960). In Figure 8, three distinct phases occur in the double-mass curve: (a) Phase 1: 1951–1968; (b) Phase 2: 1969–1992; and (c) Phase 3: 1993–2010. Figure 8 indicates a non-stationary precipitation–runoff relationship during the study period, with regime transitions occurring in 1969 and 1993.

Figure 8. Double mass curve of cumulative lake precipitation (P_L, mm) against cumulative runoff (R_i, km³) for Devils Lake for water years 1951–2010. Critical thresholds occur in 1969 and 1993.

The sub-totals for each set of variables given in Tables 1, 2, and 3 were determined for the three periods of precipitation–runoff relationship identified in Figure 8. The annual water budget variables for standard lake measurements are summarized in Table 1: P_L shows a steady increase from Phase 1 to Phase 3. The coefficient of variation for P_L is comparable for Phases 1 and 2, but much smaller for Phase 3; E_L shows a moderate increase from Phase 1 to Phase 2, and then a small decline for Phase 3; and E_L − P_L has a small increase from Phase 1 to Phase 2, but then a more substantial decline for Phase 3.

The annual water budget variables for the volumetric analysis are summarized in Table 2. All four variables increase from Phase 1 to Phase 2, and then again from Phase 2 to Phase 3: P increases by 89% and then 168% from Phase 1 to Phase 3; E increases by 85% and then 126% from Phase 1 to Phase 3; and E − P also increases progressively from Phase 1 to Phase 3, but by smaller percentages of 80 and 44%, respectively. The largest percentage increase occurs for R_i, which shows a 417% increase from Phase 1 to Phase 2, and a 392% increase from Phase 2 to Phase 3. The CVs for P, E, and E − P all increase from Phase 1 to Phase 3, while the CV for R_i decreases continuously from Phase 1 to Phase 3.

Key closed lake variables for the three phases defined by Figure 8 are shown in Table 3. Average lake area (L_A) as a percentage of closed basin area (A) increased from 1.3% in Phase 1 to 5.2% in Phase 3, which contributed to the increase in the volume of direct precipitation onto the lake (P, Fig. 4). Although the magnitude of P increased from Phase 1 to Phase 3, the relative contribution of P to the lake water budget declined, with the average precipitation factor (PF, %) decreasing from 87.7 to 51.9% from Phase 1 to Phase 3 (Table 3). The explosive growth in R_i as a contributor to the lake water budget is shown by the increase in the inflow factor from an average of 12.3% in Phase 1 to an average of 48.1% in Phase 3. By the end of the study period, R_i was the largest contributor to the lake volumetric water budget. Average runoff ratios (RR, %) for the three phases were 0.3, 1.5, and 6.4%, respectively (Fig. 7, horizontal red lines), for a nearly 20-fold increase from Phase 1 to Phase 3.

5 Discussion

The large fluctuations in reconstructed lake level during the late Holocene (Bluemle et al. 1999), and observed lake level during the late 20th century (Fig. 2), are the result of special characteristics of the Devils Lake water balance. Devils Lake is an example of an amplifier terminal lake, a term first proposed by Street (1980) following the general classification of lake water budgets presented by Szesztay (1974). The water budgets of amplifier terminal lakes are dominated by R_i and E. Their lake levels exhibit a large-amplitude response to long-term climatic fluctuations because numerous positive feedbacks amplify or attenuate runoff generation from their contributing areas (Street 1980). Two additional factors are also operative in the Devils Lake Basin that create temporal lags in the lake system response to the hydroclimatic drivers, and affect the rate of amplification/attenuation to the forcings. These concern the location of the basin along a dry sub-humid/moist sub-humid hydroclimatic boundary, and the principles of fill−spill hydrology that regulate precipitation–runoff relationships at both the basin and sub-basin scales.

Figure 9 shows a conceptual model of the positive feedback relationships in an amplifier lake associated with a long-term increase in precipitation (+ΔP_L). The same positive feedback relationships hold for a long-term decrease in precipitation (−ΔP_L), except that all of the terms in the figure would change sign. Direct precipitation on the lake (P_L) is the single greatest input to lake volume (L_V), contributing 53% of the input to the lake over the period of record. Annual +ΔP_L inputs, however, can only produce a small impact upon annual +ΔL_V. Long-term climate change producing an increase in annual P_L leads to an increase in L_V for two reasons. First, the increase in P_L contributes directly to an increase in lake precipitation (P, km³). For Devils Lake, average P_L increased by 113 mm (24%) from climate Mode 1 to Mode 2 (Todhunter 2016). Second, as lake area (L_A) expands it produces an additional linear increase in P due to the expanding lake area. The increase in L_A/A from 1 to 9% over the period of record (Fig. 5) creates a steadily increasing trend in the P time series (Fig. 4,) which is in marked contrast to the inter-annual variability in the P_L time series shown in Todhunter (2016). A long-term increase in P_L, however, leads to an absolute increase in the contribution of P to L_V, but a relative decrease in the contribution of P to L_V. Average annual P_L increases from 471 to 585 mm from Phase 1 to Phase 3 (Table 1), and average annual P increases from 0.054 to 0.274 km³ for the same period (Table 2). Table 3, however, shows a decrease in the precipitation factor of 87.7 to 51.9% from Phase 1 to Phase 3.

Figure 9. Conceptual model of the positive feedback relationship in an amplifier lake. The figure shows the positive feedback associated with a long-term increase in precipitation (+ΔP_L). All of the signs change for the positive feedback associate with a long-term decrease in precipitation (–ΔP_L). The units for all variables are depth (mm).

A long-term increase in P_L also contributes to an increase in lake volume (+ΔL_V) through its direct effect on lake evaporation (E_L) and net lake evaporation (E_L − P_L). Mean annual air temperature warmed by an average of 0.63°C from climate Mode 1 to Mode 2 (Todhunter 2016), but lake evaporation remained mostly unchanged, increasing from 752 to 774 mm for the two climate modes, or from 736 to 750 mm from Phase 1 to Phase 3 (Table 1). Increasing cloud cover and atmospheric humidity are effective at attenuating most of the increased evaporative demand associated with atmospheric warming (Schertzer 1997, Lenters et al. 2005).

More importantly, net lake evaporation, which is the only natural source of lake volume loss except for those occasions when the lake overflows through the ordinary high water mark and the lake transitions from a closed lake to a flowing lake (Szesztay 1974), decreased as the long-term climate transitioned from a drier to a wetter mode. Table 1 shows a decrease in E_L − P_L from 266 to 165 mm from Phase 1 to Phase 3. This reduces the ability of the lake to remove water, on a unit depth basis, from the lake by evaporation.

As the hydroclimate regime shifts from a drier to a wetter climate mode, Devils Lake transitions from a precipitation-dominated to a runoff-dominated water budget, or from a P-E to an IP-E lake water budget in the nomenclature of Szesztay (1974) and Street (1980). The slow, progressive, and cumulative effects of the prolonged wetter climate mode produce important changes in the water budget and runoff characteristics of the upper basin. The increasing trend in P_L and decreasing trend in E_L – P_L have a direct positive impact upon ΔL_V. These direct effects, however, are small in comparison to their indirect effects upon +ΔL_V through surface inflow (R_i). Although the increase in P_L from Mode 1 to Mode 2, from 466 to 579 mm, is substantial (+24%), it can only have a small direct effect upon +ΔL_V. The continuation of the wet cycle has a large effect upon the moisture status of the upper basin that dynamically changes R_i. Average annual soil moisture (SM) is 57 mm for climate Mode 1 and 97 mm for climate Mode 2 (Todhunter 2016), while mean annual SM increases from 67 to 100 mm from Phase 1 to Phase 3.

During Phase 1 the runoff ratio (RR) averages only 0.3%, and the average precipitation factor (PF) is seven times greater than the average inflow factor (IF), even though the basin area is 99% of the closed basin area (Fig. 5). As climate Mode 2 becomes established, SM slowly begins to increase, leading to a slow but progressive increase in RR and R_i. Once Phase 3 becomes established, basin area is only 91% of the closed basin, but RR now averages 6.4%, for a 20-fold increase from Phase 1; R_i has increased from 0.012 to 0.305 km³ from Phase 1 to Phase 3 (Table 2), for a 25-fold increase; and PF and IF are of comparable magnitude during Phase 3 (Table 3). The establishment of the wetter climate regime causes both the absolute contribution (R_i, km³) and the relative contribution (PF, %) of surface inflow to the volumetric lake water budget to increase.

The volumetric lake water budget time series illustrates these effects and transitions (Fig. 4). During Phase 1 (1951–1968) Devils Lake exhibits a clear pattern of a precipitation-dominated amplifier terminal lake (P-E). Inflow factors are always less than 35%, both P and E are large relative to R_i, and E controls the lake water budget, since net lake evaporation allows only small annual fluctuations in L_V (Szesztay 1974, Street 1980). Phase 3 (1993–2010) shows the pattern for an inflow-dominated amplifier terminal lake (IP-E). Inflow factors now range between 35 and 65%, and R_i is large relative to P and E. The lake basin now amplifies the effects of the climatic cycle fluctuations on the lake water budget leading to a new equilibrium hydrological state. Lake level and lake volume are highly unstable, rising quickly when high annual R_i overwhelms net lake evaporation, but declining slowly in years when low annual R_i allows net lake evaporation to remove water from the lake. Phase 2 (1969–1992) is a transitional period during which the direct effects of P_L fluctuations on the lake water budget become less important, and the indirect effects of P_L fluctuations on the lake water budget become dominant through their effects upon inflow. Results from Todhunter (2016) show that there is about a 15-year lag between the establishment of the wetter P_L mode and the initiation of the upward trend in L_V due to the substantial hydrological inertia present in the closed lake system (Street 1980, Mason et al. 1994).

Equilibrium is restored to the positive feedback relationship through the relative changes in P and E shown in Figure 10. Because E_L is normally greater than P_L, net lake evaporation is normally positive, producing an annual net loss in lake volume. As L_A expands, the rate of growth of E exceeds the rate of growth of P due to the surface area dependence of both terms. Figure 4 shows that both P and E show large increases during Phase 3, but the continued increase in lake volume is dependent upon annual R_i being greater than annual E − P. The weak negative relationship between E_L and P_L (Fig. 3) is overridden by the positive relationship between E − P and L_V (Fig. 10), creating a negative feedback that works to restore equilibrium to the lake water budget.

Figure 10. Plot of lake evaporation (E, km³) and lake precipitation (P, km³) against lake volume (L_V, km³) for Devils Lake for water years 1951–2010.

The precipitous rise in lake volume results from several probable causes (Fig. 2). First, the Devils Lake Basin lies along a hydrotone or a region of sharp hydroclimatical transition (Todhunter 2016). As the climate mode shifts from climate Mode 1 to Mode 2 the position of the moist–sub-humid/dry–sub-humid boundary shifts westward, producing changes in the regional moisture status that are small on an absolute basis, but much larger on a relative basis.

Second, basin antecedent wetness conditions play a critical role in precipitation–runoff partitioning. Most hydrological investigations assume stationarity in hydrological processes across the range of climatic variations (Saft et al. 2015). They also assume that the historical record of precipitation–runoff dynamics captures the full range of basin response (Saft et al. 2016b). Basins that exhibit a linear precipitation–runoff relationship for most years may exhibit a nonlinear relationship during extremes in the precipitation distribution (Risbey and Entekhabi 1996). Saft et al. (2016a, 2016b)) show that prolonged drought of decadal-scale or longer causes a reduction in the normal precipitation–runoff relationship, and that such a reduction is most likely in drier, flatter, less forested basins such as characterize the Devils Lake Basin. The prolonged drought to deluge climate change experienced in the study area would likely produce a nonlinear increase in the precipitation–runoff partitioning.

Third, the fill–spill principles of hydrology experienced at the basin scale and among the innumerable wetland complexes within each sub-basin would produce nonlinear relationship precipitation–runoff relationships across all levels of precipitation. This has been demonstrated for watersheds within the glacial drift prairie region of Canada (Shook and Pomeroy 2011, Ehsanzadeh et al. 2012, 2016, Shaw et al. 2012, Dumanski et al. 2015, Shook et al. 2015). At the basin scale, the filling and spilling of the upper basin chain of lakes would lead to an increase in the hydrological connectivity of the sub-basins contributing inflow to Devils Lake, which would produce a pattern of step function increases in the effective contributing area to the lake system. At the sub-basin scale, the innumerable wetland complexes would also experience a filling and spilling pattern leading to an increase in hydrological connectivity, but the pattern of increase in effective contributing area would be much more continuous in nature due to the vast number of wetland complexes.

Figure 8 shows three phases in the precipitation–runoff relationship, with distinct thresholds occurring in 1969 and 1993. Phase 1 illustrates the precipitation–runoff relationship of climate Mode 1, the cooler and drier climate mode, marked by sub-basin lakes and wetland complexes with storage volumes below capacity and not spilling to lower elevation areas, minimal hydrological connectivity between the lakes and wetland complexes, and low basin-wide antecedent wetness conditions. Phase 3 illustrates the precipitation–runoff relationships of climate Mode 2, the warmer and wetter climate mode. Long-term establishment of this climate mode leads to increased basin-wide antecedent wetness conditions, increased surface runoff, increased sub-basin lake and wetland complex water storage, and increased hydrological connectivity at the basin and sub-basin scales. Phase 2 is a transitional period marked by a progressive increase in basin-wide antecedent wetness conditions, increased sub-basin lake water levels, increased maximum wetland complex water storage, and improving hydrological connections between individual wetland complexes and sub-basins.

Several topics for further investigation arise. First, the explanation provided for the rapid rise in lake level at Devils Lake was based upon a bulk or basin-wide analysis of hydrological processes. A study employing a physically-based, distributed hydrological model would confirm and elaborate upon the explanation offered in this paper (Pomeroy et al. 2007). Second, following Wiche and Pusc (1994), it has been assumed that groundwater inflow was a negligible contributing factor to the lake rise. The Spiritwood Aquifer underlies the chain of lakes and bays that comprise Devils Lake, so this assumption may not be valid as the wet phase becomes fully established. Rising groundwater levels close to the lake could lead to both increasing groundwater inflow to the lake, and higher surface runoff coefficients in the adjacent area. Finally, Devils Lake acts as a low-pass filter on low-frequency climate variability (Mason et al. 1994). Climate variation is not random, but exhibits preferred time scales of residence within each transient state. Each closed basin lake, acting as a low-pass filter with specific internal system dynamics, will translate that low-frequency climate signature into a unique land surface hydrological response. Understanding this nonlinear behaviour is essential to providing meaningful assessments of hydrological response to simulations of future climate change.

6 Conclusions

The hydrology of Devils Lake exhibits features common to all amplifier terminal lakes, as well as features specific to its glacial formation. These factors combine to produce the large oscillations in lake level that characterize the late Holocene. During the drier drought phase the lake has a precipitation-dominated water budget. As the long-term climate transitions to the wetter deluge phase the lake takes on a runoff-dominated water budget, following a time lag due to basin memory effects. Fill–spill hydrology effects further accentuate the natural amplifications and dampenings of the positive feedback relationships, producing a non-stationary precipitation–runoff relationship in surface inflow to the lake.

These nonlinear precipitation–runoff relationships may explain the weak correlations obtained between measured streamflow variables and a range of standard climate indices (McCabe and Wolock 2014). The low-frequency variations in regional hydrology must be related to persistent, large-scale, low-frequency circulation anomalies, but their identification has remained elusive. Efforts to investigate the direct relationship between climate indices and precipitation (rain and/or snowfall) may be more productive than the more indirect relationship between climate indices and streamflow. The fill–spill hydrological principles and nonlinear precipitation–runoff relationships present in the northern glaciated plains also have important implications for the investigation of the impact of anthropogenic climate warming upon trends in regional hydrology.

Acknowledgements

I thank Rick Thalacker for his work in making Figure 1. Aldo Vecchia, Jr. of the USGS North Dakota Water Science Center in Bismarck, ND, provided the monthly lake water budget data for Devils Lake. Bill Schuh of the North Dakota State Water Commission provided the discharge data for the outlets on Devils Lake. Two reviewers provided several comments that helped improve the final manuscript.

Disclosure statement

No potential conflict of interest was reported by the author.

Notes

¹ http://www.swc.nd.gov/info_edu/reports_and_publications/ pdfs/wr_investigations/wr22_report.pdf.