2. Methods
To characterize changes in lake level variability, we analyzed lake
level data compiled by the Great Lakes Coordinating Committee
(Coordinating Committee on Great Lakes Basic Hydraulic Hydrologic Data,
n.d.). Monthly mean lake level data are calculated by combining
observations from a suite of gauges around the Great Lakes basin
maintained by the U.S. National Oceanic and Atmospheric Administration
(NOAA) and Canadian Hydrographic Service (CHS) (Fry et al., 2022). These
data span from 1918 to the present. In addition to lake level data, our
analysis uses Great Lakes monthly hydrologic data curated by the NOAA
Great Lakes Environmental Research Laboratory (GLERL) for net basin
supply variables: overlake precipitation, overlake evaporation, and
runoff.
The NOAA-GLERL monthly hydrometeorological database (GLM-HMD) is the
first dataset to assimilate hydrometeorological measurements across both
the U.S. and Canadian portions of the Great Lakes basin dating back to
the early and mid-1900s (Hunter et al., 2015). The GLM-HMD uses a
Thiessen weighting approach to compute overlake precipitation from daily
monitoring stations (Hunter et al., 2015), including stations up to 50
km from the lake basin depending on gauge density, though stations over
the lake surface or on islands are sparse. It should also be noted that
station availability expanded rapidly from the early 1900s before
peaking from the late 1940s to the present (see fig. 2, Hunter et al.,
2015) and GLM-HMD overlake precipitation calculations are thus available
beginning in 1940. Overlake evaporation is computed by the Large Lake
Thermodynamics Model (LLTM) as described in Croley (1989), and is
available from 1950 to present. The LLTM is a 1-dimensional
thermodynamics model that computes overlake evaporation by simulating
the energy balance above the lake surface, vertically-summed heat
storage through the lake columns, and aerodynamic evaporation (Fry et
al., 2022). GLM-HMD estimates of runoff are calculated from streamflow
measurements from the United States Geological Survey and Water Survey
of Canada gauges based on subbasin area with extrapolations for
unrepresented subbasin area and ungauged subbasins (Hunter et al.,
2015). Fry et al. (2013) determine that this approach provides reliable
estimates of total discharge to the lakes. In preparation to examine
trends in interannual variability, we transformed monthly data to annual
data by taking an annual mean (lake levels) or annual sum
(precipitation, evaporation, runoff) on the component time series as
appropriate. Finally, we note that we considered alternative data
sources, notably output from the Large Lake Statistical Water Balance
Model (L2SWBM) (Gronewold et al., 2020), and a comparison of annual
values from GLM-HMD and L2SWBM produced correlations of at least 0.97
for all lake-component combinations.
Here, we define interannual variability as the standard deviation in a
time series over a 13-year moving window. We select a 13-year window
width for two reasons. First, this width is longer than the cycles of
known 8- and 12-year modes of lake level variability (Hanrahan et al.,
2014; Cheng et al., 2021). Interannual lake level fluctuations are
driven both by modes of climate variability at various periodicities,
with notable periodicities of 8, 12, and 30 years (Hanrahan et al.,
2014; Cheng et al., 2021), as well as by long-term trends in net basin
supply components. Relevant modes of climate variability include the
Pacific Decadal Oscillation (PDO), Atlantic Multidecadal Oscillation
(AMO), North Atlantic Oscillation (NAO), El Nino-Southern Oscillation
(ENSO), and the Pacific-North American Oscillation (PNA) (Wang et al.,
2018; Ghanbari and Bravo, 2008). Second, we examined the dependence of
trends through a sensitivity analysis using window widths ranging from 5
to 25 years. This revealed a consistency in trends in interannual
variability for window widths greater than 12 years despite
inconsistencies in trends across shorter window widths. We focus our
analysis on the period from 1970 to the present. While this truncates
the length of data available for time series analysis, this is a period
of increasing influence of anthropogenic climate change. Notably, there
is evidence of an emerging intensification of the hydrologic cycle early
in this five-decade span or just prior to it (IPCC, 2021), an
intensification which is expected to continue (Seager et al., 2014).
However, we include the full extent of available data in Figure 1 to
provide historical context for recent trends.
To calculate the time dependence of lake levels and NBS components, we
use two non-parametric measures: the Kendall rank correlation
coefficient and the Trend-Free Pre-Whitened Mann-Kendall trend test. We
selected these tests due to the appropriateness of their application on
time series with high levels of autocorrelation; this restriction ruled
out more traditional linear regression techniques (e.g., the Ordinary
Least Squares and conventional Mann-Kendall trend tests) whose
statistical assumptions are not satisfied by the time series under
analysis. The Kendall rank correlation coefficient is a rank-based test
that determines if two variables are related regardless of
distributions. The Trend-Free Pre-Whitened Mann-Kendall trend test is a
modified version of the conventional Mann-Kendall trend test designed to
overcome inaccuracies in analysis due to autocorrelation (Yue et al.,
2002); both the modified and conventional Mann-Kendall trend tests are
valid independent of the distribution of the underlying data. For ease
of reading, we will refer to the Trend-Free Pre-Whitened Mann-Kendall
trend test simply as the Mann-Kendall trend test throughout this study.
In summary, to assess the robustness in trends of interannual
variability, we apply a Mann-Kendall trend test to the 13-year moving
window of standard deviation across lake levels and NBS components over
the past fifty years.