Climate Change Causes River Network Contraction and Disconnection in the H.J. Andrews Experimental Forest, Oregon, USA

Headwater streams account for more than 89% of global river networks and provide numerous ecosystem services that benefit downstream ecosystems and human water uses. It has been established that changes in climate have shifted the timing and magnitude of observed precipitation, which, at specific gages, have been directly linked to long-term reductions in large river discharge. However, climate impacts on ungaged headwater streams, where ecosystem function is tightly coupled to flow permanence along the river corridor, remain unknown due to the lack of data sets and ability to model and predict flow permanence. We analyzed a network of 10 gages with 38–69 years of records across a 5th-order river basin in the U.S. Pacific Northwest, finding increasing frequency of lower low-flow conditions across the basin. Next, we simulated river network expansion and contraction for a 65-year period of record, revealing 24% and 9% declines in flowing and contiguous network length, respectively, during the driest months of the year. This study is the first to mechanistically simulate network expansion and contraction at the scale of a large river basin, informing if and how climate change is altering connectivity along river networks. While the heuristic model presented here yields basin-specific conclusions, this approach is generalizable and transferable to the study of other large river basins. Finally, we interpret our model results in the context of regulations based on flow permanence, demonstrating the complications of static regulatory definitions in the face of non-stationary climate.


INTRODUCTION
More than 89% of the global river network is headwaters (Downing et al., 2012;Allen et al., 2018), supporting ecosystem services and the health of downstream waters (Alexander et al., 2007;US EPA, 2015). These services are associated with the frequency with which streams have surface flow (hereafter "flow permanence"), and any declines in flow permanence will effectively disconnect larger rivers from their headwaters and their functions. Flow generated in headwater streams is highly sensitive to changes in precipitation timing, magnitude, and duration based on a small number of empirical studies over short timescales (Godsey and Kirchner, 2014;Jensen et al., 2017;Zimmer and McGlynn, 2017;Prancevic and Kirchner, 2019). However, no observational studies have covered a sufficient period of record to evaluate if and how changing climate has altered flow permanence across river networks. Consequently, numerical simulations parameterized with readily available data are needed to fill this knowledge gap (Gallart et al., 2016;Ward et al., 2018a).
Changes in flow permanence can alter the transport and transformation of water, energy, dissolved and suspended materials, and organisms throughout the river network (Larned et al., 2010;Gallart et al., 2012;Steward et al., 2012;Datry et al., 2016Datry et al., , 2017Raymond et al., 2016). Evaluating how flow permanence has changed requires quantification of both the temporal variation (i.e., the frequency a given segment has surface flow) and spatial variation (i.e., the spatial connectivity of surface flow) (Covino, 2017;Wohl, 2017;Wohl et al., 2017). In headwater streams, flow permanence is controlled by the dynamic interaction of geologic setting with hydrologic forcing (Costigan et al., 2016;Prancevic and Kirchner, 2019). Climate change is primarily associated with changes to hydrologic forcing, such as altering the spatial distribution and within-year timing of precipitation. Geologic setting-such as valley width and slope, sinuosity, and hydraulic conductivity-will remain relatively static compared to the pace of climate change.
Changes in flow permanence complicate management and protection of headwater streams. Regulatory protections in the U.S. and E.U. are traditionally focused on perennially flowing waters, with emerging attention paid to temporarily flowing waters (US DoD, 1986;Nikolaidis et al., 2013;US DoD EPA, 2015;Fritz et al., 2018;US DoD US EPA, 2018;Walsh and Ward, 2019). Further complicating management, data on headwater streams, and particularly intermittent and ephemeral streams, are lacking. For example, only 3% of the rivers gaged in the U.S. are headwater streams, as gages are heavily biased toward larger rivers (Poff et al., 2006;Eng and Milly, 2007). A proposed rule would revise protected status to waters with contiguous surface flow in a "typical" 30-year period in the U.S., but does not address systematic changes in flow permanence (US DoD US EPA, 2018). The time-variable definition of the 30-year window does not consider the role of climate change and shifting norms, despite clear evidence that non-stationarity is prevalent in hydrologic systems (Milly et al., 2008). For example, systematic declines in streamflow, and particularly lower low-flows, have been observed across the U.S. Pacific Northwest (Luce and Holden, 2009). Subsequent study revealed that slowing westerlies during winter months reduced the orographic enhancement of precipitation in the region, changing both the amount of precipitation reaching the landscape and the timing for storage vs. export from catchments (Luce et al., 2013). A more recent example paints another dire picture for the future of streamflow in the Southwestern U.S. in response to shifting precipitation and temperature (Milly and Dunne, 2020).
Here we assess whether flow permanence in headwater streams has shifted over the past 65 years from the mid twentieth century baseline in response to observed changes in climate-driven hydrologic forcing. We investigate how timing and magnitude of discharge have shifted over a 65-yr period of record and yielded changes in flow permanence along mountain stream networks. Finally, we consider how our findings may inform current and future protections for streams under the proposed Waters of the U.S. (WOTUS) Rule (US DoD US EPA, 2018) and subsequently finalized in the Navigable Waters Protection Rule (US DoD US EPA, 2019). We selected the 5 th -order Lookout Creek basin (Western Cascade Mountains, Oregon, USA) because of the extensive and long-term network of gages on low-order streams (Table S2). Furthermore, this basin is representative of the broader Pacific Northwest where climate change impacts on the timing and magnitude of moisture delivery to high elevation watersheds are known to cause declines in large rivers (Luce and Holden, 2009;Luce et al., 2013). Thus, reduced orographic enhancement of precipitation due to climate change is projected at the field site. This study considers the cascading impact of this change on stream discharge, and how discharge changes in headwaters may change flow permanence and connectivity in a river network.

Site Description and Available Data
The study was conducted at the H.J. Andrews Experimental Forest (HJA), a 5th-order basin in the Western Cascades, Oregon, USA (site map in Figure S1). The basin drains about 6,400 ha of forested landscape, with elevations ranging from about 410-1,630 m a.m.s.l, making it an ideal place to evaluate the impact of a changing climate on river networks of the broader Pacific Northwest. The basin has been a long-term study site for ecological and forest management research for more than 70 years and is relatively pristine with no urban land use, no dams or reservoirs, and minimal logging during the period of this study. At the longest currently operating meteorological station (CS2MET, elev. 485 m a.m.s.l.) annual precipitation averages 2,345 mm and average annual air temperature is 9.2 deg. C. Additional summaries of temperature and precipitation including trends for the period of record for each station are summarized in Table S1 and Figures S2-S15). In general, significant trends in monthly precipitation and air temperature are infrequently detected, due largely to the short observational records at the local meteorological network (Table S1). Further details about the local climate, morphology, geology, and ecology are comprehensively described elsewhere (Dyrness, 1969;Swanson and James, 1975;Swanson and Jones, 2001;Jefferson et al., 2004;Cashman et al., 2009;Deligne et al., 2017).
The HJA includes a network of 10 stream gages with drainage areas ranging from 8.5 to 6,241.9 ha, with records of 38-69 years of data across the gage network (Table S2). Additionally, a high quality digital elevation model derived from an airborne LiDAR survey is available for the entire basin, which has been reliably processed to extract topographic metrics including valley width, valley slope, and stream sinuosity (after Corson-Rikert et al., 2016;Schmadel et al., 2017;Ward et al., 2018a,b). For each gage we also calculated annual summary metrics of discharge including annual minimum discharge, mean discharge, maximum discharge, exceedance discharges (1, 5, 10, 25, 50, 75, 90, 95, 99th percentiles), total annual discharge, and the days elapsed to various cumulative fractions of discharge (1, 5, 10, 25, 50, 75, 90, 95, 99%).
Finally, an extensive data collection effort spanning the stream orders and lithologic regions of the basin was completed in 2015, providing a database of stream and valley morphologies and hydraulic conductivities to inform network-scale model parameterization (Ward et al., 2019b).

Simulation of the River Network
Simulation of network expansion and contraction followed the methods, data sources, and implementation described by Ward et al. (2018a). Briefly, the approach conceptualizes the river corridor  Ward et al. (2018a) validated the model in a 2nd order catchment within our study basin, concluding the model was appropriate to represent network expansion and contraction based on correct prediction of flowing or dry streambed conditions for more than 95% of over 3.2 million observations. We proceed with implementation of this model across a 5th order basin on the basis of Ward's et al. (2018a) success within our study basin, particularly given the accuracy in representing network expansion and contraction in response to diurnal fluctuations driven by evapotranspiration, storms, and seasonal baseflow.
The model was intentionally derived and constructed to require geomorphic and hydrologic data that are readily estimated for unstudied catchments (Ward et al., 2018a), consistent with our application in this study. Valley width, valley slope, along-stream slope, sinuosity, and the lateral contributing area for each 10-m segment of the valley bottom were extracted from the LiDAR data using a modified version of the TopoToolbox (Schwanghart and Kuhn, 2010;Schwanghart and Scherler, 2014;Schmadel et al., 2017). Stream width at each location was estimated using a power-law regression of 62 observations collected in August 2015 (Ward et al., 2019b) where b is the width of the channel (m) and UAA is the total drainage area (ha), and the best-fit relationship had a coefficient of determination of 0.84. We assigned a uniform Manning's n of 0.045 along the entire network based on visual inspection while working in the basin. Ward et al. (2018a) identified hydraulic conductivity as the largest source of uncertainty in the model. In response, we measured hydraulic conductivity of the streambed at 57 sites in August 2015 (Ward et al., 2019b) and assigned the geometric mean, 1.53 × 10 −4 m s −1 , across the network. Porosity was assigned as 0.3 at all locations, the midpoint of past studies (Dyrness, 1969;Kasahara and Wondzell, 2003;Wondzell et al., 2009;Ward et al., 2017) and the same value used in model validation (Ward et al., 2018a). We set valley colluvium depth to a minimum value of 1-m (Gooseff et al., 2006;Ward et al., 2018a), increasing as: where h is colluvium depth (m) and w is valley width (m). This functional form was selected to reflect the limited measurements of subsurface colluvium depth that are available, including geophysical observations at several headwater locations (Crook et al., 2008;Ward et al., 2012) and along the 5th order reach of Lookout Creek (Wondzell, personal communication, and unpublished data).
To parameterize the total down-valley discharge at the upstream end of each 1st order segment, we calculated a unique power-law regression for the gage discharge and drainage area for each 15-min timestep simulated, and defined the discharge based on UAA. Thus, all available gage data, and their time variation, informed the upstream boundaries for model headwaters. Lateral inflows for each segment were estimated using the same powerlaw regression, where the change in UAA between the upand downstream end of each segment was used to calculate the associated change in discharge attributed to the lateral area Ward et al., 2018a,b). Finally, we used the threshold of at least 2.21 × 10 −4 m 3 s −1 to differentiate surface flow from dry streambeds after past studies in the basin (Ward et al., 2018a).
Consistent with Ward et al. (2018a), we underscore that reduced complexity models are intended to represent the dominant mechanisms and interactions in a system of interest. This necessarily comes at the expense of representation of complexity and heterogeneity within the system. While our model has been derived and validated for headwaters within the study basin, the parameterization detailed above requires simplifications. To that end, this model is most appropriately viewed as heuristic, consistent with common practice in the study of river corridors (e.g., Gooseff et al., 2006;Cardenas and Wilson, 2007;Trauth et al., 2013;Irvine and Lautz, 2015;Schmadel et al., 2016Schmadel et al., , 2017. At the scale of river networks, comparable models have been applied to study patterns and trends at large spatial scales at the expense of site-specific localized predictions (e.g., Gomez-Velez and Harvey, 2014;Kiel and Cardenas, 2014;Gomez-Velez et al., 2015;Schmadel et al., 2018). Thus, given the model's strong performance at the reach-scale within our study basin (Ward et al., 2018a), explicit design as a heuristic that can be implemented at minimally studied sites, the wealth of data available across our network, and the tradition of heuristic models to test hypotheses in river corridor science, we proceed with this approach.

Statistical Tests
All between-group differences were tested using one-way ANOVA, Kruskal-Wallis, and Mann-Whitney-Wilcoxson U tests. We report differences as significant only if p < 0.05 for all three tests.
For all trends in discharge metrics, flowing frequency, contiguous frequency, flowing length, and contiguous length, we used Mann-Kendall tests and Sen's slopes to define the significance and direction of trends. Decreasing trends are reported for p < 0.05 for the Mann-Kendall tests and a Sen's slope of < 0. Increasing trends are reported for p < 0.05 for the Mann-Kendall tests and a Sen's slope > 0. We report no significant trend for a Mann-Kendall test with p > 0.05 or if Sen's slope is zero.
Analysis of trends may be sensitive to the length of the data set and which years or timesteps are included (i.e., different starting or ending dates, or different trend lengths). Consequently, we analyzed trends for every moving window of 10 or more years for every metric considered in the study, including those related to discharge at gages as well as flowing and contiguous lengths. In the body of manuscript we report significance and direction based on overall trends for each analysis. In the supplemental information we also tabulate how many of the moving windows agree with the overall trends, the length of the single longest trend that opposes the overall trend, and the length of the single longest period with no significant trend. We also tabulated the mean, median, maximum, and minimum Sen's Slope for every analysis, and the number of trends that are increasing, decreasing, or exhibit no significant trend (see Table S4 for robustness of discharge trends, Table S5 for robustness of flowing length and connected length trends, and Figures S29, S30 for visualization of the annual flowing length and connected length trends).

Headwater Stream Discharge Is Declining During the Dry Season
Discharge is predominantly decreasing across all gages in the basin over a 65-year period of record (Figure 2). For example, the Lookout Creek gage at the basin outlet has decreasing discharge for 81% of the year (about 300 days), steady discharge for 4% of the year (about 15 days), and increasing discharge for 15% of the year (about 55 days) ( Figure 2C; remaining sites in Figures S16-S24; Tables S2, S3). The largest and most consistently decreasing trends are during the summer season when discharge is lowest. We found no increasing discharge for the driest 7-months of the year (April through October; Figure 2D).
Across the gage network, we find significant inter-and intraannual changes in the timing and magnitude of discharge. Annual mean, median, and total discharge are all decreasing for 9 of 10 gages across their periods of record. We also found decreasing annual minimum and maximum discharges for 7 of 10 stream gages, declining annual low-flows (75-99% exceedance flows) at all gages, and declining annual high-flows (1-25% exceedance flows) at 7 of 10 gages (Table S3). Conceptually, the changes in moisture delivery is causing an increased export of water during winter months (Luce et al., 2013; Table S3), as evidenced by the more rapid time to export the first 10% of streamflow each year. Consequently, less water is stored during the rainy season, resulting in decreased dry-season baseflow, and extended times to export the last 10% of annual discharge.

Decreased Flow Permanence Has Reduced River Network Connectivity
Using the stream gage data, topographic analysis, and published data collected in the basin, we simulated dynamic expansion and contraction of the network (Figure 3 and Figure S25) (Ward Frontiers in Water | www.frontiersin.org et al., 2018a). For the 65-year simulation period, declining discharge, and increasing early season export of water within the basin result in an overall contraction of the flowing network (Figure 3). We found the flowing network reaches a maximum length of about 40-km during the wet winter months and contracts to as short as 15-km during the driest periods of record ( Figure 3B and Figure S25). Flowing network length is a useful proxy for connectivity along the river corridor, where longer flowing lengths allow more rapid connection of hillslopes to downstream water and promote rapid export of energy and materials rather than internal transformation (van Meerveld et al., 2019).
Next, we found the connected network length plateaus at about 21-km during the wet winter months and contracts to as small as about 5-km under the lowest discharge conditions ( Figure S10). The connected network represents an average of 57.4% of the flowing network across the entire simulation (median 55.9%; range 8.9-79.8%). The connected network defines the migration corridor through which aquatic organisms may travel upstream from the basin outlet without encountering a dry streambed location.
We found significant declines in flowing length for 75.7% of the year (about 276 days) compared to 23.6% of the year with no-trend (about 86 days), and < 1% of the year (about 3 days) with increasing flowing length ( Figure 3C). The decreasing trends are common throughout much of the year except for the highest discharge conditions associated with spring storms and snowmelt runoff (April through June) when the network length is more steady. Connected length exhibits similar trends, declining for 66.7% of the year (about 243 days), no trend for 33.1% of the year (about 121 days) and increasing trends for < 1% of the year (about 2 days; Figure S25).
Decreasing Network expansion and contraction exhibit threshold behavior, generally consistent with past studies (Ward et al., 2018a;Prancevic and Kirchner, 2019). When discharge at the Lookout Creek gage is greater than about 1 m 3 s −1 , the flowing and connected lengths are nearly constant at their plateau values ( Figure 3B and Figure S27). Under these wet-season, high discharge conditions, the flowing length maximum reflects a  Figure S25.
geologic limitation on network expansion where the drainage network is sufficiently dense to drain additional precipitation from the landscape without developing additional channels. As discharge drops below 1 m 3 s −1 , the network dynamically expands and contracts in response to precipitation. Under these dynamic conditions, the capacity of the valley bottom is comparable to the down-valley discharge, resulting in the large variation in flowing length in response to minor fluctuations in discharge (Ward et al., 2018a).

Headwater Streams In Steeper and/or Wider Valleys Are the Most Sensitive to Climate Change
About 41% of the headwater stream network exhibits a decreasing surface flow frequency, with the remaining 59% exhibiting no change ( Figure 4A). No location had increasing frequency of surface flow. Similarly, 27% of locations decrease in frequency of connected flow, 73% have no change, and no sites are more frequently connected across the period of record ( Figure S11).
Declining trends in flowing and connected frequency are not evenly distributed through the year. Instead, we found few significant trends for any segment during the wet season (November through June) because maximum network extent is controlled by basin morphology and drainage density ( Figure 4B). During the dry season (July-October) we found declining frequency of surface flow and contiguous flow in many network segments, due to declining discharge during this period. Similarly, trends in flowing and connected frequency are not evenly distributed in space. The reaches with the largest declines in flowing and connected frequency have significantly smaller drainage areas, steeper valley slopes, and/or wider valleys compared to locations with no trend (Figure S28), consistent with past findings in smaller catchments and the conceptual model (Costigan et al., 2015;Ward et al., 2018a;Prancevic and Kirchner, 2019). Decreasing trends in both flowing and connected frequency are most prevalent at the most upstream extents of the network, making the lowest-order streams "canaries in the coal mine" to first detect the impacts of climate change on flow permanence.

Changing Flow Permanence Challenges Current Regulatory Strategies
Non-stationarity is now the dominant paradigm in water resources (Milly et al., 2008). In our study system, the peak and average connected lengths are significantly larger in the first 30 years than the last 30 years. From a practical perspective, some waters that would have been federally jurisdictional (herafter "jurisdictional") in 1982 (based on the period 1953-1982) may not be jurisdictional in 2018 (based on the period 1989-2018).
The revised definition for Waters of the United States, which defines the basis for a water receiving federal protections under  Figure 5). As the minimum number of flowing years for regulatory protection decreases, changes due to climate become negligible (e.g., horizontal range at Y = 15-year; Figure 5). In contrast, if the Navigable Waters Protection Rule was intended to provide a constant-in-time determination, it must explicitly adjust the definition of "typical year" with climate. The systematic contraction in our study system, and thereby loss of federal protections for streams and their nearby wetlands, is only one response to changes in climate. In a landscape where flow permanence increases due to changing climate, federal jurisdictional scope could increase. Our critique here is consistent with draft comments from the US EPA's Science Advisory Board (Honeycutt and Board, 2019).

CONCLUSIONS
While past studies have explored the reduction in discharge at downstream gages on large rivers (Luce and Holden, 2009), this study is the first to examine how known changes in precipitation (Luce et al., 2013) and discharge translate into changes in connectivity between mountain hillslopes and their headwaters. Compared to their 1953-1962 averages, the 2009-2018 network has contracted by 24.1% and 9.2% in flowing and connected length, respectively, during the driest months. The dynamic connections along the network underpin a host of ecosystem services that we expect to also vary with flowing frequency. The loss of ecological function of such streams could be irreversible, and time-variable jurisdictional protections complicate the protection of these important resources. These losses are relative to a mid twentieth century baseline, and while some function will be lost as flow permeance decreases, other functions could be amplified as a result of increased duration or frequency of non-flowing conditions. Decreases in streamflow during periods when water resources are in highest demand, as was recently observed across the western U.S. (e.g., Milly and Dunne, 2020), further highlights the need for the extending the approaches presented here to more river basins.
Simulations predict that reaches with smaller drainage areas and larger subsurface flow capacity are the most likely to change in their flowing and connected frequencies in response to climate change. Thus headwater locations with steep valley gradients, larger valley widths, and/or disproportionately high hydraulic conductivity (Wondzell, 2006;Ward et al., 2019a) should be closely monitored to assess catchment response to climate change. Importantly, there are a small number of critical locations within a valley that can cutoff entire upstream reaches from the contiguous network-one location with a wideror steeper-than-average morphology can transition to entirely subsurface flow. We observed this threshold disconnection when the Lookout Creek gage discharge dropped below about 1 m 3 s −1 . While this threshold is the result of local geologic setting, we expect other systems will exhibit similar threshold behavior as a function of subsurface flow capacity and discharge. Finally, we underscore that current regulations are not designed with climate change-induced shifts in flowing and connected frequency, which will complicate policy enforcement for protection of headwater streams. The conclusions presented here are specific to one river basin in the Pacific Northwest, but the modeling approach and interpretation were intentionally designed to be transferable to other river networks, enabling extended analysis with modest, commonly-available data.

DATA AVAILABILITY STATEMENT
Topographic data and stream discharge data are available from the H.J. Andrews Experimental Forrest LTER database (https://andrewsforest.oregonstate.edu/data) as data sets GI010 and HF004, respectively. Valley and stream geometries derived from the LiDAR are publicly available (Ward et al., 2019b). Model outputs including timeseries of flowing length, contiguous length, and flowing frequency for each modeled segment are archived in CUAHSI's HydroShare at: http://www.hydroshare. org/resource/2f9643bb5d85436ba997a466ba8ed653.

AUTHOR CONTRIBUTIONS
AW conceived of the study and led the modeling, analysis, and writing. NS completed the topographic analyses and contributed to conceptual and numerical model development. AW and SW secured funding to support the research. AW, SH, SW, and NS participated in the data analysis and writing of the study.