Loading...
HomeMy WebLinkAboutReservoir-Model-Jordan-NCSU_ObenourJordan Lake Reservoir Model Report Prepared for North Carolina Policy Collaboratory Prepared by Dario Del Giudice, Matthew Aupperle, Sankar Arumugam, & Daniel R. Obenour Environmental Modeling to Support Management and Forecasting Group Department of Civil, Construction, & Environmental Engineering North Carolina State University December 2019 Jordan Lake Reservoir Model December 2019 2 Table of Contents Executive Summary ........................................................................................................................ 3 1. Methods....................................................................................................................................... 6 Study Area ................................................................................................................................................. 6 Data ........................................................................................................................................................... 8 Watershed Load Estimation .................................................................................................................... 10 Reservoir Routing .................................................................................................................................... 11 Reservoir Temperature Estimation ......................................................................................................... 11 Nutrient Model Formulation ................................................................................................................... 11 Nutrient Model Calibration ..................................................................................................................... 14 Prior Parameter Information .................................................................................................................. 14 Chlorophyll Modeling .............................................................................................................................. 16 Scenario Analysis ..................................................................................................................................... 17 2. Results ....................................................................................................................................... 19 Total Phosphorus Model ......................................................................................................................... 19 Total Nitrogen Model .............................................................................................................................. 25 Chlorophyll Model .................................................................................................................................. 30 Scenario Analysis ..................................................................................................................................... 34 Comparison to Previous Results ............................................................................................................. 42 3. References ................................................................................................................................. 43 Jordan Lake Reservoir Model December 2019 3 Executive Summary Jordan Lake is a major water supply, flood control, and recreational reservoir located in Chatham County, North Carolina. The reservoir is highly eutrophic based on algal (i.e., chlorophyll a) levels that regularly exceed the state criterion of 40 μg/l. The lake can be separated longitudinally into 4 segments with unique water quality, based on constrictions due to road causeways and natural features. The most upstream (northern) segment receives flows and nutrient loads from watersheds that include the cities of Durham and Chapel Hill. The most downstream (southern) segment receives input from the Haw River watershed, which includes the city of Greensboro. Chlorophyll a levels are particularly elevated where these major tributaries enter the reservoir. There is a general consensus in the scientific and management community that reducing watershed nutrient (nitrogen and/or phosphorus) loads will improve water quality by reducing algal levels over time. However, it is also acknowledged in the scientific literature that internal loading from reservoir bottom sediments can continue to supply nutrients for algal growth even after watershed nutrient loading has been reduced. This phenomenon has been studied in some natural lakes, but has received less attention in man-made reservoirs. Furthermore, the degree of internal nutrient loading is likely to vary substantially among different lakes and reservoirs, considering their unique nutrient loading patterns, sediment characteristics, geometry, and climate. In this study, we develop and apply a water quality model to infer and simulate reservoir nutrient (total phosphorus and total nitrogen) dynamics over a multi-decadal time period (1983-2018). We use a parsimonious mechanistic formulation based on mass balances for the sediments and waters of the four main lake sections. The mechanistic formulation builds on previous modeling studies exploring long-term phosphorus dynamics in natural lakes (Chapra & Canale, 1991; Jensen et al., 2006). The model is calibrated in a Bayesian framework where prior knowledge of biophysical rates from relevant scientific literature is systematically updated based on the long- term calibration datasets for nitrogen and phosphorus in Jordan Lake. Furthermore, empirical relationships are used to relate seasonal nitrogen and phosphorus levels to chlorophyll a. The combined calibrated model is then used to make probabilistic scenario predictions of how the reservoir’s internal nutrient cycling and water quality will respond to potential future reductions in watershed nitrogen and phosphorus loading over time. Modeling results explain 58% and 41% of monthly phosphorus and nitrogen variability, respectively. This level of performance compares well with previous water quality modeling studies (Arhonditsis and Brett, 2004), and higher performance would not necessarily be expected given the stochasticity in sampling results within individual months and reservoir segments. Overall, these results suggest the model is well formulated to address major drivers of nutrient variability. Jordan Lake Reservoir Model December 2019 4 Phosphorus modeling results indicate that there has been a gradual (11%) decrease in phosphorus storage in the reservoir sediment from 1992 till 2018. Phosphorus storage increased in the first decade of the period of record (1983-1992) by 5.5%; a period when watershed nutrient loading was particularly elevated. In the first decade, internal phosphorus loading accounted for just 35.7% of total (internal plus watershed) nutrient loading to the water column. After 1992, internal phosphorus storage declines at a slower rate than watershed loading. As a result of these different reduction rates, in the last decade (2009-2018), internal phosphorus loading accounted for 51% of total loading. Here, total loads do not include the portion of phosphorus load lost in the most upstream portions of the reservoir due to rapid settling and burial of incoming particulate material, which is estimated to be about 46%. Overall, these results suggest that internal nutrient loading currently plays a major role in reservoir eutrophication dynamics and will mute the impact of short-term nutrient watershed loading reductions. Nitrogen modeling results, contrary to phosphorus, demonstrate a 30% increase in nitrogen storage in the reservoir sediment during the study period. In the last decade (2009-2018) average external loading was 7.4% lower than in the first decade (1983-1992), but internal loading increased by 20% during the same interval. Internal nitrogen loading accounted for 71% of total nitrogen loading to the water column from 2009 to 2018, compared to 65% from 1983 to 1992. These total loads do not account for nitrogen lost in upstream portions of the reservoir due to rapid settling of particulate material, which is estimated to be up to 24% of the external load. These results suggest internal nitrogen loading is an important contributor to the reservoir, and will mute the impacts of short-term watershed nutrient loading reductions. The empirical (multiple linear regression) model linking nutrients and chlorophyll explains about 60% of the variability in the chlorophyll data. This performance is satisfactory, considering that the model operates at daily scales, whereas similar models applied at much coarser scales have shown comparable performances (e.g., Dolman and Wiedner, 2015). Calibrated model coefficients generally indicate higher algal concentrations when nutrients (nitrogen and phosphorus) and temperature are high and when flushing is low. The influence of nutrient concentrations appears to be highest in summer (June-September), when water residence time (the inverse of flushing rate) and temperature are high and thus less likely to limit algal growth. The influence of flushing, on the other hand, seems to be highest in winter, when watershed inflows are highest and most variable. The chlorophyll model also provides information on rnp, the total nitrogen to total phosphorus ratio (TN:TP) at which algae switch between nitrogen and phosphorus limitation. We find that, on average, rnp ≈ 16 provides the best fit to the data. This value is substantially higher than the Redfield ratio of 7.2, which indicates that a portion of the nitrogen pool is not easily usable by the algae to foster their growth. Given that observed TN:TP in the reservoir typically ranges from 5 to 30, rnp ≈ 16 suggests that nitrogen and phosphorus are limiting a similar fraction of the time. Interestingly, results also show phosphorus is more limiting than nitrogen in the summer. In the upper (northern) portion of the lake, about 90% of sampled summer days show phosphorus limitation, though the frequency of phosphorus Jordan Lake Reservoir Model December 2019 5 limitation declines in lower sections of the reservoir. In general, the chlorophyll model highlights the importance of reducing both nitrogen and phosphorus in order to reduce algal biomass. The combined models (nitrogen, phosphorus, chlorophyll) are applied to simulate the reservoir’s likely response to potential changes in watershed nutrient loading over time. Because internal nutrient storage and loading may respond slowly to changes in external inputs, we perform these simulations over a 4-decade period. In these simulations, we sample from the variability in historical hydrologic conditions and from the uncertainties in the model itself. Nutrient load reductions are made relative to a 1999-2018 baseline level for hydrology and nonpoint source loading, and present-day point source (wastewater treatment plant) loading rates. If nutrient loading persists at current levels, our results indicate that lake-wide concentrations will continue to change modestly over the next 20 years (+9% w.r.t. 948 μg/l for nitrogen, -5% w.r.t. 59 μg/l for phosphorus, and -2% w.r.t. 30 μg/l for chlorophyll) due to gradual changes in sediment nutrient storage. Results indicate that a 50% external nutrient loading reduction will produce approximately 9% and 5% reductions in lake-wide phosphorus and nitrogen concentrations (respectively) after one year, 25% and 12% after 10 years, and 38% and 17% after 40 years. Further, for the highly eutrophic northern section of the reservoir (above Farrington Road), results suggest it will take about 30 years for a 75% reduction in nutrient loading to reduce the probability of exceeding 40 μg/l chlorophyll to 20% (as an April-October mean, for a given year). For the same scenario and time horizon, there is about 80% probability that concentrations in lower portions of the lake will average below 25 μg/l. When analyzing loading changes to different arms of the lake (Haw River vs. New Hope Creek), we find that nutrient loading reductions to the New Hope Arm of the lake will be most impactful in improving water quality. For instance, 2%, 10%, and 11% reductions in lake-wide mean chlorophyll concentration are expected after 10 years if loads are reduced by 50% to the Haw River Arm, New Hope Arm, and entire lake, respectively. Jordan Lake Reservoir Model December 2019 6 1. Methods Study Area This study focuses on Jordan Lake and its tributaries, including Morgan Creek, New Hope Creek, Northeast Creek, White Oak Creek, and The Haw River. For modeling and analysis purposes, Jordan Lake is divided into four segments based on the locations of flow constrictions, including the Highway 64 causeway, Farrington Road causeway, and the narrows located northeast of the dam (Figure 1). The geometry of lake segments varies substantially, with segment 1 being relatively shallow, and segment 4 having the smallest surface area but deepest mean depth (Table 1). Table 1: Jordan Lake Dimensions at USACE Normal Pool Depth of 65.8 m (216 ft) above sea level. Segments are numbered from north to south. Portion of Lake Normal Pool Area: km2 Normal Pool Volume: km3 Normal Pool Mean Depth: m whole lake 56.44 0.2653 4.70 segment 1 13.97 0.0361 2.59 segment 2 14.43 0.0639 4.43 segment 3 20.25 0.1086 5.36 segment 4 7.78 0.0567 7.28 Jordan Lake Reservoir Model December 2019 7 Figure 1: Map of Jordan Lake displaying the four segments and where they are separated. Black diamonds indicate in-lake sampling points used in this study. Blue squares represent watershed load monitoring sites. Jordan Lake Reservoir Model December 2019 8 Data We compile information from multiple institutions to provide model inputs and calibration data. Reservoir tributary flow data (Table 2, Figure 2) are retrieved from USGS (2019). Nutrient concentration data for the lake and gauged tributaries are compiled from the Water Quality Portal (2019), which provides a compilation of USGS and NC DEQ sampling data. Water quality monitoring sites corresponding to USGS flow gauging sites are shown in Table 2 and Figure 2. Water quality monitoring sites located in the middle of lake segments (used for model calibration) are shown in Table 3 and Figure 1. We developed nutrient loading estimates for wastewater treatment plants from self-reported nutrient concentration and flow data (NCDEQ, Personal Communication, 2019). Table 2: Watershed load monitoring sites directly upstream of Jordan Lake, along with the associated reservoir segment. Flow Site ID Water Quality Site ID Stream Watershed Area (km2) Years Available Reservoir Segment USGS-02096960 B2100000 Haw River 3302 1980-2018 4 USGS-02097517 B3900000 Morgan Creek 106 1983-2018 1 USGS-02097314 B3040000 New Hope Creek 197 1983-2018 1 USGS-0209741955 B3660000 Northeast Creek 55 1983-2018 1 USGS-0209782609 0209782609 White Oak Creek 31 2000-2018 2 Table 3: Lake water quality sampling locations used for analysis. segment 1 segment 2 segment 3 segment 4 B3680000 USGS-0209771550 B4010000 B2453000 B3680020 USGS-0209781125 B4010020 B2453010 B3950000 B3967000 CPF0880A B2453020 B3950020 B3967020 CPF0880Aa CPF055D CPF081A1C B4030000 CPF0880Ab CPF055E CPF086C CPF087B CPF0880Ac CPFMC03 CPF087B3 CPFMC04 Jordan Lake Reservoir Model December 2019 9 Figure 2: Map of Jordan Lake watersheds by color, with ungauged tributaries in hatched zones. Blue squares represent load monitoring sites. We obtained reservoir level and dam release data from the USACE (2019). Daily air temperature, precipitation, and open-water evaporation data are from stations at Raleigh-Durham Airport (ID: KRDU), and Chapel Hill (ID: 311677). Additionally, daily air temperature data come from stations at Raleigh St. University (ID: 317079), Raleigh Apartments (ID: 317069, NCCO, 2019). Data from Chapel Hill, Raleigh Apartments, and Raleigh St. University are averaged to create a more complete and spatially representative record of air temperature over the lake. Reservoir stage-storage relationships (Figure 3) are based on digital elevation data. We gathered elevation data for above the normal pool elevation from the USGS National Elevation Dataset (USGS, 2019) at approximately 30-m resolution. Reservoir bathymetry for below the normal pool comes from a recent University of North Carolina bathymetry survey (A. Rodriguez, personal communication, May 2019) at approximately 25-m resolution. Jordan Lake Reservoir Model December 2019 10 Figure 3: Stage-area (top) and stage-volume (bottom) relationships for the four segments of Jordan Lake. Watershed Load Estimation For each tributary monitoring station, the USGS Weighted Regressions in Time, Discharge, and Season (WRTDS) program (Hirsch & De Cicco 2015) is used to compute nutrient concentration and load on a daily timescale. WRTDS results are aggregated by lake segment and month for input to the reservoir water quality model. Loadings from wastewater treatment plants outside of the gauged area are added directly to the WRTDS segment loading estimates, assuming negligible instream nutrient removal due to their close proximity to the lake. Nonpoint source loading from ungauged watershed areas (Figure 2) are estimated using areal loading rates (kg km-2 month-1) derived from White Oak Creek, as it is the only gauged watershed without WWTP loadings. Missing monthly loads for White Oak Creek and Northeast Creek (due to gaps in the monitoring record) are filled using linear regressions developed between these sites and gauged flow in Northeast Creek and Morgan Creek, respectively. Jordan Lake Reservoir Model December 2019 11 Reservoir Routing A routing model is needed to estimate flows between the four reservoir segments. Routing is performed at a monthly time step, considering flow balances among tributary inflows, over-lake evaporation and precipitation, reservoir discharge, and changes in reservoir storage. Flows from ungauged areas, representing 14.5% of the Jordan Lake watershed, are estimated using drainage area ratios with nearby gauged streams. The overall reservoir flow balance is formulated as: 𝑃𝑜𝑟𝑟=𝑃𝑟+𝑃−𝐸+Δ𝑉+ 𝜀 Eqn 1. where 𝑃𝑜𝑟𝑟 is reservoir outflow, 𝑃𝑟 is the area-adjusted tributary inflow, E over-lake evaporation, P is over-lake precipitation, and Δ𝑉 is the change in lake storage. All of these variables are known, such that the error term (ε) results from inaccuracies in the measured values. Errors are generally small (averaging 3.8 m3/s as absolute values, compared to an average reservoir inflow of 45 m3/s). As tributary gauging is expected to be the most substantial source of uncertainty in Eqn 1, tributary inflows are adjusted up or down slightly to remove the error and close the flow balance. Flows between segments are determined by applying Eqn 1 for each segment with additional term, 𝑃𝑟, to account for flow from the upstream segment. Eqn 1 is then solved for segment outflow, 𝑃𝑜𝑟𝑟, which now represents the flow to the next segment. We note that segment flows are sometimes negative, particularly when there is a large inflow from the Haw River that fills the reservoir from its downstream end. Reservoir Temperature Estimation To provide continuous temperature inputs to the water quality model, daily water temperature is reconstructed using a linear regression between observations of air temperature (available daily) and water temperature (available more sporadically). Daily air temperature readings are used to create two-week moving averages to represent the delayed reaction of a large body of water. The linear regression is formed between these moving averages and water temperature measurements. For this purpose, the water temperature data used for the regression were collected at a depth of 4 meters or greater to better represent the temperature effects occurring deeper in the lake, as the temperature data are used to influence mass transfer rates between the sediment and water layers. Nutrient Model Formulation The mechanistic nutrient model is developed from 8 differential equations, representing nutrient mass balances in the water column and sediment layer of each of the 4 reservoir segments (Figure 4). These equations are comparable to those developed by Chapra (1990) and Jensen et al. (2006) for studying long-term phosphorus dynamics, though these previous studies represented natural lakes as a single well-mixed reactor (rather than multiple segments). The mass balance differential equation for nutrients in the water column is as follows: Jordan Lake Reservoir Model December 2019 12 𝑑𝑀𝑖 𝑑𝑟=(𝑃𝑖∗𝑃𝑟)𝜃𝑅𝑅−20 +(𝑃𝑖𝑛𝑖∗𝐵𝑖𝑛𝑖)(1 −𝜓)−𝑀𝑖(𝑘+𝑟∗𝐴𝑖 𝑉𝑖 )𝜃𝑉𝑅−20 +𝑃𝑖−1,𝑖∗𝑀𝑖−1 𝑉𝑖−1 − 𝑃𝑖,𝑖+1 ∗𝑀𝑖 𝑉𝑖 Eqn. 2.1 In the case with reversed flow (south to north) the last two terms are replaced with the following: 𝑃𝑖+1,𝑖∗𝑀𝑖+1 𝑉𝑖+1 −𝑃𝑖,𝑖−1 ∗𝑀𝑖 𝑉𝑖 Eqn. 2.2 The mass balance differential equation for nutrients in the sediment layer is as follows: 𝑑𝑅𝑖 𝑑𝑟=(𝑀𝑖∗(𝐾𝑟+𝑉𝑟∗𝐴𝑖 𝑉𝑖 ))(𝜃𝑉𝑅−20)−(𝑃𝑖∗𝑃𝑟)(𝜃𝑅𝑅−20)−(𝑃𝑖∗𝐵𝑟) Eqn. 3 The terms in the preceding equations are described as follows: Mi = Mass of phosphorus in segment i water layer. [kg] Si = Mass of phosphorus in segment i sediment layer. [kg] v = Transfer rate of nutrients to the sediment layer, as an effective settling velocity. [m•month-1] k = Transfer rate of nutrients to the sediment layer, as a first order removal rate. [month-1] ψ = Watershed nutrient load adjustment factor. [n/a] R = Recycling rate of nutrients from the sediment back into the water layer. [month-1] B = Removal rate from the sediment to permanent burial and/or denitrification. [month-1] 𝜃𝑟 = Temperature adjustment parameter for the transfer of nutrients from water column to sediments. [n/a] 𝜃𝑅 = Temperature adjustment parameter for the sediment nutrient recycling rate. [n/a] Ai = Area of segment i water layer at the surface (varies with time). [106 m2] Vi = Volume of segment i water layer (varies with time). [106 m3] Qini = Watershed inflow to segment i (varies with time). [106 m3•month-1] Cini = Concentration of nutrient in watershed inflow (varies with time) [mg•m-3] 𝑃𝑖−1,𝑖 = Flow from upstream segment to segment i (varies with time) [106 m3•month-1] 𝑃𝑖+1,𝑖 = Flow from segment i to downstream segment (varies with time) [106 m3•month-1] These differential equations were solved numerically on a monthly time scale using the ODIN package (FitzJohn, 2019) in R (R Core Team, 2018). The model was run from 1983 through the end of 2018. Jordan Lake Reservoir Model December 2019 13 Figure 4: Diagram of the nutrient model used to track flow through Jordan Lake. During the months of May through September, segments 2-4 demonstrate consistent stratification (R. Luttich, personal communication, 2019). To account for such stratification, the model output is adjusted to predict surface concentration (upper 3 m of water column). Surface concentration is calculated from the predicted overall water-column nutrient mass adjusted for the observed difference between bottom and surface water column concentrations: 𝐵𝑟𝑟𝑟=𝑀𝑖[𝑉𝑖+𝑉𝑖,𝑏(𝑃𝑖,𝑏:𝑟−1)]−1 Eqn. 4 where 𝑉𝑖,𝑏 is the volume of segment i below a depth of 3 meters, 𝑃𝑖,𝑏:𝑟 is the nutrient concentration ratio for observations below/above 3 m, for months where consistent differences were observed in the historical data (Mar-Sept for TP, Apr-Sept for TN). These ratios are aggregated as medians to avoid the influence of extreme values. Also, because of the paucity of paired bottom/surface observations, we grouped months together based on an inspection of the typical seasonality of these ratios, in order to achieve more robust and realistic results (Table 4). We note that segment 1 is shallow and does not persistently stratify, so that no adjustment was required. Jordan Lake Reservoir Model December 2019 14 Table 4: Ratios of bottom to surface nutrient concentration, 𝑃𝑖,𝑏:𝑟, by month segment 2 ratios segment 3 ratios segment 4 ratios Month TP TN TP TN TP TN 3 1.1 1.6 1.1 4 1.1 1.0 1.6 1.4 1.1 1.2 5 1.1 1.0 1.6 1.4 1.1 1.2 6 1.1 1.0 1.6 2.2 1.1 1.7 7 1.2 1.0 1.6 2.2 1.7 1.7 8 1.2 1.0 1.6 2.2 1.7 1.7 9 1.2 1.0 1.6 1.4 1.7 1.2 Nutrient Model Calibration The mechanistic model is calibrated in a Bayesian framework to produce probabilistic estimates of key model parameters. Calibration is achieved by updating prior rate and nutrient information from previous research using the observed surface concentration data gathered from Jordan Lake over the past 30 years. The parameters calibrated in the model are listed below in Table 5. Bayesian model calibration is implemented using an adaptive Markov Chain Monte-Carlo approach as implemented in the “adaptMCMC” package in R (Scheidegger 2018). The model is calibrated in a log10 transformed scale to help account for the right skew typical of pollutant concentration data (Ott, 1990). Parameters calibrated with this approach are centered at their optimal value and additionally incorporate a measure of their uncertainty. Prior Parameter Information Before model calibration, prior distributions for model parameters are formulated from existing literature and knowledge of the system. Most priors are Gaussian (normal) and truncated at zero to avoid unrealistic negative values (Table 5). Priors for nutrient removal rates are developed based on previous literature exploring internal phosphorus cycling in lakes. Effective phosphorus settling rates have been estimated to be within 1 to 4 m/mo (Chapra, 1975; Chapra & Canale, 1981; Nürnberg, 1984; Jensen et al., 2006). Converting these settling rates to 1st-order removal rates based on typical lake depths indicates a range of 0.1 to 0.7 mo-1. Because we include both options within the model, our prior means are half of the center of these ranges: 1.3 m/mo and 0.20 mo-1 for v and k, respectively. Prior standard deviations are set to 1.4 m/mo and 0.26 mo-1, respectively, so that that the upper 0.975 quantile of the prior distributions align with the upper end of the literature ranges. In addition, the load adjustment factor ψ is assigned a mean of 0.1 and standard deviation of 0.1, considering the potential for bias in loading estimates (Hirsch et al., 2014), and the potential for removal of particulate P in the upstream reaches of the reservoir (Duan et al., 2014; River & Richardson, 2018). For eutrophic lakes, reported phosphorus fluxes typically range from around 0.2 to 3.0 g/m2/yr (Nürnberg, 1988; Chapra & Canale, 1991; Moore et al., 1998; Welch & Jacoby, 2001; Haggard Jordan Lake Reservoir Model December 2019 15 et al., 2005; Nürnberg, 2009; Matisoff et al., 2016). Converting these fluxes to release rates requires information on sediment phosphorus content, reported at typical ranges of 20 to 80 g/m2 for eutrophic lakes (Larsen et al., 1981; James et al., 2005; Jensen et al., 2006). Based on these values, the prior recycle rate is determined to have a mean of 0.0030 mo-1 with standard deviation of 0.0018 mo-1. There is little prior information available for the phosphorus burial rate, but it is expected to be smaller than the recycle rate, as it is not included in all models (Jensen et al., 2006). Here, we use a prior with mean and standard deviation of 0.001 mo-1, based loosely on Chapra & Canale (1991). Table 5: Priors Applied to Parameters Subject to Calibration. Parameter Distribution Lower Upper Description Units Phosphorus Model Priors ψ N(0.1,0.1) -1 1 watershed load reduction factor n/a R N(0.003,.0018) 0 Inf 1st order sediment recycling rate mo-1 B N(0.001,0.001) 0 Inf 1st order burial rate mo-1 Sinit N(0.0456,0.025) 0 Inf initial sediment concentration kg•m-2 Nitrogen Model Priors ψ N(0.05,0.05) -1 1 watershed load reduction factor n/a R N(0.033,.017) 0 Inf 1st order sediment recycling rate mo-1 B N(0.15,0.1) 0 Inf 1st order burial/denitrification rate mo-1 Sinit N(.25,0.075) 0 Inf initial sediment concentration kg•m-2 General Priors v N(1.3,1.4) 0 Inf effective settling velocity to sediment m•mo- 1 k N(0.2,0.26) 0 Inf 1st order transfer rate to sediment mo-1 θv N(1.05,0.03) 1 Inf settling/transfer temperature adjustment n/a θR N(1.05,0.03) 1 Inf recycling temperature adjustment n/a Jordan Lake Reservoir Model December 2019 16 Prior information for nitrogen water column removal rates is less available than for phosphorus. However, because nitrogen removal from the water column is mediated by largely the same processes to phosphorus removal (e.g., phytoplankton growth, consumption, settling of detritus), we apply the same priors for v and k used in the phosphorus model. We note that 1.1 m/mo has been reported as an average apparent settling rate for reservoirs (Harrison et al., 2009) and is generally consistent with our prior, though it is a net rate reflecting both nitrogen settling and recycling. For the load adjustment factor, we assign a prior mean and standard deviation of 0.05 to account for potential settling of particulate organic matter in the upper reaches of the reservoir as well as potential biases in load estimation. This prior is notably lower than for phosphorus, as nitrogen does not adsorb to settling inorganic particulate matter. Due to denitrification, many long-term studies assume little or no buildup of nitrogen within the sediment layer over time (e.g., Jensen, 1992; David et al., 2006). Lake and reservoir denitrification rate estimates tend to vary widely, from about 10 to 800 g/m2/yr in literature (Windolf et al., 1996; Tomaszek & Czerwieniec, 2000; David et al., 2006). Based on typical areal TN concentrations from 100 to 400 g/m2 (Lane & Koelzer, 1943; Tomaszek & Czerwieniec, 2000; Fisher et al., 2001; James et al., 2005), these values suggest a denitrification rate with a prior mean of 0.15 mo-1 with standard deviation of 0.10 mo-1. In the nitrogen model, this denitrification rate replaces the phosphorus burial rate and is mathematically equivalent to it. As an alternative to denitrification, sediment nitrogen may also be returned to the water column as ammonia (Di Toro, 2001). Reported summer ammonia flux rates typically vary from 5-25 g/m2/yr for hypoxic summer conditions (Graetz et al., 1973; Beutel, 2001). On a yearly basis, these values suggest a recycle rate prior with a mean of 0.033 mo-1 with standard deviation of 0.017 mo-1. Chlorophyll Modeling In-lake concentrations of total nutrients (TN and TP) are used to predict chlorophyll via regression. The model is composed of three terms that are linearly combined to predict chlorophyll concentration (c), namely, the nutrient concentration term, the flushing rate (F, the inverse of water residence time) and water temperature (T): 0log( ) log min , log( ) log( )np f t np TNc TP F Tr    = + + + , Eqn. 5 where ’s and rnp are parameters to be calibrated using ‘penalized least squared’ optimization implemented using the ‘penalized’ package (Goeman, 2010). This approach is equivalent to maximum likelihood estimation with np and t being constrained to be nonnegative, as negative coefficients for nutrients and temperature would be physically implausible. The nutrient concentration term of the chlorophyll model is based on the ‘equivalent nutrient approach’ proposed by Dolman and Wiedner (2015). This approach represents chlorophyll as being controlled by the nutrient in shortest supply, either TN or TP. To account for the fact that algae require more TN than TP to grow, TN is divided by rnp, a calibration parameter Jordan Lake Reservoir Model December 2019 17 representing the TN:TP at which TN and TP are equally limiting. rnp is constrained between 5 and 30, a range consistent with the observed TN:TP in Jordan Lake and with previous literature (Dolman and Wiedner 2015). Calibration is conducted using daily input and output data derived from lacustrine measurements, except for the flushing rate that was derived from monthly outputs of the routing model (see Reservoir Routing section). Each predictor term as well as the dependent variable are log10 transformed, a solution often adopted in regression modeling of chlorophyll to mitigate the impact of extremely high values (Dolman and Wiedner 2015; Filstrup and Downing, 2017; Prairie et al. 1989). Additionally, the most extreme outliers are removed from model calibration to avoid corrupting parameter estimates. Overall, 12 (out of 1020) daily data points are dropped because (transformed) chlorophyll, TN or TP were above their upper quartile plus twice the interquartile range or below their lower quartile minus twice the interquartile range. The model is allowed to calibrate different parameters according to the season and segment, as preliminary analyses show non-stationarity in the relationship between predictors and dependent variable based on the time of the year and the portion of the lake. Scenario Analysis The three developed models (for TP, TN and chlorophyll) are combined to predict future water column and sediment concentrations at monthly scales. Specifically, the combined model is used to simulate the reservoir’s response in terms of TP, TN, and chlorophyll concentrations to past and future nutrient loading. While model calibration is performed over the monitoring period of 1983-2018, simulations for scenario analysis are conducted over the 2019-2058 period. To realistically account for seasonality in hydrology over this latter period, for each simulation, an input time series is created by randomly rearranging the historical (1999-2018) inputs by year. This historical period is chosen as representative because around the late 1990s nutrient concentrations in the lake became approximately stationary. To account for two major and sudden reductions in nitrogen loading due to improvements in sewage treatment, we adjust the concentration of segment 1 riverine inputs used for scenario analysis. Specifically, during preliminary data analysis, we find that, after December 2004, Northeast Creek average TN loading decreases by 58%, whereas after March 2010, Morgan Creek average TN loading decreases by 28%. Therefore, for scenario analysis, riverine inputs from 1999 to 2010 are reduced accordingly, in order to not project forward loading conditions that irreversibly changed. The main nutrient loading scenarios are calculated by varying tributary concentrations from - 100% to +100% of historical values. Ancillary scenarios are also computed by varying of ±50% loadings from either the New Hope Arm or the Haw River Arm. To account for model parameter uncertainty for a given loading scenario, each model is run 1000 times, each time with a different sample from the calibrated posterior parameter distribution. Additionally, for each simulation, realizations of the residual error distributions are added to TP, TN, and chlorophyll to account Jordan Lake Reservoir Model December 2019 18 for model structure and observation uncertainty. Propagating parameter and residual uncertainties through the three models is paramount to quantify the probability of exceeding a given concentration in the future. Jordan Lake Reservoir Model December 2019 19 2. Results Total Phosphorus Model Overall, 58% of the variability of the data is explained by the model (Figure 5). The coefficient of determination (R2 or fraction of variability in the data explained by the model) is equal to 0.42, 0.25, 0.30, and 0.23 for segments 1, 2, 3, and 4, respectively. The TP model effectively captures the higher levels of recorded phosphorus from 1983 to 1990, and the generally lower TP concentrations thereafter (Figure 6). In addition, the model captures intra-annual variability, with high concentrations generally occurring in the fall-winter. Figure 5: In-lake surface TP observations versus predictions (R2 = 0.58). Jordan Lake Reservoir Model December 2019 20 Figure 6: Monthly time series showing the TP model predictions along with monthly observations for segment 1 (above Farrington Rd). Posterior distributions generally lie within the prior distributions suggested by previous literature (Figure 7). Posterior distributions are generally narrower than prior distributions, indicating substantial reductions in uncertainty after assimilating calibration data. We note that the posterior distribution for ψ is substantially shifted from the prior, indicating more phosphorus loss (i.e., settling and permanent burial) in upstream portions of the reservoir than expected and/or a bias in WRTDS loading estimates. Additionally, many tributaries (e.g., Morgan Creek, New Hope Creek, and Northeast Creek) also pass through wetlands and impoundments before entering the lake, which can contribute to TP removal. On the other hand, the burial rate, B, is quite small, indicating less permanent phosphorus removal within the deeper main body of the lake. In the main body of the lake, phosphorus appears to be efficiently recycled, as indicated by the relatively high R. Jordan Lake Reservoir Model December 2019 21 Figure 7: Comparison between prior (dashed) and posterior (solid) distributions for the TP model. The y-axis represents probability density. Prior distributions reflect knowledge gathered from previous studies. Posterior distributions represent calibrated parameters. The total phosphorus model shows an overall reduction in the water column concentration of phosphorus over the course of the 35 year study period (Figure 8). Water column concentrations drop by 50%, 29%, 26%, and 40% in segments 1 through 4, respectively, when comparing the first (1983-1992) and last (2009-2018) decades of the study period. Sediment concentrations increase during the first decade of study and then, starting from the early 1990s, they gradually decrease. Sediment concentrations drop by 8% when comparing the first and last decades. This concentration reduction is more prominent in segments 1 and 4, which are most exposed to the inflow of TP from the watershed than the innermost segments. The long term trend of sediment TP is important to note, as it will influence the rates of internal loading in the future (see section on scenario analysis). The majority of the phosphorus enters the lake from the watershed during the beginning of the period, but over time, the majority of the loading to the lake (i.e., water column) starts to come from internal recycling from the sediment layer in most years (Figure 8, top). This change in proportion is mostly due to reduced watershed TP loading and, to some extent, to increased sediment loading caused by an increase in sediment TP content in the beginning of the period. Years with higher inflow generally have higher loadings, especially in segment 4 (Figures 8 and 9). Segment 4 has the strongest correlation between flushing rate and inflow, as expected due to the dominance of the Haw River. Jordan Lake Reservoir Model December 2019 22 Figure 8: In the Top plot, mass flow of Total Phosphorus for Jordan Lake overall is represented. In the bottom plot, total water column concentration of TP over the course of the entire period, as well as average areal concentration of TP in the sediments are shown. Jordan Lake Reservoir Model December 2019 23 Figure 9: Average yearly flushing rate for individual segments using (left y-axis), calculated in terms of how many times that segment’s volume flowed through during the period. Also plotted is the total inflow to the reservoir (right y-axis). The reservoir demonstrates substantial intra-annual variability in TP gains and losses (Figure 10). Internal loading is the primary contributor of TP during summer months, whereas watershed loading dominates during the winter. Higher settling rates occur during the summer, due to the temperature response parameters. External watershed loads and releases of nutrients through the dam follow the same pattern as inflow and flushing rates (Figure 11). Jordan Lake Reservoir Model December 2019 24 Figure 10: Top: Average monthly mass transfer rates for the lake as a whole. Bottom: Average monthly segment concentrations; dotted lines and empty symbols represent surface concentrations during stratified periods for segments 2-4. Figure 11: Average monthly flushing rate for individual segments, calculated in terms of how many times that segment’s volume flowed through during the period. Jordan Lake Reservoir Model December 2019 25 Total Nitrogen Model Overall, 41% of the variability of the TN data is explained by the model (Figure 12). For individual segments of the lake the coefficient of determination is equal to 0.21, 0.40, 0.34, and 0.53 for segments 1, 2, 3, and 4, respectively. The model generally captures temporal variability in TN concentrations well. However, for segment 1 (lowest R2), the model tends to over-predict TN concentrations in the 1990s (Figure 13). Figure 12: In-lake surface TN observations versus predictions (R2 = 0.41). Jordan Lake Reservoir Model December 2019 26 Figure 13: Monthly time series showing the TN model predictions along with monthly observations for segment 1. Estimated model parameters (Figure 14) generally fall within the ranges indicated by the prior distributions. The relatively high effective settling rate and low first order transfer rate generally suggest that nutrient transfer to the sediment is better represented as a settling term. Parameter B represents permanent burial and/or denitrification, which is estimated to be nearly zero, indicating denitrification is not a significant process within the main body (i.e., center) of the lake. The posterior distribution for ψ is substantially larger compared to the prior distribution, indicating nitrogen loss occurs before entering the main body of the lake and/or a bias in WRTDS loading estimates. The high ψ posterior distribution is likely explained as an initial settling of particulate matter due to water slowing down before it enters the body of the lake, a trend documented before (Jiao et. al. 2018). Also, as noted for phosphorus, many tributaries (e.g., Morgan Creek, New Hope Creek, and Northeast Creek) pass through wetlands and impoundments before entering the lake, which can contribute to burial and denitrification. At the same time, the ψ for the TN model is substantially lower than for the TP model, and together with the low B, partially explains why there is more accumulation of TN than TP in reservoir sediments over the study period. Jordan Lake Reservoir Model December 2019 27 In general, inferred TN and TP sediment rates (Figures 7 and 14) are at least qualitatively consistent with recent nutrient flux measurements (M. Piehler, personal communication, 2019). These flux measurements indicated the sediment nitrogen releases (primarily as ammonia) were about 40 times higher than phosphate releases, indicating that nitrogen is more efficiently recycled from the sediment than phosphorus. Also, measured releases of nitrogen gas, which indicate denitrification, were minimal except for at a location in the upstream Haw River portion of segment 4. This is consistent with the minimal estimate for B and substantially negative estimate of ψ, which also indicate that permanent nitrogen losses are largely limited to around the mouths of the tributaries. Further comparisons with measured sediment nutrient concentrations and fluxes will be explored in the future. Figure 14: Comparison between prior (dashed) and posterior (solid) distributions for TN model. The y-axis represents probability density. Prior distributions reflect knowledge gathered during the literature review process. Posterior distributions represent calibrated parameters. There is a higher proportion of internal loading for TN (Figure 15) than for TP (Figure 8). Over 1983-2018, internal loading contributes 71% and 53% of total loading for these two nutrients, respectively. There is also an accumulation of sediment TN over the course of the entire study period, while sediment TP begins to decline after the first decade. Also, where the TP model shows notable decreases in average water column concentration, the TN model shows minor increases of 1.7%, 6.8%, 1.4%, and 6.0% in segments 1 through 4, respectively, when comparing the last decade to the first decade of the study period (2009-18 vs. 1983-92). It is important to note that there was no substantial external watershed loading reductions in TN, unlike TP, which Jordan Lake Reservoir Model December 2019 28 explains part of the difference in water column concentration changes over the study period. Similarly, sediment TN concentrations increase by 20% over the same period. Segment 4 again shows the strongest relationship between higher loading on a yearly scale and higher inflow (Figures 15 and 9). Figure 15: In the top plot, mass flow of Total Nitrogen for Jordan Lake overall. In the bottom plot, total water column concentration of TN over the course of the entire period, as well as average areal concentration of TN in the sediments. Internal loading of TN is a substantial contributor during all months (Figure 16). Internal loading is only overtaken by watershed loading during the early months of the year when watershed loading is the highest, aligning with high inflow and flushing rates (Figure 11). There is a much greater difference in surface concentration relative to total water column concentration during the summer months in the segments 3 and 4, which reflect stratification behavior (Figure 16, bottom). This is consistent with the higher rates of internal loading of TN, as this will create greater differences between concentrations in the epilimnion versus the hypoliminion. Jordan Lake Reservoir Model December 2019 29 Figure 16: Top: Average monthly mass transfer rates for the lake as a whole. Bottom: Average monthly segment concentrations; dotted lines and empty symbols represent surface concentrations during stratified periods for segments 2-4. Jordan Lake Reservoir Model December 2019 30 Chlorophyll Model The chlorophyll regression model with three predictor variables explains about one third of the variability in daily chlorophyll for a given segment and season (see Figure 17 for examples of summer). However, when considering the entire dataset, the model is able to explain more than half of the observed variability (R2=0.59). For most of the 16 cases considered (4 seasons × 4 segments) we find that chlorophyll is positively related to in-lake nutrients (nitrogen or phosphorus, depending on which one is limiting) and temperature, whereas a mostly negative relationship is found with segment flushing rate (Table 6). This is consistent with previous literature, as higher temperatures and nutrients and low flushing are known to favor algal growth (Paerl and Otten, 2013). Interestingly, in some cases such as wintertime in segment 4, nutrients do not show significant positive relationships with nutrients. This implies that, in those conditions, variability in algal biomass is mainly controlled by hydrometeorology. The chlorophyll model also provides information on rnp, the total nitrogen to total phosphorus ratio (TN:TP) at which algae switch between nitrogen and phosphorus limitation. We find that, on average, rnp ≈ 16 provides the best fit to the data, although the exact best value varies by segment and season. In all cases except one (Table 6), the value of rnp is substantially higher than the Redfield ratio of 7.2, which indicates that a portion of the nitrogen pool is not easily usable by algae to stimulate growth. Interestingly, rnp = 16 is approximately equal to the median of the observed TN:TP, which suggests that in Jordan Lake nitrogen and phosphorus are limiting a similar fraction of the time. Additionally, for all segments, rnp reaches a minimum in summer, which suggests that in summer phosphorus is more limiting than nitrogen. Specifically, the model results indicate that about 80% of summer days experience some degree of phosphorus limitation. Consequently, summer is the season when the importance of phosphorus for controlling algal biomass is high, relative to other seasons and to nitrogen. Additionally, segment 1 appears to be limited by phosphorus most frequently (about 90% of sampled summer days) whereas segment 4 is phosphorus limited the least (about 50% of sampled summer days). In general, however, results of chlorophyll model calibration highlight the importance of reducing both nitrogen and phosphorus in order to reduce algal biomass throughout the year and throughout the entire lake. Jordan Lake Reservoir Model December 2019 31 Figure 17: Observed vs predicted daily log-transformed chlorophyll concentration during summer (June-September). The diagonal lines represent the 1:1 line. The model captures reasonably well the daily variability in different portions of the lake using only three predictors or less. Jordan Lake Reservoir Model December 2019 32 Table 6: maximum likelihood estimates of the regression model parameters calibrated to different segments (1 to 4) and seasons. Note that in cases where the nutrient term is not significantly positive the ratio rnp is not applicable. The units of the β coefficients are (from the first to the fourth row) log(μgchlorophyll/l) divided by: 1, log(μgphoshorus/l), log(1/mo), and log(°C). 1,winter 1,spring 1,summer 1,autumn 2,winter 2,spring 2,summer 2,autumn 0 -0.72 0.58 -0.58 0.35 -0.14 -0.44 -3.69 0.70 np 0.86 0.53 1.01 0.61 0.65 1.14 0.88 0.15 f -0.23 -0.11 0 -0.11 -0.06 0.01 0.05 -0.11 t 0.72 0.06 0.22 0.18 0.54 0 2.60 0.43 npr 16 14 5 12 24 18 11 21 3,winter 3,spring 3,summer 3,autumn 4,winter 4,spring 4,summer 4,autumn 0 0.88 -0.12 -2.69 -1.22 1.05 1.23 -0.13 -0.41 np 0 0.9 0.92 0.92 0 0 0.62 0.20 f -0.04 -0.10 -0.04 -0.06 -0.24 -0.28 -0.14 -0.11 t 0.47 0 1.83 0.88 0.22 0.19 0.41 1.15 npr NA 20 13 17 NA NA 15 25 Combining the chlorophyll model with the TP and TN models allows us to predict chlorophyll for each month rather than only when nutrient measurements are available (see Figure 18 for the example of segment 1). The combined model also enables us to analyze effects of changing riverine inputs (see next section). The model captures well the high chlorophyll events such as the one observed in July 1986 and enables us to reconstruct other probable algal blooms, such as the one of September 2002, which had not been monitored (Figure 18). The combined models provide robust chlorophyll prediction. We note that the stand-alone chlorophyll regression model was developed using actual nutrient observations, and that substitution of modeled nutrient values would be expected to degrade predicative performance. Still, for the entire dataset, the combined model explains a large fraction of the observed variance (R2=0.44). This testifies to the robustness of the underlying nutrient models and the overall modeling chain required to predict chlorophyll. Jordan Lake Reservoir Model December 2019 33 Figure 18: Monthly time series showing the output combined TN-TP-chlorophyll model along with monthly observations for segment 1. Jordan Lake Reservoir Model December 2019 34 Scenario Analysis The combined model is used to assess future response of the reservoir to a broad range of changes in riverine nutrient (TN and TP) inputs. Considering the case with nutrient loading persisting at approximately current levels, model results show that the decreasing trend of TP observable in Figure 8 is going to continue. For instance, after 20 years, average lacustrine concentrations are expected to be 5% lower than the recent historical average (59 μg/l for 1999- 2018). However, halving phosphorus loads would facilitate larger lake-wide TP concentration reductions both immediately (due to halved watershed contributions) and gradually (due to a progressive decrease in sediment loading), as shown in Figure 19 for segment 1. Considering the entire lake, a 50% reduction in incoming TP after 20 years would lead to lacustrine concentrations about 30% lower than in the historical period. Nitrogen is expected to behave somewhat differently than phosphorus. Under a scenario of no loading changes, in 2038 TN concentrations are expected to be 9% higher than the recent historical average of 948 μg/l. This increase in lacustrine concentrations is explainable with the continued gradual accumulation of nitrogen in the sediments (Figure 15), which also implies a gradual increase in internal N loading. Given this increasing internal input, halving external (watershed) inputs would lead to a less steep reduction in TN concentrations when compared to TP (cf. Figures 19 and 20). Specifically, for the entire lake, reducing external loads of 50% would lead, after 20 years, to TN concentrations only 15% lower than in the period 1999-2018. Jordan Lake Reservoir Model December 2019 35 Figure 19: Probabilistic time series of predicted phosphorus for the recent historical period and, separated by the dashed line, the first two decades of the projection period. In this example, riverine nutrient inputs from all watersheds are reduced of 50%. Prediction intervals incorporate the effect of model/data and hydrologic uncertainty. Jordan Lake Reservoir Model December 2019 36 Figure 20: Probabilistic time series of predicted nitrogen for the recent historical period and, separated by the dashed line, the first two decades of the projection period. In this example, riverine nutrient inputs from all watersheds are reduced of 50%. Prediction intervals incorporate the effect of model/data and hydrologic uncertainty. The interplay of opposite trends in nitrogen and phosphorus internal loads is largely going to balance out in terms of chlorophyll concentration. Specifically, if future riverine nutrient loads are maintained at recent historical levels, chlorophyll is expected to approximately remain at its historical mean (30 μg/l, annual) for the first decade. However, over a longer time horizon, chlorophyll is expected to decrease slightly (Table 7), due to a gradual depletion in sediment TP. Given this mild decreasing trend, changes in external loads of the same magnitude but opposite sign are going to lead to slightly asymmetric changes in chlorophyll concentration. For instance, a load change of +50% is only going to lead to a chlorophyll change of +11% after 20 years, compared to -15% with a -50% load change. Jordan Lake Reservoir Model December 2019 37 Table 7: Percentage of average changes in lake-wide chlorophyll concentration. The rows represent the loading scenarios, while the columns indicate time periods. The values in column 1 indicate the percentage changes in loading for all watersheds, except for the last four rows which show the effect of modifying the loads only for the New Hope (NH) Arm or Haw River (HR) Arm. Loading adjustment Historical 2019 2028 2038 2058 -100 0 -14 -23 -30 -41 -75 0 -9 -17 -22 -30 -50 0 -6 -12 -15 -21 -25 0 -3 -6 -8 -12 0 0 1 -1 -2 -4 25 0 4 4 4 4 50 0 7 9 11 12 75 0 10 14 16 19 100 0 13 18 22 27 -50NH 0 -5 -10 -13 -19 50NH 0 6 7 8 10 -50HR 0 0 -2 -4 -7 50HR 0 2 1 0 -2 Changing loads in the New Hope (NH) Arm would have a much larger impact on lake-wide chlorophyll than changing loads in the Haw River (HR) Arm (Table 7), even though NH is only 20% the size of the total (NH + HR) watershed and it only contributes roughly 20% of total incoming nutrients (see Watershed Modeling Report). Taking a 50% load reduction and a 20 year horizon as an example, changes in NH loads are predicted to cause a chlorophyll change of - 13% whereas changes in HR loads would only lead to -4%. The higher impact of changes in NH loads can be explained by the following considerations. First, the main pool of segment 4 has lower chlorophyll (22 μg/l) compared to the lake average (30 μg/l) and segment 1 (46 μg/l). Therefore, chlorophyll changes in segment 4 will have a lower impact on overall lake chlorophyll compared to the same percentage change for segment 1. Additionally, in Tables 8 and 9, it is evident that meaningful changes in segment-specific chlorophyll only occur when loads from directly contributing tributaries are reduced. For instance, changing NH loads by ±50%, which have direct impacts on segments 1-3, will produce no significant change in segment 4. Finally, regression results (Table 6) show that chlorophyll in segment 4 is less sensitive to nutrient changes than chlorophyll in segment 1. As a result, for a -50% change in overall loading, after 20 years segment 1 will experience 21% lower chlorophyll (Table 8), whereas segment 4 chlorophyll will only be reduced by 9% (Table 9). Time series of chlorophyll predictions further illustrate these considerations. Specifically, Figure 21 shows that segment 1 has both high chlorophyll and marked responsiveness to load reduction. Segment 4 (Figure 22), Jordan Lake Reservoir Model December 2019 38 on the other hand, has lower chlorophyll and a barely perceptible reduction in chlorophyll even after a 50% reduction in riverine loads. Table 8: Percentage of average changes in chlorophyll for segment 1. The rows represent the loading scenarios, while the columns indicate time periods. The values in column 1 indicate the percentage changes in loading for all watersheds, except for the last four rows which show the effect of modifying the loads only for the New Hope (NH) Arm or Haw River (HR) Arm. Loading adjustment Historical 2019 2028 2038 2058 -100 0 -18 -33 -42 -55 -75 0 -12 -24 -31 -41 -50 0 -7 -16 -21 -28 -25 0 -2 -8 -11 -16 0 0 4 -1 -2 -5 25 0 9 7 6 6 50 0 13 13 15 16 75 0 18 20 22 26 100 0 22 26 31 36 -50NH 0 -6 -16 -21 -28 50NH 0 13 13 14 17 -50HR 0 4 0 -3 -5 50HR 0 4 -1 -2 -5 Table 9: Percentage of average changes in chlorophyll for segment 4. The rows represent the loading scenarios, while the columns indicate time periods. The values in column 1 indicate the percentage changes in loading for all watersheds, except for the last four rows which show the effect of modifying the loads only for the New Hope (NH) Arm or Haw River (HR) Arm. Loading Adjustment Historical 2019 2028 2038 2058 -100 0 -12 -18 -23 -30 -75 0 -9 -13 -16 -19 -50 0 -5 -8 -9 -11 -25 0 -3 -3 -4 -5 0 0 1 0 1 0 25 0 3 4 4 5 50 0 5 7 9 10 75 0 8 10 12 14 100 0 9 13 16 18 -50NH 0 0 0 0 0 50NH 0 0 0 1 1 -50HR 0 -5 -7 -9 -11 50HR 0 5 7 8 9 Jordan Lake Reservoir Model December 2019 39 Figure 21: Probabilistic time series of predicted chlorophyll a in segment 1 for the recent historical period and, separated by the dashed line, the first two decades of the projection period. In this example, riverine nutrient inputs from all watersheds are reduced of 50%. Prediction intervals incorporate the effect of model/data and hydrologic uncertainty. Jordan Lake Reservoir Model December 2019 40 Figure 22: Probabilistic time series of predicted chlorophyll a in segment 4 for the recent historical period and, separated by the dashed line, the first two decades of the projection period. In this example, riverine nutrient inputs from all watersheds are reduced of 50%. Prediction intervals incorporate the effect of model/data and hydrologic uncertainty. Besides elucidating mean trends, the probabilistic model projections in this study take into account parameter, hydrologic, model structure, and observation uncertainty. Consequently, these projections enable us to quantify the probability of meeting the chlorophyll state criterion of 40 μg/l (NCDENR, 2007) for a variety of loading scenarios and temporal horizons. Results of these probabilistic analyses, reported in Table 10, are centered on segment 1, which is the portion of the lake with the highest chlorophyll concentrations (46 μg/l during the historical period). In general, Table 10 shows that substantial load reduction and multiple decades are necessary for segment 1 to consistently achieve chlorophyll concentrations below the criterion with high probability. For instance, with a 75% load reduction, it would take four decades to attain April- October chlorophyll of 40 μg/l or lower with 86% confidence. For the same temporal horizon but only a 25% reduction, the state criterion would instead be exceeded with about 72% probability. Jordan Lake Reservoir Model December 2019 41 Table 10: Probability that mean April-October chlorophyll in segment 1 does not exceed the state criterion of 40 μg/l in any given year. The rows represent the loading scenarios, while the columns indicate time periods. The values in column 1 indicate the percentage changes in loading for all watersheds, except for the last four rows which show the effect of modifying the loads only for the New Hope Arm or the Haw River Arm. Loading Adjustment Historical 2019 2028 2038 2058 -100 0.05 0.47 0.69 0.84 0.98 -75 0.05 0.3 0.49 0.66 0.86 -50 0.05 0.17 0.29 0.40 0.62 -25 0.05 0.09 0.13 0.18 0.28 0 0.05 0.04 0.06 0.09 0.12 25 0.05 0.02 0.02 0.02 0.02 50 0.05 0.00 0.01 0.01 0.00 75 0.05 0.00 0.00 0.00 0.00 100 0.05 0.00 0.00 0.00 0.00 -50NH 0.05 0.17 0.30 0.39 0.64 50NH 0.05 0.01 0.01 0.00 0.00 -50HR 0.05 0.05 0.05 0.06 0.10 50HR 0.05 0.04 0.04 0.09 0.10 Jordan Lake Reservoir Model December 2019 42 Comparison to Previous Results The last major policy informing modeling project for Jordan Lake was performed by TetraTech (2003), based on a study period of 1997 to 2001. This study used the Water Analysis Simulation Program (WASP) for modeling water quality and the Environmental Fluid Dynamics Code (EFDC) for hydrodynamics. Although the model was calibrated to a shorter time period, it had higher spatio-temporal resolution, and represents eutrophication processes with increased mechanistic detail. However, there was no focus on long-term sediment nutrient dynamics or on probabilistic model calibration. Additionally, in our work we have done a more complete analysis of error and model performance relative to observed quantities. To compare between this study and our results, TetraTech model segments 1-4, 5-8, 9-13, and 14-15 map to segments 1, 2, 3, and 4 of this study, respectively. TetraTech (2003) reported the required watershed reductions in TP and TN to reach a 10% or lower frequency of chlorophyll a concentrations above 40 μg/l during the “growing season” of May-September. Both our results (Table 8) and TetraTech confirm that Haw River reductions will have minimal influence on water quality in the highly eutrophic segment 1. Based on the TetraTech report, a 50% reduction of both TN and TP would instantaneously reduce the frequency of days with chlorophyll a above 40 μg/l to 10% or less in segment 1. From our own analysis, a 50% reduction in loading is expected to produce only a 17% chance of reducing the average April-October chlorophyll a concentration to below 40 μg/l at first. However, after 40 years of sustained 50% reductions, the probability would increase to 62% (Table 10). While the two studies assessed somewhat different improvement goals, both studies indicate that reducing external nutrient loading reductions can substantially improve water quality in the reservoir. The primary difference is that our study indicates such improvements will likely take decades to fully realize. Jordan Lake Reservoir Model December 2019 43 3. References Beck, M. B. (1987). Water quality modeling: A review of the analysis of uncertainty. Water Resources Research, 23(8), 1393–1442. https://doi.org/10.1029/WR023i008p01393 . Beutel, M. W. (2001). Oxygen consumption and ammonia accumulation in the hypolimnion of Walker Lake, Nevada. In Saline Lakes (pp. 107-117). Springer, Dordrecht. Billen, G., Beusen, A., Bouwman, L., & Garnier, J. (2010). Anthropogenic nitrogen autotrophy and heterotrophy of the world’s watersheds: Past, present, and future trends. Global Biogeochemical Cycles, 24(2), 1–12. https://doi.org/10.1029/2009GB003702 . Chapra, S. C. (1997). Surface Water-Quality Modeling. Waveland Press, Inc. Chapra, S. C., & Canale, R. P. (1991). Long-term phenomenological model of phosphorus and oxygen for stratified lakes. Water Research, 25(6), 707–715. https://doi.org/10.1016/0043- 1354(91)90046-S Chapra, S. C., & Canale, R. P. (1991). Long-term phenomenological model of phosphorus and oxygen for stratified lakes. Water Research, 25(6), 707–715. https://doi.org/10.1016/0043- 1354(91)90046-S Di Toro, D. M. (2001). Sediment flux modeling (Vol. 116). New York: Wiley-Interscience. Dolman A.M. & Wiedner C. (2015) Predicting phytoplankton biomass and estimating critical N : P ratios with piecewise models that conform to Liebig’s law of the minimum. Freshwater Biology, 60, 686–697 Duan, S., Kaushal, S. S., Groffman, P. M., Band, L. E., & Belt, K. T. (2012). Phosphorus export across an urban to rural gradient in the Chesapeake Bay watershed. Journal of Geophysical Research: Biogeosciences, 117(G1). Filstrup, C.T., Downing, J.A., 2017. Relationship of chlorophyll to phosphorus and nitrogen in nutrient-rich lakes. Inland Waters 7, 385–400. Fisher, M. M., Reddy, K. R., & James, R. T. (2001). Long-term changes in the sediment chemistry of a large shallow subtropical lake. Lake and Reservoir Management, 17(3), 217- 232. FitzJohn, R. (2019). odin: Ode Generation and Integration. R package version 0.2.4. https://github.com/mrc-ide/odin Goeman, J. J. (2010). L1 penalized estimation in the Cox proportional hazards model. Biometrical Journal, 52(1), 70–84. https://doi.org/10.1002/bimj.200900028 Jordan Lake Reservoir Model December 2019 44 Graetz, D. A., Keeney, D. R., & Aspiras, R. B. (1973). Eh status of lake sediment‐water systems in relation to nitrogen transformations 1. Limnology and oceanography, 18(6), 908-917. https://doi.org/10.4319/lo.1973.18.6.0908. Hadley Wickham (2007). Reshaping Data with the reshape Package. Journal of Statistical Software, 21(12), 1-20. URL http://www.jstatsoft.org/v21/i12/ . Haggard, B. E., Moore, P. A., & DeLaune, P. B. (2005). Phosphorus flux from bottom sediments in Lake Eucha, Oklahoma. Journal of Environmental Quality, 34(2), 724–728. https://doi.org/10.2134/jeq2005.0724 Harrison, J. A., Maranger, R. J., Alexander, R. B., et. al. (2009). The regional and global significance of nitrogen removal in lakes and reservoirs. Biogeochemistry, 93(1-2), 143-157. Hirsch, R. M. (2014). Large biases in regression‐based constituent flux estimates: Causes and diagnostic tools. JAWRA Journal of the American Water Resources Association, 50(6), 1401-1424. Hirsch, R.M., and De Cicco, L.A., 2015, User guide to Exploration and Graphics for RivEr trends (EGRET) and dataRetrieval: R packages for hydrologic data (version 2.0, February 2015): U.S. Geological Survey Techniques and Methods book 4, chap. A10, 93 p., https://doi.org/0.3133/tm4A10 . James, R. T., Bierman Jr, V. J., Erickson, M. J., & Hinz, S. C. (2005). The Lake Okeechobee water quality model (LOWQM) enhancements, calibration, validation and analysis. Lake and Reservoir Management, 21(3), 231-260. Jensen, J. P., Jeppesen, E., Kristensen, P., Christensen, P. B., & Søndergaard, M. (1992). Nitrogen loss and denitrification as studied in relation to reductions in nitrogen loading in a shallow, hypertrophic lake (Lake Søbygård, Denmark). Internationale Revue der gesamten Hydrobiologie und Hydrographie, 77(1), 29-42. Jensen, J. P., Pedersen, A. R., Jeppesen, E., & Søndergaard, M. (2006). An empirical model describing the seasonal dynamics of phosphorus in 16 shallow eutrophic lakes after external loading reduction. Limnology and Oceanography, 51(1 II), 791–800. https://doi.org/10.4319/lo.2006.51.1_part_2.0791 . Jeppesen, E., Søndergaard, M., Meerhoff, M., Lauridsen, T. L., & Jensen, J. P. (2007). Shallow lake restoration by nutrient loading reduction - Some recent findings and challenges ahead. Hydrobiologia, 584(1), 239–252. https://doi.org/10.1007/s10750-007-0596-7. Jiao, Y., Yang, C., He, W., Liu, W. X., & Xu, F. L. (2018). The spatial distribution of phosphorus and their correlations in surface sediments and pore water in Lake Chaohu, China. Environmental Science and Pollution Research, 25(26), 25906–25915. https://doi.org/10.1007/s11356-018-2606-x. Lane, E. W., & Koelzer, V. A. (1943). Density of sediments deposited in reservoirs (Report No. 9). St Paul US Engineer District Sub-Office, Hydraulic Laboratory, University of Iowa. Jordan Lake Reservoir Model December 2019 45 Larsen, D. P., Schults, D. W., & Malueg, K. W. (1981). Summer internal phosphorus supplies in Shagawa Lake, Minnesota. Limnology and Oceanography, 26(4), 740-753. Lin, J., & Li, J. (2011). Nutrient response modeling in falls of the Neuse Reservoir. Environmental Management, 47(3), 398–409. https://doi.org/10.1007/s00267-011-9617-4 . Lürling, M., & van Oosterhout, F. (2013). Case study on the efficacy of a lanthanum-enriched clay (Phoslock®) in controlling eutrophication in Lake Het Groene Eiland (The Netherlands). Hydrobiologia, 710(1), 253–263. https://doi.org/10.1007/s10750-012-1141-x . Matisoff, G., Kaltenberg, E. M., Steely, R. L. et. al. (2016). Internal loading of phosphorus in western Lake Erie. Journal of Great Lakes Research, 42(4), 775–788. https://doi.org/10.1016/j.jglr.2016.04.004 . Mersmann, O., Trautmann, H., Steuer, D., and Bornkamp, B. (2018). truncnorm: Truncated Normal Distribution. R package version 1.0-8. https://CRAN.R- project.org/package=truncnorm.. Moore, P. A., Reddy, K. R., & Fisher, M. M. (1998). Phosphorus flux between sediment and overlying water in Lake Okeechobee, Florida: Spatial and temporal variations. Journal of Environmental Quality, 27(6), 1428–1439. https://doi.org/10.2134/jeq1998.00472425002700060020x. NCCO (2019). North Caroline Climate Office. Retrieved August, 1, 2019, from: https://climate.ncsu.edu/cronos. NCDENR (2007). Classifications and Water Quality Standards Applicable to Surface Waters and Wetlands of N.C. 15A NCAC 02b .0211. from: http://h2o.enr.state.nc.us/admin/rules/codes_statutes.htm . Nürnberg, G. (1988). Prediction of Phosphorus Release Rates from Total and Reductant-Soluble Phosphorus in Anoxic Lake Sediments. Canadian Journal of Fisheries and Aquatic Sciences, 45(3), 453–462. https://doi.org/https://doi.org/10.1139/f88-054 . Nürnberg, G. K. (2009). Assessing internal phosphorus load - Problems to be solved. Lake and Reservoir Management, 25(4), 419–432. https://doi.org/10.1080/00357520903458848 . Ott, W.R. (1990). A physical explanation of the lognormality of pollutant concentrations. United States. doi:10.1080/10473289.1990.10466789. Paerl, H., Otten, T. (2013). Harmful cyanobacterial blooms: causes, consequences, and controls. Microb. Ecol., 65 (4) (2013), pp. 995-1010 Paerl, H. (2014). Mitigating Harmful Cyanobacterial Blooms in a Human- and Climatically- Impacted World. Life, 4(4), 988–1012. https://doi.org/10.3390/life4040988 . Peder Jensen, J., SØndergaard, M., Landkildehus, F., Jeppesen, E., & Lauridsen, T. (2000). Trophic structure, species richness and biodiversity in Danish lakes: changes along a Jordan Lake Reservoir Model December 2019 46 phosphorus gradient. Freshwater Biology, 45(2), 201–218. https://doi.org/10.1046/j.1365- 2427.2000.00675.x . Poor, N. D. (2010). Effect of lake management efforts on the trophic state of a subtropical shallow lake in Lakeland, Florida, USA. Water, Air, and Soil Pollution, 207(1–4), 333–347. https://doi.org/10.1007/s11270-009-0140-7 . Prairie, Y.T., Duarte, C.M., and Kalff, J. 1989. Unifying nutrient–chlorophyll relationships in lakes. Can. J. Fish. Aquat. Sci. 46: 1176–1182 R Core Team (2018). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. URL: https://www.R-project.org/. River, M., & Richardson, C. J. (2018). Particle size distribution predicts particulate phosphorus removal. Ambio, 47(1), 124-133. Romo, S., Villena, M. J., Sahuquillo, M., Soria, J. M., Giménez, M., Alfonso, T., … Miracle, M. R. (2005). Response of a shallow Mediterranean lake to nutrient diversion: Does it follow similar patterns as in northern shallow lakes? Freshwater Biology, 50(10), 1706–1717. https://doi.org/10.1111/j.1365-2427.2005.01432.x . Sastrasinh, M., Weinberg, J. M., & Humes, H. D. (1982). Gentamicin (G) inhibits mitochondrial (M) calcium transport. Kidney International, 21(1), 224. Scheidegger, Andreas (2018). adaptMCMC: Implementation of a Generic Adaptive Monte Carlo Markov Chain Sampler. R package version 1.3. https://CRAN.R- project.org/package=adaptMCMC . Soetaert, K., Petzoldt, T., Setzer, R. Woodrow. (2010). Solving Differential Equations in R: Package deSolve. Journal of Statistical Software, 33(9), 1--25. https://doi.org/10.18637/jss.v033.i09 . Testa, J. M., Brady, D. C., Di Toro, D. M., Boynton, W. R., Cornwell, J. C., & Kemp, W. M. (2013). Sediment flux modeling: Simulating nitrogen, phosphorus, and silica cycles. Estuarine, Coastal and Shelf Science, 131(2013), 245–263. https://doi.org/10.1016/j.ecss.2013.06.014 . Tetratech. (2003). B . Everett Jordan Lake Nutrient Response Model Enhancement. Tomaszek, J. A., & Czerwieniec, E. (2000). In situ chamber denitrification measurements in reservoir sediments: an example from southeast Poland. Ecological engineering, 16(1), 61- 71. USACE (2019). B. Everett Jordan - Cape Fear River Basin. Retrieved from https://epec.saw.usace.army.mil/jord.htm USGS (2019). National Water Information System. Retrieved August 29, 2019, from https://maps.waterdata.usgs.gov/mapper . Jordan Lake Reservoir Model December 2019 47 USGS National Geospatial Program. (2019). The National Map. Retrieved from https://www.usgs.gov/core-science-systems/national-geospatial-program/national-map.Water Quality Portal (2019). National Water Quality Monitoring Council; Working together for clean water. Retrieved August 29, 2019, from https://www.waterqualitydata.us/ . Welch, E. B., & Jacoby, J. M. (2001). On determining the principal source of phosphorus causing summer algal blooms in Western Washington Lakes. Lake and Reservoir Management, 17(1), 55–65. https://doi.org/10.1080/07438140109353973 . Withers, P. J. A., Neal, C., Jarvie, H. P., & Doody, D. G. (2014). Agriculture and eutrophication: Where do we go from here? Sustainability (Switzerland), 6(9), 5853–5875. https://doi.org/10.3390/su6095853 . Wu, Z., Liu, Y., Liang, Z., Wu, S., & Guo, H. (2017). Internal cycling, not external loading, decides the nutrient limitation in eutrophic lake: A dynamic model with temporal Bayesian hierarchical inference. Water Research, 116, 231–240. https://doi.org/10.1016/j.watres.2017.03.039 . Zhang, M., Dolatshah, A., Zhu, W., & Yu, G. (2018). Case study on water quality improvement in Xihu lake through diversion and water distribution. Water (Switzerland), 10(3). https://doi.org/10.3390/w10030333 .