Loading...
HomeMy WebLinkAboutWatershed-Model-Jordan-NCSU_MillerJordan Lake Watershed Model Report Prepared for North Carolina Policy Collaboratory Prepared by Jonathan Miller, Kimia Karimi, Sankar Arumugam, & Daniel R. Obenour Environmental Modeling to Support Management and Forecasting Group Department of Civil, Construction, & Environmental Engineering and Center for Geospatial Analytics North Carolina State University December 2019 Jordan Lake Watershed Model December 2019 2 Table of Contents List of Tables .................................................................................................................................. 4 List of Figures ................................................................................................................................. 5 Executive Summary ........................................................................................................................ 6 1. Introduction ................................................................................................................................. 8 2. Methods....................................................................................................................................... 8 2.1 Study Area ............................................................................................................................. 8 2.2 Loading Monitoring Sites (LMSs) ........................................................................................ 9 2.2.1 Identification of LMSs .................................................................................................... 9 2.2.2 Incremental LMS watersheds ....................................................................................... 12 2.3 Model Segmentation - Subwatersheds ................................................................................ 12 2.4 Anthropogenic Factors ........................................................................................................ 12 2.4.1 Land covers................................................................................................................... 12 2.4.2 Point source dischargers ............................................................................................... 12 2.4.3 Livestock ...................................................................................................................... 14 2.5 Natural Factors .................................................................................................................... 14 2.5.1 Meteorological data ...................................................................................................... 14 2.5.2 Soil types ...................................................................................................................... 15 2.6 Nutrient Load Calculations ................................................................................................. 15 2.6.1 TN and TP loading estimates ........................................................................................ 15 2.6.2 Study period boundaries for yearly WRTDS estimates ................................................ 16 2.6.3 Uncertainty in WRTDS loading estimates ................................................................... 16 2.6.4 Incremental nutrient loading and error variance ........................................................... 16 2.7 Model Construction ............................................................................................................. 17 2.8 Nutrient Retention ............................................................................................................... 20 Jordan Lake Watershed Model December 2019 3 2.9 Model Assessment............................................................................................................... 20 3. Results ....................................................................................................................................... 21 3.1 In-stream Nutrient Loading ................................................................................................. 21 3.2 Model Posterior Parameter Estimates ................................................................................. 22 3.3 Spatial Variations in Nutrient Export and Retention .......................................................... 24 3.4 Nutrient Source Allocations over Time .............................................................................. 24 3.5 Nutrient Load Reduction Scenarios .................................................................................... 27 3.6 Watershed-Level Random Effects ...................................................................................... 29 3.7 Model Skill Assessment ...................................................................................................... 29 4. Discussion ................................................................................................................................. 32 4.1 Comparison to previous JL Watershed Model .................................................................... 32 4.2 Recommendations for Potential Nutrient Reductions ......................................................... 33 References ..................................................................................................................................... 35 Supplementary Material- Figures ................................................................................................. 38 Supplementary Material- Tables ................................................................................................... 41 Jordan Lake Watershed Model December 2019 4 List of Tables Table 2.2.1 Nutrient loading stations in Jordan and Falls lake watersheds .................................. 11 Table 2.7 Prior distributions of model parameters........................................................................ 19 Table 3.2 Mean parameter estimates for TN and TP models ....................................................... 22 Table 3.4 Percent of nutrient source contributions to Haw and New Hope Creek watersheds .... 27 Table 4.1 Summary of export coefficients for previous JL watershed model ............................. 32 Table S1 Yearly TN and TP load for major and minor WWTP ................................................... 41 Table S2 LMS sites with breaks (walls) in WRTDS loading estimates ...................................... 43 Jordan Lake Watershed Model December 2019 5 List of Figures Figure 2.1 Jordan and Falls Lake nutrient load monitoring sites (LMSs) and watersheds ........... 9 Figure 2.2.1 LMSs and their incremental watersheds................................................................... 10 Figure 2.4.1 Land cover, discharger, and livestock trends from 1994-2017 ............................... 13 Figure 2.5.1 Box and whisker plots of yearly precipitation .......................................................... 15 Figure 3.1 WRTDS in-stream and flow-normalized nutrient estimates ....................................... 21 Figure 3.2 TN and TP model posterior values. ............................................................................. 23 Figure 3.3a Model estimates for TN and TP export by subwatershed .......................................... 25 Figure 3.3b Expected stream and reservoir retention by subwatershed ........................................ 26 Figure 3.4 Nutrient source allocation by basin from 1994-2017. ................................................. 28 Figure 3.6 Watershed-level random effects for TN and TP. ......................................................... 30 Figure 3.7 Predicted vs. Observed plots for TN and TP ............................................................... 31 Figure S1 TN and TP yearly loadings from major WWTPs close to JL and FL .......................... 38 Figure S2 Major geologic regions in JL and FL region ................................................................ 39 Figure S3 Power relationship quantifying the uncertainty of TN and TP yearly loadings .......... 40 Jordan Lake Watershed Model December 2019 6 Executive Summary Jordan Lake is a major water supply, flood control, and recreational reservoir located in Chatham County, North Carolina. Watershed nutrient loading promotes excessive algal growth in the reservoir, resulting in chlorophyll concentrations that regularly exceed the state criterion of 40 ug/L. The Jordan Lake watershed can be divided into two major sections: the Haw River watershed and the New Hope Creek watershed. The Haw watershed, which includes the Cities of Greensboro and Burlington, discharges into the reservoir near its downstream (southern) end. The New Hope Creek watershed includes major portions of the Cities of Durham and Chapel Hill, and is an important contributor of flow and load to the reservoir’s upstream (northern) end. To efficiently manage watershed nutrient loading, it is critical to identify the major sources and locations of nutrient loading within these watersheds. In this study, we develop and apply a “hybrid” watershed modeling approach to characterize loading rates from point (i.e., wastewater treatment plant discharges) and nonpoint (i.e., diffuse washoff from the land surface) sources of nitrogen and phosphorus. The hybrid modeling approach combines the benefits of parsimonious mechanistic modeling with a rigorous statistical framework for data-driven inference and uncertainty quantification (Strickling & Obenour, 2018). The approach is comparable to the well-established USGS SPARROW model, but it is enhanced to account for interannual variability and to systematically incorporate and update prior information from previous studies through Bayesian inference. By modeling interannual variability, the model provides an assessment of how land use change and hydroclimatological variations have affected nutrient loading over time. The spatial scope of this study includes the Upper Falls Lake watershed in addition to the Jordan Lake (Haw and New Hope) watershed. The Upper Falls Lake watershed is intensively monitored, and thus provided substantial additional data for discerning differences in loading rates from various source types within the hybrid model. The primary source categories considered were point source discharges, undeveloped land, urban land, agricultural land, and livestock. Urban land was further divided into pre versus post-1980 development, to assess whether the age of the development is related to nutrient loading rates. The model was calibrated to data collected from 1982 to 2017. Model results show that urban land contributes the greatest total nitrogen (TN) load (i.e., export) on a per unit area basis. In particular, pre-1980 development contributes 9.5 kilograms per hectare per year (kg/ha/yr) of TN, while post-1980 development contributes 3.9 kg/ha/yr. Agriculture also contributes a substantial 4.0 kg/ha/yr of TN, while undeveloped land contributes a relatively low 0.7 kg/ha/yr. Nutrient removal within the watershed stream network is generally low (13%), except where large reservoir impoundments allow for greater removal rates of up to 75%, due to their relatively long residence times. Currently, the New Hope Creek watershed contributes 19% of the total nitrogen load to Jordan Lake, while the Haw River watershed makes up the rest. Point sources currently make up 48% of the load to the lake, while nonpoint sources make up the rest. Excess livestock export (beyond regular agricultural land export), which is considered part of the nonpoint source category, contributes around 2% of the TN load to the lake. Since 1994, there has been a slight reduction in TN loading to the New Hope Arm of Jordan Lake. Though point sources in the New Hope Creek watershed have decreased, TN loading from Jordan Lake Watershed Model December 2019 7 post-1980 urbanization has been increasing. TN loading to the Haw River watershed varies greatly among years due to the responsiveness of agricultural land export to changes in annual precipitation. In addition, TN loading due to point sources and urbanization are both increasing in the Haw River watershed. Some watersheds were found to produce more or less nutrients than would be expected based on the estimated nutrient export rates. Of note, two highly urbanized watersheds, Third Fork Creek, and Sandy Creek in Durham, exported less TN (and TP) than expected, while Cane Creek, Morgan Creek, and Ellerbe Creek all exported more TN (and TP) than mean model parameters would suggest. Modeling results suggest potential opportunities for reducing nitrogen loads to Jordan Lake. In general, nitrogen loads from undeveloped areas were found to be relatively low, indicating that most of the nitrogen loading in the watershed is due to human activities. For example, if point source discharges of nitrogen could be reduced by 25%, this would reduce total loading to the lake by 12%. Moreover, if pre-1980 development export rates could be reduced to the level of post-1980 development, this would result in an additional 13% load reduction to the lake. Also, if agricultural and livestock sources could be reduced by 25%, this would result in a 5% load reduction. Taken together, these example improvements would result in an overall lake nitrogen load reduction of 30% to Jordan Lake. The total phosphorus (TP) model showed similar trends to the TN model. Pre-1980 development contributes a substantial 1.5 kg/ha/yr of TP, while post-1980 development contributes 0.6 kg/ha/yr. Agriculture also contributes 0.6 kg/ha/yr of TP, while undeveloped land contributes a relatively low 0.05 kg/ha/yr. Mean nutrient removal throughout the stream network for TP (17%) was higher than for TN, mainly attributable to increased removal in reservoirs. For average precipitation, the New Hope Creek watershed contributes 20% of the TP load to JL, while the Haw River watershed contributes 80%. Point sources contribute a significantly lower percentage of the total TP loadings to JL (24%) as compared to TN (48%), with the rest of TP loading coming from nonpoint sources. Excess livestock export contributes around 6% of the TP load to the lake. Point source reductions of 25% would lower TP loading to JL by just 6%. If pre-1980 development loads could be reduced to the level of post-1980 development, loadings could be reduced 21%. A reduction in 25% of agricultural export would reduce TP loadings in JL 9%. If all three example improvements were implemented together, the overall lake phosphorus load reduction would be 35% to Jordan Lake. Finally, we note that hydroclimatological variability plays an important role in the interannual variability in nitrogen and phosphorus loading. Wet years (i.e., upper 33% of years) currently produce 60% and 82% more load than dry years (lower 33%) for nitrogen and phosphorus, respectively. However, since wet years also result in 65% more flow than dry years, TN concentrations are 6% lower while TP concentrations are 10% higher. Jordan Lake Watershed Model December 2019 8 1. Introduction Jordan Lake (JL; i.e., B. Everett Jordan Reservoir) and Falls Lake (FL) have both been identified as impaired due to high levels of chlorophyll-a caused by excessive nutrient loading. To address this, watershed-level management strategies were developed for both reservoirs by the North Carolina Division of Water Resources (DWR) in the early 2000s (Tetra Tech 2007; NC DWR 2009). In order to understand the sources of nutrient loading to the lakes, various process based models have been developed. The Watershed Analysis Risk Management Framework (WARMF) was used to model FL in 2009 (NC DWR 2009), whereas a Generalized Watershed Loading Function based model (GWLF) and Loading Simulation Program in C++ (LSPC) were developed for JL in 2002 and 2014 (Tetra Tech 2002; 2014). These models determined baseline nutrient loads entering JL and FL, separated loadings by source type and watersheds, and became the basis for watershed planning regulations and recommendations. The primary objective of this study is to improve our understanding of nitrogen and phosphorus export rates from various land uses, livestock, and point sources, as well as instream nutrient retention rates, so that we can determine how various source types contribute loading to JL. Compared to previous nutrient modeling efforts in the region, this study benefits from an extended calibration dataset and focuses on data-driven rate (i.e., parameter) estimation. Our modeling effort is an extension of the Strickling and Obenour (2018) hybrid watershed model, adapted and downscaled for the JL and FL watersheds. The hybrid modeling approach combines the benefits of parsimonious mechanistic modeling with a rigorous statistical framework for data-driven inference and uncertainty quantification, the latter being important for risk-based watershed management (Reckhow 1994; NRC 2001). The approach is comparable to the well- established USGS SPARROW model (Smith et al. 1997, Hoos and McMahon 2009; Garcia et al. 2011; Gurley et al. 2019), but it is enhanced to account for interannual variability and to systematically incorporate and update prior information from previous studies through Bayesian inference (Strickling & Obenour, 2018). By modeling interannual variability over multiple decades, the model benefits from a richer calibration dataset, and provides an assessment of how land use change and hydroclimatological variations have affected nutrient loading over time. While our primary focus is JL, including the FL watershed expands our calibration dataset, enhancing the ability of the model to characterize different nutrient loading sources. 2. Methods 2.1 Study Area JL and FL are located in the NC Piedmont (Figure 2.1). Both reservoirs were completed in the early 1980s and immediately declared nutrient sensitive by the state. The upper New Hope (NH) Arm of JL was determined to be impaired and placed on the 2002 303d list for high levels of chlorophyll-a, and the lower NH Creek and Haw River Arms (HR) were added in 2006 (Tetra Tech 2007). In addition, the Upper HR Arm of JL had elevated pH levels. JL watershed planning has been ongoing since the early 2000s and is still in the process of being formalized. Initial nutrient load reductions were set as 35% TN and 5% TP for the Upper NH Creek Arm and 8% TN and 5% TP for the HR Arm (Tetra Tech 2007). Intensive sampling in FL began in the early 2000s which led to nutrient watershed planning mandated by the state (Session Law 2009-486). Phase I goals of the Falls Lake Rules included a reduction of 40% TN and 77% TP from major sources in the watershed (http://portal.ncdenr.org/web/fallslake/). JL and FL belong to different Jordan Lake Watershed Model December 2019 9 watersheds (Cape Fear and Neuse, respectively), but both share similar underlying soil conditions and relative levels of urbanization and development such that nutrient export dynamics were expected to be similar. 2.2 Loading Monitoring Sites (LMSs) 2.2.1 Identification of LMSs Nutrient load monitoring sites (LMSs) were identified based on current and historical sampling locations that had sufficient flow and nutrient sampling data to calculate yearly and seasonal TN and/or TP loads. To be included as an LMS, a site needed a minimum of five years of daily flow records and at least 50 water quality samples during that period of record. These minimum conditions are consistent with previous studies using WRTDS load estimates and model defaults (Hirsch and De Cicco 2015; Chanat et al. 2016). All flow data were obtained from the United States Geological Survey (USGS), whereas nutrient data were obtained from the Water Quality Portal (WQP; Read et al. 2017) as well as local city managers (e.g., city of Durham). The two largest sources of nutrient data (downloaded from the WQP) came from the USGS and the North Carolina Department of Environmental Quality (NCDEQ). Sites from these different entities were often located in close proximity. Data from water quality sites with less than 5% deviations in watershed area and no intervening point sources were compiled together (Table 2.2.1). Figure 2.1: Jordan and Falls Lake watersheds. Jordan Lake is further split between the Haw River and New Hope creek watersheds. Nutrient load monitoring sites (LMSs) used to collect nutrient data for calibrating the model, and major and minor wastewater treatment plants (WWTPs) are also shown. Shading represents the entire Jordan and Falls Lake watersheds. Jordan Lake Watershed Model December 2019 10 In many cases, ample water quality data were available at the location of the USGS flow monitoring station. However, if little or no water quality data were located at the flow station, nearby water quality stations were used instead, assuming there was less than a 20% change in watershed area between the flow and water quality monitoring stations. If multiple water quality sites were located close to the flow station, only the site with the longest record was chosen. In one exception, two water quality sites shared the same flow monitoring station (NH1 and NH6; Table 2.2.1), which was done to include two substantial data records collected above and below a major wastewater discharger on Morgan Creek. In this case, the LMSs were represented at the location of the water quality monitoring sites, and flows are adjusted based on the drainage area ratio between the two sites discounting any intervening wastewater flow. There were 26 LMSs located upstream of JL and FL (Figure 2.2.1). Stations were split into three major basins for classification purposes: the Haw River (HR) watershed of Jordan Lake, the New Hope (NH) Creek watershed of Jordan Lake, and Falls Lake (FL). Five LMSs were located on tributaries that drain directly into JL (HR1, NH 1-4; Table 2.2.1): the Haw River at Bynum, Morgan Creek, New Hope Creek, Northeast Creek and White Oak Creek; and they represent 85% of the HR watershed and 49% of the NH Creek watersheds. Seven LMSs (HR 2-8) were located upstream of the HR at Bynum LMS (HR1), five of which were in or near Greensboro, NC. Ten LMSs were upstream of Falls Lake (FL1-FL10) representing 62% of the FL watershed. Three sites were located directly downstream of large reservoirs (HR4, FL6, and FL 9). Figure 2.2.1: Load monitoring stations shown with their incremental watersheds. 79 subwatersheds used to aggregate nutrient export to LMSs are also shown. Jordan Lake was separated between the Haw River and New Hope Arms of Jordan Lake. Jordan Lake Watershed Model December 2019 11 Table 2.2.1: Loading stations located in the Jordan and Falls watersheds along with their drainage areas. Years of record corresponds to time that loadings could be estimated (i.e., when daily flow and monthly water quality sampling was performed). The number of water quality samples available is also shown. TN TP LMS Name Res Drainage area (km2) Years of record # of samples Years of record # of samples NH1 Morgan Creek, Jordan Lake JL 121.4 1994-2017 578 1994-2017 588 NH2 New Hope Creek JL 203.9 1994-2017 575 1994-2017 487 NH3 Northeast Creek JL 53.6 1996-2017 430 1996-2017 434 NH4 White Oak Creek JL 31.1 2000-2017 106 2004-2017 104 NH5 Morgan Creek, White Cross JL 21.4 2000-2017 116 1999-2017 128 NH6 Morgan Creek, Chapel Hill JL 103.2 2001-2013 141 2000-2013 159 NH7 Sandy Creek, Cornwallis JL 12.1 2009-2017 133 2009-2017 142 NH8 Third Fork Creek JL 41.2 2009-2017 107 2009-2017 106 HR1 Haw River, Bynum JL 3296.4 1994-2017 590 1994-2017 588 HR2 Cane Creek JL 19.6 1989-2017 227 1989-2017 235 HR3 Haw River, Burlington JL 1562.1 1994-2017 268 1994-2017 268 HR4 Reedy Fork , Gibsonville JL 316.6 1981-1986 2001-2017 341 1982-1986 2001-2017 328 HR5 N. Buffalo Creek JL 96.2 1999-2017 394 2001-2017 300 HR6 S. Buffalo Creek JL 88.6 2000-2017 343 2000-2017 363 HR7 Reedy Fork, Oak Ridge JL 53.4 2001-2017 255 2001-2017 213 FL1 Ellerbe Creek, Gorman FL 54.8 2006-2017 280 2006-2017 283 FL2 Ellerbe Creek, Murray FL 11.2 2009-2013 100 2009-2013 113 FL3 Eno River, Durham FL 367.2 1994-2000 2004-2017 375 1994-2000 2004-2017 319 FL4 Eno River, Hillsborough FL 171.0 1990-2017 223 1990-2017 205 FL5 Little River, Orange Factory FL 202.7 1988-2000 2005-2017 381 2006-2017 236 FL6 Little River, Fairntosh FL 246.4 1996-2011 196 1996-2011 193 FL7 Mountain Creek FL 20.8 1995-2011 156 1995-2011 140 FL8 Flat River, Bahama FL 385.9 1981-2011 472 1981-2011 471 FL9 Flat River, Dam FL 434.4 1983-1990 2003-2017 225 2003-2017 216 FL10 Knap of Reeds Creek FL 111.4 2006-2017 142 2006-2017 156 Jordan Lake Watershed Model December 2019 12 2.2.2 Incremental LMS watersheds LMS watersheds were delineated using Spatial Analyst tools in ArcMap 10.6.1 (ESRI 2018). Watershed drainage areas ranged from 11 km2 to over 3000 km2 (Table 2.2.1) with a median value of 106 km2. Often, LMS watersheds had one or more LMS contained within their upstream watershed (Figure 2.2.1). In order to use multiple LMSs that drained to each other, our model predicted nutrient loadings of incremental LMS watersheds. We determined incremental LMS watersheds by subtracting out any upstream LMS watersheds that were contained in a larger LMS watershed. If a LMS did not have an upstream LMS in its watershed, its incremental watershed was equal to its total watershed. For periods when loading data were unavailable for a given LMS, the incremental watershed for that LMS was combined with the next downstream incremental watershed. 2.3 Model Segmentation - Subwatersheds To more accurately account for nutrient transport and retention, incremental LMS watersheds were divided into subwatersheds (Figure 2.2.1). All predictor data (e.g., land covers, precipitation, livestock) were compiled at the subwatershed level. The largest possible subwatershed corresponded to a USGS 12-digit hydraulic unit code (HUC; https://water.usgs.gov/GIS/huc.html). If a LMS was located in the middle of a HUC, the HUC was split into two parts: an upstream portion that drained to the LMS, and a downstream portion that drained to a subsequent downstream LMS. Eleven LMSs had watersheds that were smaller or equal to one HUC 12 (Figure 2.2.1). Seventy-nine subwatersheds were located within the study area, with a mean drainage area of 63.0 km2, minimum of 11.2 km2, and maximum of 146.4 km2. 2.4 Anthropogenic Factors 2.4.1 Land covers Land cover variables were derived from the U.S. Conterminous Wall-to-wall Anthropogenic Land use Trends (NWALT) dataset (Falcone 2015). We aggregated NWALT land cover designations into three major categories: Urban (low, medium, and high density residential areas, transportation, industrial, and commercial development), Agriculture (pasture and crop), and Undeveloped (semi-developed and low use and wetlands). Semi-developed land was included with undeveloped because it is mostly comprised of forested land (Miller et al. 2019). We further split urbanization into two additional categories: urban development before 1980 and development after 1980. In order to determine when urbanization occurred in the region, we interpolated available NWALT data (1974, 1982, 1992, 2002, and 2012) to obtain year-specific land cover values for each subwatershed. Since our study extended beyond 2012, we also used linear extrapolation for years 2013-2017 based on 2002 and 2012 values. Land cover trends throughout the study period were generally gradual, such that linear extrapolation was considered reasonable for post 2012 years. Land cover trends for the three main basins are shown in Figure 2.4.1. 2.4.2 Point source dischargers Point source dischargers included major (> 1 million gallons per day) and minor WWTPs in the study area (Figure 2.1). Discharge data were collected starting in 1994 from the NC DEQ (personal communication 2019) and included monthly TN, TP, and flow values. However, many Jordan Lake Watershed Model December 2019 13 WWTPs had numerous missing months, so we determined annual loads by multiplying the yearly median concentration and flow values together for each WWTP. Yearly TN and TP loadings from major WWTPs located directly upstream of JL and FL are shown in Figure S1, and a full list of WWTPs used in this study, their location, and their nearest loading station are in Table S1. LMSs with major WWTPs in their watersheds were only modeled starting in 1994 due to a lack of discharge data before that year (i.e., HR1,3,5, NH1-3, and FL1,3,10), while LMSs without major WWTPs were potentially modeled to 1982 depending on data availability. Only one LMS (HR4) had both monitoring data and a minor WWTP that discharged before 1994. Since the minor WWTP only represented < 3% of its mean yearly load, we assumed pre-1994 discharge from the WWTP was equal to its post 1994 discharge. TN and TP trends of point dischargers aggregated by basin are shown in Figure 2.4.1. Figure 2.4.1: Land cover, discharger, and livestock trends from 1994-2017 in the Haw River, New Hope Creek, and Falls Lake watersheds. Jordan Lake Watershed Model December 2019 14 2.4.3 Livestock The numbers of cows, chickens, and hogs in subwatersheds were calculated by using county- level United States Department of Agriculture (USDA) census and survey reports (https://www.nass.usda.gov/). Cow and hog census and survey results covered our entire study period (1982-2017), while chickens had census data collected every five years beginning in 1997. For missing years between census dates, chicken counts were interpolated, whereas chicken counts before 1997 were assumed to be equal to 1997 values. Only two incremental watersheds (HR1, 3) had large chicken counts (> 1,000,000 and > 150,000, respectively) and these watersheds were not modeled before 1994 (as they were also missing major WWTP discharge data). Because USDA livestock counts were only reported by county, there was uncertainty as to where actual farms were located in our study area. In order to accurately represent the spatial locations of livestock throughout the region, county-level data were assigned to incremental watersheds based on an area ratio. Major urban areas were excluded when calculating these proportions, as livestock were assumed to be located outside of city areas. Livestock counts were then further divided into the subwatersheds (using area ratios again) that comprised the LMS (see Table 2.2.1). However, chickens in Chatham County were accounted for differently because a large majority of Chatham’s chicken farms (>90%) are located outside of the JL watershed, and the county has an extremely high chicken count (>3,000,000; USDA). Chatham Co. records of current chicken farm locations were available (opendata-chathamncgis.opendata.arcgis.com) and used to find the proportion of chicken farms in Chatham Co. that were in the JL watershed (i.e., 8.2%). This ratio was then used, instead of an area ratio, to more accurately represent the chicken count from Chatham County. Basin level trends of livestock are shown in Figure 2.4.1. 2.5 Natural Factors Natural variations can greatly influence the export of nutrients from watersheds. Fluctuations in precipitation across multiple scales have been shown to be a controlling mechanism on downstream loads (Howarth et al., 2012; Sinha and Michalak, 2016; Strickling and Obenour, 2018). In addition, soil and geologic properties have been linked to both nutrient sources as well as nutrient transport in watersheds (Preston et al. 2011). Therefore, we investigated variation in both meteorological and geologic conditions across the region for consideration in model construction. 2.5.1 Meteorological data Precipitation was obtained from the PRISM Climate Group (http://www.prism.oregonstate.edu/). We downloaded monthly precipitation rasters from 1982 to the present and used spatial analysis in the R package “raster” (Hijmans et al. 2015) to determine mean precipitation across the study area. There was substantial variation across years as well as across subwatersheds (Figure 2.5.1). Jordan Lake Watershed Model December 2019 15 Figure 2.5.1: Box and whisker plots of yearly precipitation. The dashed line represents mean yearly precipitation for the JL and FL watersheds. 2.5.2 Soil types Different soil characteristics were analyzed for possible inclusion in the model formulation. The largest variation in soil conditions in the JL and FL watershed are associated with Triassic basin soils, which have a distinct geologic history. Triassic soils generally have lower infiltration rates, resulting in higher erosion potential and lower baseflows in streams, such that export and retention rates have been parameterized differently for Triassic soils in previous mechanistic JL models (Tetra Tech 2002, 2014). We did not explicitly model Triassic soils in our model, but watershed-level random effects (see section 2.7) were analyzed to determine if there were any trends (like Triassic regions) in nutrient export. Six LMSs were located in predominantly Triassic soils (Figure S2; NH3,4,7,8 and FL1,2). 2.6 Nutrient Load Calculations 2.6.1 TN and TP loading estimates Our model required yearly TN and TP loadings, which is the product of concentration and discharge. Most riverine monitoring programs measure streamflow daily, whereas water quality variables are sampled less frequently (e.g., monthly). Regression-based methods can be used to estimate missing daily concentration values, thus enabling the estimation of loadings across time. In this study, daily TN and TP concentrations and loads were estimated using the USGS Weighted Regressions on Time, Discharge, and Season (WRTDS; Hirsch et al., 2010). This method estimates daily concentration using time, discharge, and the time of the year (seasonality) as explanatory variables. WRTDS develops a unique regression model for each day in the estimation period. It uses a semi-parametric regression where observations that are collected under similar conditions to the estimation date (in terms of time, discharge, and season) are more heavily weighted. Jordan Lake Watershed Model December 2019 16 Several major WWTPs had large, sudden increases or decreases in their yearly loadings presumably due to plant upgrades. In order to not bias WRTDS loading estimates near these events, the period of record was split into two (pre-change, post-change) with a minimum of 50 observations in each (Table S2). 2.6.2 Study period boundaries for yearly WRTDS estimates LMSs needed to have a minimum of 5 consecutive years of daily flow data with at least 50 observations within that time period to run yearly WRTDS estimates (Table 2.2.1). In addition, water quality samples in the beginning and ending year of record had to span more than half a year (i.e., > 6 samples). Water quality sampling did not always start or end concurrent with a calendar year though, and this meant that years at the beginning and end of water quality records were often incomplete (< 6 samples). Though, yearly loading estimates for that year were not included in full model construction, sampling observations from those partial years (i.e., years with < 6 samples) were included to parameterize WRTDS estimates for other years. In addition, several LMSs had significant gaps in their water quality or flow monitoring data during their period of record (Table 2.2.1). For years with missing daily flow data, no nutrient loading estimates could be made. Gaps in water quality samples of one year were considered acceptable, as preliminary analysis (using simulations) showed that loading estimates during a one year gap varied < 1% from estimates using the entire time series. If a longer gap occurred for a given LMS, the above rule (> 6 samples) for identifying the beginning and ending of the loading estimation period was applied. 2.6.3 Uncertainty in WRTDS loading estimates Uncertainty in loading estimates were determined through subsampling of three NC stations that had nearly daily TN and TP observations for at least 7 consecutive years. The method is adapted from Strickling & Obenour (2018). Initially, WRTDS was run with the full record of nutrient data for the three sites, and those estimates were considered to be the “true” nutrient load. Then, subsamples of the full dataset were taken randomly at frequencies of 6, 12, 25, 50, 90, 120 and 150 samples-per-year. For each sampling frequency, 300 different random subsamples were created. Then, for each subsample, yearly TN and TP loads were estimated using WRTDS. These yearly loads (based on the subsamples), were then compared to the “true” loads, and their residuals were evaluated to determine a coefficients of variation (CV) for each site. In order to determine uncertainty, the standard deviation of the residuals (𝜎𝜎𝑡𝑡) of the 300 subsampled loading estimates were evaluated and divided by the mean “true” load (𝜇𝜇𝑡𝑡) for each year (t): 𝐶𝐶𝐶𝐶𝑡𝑡=𝜎𝜎𝑡𝑡𝜇𝜇𝑡𝑡 Eq. 1 Finally, a power relationship was used to model 𝐶𝐶𝐶𝐶 based on annual sampling frequency, considering subsampling results for all frequencies and subsampled sites (Figure S3). This allowed us to determine the CV (and consequently the standard deviation) of WRTDS estimates based on the number of water quality samples available for a given year. 2.6.4 Incremental nutrient loading and error variance The response variable in our model was the change in nutrient load across an incremental watershed. This is defined as the difference between the load at an incremental watershed’s downstream LMS and the cumulative load from any upstream LMSs. For sites with no upstream Jordan Lake Watershed Model December 2019 17 LMSs, the incremental load is equal to its total load. The uncertainties of incremental loads were calculated based on the relationship between correlated random variables (Eq. 2; Kottegoda & Rosso, 2008): 𝜎𝜎�𝑥𝑥,𝑡𝑡2 = 𝜎𝜎𝑥𝑥,𝑡𝑡2 −2 ∑𝜌𝜌𝑥𝑥,𝑖𝑖𝜎𝜎𝑥𝑥,𝑡𝑡𝜎𝜎𝑖𝑖,𝑡𝑡𝑛𝑛𝑖𝑖=1 +∑∑𝜌𝜌𝑖𝑖,𝑗𝑗𝜎𝜎𝑖𝑖,𝑡𝑡𝜎𝜎𝑗𝑗,𝑡𝑡𝑛𝑛𝑗𝑗=1𝑛𝑛𝑖𝑖=1 Eq. 2 Where 𝜎𝜎�𝑥𝑥,𝑡𝑡2 is the incremental load variance for a given incremental LMS watershed (x) in year (t), where (i) and (j) represent upstream LMSs. Here, 𝜎𝜎𝑥𝑥,𝑡𝑡2 is the WRTDS error variance for the downstream LMS, 𝜎𝜎𝑖𝑖,𝑡𝑡 and 𝜎𝜎𝑗𝑗,𝑡𝑡 are the WRTDS standard deviations for upstream LMSs, and 𝜌𝜌𝑥𝑥,𝑗𝑗 and 𝜌𝜌𝑖𝑖,𝑗𝑗 are correlation coefficients between the LMS loadings. This relationship was extended to up to n=3 upstream LMS sites (Figure 2.2.1). 2.7 Model Construction Our model was formulated similar to Strickling & Obenour (2018). Within a Bayesian framework, we related deterministically predicted incremental nutrient loads (𝑦𝑦�i,t; Eq. 3) to an inferred incremental load (yi,t) for every incremental watershed (i) and year (t). The watershed-level random effect (αi; Gelman et al. 2014) accounted for spatial variability across incremental watersheds not explained by the deterministic prediction (𝑦𝑦�i,t), and the residual error (with standard deviation of σε) primarily accounts for temporal variability unexplained by the deterministic prediction. The hyperdistribution of the normally distributed watershed-level random effect was centered on zero, with variance σLMS2. L(yi,t) ~ N( L(𝑦𝑦�i,t + αi), σε) Eq. 3 αi ~ N(0, 𝜎𝜎𝐿𝐿𝐿𝐿𝐿𝐿) L(y) is the natural log transformation of y + 105 (kg/year of TN; y + 104 for TP). This transformation addresses non-normality of the loading residuals, while the offset accounts for any negative incremental loads that would produce non-real values when log transformed. Negative incremental loads were common for incremental watersheds that included large reservoirs that retained a substantial portion of the upstream nutrient load. The inferred incremental load (yi,t; Eq. 4) is related to the WRTDS incremental estimates (𝑦𝑦�i,t) by taking into account the uncertainty of those loading estimates (𝜎𝜎�i,t; as determined by the number of samples per year, see section 2.6.3). By accounting for the uncertainty of WRTDS predictions, loading estimates at sites with more observed data were given more weight in the Bayesian calibration than sites with less observed data. 𝑦𝑦�i,t ~ N(yi,t, 𝜎𝜎�i,t) Eq. 4 Within the model, the deterministic prediction of incremental nutrient load (𝑦𝑦�i,t) was calculated by aggregating watershed source contributions (i.e., land use, dischargers, livestock) and subtracting in-stream losses from upstream nutrient loads (Eq.5). 𝑦𝑦�i,t = Li,t,ur1 + Li,t,ur2 + Li,t,ag + Li,t,und + Li,t,ps + Li,t,ch + Li,t,h + Li,t,cw – Ui,t * ri,z + ԑi,t Eq. 5 Contributions were calculated for pre- and post-1980 urban (Li,t,ur1; Li,t,ur2), agricultural (Li,t,ag), and undeveloped (Li,t,und) lands, point sources (i.e., WWTPs; Li,t,ps), chickens (Li,t,ch), hogs (Li,t,h), and cows (Li,t,cw). Loads from upstream incremental watersheds (Ui,t) were reduced by their Jordan Lake Watershed Model December 2019 18 expected in-stream and reservoir losses (ri,z; see section 2.8; Eq.7). Each source specific load was calculated as follows: Li,t,x = β,x ( 𝑝𝑝�𝑖𝑖,𝑡𝑡 𝛾𝛾,𝑥𝑥) * aTi,t,x * (1 - ri,x) Eq. 6a (land covers) Li,t,x = β,x (𝑝𝑝�𝑖𝑖,𝑡𝑡 𝛾𝛾,𝑥𝑥) * hTi,t,x * (1 - ri,x) Eq. 6b (livestock) Li,t,x = β,ps * wTi,t * (1 - ri,x) Eq. 6c (dischargers) where Li,t,x is calculated in kg/year and represents the total contributed load from a given source type. Parameter β,x represents a land cover’s export coefficient (kg/ha/yr) or a livestock or WWTP delivery coefficient (unitless, 0-1). Parameter γ,x is the precipitation impact coefficient (PIC, unitless) for a given nonpoint source which is parameterized as a power relationship with the export coefficient. 𝑝𝑝�𝑖𝑖,𝑡𝑡 is a scaled precipitation value (actual precipitation for incremental watershed divided by mean precipitation) specific to each incremental watershed and year. Often, there were multiple sources within incremental watersheds that had to be aggregated to obtain a total source contribution. For land covers and livestock, this occurred when an incremental watershed was made up of more than one subwatershed. For point sources, this meant there was more than one WWTP in the incremental watershed. To account for multiple contributions to a source, aTi,t,x, hTi,t,x, and wTi,t,x are transposed vectors of individual source locations (i.e., ha of land covers, livestock counts, and kg/yr nutrient output from WWTPs, respectively) that were multiplied by a vector (ri,x) of location-specific stream and reservoir retention losses. Export coefficients (β,x) accounted for the long-term median expected flux of nutrients off land surfaces in the region, while discharge coefficients represented the fraction of nutrients from point source or livestock that reached stream networks. PICs (γ,x) accounted for intra-annual variability of nutrient loads based on annual precipitation. For mean annual precipitation (i.e., 𝑝𝑝�𝑖𝑖,𝑡𝑡=1), PIC values do not affect loading estimates. PICs differ by source, but are related to each other through a common hyperdistribution with mean µ𝛾𝛾, and standard deviation 𝜎𝜎𝛾𝛾 (Table 2.7). Point sources do not have a PIC term (Eq.6c) as yearly WWTP discharge values were reported and accounted for potential yearly variation. All model parameters were assigned a prior distribution within the model (Table 2.7). Informative priors were used when previous studies reporting similar parameters were available. Prior distributions for land export rates were taken from Dodd (1992), while stream retention rates were adapted from previous SPARROW models (Hoos and McHahon 2009; Garcia et al. 2011). Prior distributions for chicken and hog TN delivery coefficients were adapted from Strickling and Obenour (2018) to represent kilograms of TN per animal per year. Priors for TP delivery coefficients were adapted from the TN value based on common TN:TP ratios (2:1) in chicken and hog manure (NCSU 2019). Uninformative priors (i.e., uniform priors) were used for parameters without any relevant prior information. Models were parameterized within a Bayesian framework using RStan software in R (R. Core Team, 2018; Stan Development Team, 2018). RStan uses a Hamiltonian Monte Carlo sampling method for determining posterior distributions and is considered to be faster than other Bayesian samplers (Gelman et al. 2015). 20,000 iterations were run in three parallel chains with a burn in period of 5,000 iterations (that were discarded), creating 9,000 posterior samples after thinning Jordan Lake Watershed Model December 2019 19 (accepting every fifth posterior value). Parameters were considered to have converged if their scale reduction coefficient (𝑅𝑅�) was approximately equal to one (Gelman and Rubin 1992). Table 2.7: Prior distributions (or hierarchical distributions for PIC) of watershed parameters. Note, SD is standard deviation, PIC is precipitation impact coefficient, EC is export coefficient and DC is delivery coefficient. EC units are kilograms per hectare per year (kg/ha/yr) and DC units are kilograms per animal per year (kg/an/yr). Para-meter Name Units Prior distribution TN model Prior distribution TP model Source β,ag Agriculture EC kg/ha/yr N(9,7) N(1,0.7) Dodd et al. 1992 β,ur1 Pre-1980 Urban EC kg/ha/yr N(8,3) N(1,0.9) Dodd et al. 1992 β,ur2 Post-1980 Urban EC kg/ha/yr N(8,3) N(1,0.9) Dodd et al. 1992 β,und Undeveloped EC kg/ha/yr N(2,2) N(0.2,0.1) Dodd et al. 1992 β,ch Chicken DC kg/an/yr N(0.001,0.0003) N(0.0005,0.0002) Strickling and Obenour 2018 β,h Hog DC kg/an/yr N(0.04,0.02) N(0.02,0.01) Strickling and Obenour 2018 β,cw Cow DC kg/an/yr U(0,5) U(0,5) - β,ps Point source DC - N(1,0.10) N(1,0.10) - γ,ag Agriculture PIC - N(µγ, σγ) N(µγ, σγ) - γ,ur1 Pre-1980 Developed PIC - N(µγ, σγ) N(µγ, σγ) - γ,ur2 Post-1980 Developed PIC - N(µγ, σγ) N(µγ, σγ) - γ,und Undeveloped PIC - N(µγ, σγ) N(µγ, σγ) - γ,ch Chicken PIC - N(µγ, σγ) N(µγ, σγ) - γ,h Hog PIC - N(µγ, σγ) N(µγ, σγ) - γ,cw Cow PIC - N(µγ, σγ) N(µγ, σγ) - γ,ret Retention rate PIC - N(0,1) N(0,1) - ds Stream decay rate d-1 N(0.14,0.05) N(0.20,0.08) Hoos and McHahon 2009, Garcia et al. 2011 ρr Reservoir loss rate m/yr N(11,2) N(30,8) Hoos and McHahon 2009, Garcia et al. 2011 σε Model residual SD kg/yr U(0,1x106) U(0,1x106) - σLMS LMS random effect SD kg/yr U(0,1x106) U(0,1x106) - µγ PIC mean - U(1,1) U(1,1) - σγ PIC SD - U(0,1) U(0,1) - Jordan Lake Watershed Model December 2019 20 2.8 Nutrient Retention Nitrogen retention was incorporated into our model to account for the effects of settling and denitrification or burial that occurs in streams and impoundments. Nutrient retention in streams was modeled as a first-order decay (−𝑑𝑑𝑠𝑠, days-1) related to the mean stream residence time in days (𝜏𝜏𝑖𝑖,𝑧𝑧) for each incremental watershed (i) and path (z) from a given source (subwatershed or discharger) to their downstream LMS. Nutrient retention in reservoirs was modeled as a function of their hydraulic loading rate (ratio of flow to surface area) and mass transfer coefficient (−𝑝𝑝𝑟𝑟; Kelly 1987). Often, several waterbodies (w) were located along the flow path from a source such that multiple hydraulic loading rates (𝑘𝑘𝑖𝑖,𝑧𝑧,𝑤𝑤 m/yr) were included. Stream characteristics (flow and length) were obtained from NHD+ flowlines (Moore and Dewald 2016) and used to determine stream residence times from all possible sources to their downstream LMSs. Retention rates for point discharges (wTi,t,x) and upstream LMS loadings (Ui,t; Eq. 5) were calculated from their locations. Nonpoint sources (aTi,t,x, hTi,t,x) distributed within a given subwatersheds were assumed to have the same travel time, calculated using one-half the distance of the longest flow path within that subwatershed plus the distance from the outlet of the subwatershed to the downstream LMS. An overall retention rate that combined both stream and reservoir retention (ri,z; Eq. 7), was then calculated for each nonpoint source, point source discharger, and upstream LMS loading: ri,z = 1 – exp�−𝑑𝑑𝑠𝑠∗𝜏𝜏𝑖𝑖,𝑧𝑧�∗∏exp �−𝑝𝑝𝑟𝑟𝑘𝑘𝑖𝑖,𝑧𝑧,𝑤𝑤�𝑤𝑤 Eq. 7 where stream retention times (𝜏𝜏𝑖𝑖,𝑧𝑧) and reservoir hydraulic loading rates (𝑘𝑘𝑖𝑖,𝑧𝑧,𝑤𝑤) were adjusted every year (i) using a PIC (γ,ret) that related to the normalized yearly precipitation (pi,t; yearly precipitation minus mean precipitation divided by its standard deviation) specific to each incremental watershed (i) in a given year (t) (Eq. 8a,b). This was done to account for variation in nutrient retention that might occur due to higher (or lower) than normal flows. 𝜏𝜏𝑖𝑖,𝑧𝑧= 𝜏𝜏𝑖𝑖,𝑧𝑧(1+𝛾𝛾,𝑟𝑟𝑟𝑟𝑡𝑡∗ 𝑝𝑝𝑖𝑖,𝑡𝑡) Eq. 8a 𝑘𝑘𝑖𝑖,𝑧𝑧,𝑤𝑤 = 𝑘𝑘𝑖𝑖,𝑧𝑧,𝑤𝑤 * (1+ γ,ret* pi,t) Eq. 8b 2.9 Model Assessment Model performance was summarized using the coefficient of determination (R2) calculated based on annual incremental nutrient loads. Predicted incremental loading estimates were derived using the mean posterior values from the Bayesian calibration and compared to WRTDS loading estimates (𝑦𝑦�i,t). R2 was determined for model predictions with and without the watershed-level random effect for the HR, NH, and FL watersheds. To test the ability of the model to perform out-of-sample predictions, we performed a 3-fold cross-validation (Elsner and Schmertamann, 1994). The data was split into three groups by major watershed (HR, NH, and FL), and the model was trained on 2 of the 3 watersheds in turn. Predictions were then made on the excluded group, in turn, such that we generated predictions for all observations using separate data. Jordan Lake Watershed Model December 2019 21 3. Results 3.1 In-stream Nutrient Loading WRTDS in-stream nutrient estimates and flow normalized estimates can be observed to determine long term trends in nutrient loadings in JL and FL (Figure 3.1). While there is much intra-annual variation for in-stream estimates, flow normalized loadings filter out year to year variability making them a good visualization of long-term trends (Hirsch et al. 2010). Initially, nutrient concentrations in all basins decreased substantially from 1980-2000. After 2000, though, there are inconsistent trends both within and across basins. In the HR watershed, TN loadings have steadily increased since 2000 while TP loadings remained constant. In the NH Creek watershed, substantial increases in loads were seen in the early 2000s, but then large TN reductions occurred in watersheds associated with major WWTP improvements (NH1-3; Figure 3.1). TP loadings in the NH Creek watershed have decreased only slightly. Finally, in FL, large reductions were observed in FL1 and FL10, both of which have major WWTP dischargers, while TN and TP loadings in other FL watersheds have remained constant or trended upwards. Figure 3.1: Weighted-regression on Time, Discharge, and Season (WRTDS) annual nutrient loading estimates (points) and flow-normalized estimates (lines) for (A) TN and (B) TP. Load monitoring stations (LMS) shown are the farthest downstream sites on main tributaries to Jordan and Falls Lake. LMSs shown account for 85%, 49%, and 62% of the drainage area to the HR Arm, NH Creek Arm, and FL. Jordan Lake Watershed Model December 2019 22 3.2 Model Posterior Parameter Estimates Model results show that urban land contributes the greatest nitrogen load on a per unit area basis. In particular, pre-1980 development contributes 9.5 kilograms per hectare per year (kg/ha/yr) of TN, while post-1980 development only contributes less than half of that amount, 3.9 kg/ha/yr. Agriculture also contributes a substantial 4.0 kg/ha/yr of TN, while undeveloped lands export a relatively low 0.7 kg/ha/yr (Table 3.2; Figure 3.2a). Land cover ECs (β,x) represent expected nutrient export for mean yearly precipitation (i.e., scaled precipitation (𝑝𝑝�𝑖𝑖,𝑡𝑡) = 1; Eq. 6a, b). Model posterior distributions show the uncertainty of model parameters (i.e. 95% credible intervals), and appeared to be more influenced by the data, than their widely distributed prior distributions. For TP, pre-1980 urban land exported 1.5 kg/ha/yr of TP, post-1980 urban land 0.6 kg/ha/yr, agriculture 0.6 kg/ha/yr, and undeveloped land 0.05 kg/ha/yr (Table 3.2; Figure 3.2b). Table 3.2: Mean parameter estimates for the TN and TP models along with 95% credible intervals (CI). Export coefficients and retention rates Precipitation Impact Coefficients TN TP TN TP Parameter Mean 95% CI Mean 95% CI Parameter Mean 95% CI Mean 95% CI β,ag 4.0 2.3-5.7 0.6 0.4-0.8 γ,ag 4.1 2.9-5.0 4.0 2.9-5.1 β,ur1 9.5 7.4-11.4 1.5 1.1-1.8 γ,ur1 1.2 0.7-1.7 1.8 1.1-2.5 β,ur2 3.9 0.7-7.3 0.6 0.03-1.4 γ,ur2 2.1 0.4-4.0 2.0 0.2-3.9 β,und 0.7 0.1-1.5 0.05 0-0.13 γ,und 2.8 0.6-5.2 2.4 0.5-4.5 β,ch 0.01 0-0.02 0.004 0-0.009 γ,ch 1.9 0.3-3.8 2.4 0.5-4.8 β,h 0.04 0.01-0.07 0.02 0-0.04 γ,h 1.8 0.3-3.7 2.0 0.3-4.1 β,cw 0.5 0.1-1.0 0.16 0-0.55 γ,cw 1.8 0.3-3.7 2.3 0.4-4.4 β,ps 0.83 0.75-0.91 0.87 0.70-1.03 γ,ret 0.07 0.01-0.17 0.09 0 -0.22 ds 0.04 0.01-0.07 0.03 0-0.07 µγ 1.6 1.1-2.0 1.9 1.1-2.7 ρr 11.2 8.7-13.6 25.9 17.7-34.8 σγ 1.2 0.8-1.6 1.0 0.5-1.6 σε 0.07 0.07-0.08 0.16 0.14-0.17 σLMS 1.3 0.9-1.9 1.8 0.8-3.9 Precipitation impact coefficients (PIC) identify the variability of export rates for different levels of precipitation. Mean nutrient export is adjusted by the PICs to determine export during low and high flow years. Agriculture had the largest PIC for both TN and TP implying that export from agricultural lands (crop and pasture) vary the most due to rainfall. During a high flow year (90th percentile, 𝑝𝑝�𝑖𝑖,𝑡𝑡=1.18), mean nutrient export for agriculture would almost double from 4.0 to 7.9 kg/ha/yr for TN and 0.6 to 1.2 kg/ha/yr for TP (see Table 4.1 in Discussion). For a low flow year (10th percentile, 𝑝𝑝�𝑖𝑖,𝑡𝑡=0.81), mean nutrient export for agriculture would close to half the normal export: 1.7 from 4.0 kg/ha/yr for TN and 0.3 to 0.6 kg/ha/yr for TP. Variation in nutrient export due to precipitation is discussed in more detail in section 4.1. Jordan Lake Watershed Model December 2019 23 Figure 3.2: TN (A) and TP (B) model posterior values (solid lines) shown with their prior distribution (dotted line). Jordan Lake Watershed Model December 2019 24 3.3 Spatial Variations in Nutrient Export and Retention The amount of TN and TP exported from nonpoint sources (i.e. land covers and livestock) was calculated for each subwatershed (Figure 3.3a) using 2017 land covers, mean precipitation, and mean posterior ECs (βec; Table 3.2). Since the largest nonpoint source of nutrient export came from pre-1980 urban lands, subwatersheds located in the urban cores of Burlington, Durham, and Greensboro (Figure 3.3a) had the largest expected export. On average, predominantly undeveloped watersheds exported between 1-3 kg/ha/yr of TN, and 0.2-0.4 kg/ha/yr of TP while urban cores generally exported over 6 kg/ha/yr of TN and 1 kg/ha/yr of TP. Stream retention rates for TN and TP were relatively similar (0.04; 0.03 day-1, respectively; Table 3.2), while more retention occurred in reservoirs for TP than TN (25.9 m/d vs. 11.3 m/d; Table 3.2). These rates imply that on average, 13% of TN and 17% of TP is retained within the JL stream network (Figure 3.3b). Maximum TN and TP removal rates were 75% and 95%, respectively, for several subwatersheds north of Greensboro located hydraulically behind several large reservoirs. Residence times and hydraulic loading rates were affected by the PIC for stream retention (γ,ret; 0.07, 0.09; Table 3.2). These values implied that for one standard deviation increase in yearly precipitation, expected residence times would decrease around 7-8% while the hydraulic loading rates would increase 8-10% (Eq. 8a,b). 3.4 Nutrient Source Allocations over Time Yearly TN and TP loadings from 1994-2017 were calculated based on mean model parameters, land use, livestock counts, and precipitation (Figure 3.4). Temporal loadings varied greatly in the HR watershed as compared to the NH Creek watershed due to high levels of agriculture which have substantially increased export during high flow years (Table 3.2). TN and TP loadings nearly triple to the HR Arm, from the lowest to highest precipitation years, while loadings to the NH Creek Arm never double (Figure 3.4). Jordan Lake Watershed Model December 2019 25 Figure 3.3a: Exports of (A) TN and (B)TP export from land uses and livestock by subwatershed; point sources are shown separately as dots. Jordan Lake Watershed Model December 2019 26 Figure 3.3b: Expected nutrient retention for (A) TN and (B) TP in streams and reservoirs prior to reaching Jordan and Falls Lake. Jordan Lake Watershed Model December 2019 27 Nutrient loadings were converted to percent contributions to JL by basin (HR and NH) and summarized based on normal (33-67 percentile flow years), low (lower 33%) and high (upper 67%) flow years (Table 3.4). During normal flow years, point sources accounted for 44% of TN loadings to JL (33% from HR; 11% from NH) and 26% of TP loadings (Table 3.4). Pre-1980 urbanization accounted for 24% of TN loadings and 33% of TP loadings to JL, while post-1980 urbanization accounted for only 4% of TN and TP loadings. Agriculture nutrient sources are most prevalent in the HR watershed (18% for TN, 26% for TP), while livestock contributions are relatively small (2% for TN and 6% for TP). Table 3.4: Percent of nutrient sources that contributes to total Jordan Lake loadings from the Haw River (HR) and New Hope (NH) watersheds for normal flow years (33-67 percentile flow years). In parenthesis are the percent of each nutrient source during low flow years (lower 33%) and high flow years (upper 67%), respectively. % TN % TP Nutrient source HR NH HR NH Agriculture 17 (11,25) 1 (1,2) 24 (17,31) 2 (1,2) Urban, pre-1980 18 (18,17) 6 (6,5) 25 (23,22) 8 (7,7) Urban, post-1980 2 (1,1) 2 (1,1) 2 (2,2) 2 (2,1) Undeveloped 6 (5,7) 2 (2,2) 4 (3,4) 1 (1,1) Livestock 2 (2,2) 0 (0,0) 5 (5,6) 1 (1,1) Discharger 33 (40,28) 11 (15,10) 21 (32,20) 5 (6,3) 3.5 Nutrient Load Reduction Scenarios Simple forecasts were made to understand the ability to reduce nutrients through mitigation of different sources specifically for JL based on mean yearly loading rates from 1994-2017. For example, if point source discharges of nitrogen could be reduced by 25%, this would result in a 12% load reduction to the lake. Moreover, if pre-1980 development loads could be reduced to the level of post-1980 development, this would result in a 13% load reduction to the lake. Also, if agricultural and livestock sources could be reduced by 25%, this would result in a 5% load reduction. Taken together, these example improvements would result in an overall lake nitrogen load reduction of 30% to Jordan Lake. These potential reduction estimates are dependent on yearly precipitation as that determines the actual percentage of each nutrient source in JL (Table 3.4). For TP, point source reductions of 25%, would only lower TP reaching JL 6%. If pre-1980 development loads could be reduced to the level of post-1980 development, loadings could be reduced 21%. A reduction in 25% of agricultural impact would reduce TP loadings in JL 9%. If all three example improvements were done together, the overall lake phosphorus load reduction would be 35% to Jordan Lake. In the example nutrient loading reduction scenarios provided above, we assume equal percent reductions across both wet and dry conditions. However, some best management practices may be more effective at reducing nutrient loads under dry conditions than under wet conditions. Jordan Lake Watershed Model December 2019 28 Figure 3.4: Nutrient source allocation by basin from 1994-2017 representing (A) TN and (B) TP that reached JL. Jordan Lake Watershed Model December 2019 29 3.6 Watershed-Level Random Effects Watershed-level random effects account for spatial variation in nutrient loading among incremental watersheds that is not explained by deterministic model parameters (i.e., export rates, retention, and PICs). Random effects show watersheds that have either increased or reduced nutrient loadings as compared to expected model predictions. Three highly urban Durham watersheds, Sandy (NH7), Third Fork (NH8), and the upper portion of Ellerbe Creek (FL2), two watersheds downstream of major reservoirs, Little and Flat River (FL6,9), and the upper portion of the Eno River watershed (FL4) exported somewhat less TN than expected (Figure 3.6a). The random effect for Third Fork, Sandy, and Ellerbe Creek, which are mainly comprised of pre-1980 urban development, implies that these watersheds export less TN (-2.5, - 1.8, and -1.2 kg/ha/yr, respectively) than shown in Figure 3.3a. The Little and Flat River include reservoirs which may be particularly efficient at trapping nutrients. On the other hand, Cane (HR2), Morgan (NH1), North Buffalo (HR5), and Ellerbe Creek (FL1) all exported significantly more TN than mean model parameters would suggest. For Morgan, Ellerbe, and North Buffalo Creeks, this elevated TN watershed-level effect is possibly associated with major WWTP dischargers which are located in close proximity to their LMSs. For Cane Creek, this increased export is possibly related to nonpoint sources of nutrients. Random effects for the TP model had similar trends to TN (Figure 3.6b). Watershed-level random effects were also analyzed to determine if larger regional trends occurred across the study area. There were no clearly observable trends, either between basins or within larger stream networks with multiple LMSs. In addition, there were no clear trends in TN or TP export from Triassic soils which were prevalent in five LMS watersheds: Northeast, White Oak, Sandy, Third Fork and Ellerbe Creeks (NH3, NH4, NH7, NH8, FL1; Figure S2). Two Triassic watersheds had substantially positive random effects (Ellerbe and Northeast) and two had substantially negative effects (Sandy and Third Fork Creek), while White Oak had a relatively minor positive random effect (Figure 3.6). However, the effects of Triassic soils could potentially be obscured my more dominant drivers of loading heterogeneity (e.g., point sources, urban development). 3.7 Model Skill Assessment The full TN model, including random effects, explained 93% of the variation in TN loading estimates at LMSs. Discounting the random effects, the model still explains 92% of TN loading variability. The model (with watershed-level random effects) explained 95% of the variation in the HR, 92% in NH, and 81% in FL (Figure 3.7a). The TP model was slightly less predictive than the TP model explaining 86% of the variation in the data (84% without watershed-level random effect). It was able to explain 92% of the variation in the HR, 84% in NH, and 62% in FL (Figure 3.7b). In cross validation, the predictive ability of the model remained high. The R2 of the full TN model (without watershed-level random effects) lowered slightly from 92% to 89%. For TP, R2 lowered from 84% to 80% in cross validation. Jordan Lake Watershed Model December 2019 30 Figure 3.6: Watershed-level random effects for each incremental watershed for (A) TN and (B) TP. Jordan Lake Watershed Model December 2019 31 Figure 3.7: Predicted vs. Observed plots for (A) TN and (B) TP including watershed-level random effects Jordan Lake Watershed Model December 2019 32 4. Discussion 4.1 Comparison to previous JL Watershed Model Model results show similar trends in land cover nutrient export to the Tetra Tech (2014) JL watershed model, even though exact comparisons are not possible because land covers were aggregated differently. In the Tetra Tech model, high density residential and commercial land exported the most nutrients per hectare and showed the least variation in export due to precipitation (Table 4.1). In our model, pre-1980 urban land, similarly, had the highest nutrient export (9.5 kg/ha/yr of TN; 1.5 kg/ha/yr of TP) while having the lowest PICs (Table 3.2). There is a high correlation (> 0.8) between high-density and pre-1980 urbanization in the JL watershed, and both models suggest that the most intense nutrient export in the JL watershed comes from the urban cores of large cities (where high-density and pre-1980 urban lands occur most frequently). Further research is needed to determine if the excess nutrients are due to legacy effects (e.g., older wastewater infrastructure, scoured or buried stream networks, lack of best management practices (BMPs)), or merely due to increased imperviousness associated with high-density areas. As new developments are becoming more clustered and high-density in nature, understanding the precise source of nutrients from high-density development and the ability of BMPs to reduce them will be important for future management decisions. Table 4.1: Summary of export coefficients for previous JL watershed model (Tetra Tech, 2014) and this study. Ranges for parameters represent export rates due to variations in precipitation, not the uncertainty of model parameters. TN TP Model Nutrient source (kg/ha/yr) (kg/ha/yr) Tetra Tech (2014) High-density residential/commercial 5.7-9.2 0.9-1.6 Low/Medium-density residential 2.3-6.5 0.3-0.9 Row crops 2.4-11.4 0.2-1.4 Pasture/grassland 2.0-5.7 0.1-0.3 Forest 1.1-3.4 0.05-0.2 Current model (2019) Pre-1980 urban 7.4-11.6 1.0-2.0 Post-1980 urban 2.5-5.5 0.4-0.8 Agriculture 1.7-7.9 0.3-1.2 Undeveloped 0.4-1.1 0.03-0.1 Export from other land covers were also consistent between the two models. In both models, undeveloped land (primarily forest) was the lowest exporter of nutrients. Also, in both models, export from agricultural areas were most strongly affected by variations in precipitation. In addition, Tetra’s Tech lower urban category (low and medium-density residential) exported nutrients at a similar rate to pasturelands as well as croplands during low to normal precipitation (Table 4.1). In our model, post-1980 urban lands (our lower urban category) exported similar rates to agriculture (crop plus pasture; Table 3.2). Though both models differ in how urban lands were split and aggregated, they both imply that different urban lands can have substantially different export rates. Jordan Lake Watershed Model December 2019 33 Two differences between the JL models were the magnitude of nutrient export rates from undeveloped (i.e. forest) areas and retention rates within the stream network. In our model, more than 50% less nutrient export was expected to occur from undeveloped areas than the Tetra Tech model across a wide range of hydro-climatological conditions (Table 4.1). Also, our model implies there is less stream retention from upstream subwatersheds (Figure 3.3b) as compared to the Tetra Tech model. For example, retention from Greensboro non-point sources are predicted to be 15-20% for TN in our model, while it was closer to 30% in the Tetra Tech model (2014). One explanation for this difference might be the point discharger delivery coefficient included in our model (0.83; Table 3), which possibly accounts for additional stream retention of WWTP nutrient discharges (i.e., inorganic nutrients), that are potentially more readily attenuated than nutrients from nonpoint sources. 4.2 Recommendations for Potential Nutrient Reductions Model results strongly support the theory that the majority of nutrient inputs to Jordan Lake are from anthropogenic sources. Reducing these sources will be instrumental to lowering future nutrient loads to the reservoir. Based on identified sources of TN and TP in the watershed, four management strategies would lead to large potential nutrient load reductions in the lake: 1) reduction of point source loadings (i.e., WWTPs), which are the largest source of TN and 2nd largest of TP, 2) retrofitting or replacing infrastructure in older urban environments (i.e. pre-1980 urban), which are the largest nonpoint source of nutrients per unit area, 3) mitigating TN and TP loading from agricultural lands, especially during wet conditions, and 4) limiting, reducing, or offsetting the removal of undeveloped land which is the lowest exporter of nutrients. WWTP point sources are responsible for nearly 50% of TN and 25% of TP loadings top JL, and these percentages rise substantially during low flow years (Figure 3.4). In the NH Creek watershed, TN loadings have been reduced by nearly 50% from 1994 to 2017, while TP reductions have been slight (Figure 2.4.1). However, in the HR watershed, the opposite trends have occurred for TN and TP. 2017 TN loadings are similar in magnitude to 1994 levels, while TP loadings have been reduced nearly 50% (Table 2.4.1). In particular, significant increases in both TN and TP loadings from WWTPs in the Upper HR Arm have been observed since 2010 (Figure 2.4.1). Due to low levels of stream retention observed in the JL watershed (Figure 3.3), the majority of these loadings (> 70% for TN and TP) are expected to reach the Jordan Lake reservoir. Pre-1980 urban lands were the largest nonpoint source of nutrients (normalized by area) in the watershed, implying that the urban cores of major cities are a large source of nutrients (Figure 3.3a). Lower precipitation impact coefficients for these pre-1980 urban lands indicate they are also a more constant source of nutrients than other land covers (Table 3.2). The increased export of nutrients from urban cores could be attributable to many causes (e.g., increased pollutant washoff due to impervious cover, lack of BMPs, leaky wastewater infrastructure) and more research is necessary to determine the exact source for each subwatershed and LMS. If pre-1980 development export rates could be reduced to the level of post-1980 development, large reductions could be attained for TN (12%) and TP (21%; see section 3.5). Agricultural-related nutrient export (agricultural lands plus livestock) is dominated by land export as opposed to excess livestock waste (>85% TN, > 80% TP; Table 3.4). Agricultural land export rates were generally lower than previous research might suggest (Dodd et al. 1992; Strickling and Obenour 2018), and imply agricultural lands export at rates similar to post-1980 Jordan Lake Watershed Model December 2019 34 urban lands (Table 3.2). This is due to the fact that the majority of agriculture in the JL watershed is dedicated to pasture lands (> 90%; Miller et al. 2019), as opposed to croplands, which often have higher export rates. Though agricultural TN export only account for 18% of JL loadings during normal flow years (26% for TP), that number increased to 27% during high flow years (33% for TP; Table 3.4). Therefore, strategies aimed at mitigating loading from agricultural areas will have a larger effect on lake loadings in high flow years. Undeveloped lands export the lowest amount of nutrients of all non-point sources (Table 3.2), which is consistent with findings from targeted water quality monitoring recently done in JL (Delesantro & Riveros-Iregui, personal communication, 2019). However, the amount of undeveloped land in JL is decreasing due to the high level of development occurring in the watershed. Even though more recent development (post-1980) has significantly lower TN and TP export when compared to pre-1980 development, it still exports more TN (5.5x) and TP (12x) than undeveloped land (Table 3.2). Stormwater nutrient accounting tools such as SNAP (Stormwater Nitrogen and Phosphorus; NC DEQ 2019) are currently being used in the JL watershed for newly constructed lands. Model results can be used to assess and update those tools to ensure current strategies are accomplishing their intended goals. Finally, post-1980 urban export rates are essentially equal to agricultural rates (Table 3.2). This implies that new urban construction in agricultural areas may have limited impact on nutrient loading. Jordan Lake Watershed Model December 2019 35 References Chanat, J. G., Moyer, D. L., Blomquist, J. D., Hyer, K. E., & Langland, M. J. 2016. Application of a weighted regression model for reporting nutrient and sediment concentrations, fluxes, and trends in concentration and flux for the Chesapeake Bay Nontidal Water-Quality Monitoring Network, results through water year 2012 (No. 2015-5133). US Geological Survey. Dodd, R. C., Cunningham, P. A., Curry, R. J., Liddle, S. K., McMahon, G., Stichter, S., & J. P. Tippett. 1992. Watershed planning in the Albemarle-Pamlico estuarine system (No. 3). Albemarle-Pamlico Estuarine Study, NC Dept. of Environment, Health, and Natural Resources. Elsner, J. B., and C. P. Schmertmann. 1994. Assessing forecast skill through cross validation. Weather and Forecasting, 9(4): 619-624. Falcone, J. A. 2015. US conterminous wall-to-wall anthropogenic land use trends (NWALT), 1974–2012 (No. 948). US Geological Survey. Garcia, A. M., A. B. Hoos, & S. Terziotti. 2011. A Regional Modeling Framework of Phosphorus Sources and Transport in Streams of the Southeastern United States Journal of the American Water Resources Association 47(5): 991-1010. Gelman, A., & D. B. Rubin. 1992. Inference from iterative simulation using multiple sequences. Statistical science, 7(4): 457-472. Gelman, A., Carlin, J. B., Stern, H. S., & Rubin, D. B. 2014. Bayesian data analysis (Vol. 2). Boca Raton, FL, USA: Chapman & Hall/CRC. Gelman, A., D. Lee, & J. Guo. 2015. Stan: A probabilistic programming language for Bayesian inference and optimization. Journal of Educational and Behavioral Statistics, 40(5): 530-543. Gurley, L.N., A. M. Garcia, S. Terziotti, and A. B. Hoos. 2019. SPARROW model datasets for total nitrogen and total phosphorus in North Carolina, including simulated stream loads: U.S. Geologcial Survey data release, https://doi.org/10.5066/P9UUT74V. Hijmans, R. J., van Etten, J., Cheng, J., Mattiuzzi, M., Sumner, M., Greenberg, J. A., ... & M. R. J. Hijmans. 2015. Package ‘raster’. R package. Hirsch, R. M., D.L. Moyer, & S.A. Archfield. 2010. Weighted regressions on time, discharge, and season (WRTDS), with an application to chesapeake bay river inputs. Journal of the American Water Resources Association 46(5): 857–880. 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. Hoos, A. B. & G. McMahon. 2009. Spatial analysis of instream nitrogen loads and factors controlling nitrogen delivery to streams in the southeastern United States using spatially referenced regression on watershed attributes (SPARROW) and regional classification frameworks. Hydrological Processes: An International Journal 23(16): 2275-2294. Jordan Lake Watershed Model December 2019 36 Howarth, R., Swaney, D., Billen, G., Garnier, J., Hong, B., Humborg, C., ... & Marino, R. 2012. Nitrogen fluxes from the landscape are controlled by net anthropogenic nitrogen inputs and by climate. Frontiers in Ecology and the Environment 10(1), 37-43. Kelly, C. A., Rudd, J. W. M., Hesslein, R. H., Schindler, D. W., Dillon, P. J., Driscoll, C. T., ... & R.E. Hecky. 1987. Prediction of biological acid neutralization in acid-sensitive lakes. Biogeochemistry, 3(1-3): 129-140. Kottegoda, N. T., & R. Rosso. 2008. Applied statistics for civil and environmental engineers. Malden, MA: Blackwell. Miller, J. W., Paul, M. J., & D. R. Obenour. 2019. Assessing potential anthropogenic drivers of ecological health in Piedmont streams through hierarchical modeling. Freshwater Science, 38(4): 771-789. Moore, R. B. and T. G. Dewald. 2016. The Road to NHDPlus—Advancements in Digital Stream Networks and Associated Catchments. Journal of the American Water Resources Association 52(4): 890-900. NC DEQ (North Carolina Department of Environmental Quality). 2019. Nutrient Practices and Crediting. https://deq.nc.gov/about/divisions/water-resources/planning/nonpoint-source- management/nutrient-offset-information. Accessed November, 2019. NC DWR (North Carolina Department of Water Resources). 2009. Falls Lake Watershed Analysis Rick Management Framework (WARMF) Development. Division of Water Quality. Planning Section. Modeling/TMDL Unit. NCSU (North Carolina State University). 2019. North Carolina Agricultural Chemicals Manual. UNC Press. Retrieved from: https://content.ces.ncsu.edu/north-carolina-agricultural-chemicals-manual. NRC (National Research Council). 2001. Assessing the TMDL approach to water quality management. National Academies Press. Read, E. K., Carr, L., De Cicco, L., Dugan, H. A., Hanson, P. C., Hart, J. A., ... & Winslow, L. A. 2017. Water quality data for national‐scale aquatic research: The Water Quality Portal. Water Resources Research, 53(2), 1735-1745. Reckhow, K. H. 1994. Water quality simulation modeling and uncertainty analysis for risk assessment and decision making. Ecological Modelling, 72(1-2): 1-20. Sinha, E., & Michalak, A. M. 2016. Precipitation dominates interannual variability of riverine nitrogen loading across the continental United States. Environmental science & technology 50(23): 12874-12884. Smith, R. A., G. E. Schwarz, & R. B. Alexander, 1997. Regional interpretation of water‐quality monitoring data. Water resources research, 33(12), 2781-2798. Strickling, H.L. & Obenour, D.R. 2018. Leveraging Spatial and Temporal Variability to Probabilistically Characterize Nutrient Sources and Export Rates in a Developing Watershed. Water Resources Research 54(7): 5143-5162. Jordan Lake Watershed Model December 2019 37 Preston, S. D., R.B. Alexander, G.E. Schwarz, & C.G. Crawford. 2011. Factors Affecting Stream Nutrient Loads: A Synthesis of Regional SPARROW Model Results for the Continental United States. Journal of the American Water Resources Association 47(5): 891-915. Tetra Tech. 2002. Jordan Lake Nutrient Response Model. The Jordan Lake Project Partners. Tetra Tech. 2007. B. Everett Jordan Reservoir, North Carolina Phase I Total Maximum Daily Load. Final Report. Department of Environment and Natural Resources. Division of Water Quality. Tetra Tech. 2014. Lake B. Everett Jordan Watershed Model Report. NC Nutrient Science Advisory Board, NC Division of Water Resources, and Triange J. Council of Governments. Jordan Lake Watershed Model December 2019 38 Supplementary Material- Figures Figure S1: TN and TP yearly loadings from major WWTPs located closest to Jordan and Falls Lake. Jordan Lake Watershed Model December 2019 39 Figure S2: Major geologic regions in JL and FL region. Jordan Lake Watershed Model December 2019 40 Figure S3: Power relationship quantifying the uncertainty (i.e., coefficient of variation (CV)) of (A) TN and (B) TP yearly (Weighted Regression on Time Discharge and Season) loadings based on the number of water quality samples available. Jordan Lake Watershed Model December 2019 41 Supplementary Material- Tables Table S1: Yearly TN and TP loads, discharge and residence time to downstream loadings stations for major and minor WWTPs . NC DEQ ID# Name TN load (kg) (2017) TP load (kg) (2017) Discharge (MGD) (2017) Downstream LS Mean RT (days) to LS Mean HL (m/yr) to LS NC0026824 John Umstead Hospital WWTP 3396 192 1.610 FL10 0.05 NC0023841 North Durham WRF 25,765 1,230 8.280 FL1 0.24 NC0026433 Hillsborough WWTP 2016 206 0.919 FL3 1.10 NC0037869 Arbor Hills Mobile Home Park 18 4 0.001 FL3 0.55 312.0 NC0056731 Grande Oak Subdivision WWTP 41 7 0.002 FL3 0.20 107.4 NC0023876 Southside WWTP 43,342 1,420 5.872 HR1 1.26 NC0021211 Graham WWTP 0 0 0.000 HR1 1.30 NC0021474 Mebane WWTP 6,177 1,444 1.269 HR1 1.76 96.0 NC0035866 Bynum WWTP 122 40 0.006 HR1 0.02 NC0042285 Trails WWTP 485 105 0.020 HR1 1.19 NC0022675 Birmingham Place 0 0 0.000 HR1 2.94 NC0022098 Cedar Valley WWTP 126 25 0.006 HR1 3.29 NC0038164 Nathanael Greene Elem. School WWTP 0 0 0.000 HR1 2.64 NC0045128 Sylvan Elementary School 172 15 0.001 HR1 1.78 NC0045152 Jordan Elementary School 130 18 0.002 HR1 1.08 NC0042528 B Everett Jordan 1927 LLC 123 10 0.008 HR1 1.04 NC0025241 Mason Farm WWTP 55,655 3,387 5.878 NH1 0.20 NC0056413 Carolina Meadows WWTP 2115 39 0.148 NH1 0.01 NC0047597 South Durham WRF 78,164 3,159 8.023 NH2 0.55 NC0042803 Birchwood Mobile Home Park 102 22 0.007 NH2 1.20 NC0074446 Hilltop Mobile Home Park 335 54 0.010 NH2 1.54 NC0026051 Triangle WWTP 35,557 1,731 4.502 NH3 0.14 NC0023868 Eastside WWTP 55,515 1,810 3.883 HR3 0.04 NC0024881 Reidsville WWTP 48,902 2,429 2.280 HR3 1.50 NC0066966 Quarterstone Farm WWTP 1,438 111 0.038 HR3 1.28 Jordan Lake Watershed Model December 2019 42 NC0046019 The Summit WWTP 0 0 0.000 HR3 1.76 NC0046809 Pentecostal Holiness Church 19 2 0.000 HR3 1.59 NC0060259 Willow Oak Mobile Home Park 164 41 0.008 HR3 1.14 NC0077968 Horners Mobile Home Park 127 3 0.003 HR3 0.52 NC0045161 Altamahaw/Ossipee Elementary School 283 34 0.003 HR3 0.50 NC0045144 Western Alamance High School 478 44 0.006 HR3 0.40 NC0031607 Western Alamance Middle School 0 0 0.000 HR3 0.45 NC0055271 Shields Mobile Home Park 0 0 0.000 HR3 0.41 NC0065412 Pleasant Ridge WWTP 40 18 0.008 HR3 1.14 NC0073571 Countryside Manor WWTP 300 17 0.008 HR3 2.78 NC0046043 Oak Ridge Military Academy 255 29 0.004 HR3 2.74 NC0022691 Autumn Forest Manuf. Home Community 160 41 0.018 HR4 0.44 NC0047384 T.Z. Osborne WWTP 582,529 31,072 27.680 HR3 1.55 NC0038172 McLeansville Middle School WWTP 0 0 0.000 HR3 1.55 NC0029726 Guilford Correctional Center WWTP 0 0 0.000 HR3 1.53 NC0024325 North Buffalo Creek WWTP 161,171 4,724 6.619 HR5 0.29 Jordan Lake Watershed Model December 2019 43 Table S2: LMS sites that were run for TN and TP loading in WRTDS with a break (i.e. wall) in their record due to a large sustained change (> 25% up or down) in an upstream point source load. Years of record # of samples Point source LMS Name Nutrient Pre Post Pre Post Permit WWTP NH1 Morgan Creek, JL TN 1994-2009 2010-2017 360 228 NC0025241 Mason Farm NH3 Northeast Creek TN 1995-2004 2005-2017 108 322 NC0026051 Triangle UH3 Haw River, Burlington TN 1994-2013 2014-2017 223 56 NC0047384 T.Z. Osborne UH6 N. Buffalo Creek TN 1999-2007 2008-2017 148 241 NC0024325 North Buffalo Creek FL3 Eno River, Durham TN 1994-2006 2007-2017 158 217 NC0026433 Hillsborough UH6 N. Buffalo Creek TP 2000-2007 2008-2017 59 241 NC0024325 North Buffalo Creek