Multiple Indicators of Extreme Changes in Snow ‐ Dominated Streamflow Regimes, Yakima River Basin Region, USA

: Snow plays a major role in the hydrological cycle. Variations in snow duration and timing can have a negative impact on water resources. Excluding predicted changes in snowmelt rates and amounts could result in deleterious infrastructure, military mission, and asset impacts at military bases across the US. A change in snowpack can also lead to water shortages, which in turn can affect the availability of irrigation water. We performed trend analyses of air temperature, snow water equivalent (SWE) at 22 SNOTEL stations, and streamflow extremes for selected rivers in the snow ‐ dependent and heavily irrigated Yakima River Basin (YRB) located in the Pacific Northwest US. There was a clear trend of increasing air temperature in this study area over a 30 year period (water years 1991–2020). All stations indicated an increase in average air temperatures for December (0.97 °C/decade) and January (1.12 °C/decade). There was also an upward trend at most stations in Feb ‐ ruary (0.28 °C/decade). In December–February, the average air temperatures were 0.82 °C/decade. From these trends, we estimate that, by 2060, the average air temperatures for December–February at most (82%) stations will be above freezing. Furthermore, analysis of SWE from selected SNOTEL stations indicated a decreasing trend in historical SWE, and a shift to an earlier peak SWE was also assumed to be occurring due of the shorter snow duration. Decreasing trends in snow duration, rain ‐ on ‐ snow, and snowmelt runoff also resulted from snow modeling simulations of the YRB and the nearby area. We also observed a shift in the timing of snowmelt ‐ driven peak streamflow, as well as a statistically significant increase in winter maximum streamflow and decrease in summer max ‐ imum and minimum streamflow trends by 2099. From the streamflow trends and complementary GEV analysis, we show that the YRB basin is a system in transition with earlier peak flows, lower snow ‐ driven maximum streamflow, and higher rainfall ‐ driven summer streamflow. This study highlights the importance of looking at changes in snow across multiple indicators to develop future infrastructure and planning tools to better adapt and mitigate changes in extreme events.


Introduction
About one-sixth of the world's population is dependent on seasonal snowpacks and glaciers for water resources [1]. There is mounting evidence that the Northern Hemisphere is experiencing changes in snowpack characteristics [2][3][4][5][6][7]. More recent studies in-dicate that snowpacks are continuing to change in subarctic, Arctic, alpine, and mid-latitude regions [8][9][10]. For the period of 1980 to 2018, Pulliainen et al. [11] reported a decreasing trend of annual maximum snow mass for the Northern Hemisphere. Additionally, Coupled Model Intercomparison project 6 (CMIP6) models project strong negative snow extent and mass trends in the Northern Hemisphere [10]. In the western US, changes in snowpack could result in major consequences for the overall hydrology [3,5,[12][13][14][15][16][17]. For example, in western watersheds, snow is a critical surface water resource; 50% to 80% of total runoff is from spring/summer runoff [18][19][20], leading to a significant loss in water availability. Several studies show evidence that snowmelt is occurring earlier in the spring (e.g., [12,18,21,22]). For example, in the western US, there is an indication that snowmelt is taking place up to 26 days earlier [6]. Moreover, springtime snow water equivalent (SWE) has been declining since 1925 [12]. According to a more recent study by Mote et al. [23], a majority of the sites continue to show a decreasing trend. Not only is the timing of snowmelt changing, but snow duration is also declining [24]. Future predictions indicate a reduction in both total snow amount and snow duration for the western US (e.g., [18,25,[26][27][28][29][30][31]). In the western US mountain regions, [32] projects that snow duration will decrease by about 2 months by the mid-21st century, and there will be a 30% reduction in areal extent of wintertime snow-dominated area. Similarly, Lute et al. [33] estimated a decline in both the total number of days with snowfall and the SWE in the same period.
Spring runoff occurs earlier in snowmelt-dominated rivers, with a shift of up to 3 weeks earlier [34]. Even with no change in precipitation intensity, predictions indicate that this trend will continue, and that, by 2050, the maximum peak will take place approximately 1 month earlier [1]. Additionally, Casola et al. [35] indicated that a warming climate will result in a loss of about 20% of the 1 April snowpack in the Washington Cascade Mountains [35].
The Columbia River, located in the Pacific Northwest region of the western US, is the fourth largest river in the US, when comparing annual flow. The Columbia is mainly fed by mountain snowmelt, and approximately half of the annual flow is stored for flood control, hydropower, and irrigation [19]. A tributary of the Columbia River, the Yakima River, in the south central and eastern Washington State, houses important US military investments, such as Hanford and Yakima Training Center. To help with future land management planning in the Yakima River basin (YRB) region, Washington (WA), USA, a Strategic Environmental Research and Development Program (SERDP) study was funded to examine future impacts to snowpack at this site [36].
In this paper, we present the results of a 3 year SERDP study. This study included an assessment of our hypothesis that there are multiple indicators of a change in climate in the Yakima River region, with a specific focus on impacts on the US Army's Yakima Training Center (YTC). The paper considers trends in air temperature, snow accumulation, SWE arrival and departure dates, SWE peak, and streamflow within YRB and close-by watersheds as indicators of climate change.

Study Area
The YRB is located in the northwestern corner of the Columbia Basin, Washington ( Figure 1). Within the YRB is YTC, a 327,000 acre maneuver and live-fire training area [37]. The Yakima River is a main branch of Columbia River with 15,941 km 2 to its drainage outlet, and the upper basin above Union Gap is 9018 km 2 . The mean annual precipitation varies from 180 mm at lower elevations to 3000 mm at higher elevations [38]. Most of the runoff from the basin is generated from the snowpack [39]. The Bureau of Reclamation initiated the Yakima Project for irrigation efforts in the YRB in 1910. In 1932, the basin was also developed for hydropower but the power is presently mostly generated for agricultural use [40]. Several streams in YRB are regulated, where excess water is stored in reservoirs and released during low flow. From this process, Yakima River and its tributaries provide approximately 180,000 hectares of irrigation water to agriculturalists operating in the region [39]. Vaccaro [41] estimated that regulation of flow reduced the mean annual flow by about 20% at Union Gap and 48% near Parker.  Table S2), and SNOTEL stations (red solid dots, see Table S1).
(b) Elevation and (c) land cover from 2013 National Land Cover Data [42] for the snow modeling area.

SNOTEL Sites and Climate Data
In the US, repeated manual measurements of snowpack along staked out snow courses started in early 1900s. By the 1970s, the seasonal snowpack in the US was manually measured along snow courses at close to 2000 stations [19]. In the mid-1960s, a new fully automated and unattended snow telemetry (SNOTEL) network to measure snow depth was initiated. The network is operated by the US Department of Agriculture (USDA) Natural Resources Conversation Service (NRCS) and currently includes over 900 stations. Snow water equivalent (SWE) measurements are made using snow pillows, and other measurements include snow depths, minimum and maximum air temperatures, and precipitation.
For this work, we used SWE and air temperature from 22 SNOTEL stations, where seven stations were located within the YRB and the remaining stations were in close-by watersheds ( Figure 1). The selected stations had a record of at least 30 years and included data from water year (WY) WY1991-WY2020 (Table S1). The stations varied in elevation (el.) from 1100 m at Fish Lake, located in the Headwaters Cle Elum River watershed, to 2000 m at Harts Pass, located in the Upper Canyon Creek watershed. Air temperature and precipitation data for the Yakima Airport site were extracted from the Natural Resources Conservation Service (NRCS) database.

Streamflow Data
Streamflow gage information was compiled for eight US Geological Survey (USGS) gages within our study area (Table S2). Historical USGS streamflow gage and naturalized, modeled streamflow information was generated by the University of Washington's study of the Columbia River basin (Chegwidden,et al. [43], referred to as UW17). The UW17 dataset includes output from several hydrologic models, run with different parameterizations (P2 was selected for this work). For the purposes of this work, we selected output from the Variable Infiltration Capacity (VIC) hydrologic model that matched well with historical gage data [44]. We additionally included gages outside the Yakima River basin for comparison with the Yakima gages, which are heavily impacted by dams, diversion structures, and water withdrawals. The UW17 data also included future model projections, as described below.

Climate Scenarios
For the streamflow analysis, we utilized future projections from two different regional climate models (RCMs) and earth system models (ESMs)-Canadian Earth System Model v2, CanESM2 (CanRCM4) and Global Fluid Dynamical Lab Environmental System Model v2 Weather and Research Forecast model, GFDL-ESM2M (WRF)-that represented high (wet) and low (dry) scenarios of change as reported in the aforementioned SERDP study [36]. To project future streamflow changes, we focused on representative concentration pathway (RCP) 8.5 regional and earth system models scenarios, extracted from the Coordinated Regional Downscaling Experiment (CORDEX) [45]. Additional analysis of precipitation comparison between 2006 and 2099 for the CORDEX RCP 8.5 models that support the selection of our two future models is presented in Supplemental Figure S1.
Future streamflow predictions were extracted from the University of Washington's Columbia River Basin study and utilized the CanESM2 and GFDL-ESM2M models, based on multivariate adaptive constructed analogs (MACA) downscaled climate data and the VIC hydrology model. More details on the datasets can be found in the Supplementary Materials [46]. Future streamflow was obtained at eight USGS gages (ESMs; see Table S2). CanESM2 (CanRCM4) simulates regional climate model feedbacks and responses using the Canadian Earth System Model as boundary conditions and the Canadian Regional Climate Model for internal energy and water simulations [47]. The GFDL-ESM2M Earth System Model provides boundary conditions to the WRF model [47].

SnowModel Data
We present snow modeling simulation results derived from our SERDP study [36] to illustrate the spatial representation of the snow coverage and trends of our study area. In that study, we used SnowModel [48,49], a physically based model. We applied a 3 h time step for 36 years starting on 1 September 1979 to 2015. We used a 300 m grid increment, and the simulation domain was approximately 50,600 km 2 (Figure 1b,c). The dominant vegetation is cropland with upland shrubs and coniferous forest at higher elevations. We used topography from the National Elevation Database (NED) and land-cover data from the National Land Cover Dataset [42]. The meteorological data were from the National Aeronautics and Space Administration's (NASA's) North American Land Data Assimilation System (NLDAS, accessed at https://disc.gsfc.nasa.gov/datasets?keywords=NLDAS, accessed on: 5 December 2015). For this work, the spatial changes in snow duration, rainon-snow (ROS), and SWE runoff are presented to illustrate the general basin trends.

Air Temperature and SWE Trends
For each SNOTEL station and the Yakima Airport, we calculated both the average monthly air temperatures for each month from December through February and the average air temperature for the winter months (December-February) for WY1991-WY2020. From these averages, we performed a linear fit analysis and calculated the yearly increase in monthly air temperatures over our studied time period. We recognize that assuming a linear trend in time may underestimate temperature changes. If we were to fit a more complex model to the future project climate data, there would be considerable feedbacks that may lead to stronger increases in temperature. Thus, we believe that, by using a simple linear approximation for modeling, our extrapolation of the air temperature is an appropriate approach, particularly in combination with our other multiple (and nonlinear) methods applied in our work. We used this linear increase to estimate the future year when the average air temperatures would be above freezing. This was used to calculate the percentage of stations that would be above freezing (or already at positive temperatures) by 2020, 2040, 2060, 2080, and 2100.
We also identified maximum peak SWE and 1 April SWE for each of the years available for the SNOTEL stations. For each station, we performed a linear fit trend analysis from the time series (WY1991-WY2020) of maximum peak SWE and 1 April SWE. These results were used to determine the historical SWE change. Additionally, we investigated the trends of the SWE accumulation and snow-free dates for all years in the time period. These dates were chosen from longest period of continuous snow cover for each WY following the definition by Liston and Hiemstra [4].

Streamflow Trends and Generalized Extreme Value (GEV)
We analyzed the weekly, seasonal, and annual trends of maximum and minimum streamflow trends for the eight USGS gages using GFDL-ESM2M and CanESM2. We consider a historical time period (1991-2010) and a future time period (2070-2099), along with the entire duration of time . Both trends and the GEV approach examine the annual (YR) and seasonal (December, January, February (DJF); March, April, May (MAM); June, July, August (JJA); September, October, November (SON)) time series trends. Nonparametric trends were calculated using R's zyp package, using Sen's slope to estimate the magnitude of the trend and the Mann-Kendall to estimate the trend significance; significance was determined at a p-value ≤0.05.
We utilized a GEV distribution to further consider the type of changes projected in streamflow. GEV is based on the theory of extremal limits, which states that "a sufficiently long time series of block maxima will approach the GEV distribution asymptotically at large sample sizes" [50]. We applied a technique that allowed us to consider the best fit form of nonstationarity for the streamflow. For this component of the work, we considered only the maximum streamflow and examined the time period from 1991-2099. Additionally, for the GEV, we focused only on GFDL-ES2M streamflow projections, because of the strong similarities found between CanESM2 and GFDL-ES2M trend results (see Figures S2 and S3, Table S3).
The GEV approach in this study used the R-project's Generalized Extreme Values conditional density estimation network (GEVcdn) package [51,52], which was described in detail in Bennett et al. [53]. For this application, we followed the same criteria to select the 'best' model on the basis of the minimized Akaike information criteria, corrected for small sample sizes (AICc) [54]. We considered five candidate models, including a stationary (S), a linear nonstationary (LNS), and three nonlinear nonstationary models with changing the number of nodes from one to three. In the S model, the parameters were not allowed to vary in time, whereas location and scale parameters were allowed to vary in time in the other candidate models.

Air Temperature at SNOTEL Sites and Yakima Airport
To investigate regions vulnerable to changes in water storage and availability, as well as future increases in air temperatures, monthly and seasonal air temperature trends for each SNOTEL site and Yakima Airport were calculated for WY 1991-2020. From these trends, we projected air temperatures for each SNOTEL site and the years when temperatures above freezing would occur (Figure 2a). Such an increase in air temperature to above freezing indicates a lesser snowpack in that year, which would most likely result in a decrease in water storage. As expected, the air temperatures were lower at higher elevations. Surprisingly, the air temperatures in February increased at a lower rate than during other months. In fact, in February, three stations (Trough, Corral Pass, and Harts Pass) did not result in an increasing trend in air temperature, indicating that temperatures at these stations will continue to be below freezing during this century. We found a significant increasing trend in air temperatures that varied both temporally and spatially. The monthly average air temperature across stations was 0.97 °C/decade (0. 59   From our linear trend analyses, we also estimated the percentage of studied SNOTEL stations where positive air temperatures were encountered for every 20 years starting in 2020 through 2100, to see if similar trends continued (Figure 2b). In year 2040, several stations were already at above freezing air temperatures. In fact, more than one-third of the stations (41%) experienced above freezing air temperatures in the period from December through February. This percentage increased to 82% in 2060, and, by 2080, only one station (Harts Pass) was below freezing temperatures. In 2100, below freezing was only encountered at about half of the stations in the month of February.   Table S1). Gray lines indicate SWE for each WY and the black solid line is the average SWE for the WY1991-WY2020. The marks on the x-axis indicate the first day of each month.

SWE Analysis
From the SWE dataset (Figure 3), we calculated the 30 year linear trends of SWE accumulation ( Figure 4a) and snow-free date (Figure 4b) for all stations. Most stations were experiencing a later SWE start date and an earlier SWE end date. There was a strong negative trend in SWE starting date for stations with low to high elevation, where the longest delay in start date was at the lower elevation stations. The longest delay in start date was 18 days (Sasse Ridge). Opposite to this trend, some of the highest elevation stations were experiencing an earlier start date (Pigtail Peak, Lyman Lake, and Harts Pass). Green Lake, also a high-elevation station, was experiencing a later SWE start date of 11 days. A positive trend was seen in SWE end dates where the end date was earlier at lower elevations compared to higher elevations. The earliest SWE end date was 19 days earlier, and this occurred at Blewett Pass. This station was also the station experiencing the greatest difference in the number of total SWE days, with 29 fewer days in WY 2020 compared to WY 1991.
We also calculated the trend in peak SWE; there was a clear decreasing trend in peak SWE for the 30 year period (Figure 4c). The greatest decrease in peak SWE was 393 mm found at Stampede Pass. Among the 22 stations, this station was the third lowest elevation station. Even though it was one of the lowest elevation stations, it had a high SWE (max 1 April SWE peak of 2042 mm). In addition to investigating the trends in maximum SWE on 1 April, which is typically greatest and also commonly used by water resource planners, we also analyzed trends in maximum peak (Figure 4c). In general, the trends were similar; however, for the lower stations, the 1 April SWE was lower than the SWE peak.

Snow Model Simulations
Snow model simulations were performed to develop a spatial coverage of snow duration, ROS, and snowmelt runoff averages and trends from 1 September 1979-2015 for the studied area. Figure 5a-c shows the snow duration and trends during the core snow season. Following Liston and Hiemstra [4], we defined the core season as the longest period of continuous snow cover in each year (see the inset in Figure 5a). The average snow duration at lower elevations was less than 90 days during the 36 year period. For higher elevations, snow duration was over 300 days. Trends in snow-cover duration varied between ±40 days/decade. In general, lower elevations showed an increase in snow trends and, at intermediate and higher elevations, snow-cover duration increased. When averaging the yearly snow duration for the study area, it showed a declining trend of 20 days/decade for the 36 year period (Figure 5c). The number of ROS events was low at lower elevations and up to 50 days/year at higher elevations where snowpacks are longlived (Figure 5d). Most of the studied area indicated almost no change in ROS trends (Figure 5e). Overall, there was a decreasing trend of ROS events (Figure 5f). As expected, higher elevations corresponded to higher snowmelt runoff values (400 cm) compared to much smaller values at lower elevations ( Figure 5g). The entire modeling domain showed mostly negative snowmelt runoff trends, with −10 cm/decade at the highest elevations not uncommon. A decline in average snowmelt runoff was, therefore, no surprise ( Figure 5i).  [36]).

Streamflow Patterns and Trends
We analyzed historical streamflow and two ESMs at eight USGS gages (Table S2) for the mean, maximum, and minimum streamflow historical and future streamflow patterns ( Figure 6). The historical streamflow was lowest in August-September and increased gradually until peak flow in late spring or early summer, with several systems experiencing higher flows in October-December. In all but the Tieton River watershed, peak flow occurred close to mid-June to mid-July, with higher-elevation gages peaking later in the year. Future streamflow patterns exhibit large shifts, with streamflow peaking much earlier in the season (i.e., January through May) in the Yakima River basin, with the higherelevation systems (Wenatchee, Stehekin, Methow) peaking approximately 1 month earlier. Watersheds ranged in size from the Yakima River USGS station at Parker to the smaller American River at Nile station; hence, flows also varied on the basis of this factor.
Over the 109 year period, minimum streamflow trends are projected to decrease significantly in summer and fall at almost all stations with the exception of Tieton River and Stehekin River, where only small declines are projected. Spring minimum streamflow trends are expected to increase at all stations, with the exception of Methow River, where streamflow trends are projected to decrease in all seasons and annually. Maximum streamflow is expected to increase significantly in winter at all stations with the exception of Methow River. Furthermore, for maximum streamflow, trends for all stations are projected to significantly decrease in summer, again with the exception of Stehekin River, where summer streamflow is expected to increase (not significant). Annually, maximum streamflow patters are projected to increase over the 109 year period, while minimum streamflow trends are expected to decrease overall at most stations.

GEV Analysis
The GEV results illustrate how the maximum winter and maximum summer streamflow is projected to change in the future for the Yakima River basin region. For the GFDL-ESM2M data, the LNS model was minimized for the annual, winter, and summer streamflow maximums (Figures S6 and S7). The linear, nonstationary responses indicate that shifts in the systems will occur linearly with strong trends occurring in the systems over time, indicative of the overarching impacts of snowfall reductions and rainfall increases observed in the region. The lack of nonlinear responses in these systems reinforces that the changes are projected to be unidirectional under the snow-to-rain transition in regime. The GEV revealed strongly increasing maximum streamflow in winter (DJF) over the historical-to-future periods, while summer (JJA) maximum streamflow decreased over the period at all stations, further supporting the trend analysis results for the region. Annually, the three systems north of Yakima (Wenatchee, Stehekin, and Methow) are projected to exhibit some slight downward shift in trends by the 2050s, indicating a potential change in the systems, possibly indicative of the shift toward a loss of snowpack, transitioning toward a rainfall-dominant system by this time.

Discussion
Our study presented historical and future trends in multiple indicators (air temperature, precipitation (snow and rain), and timing of minimum, average, and maximum streamflow) for YRB, a snow-dominated region in the western US. We present three decades of air temperature trends, SWE measurements, and snow duration trends, in addition to simulated streamflow trends and GEV analysis from present to 2100. This study offers insight into how changes in climate affect this region, and how these changes impact streamflow. Even though hydrological simulations in a heavily irrigated watershed like the YRB can be difficult [55], results as presented in this study across multiple indicators represent the best estimate of how changes in climate could impact the region.
The air temperature trends for a 30 year historical period indicates that air temperatures are increasing at the majority of YRB stations and during most winter months except for February. Confirming these colder trends is the study by Kapnick and Hall [56], which reported lower warming trends for February for 12 western US states compared to other months at more isolated locations. As such, our results indicate that YRB is in a vulnerable transitional phase between precipitations falling as snow or rain. Hamlet and Lettenmaier [57] reported that, in a transient watershed, any change in temperature could have a great impact on runoff. Therefore, these increasing trends in air temperature will affect future water storage in the basin. For our study period (WY1991-WY2020), the highest increase in monthly air temperature observed was 1.82 °C/decade (January), recorded at a midelevation SNOTEL location (Park Ridge). We saw an increase in winter temperature (December-February) of 0.82 °C/decade as an average between the 22 stations. Similar warming was reported by Hu and Nolin [58], who indicated a warming trend of 0.4-1.2 °C/decade for the Interior West. Other studies also estimated warming rates of about 0.5 °C/decade [56,59,60], which are similar to the trends in yearly average air temperatures reported in this study (0.33 to 1.03 °C/decade). This extremely high rate of warming across the region that is already being observed historically is indicative of the rate of warming that has already occurred across the YRB.
As expected, associated with the increasing trend in air temperature, our analyses showed a decrease in both peak SWE and 1 April SWE over the 30 year period (WY1991-WY2020). Moreover, trends at most stations have experienced a delay in SWE start date and an earlier SWE end date, resulting in a shorter snow duration. A trend in shorter snow duration over the period of WY1980-WY2015 for the snow modeling area was also evident. Considering these changes, as recommended by other studies (e.g., [32,61,62]), a different date for peak SWE than 1 April should be considered for water resource planning.
Previous studies in the Washington area showed similar outcomes to our analyses, where major shifts in snowpack and precipitation transitions from snow to rain were identified. Klos, Link, and Abatzoglou [32] projected changes in winter precipitation (December-February) from late 20th century (1979-2012) for HUC4 watersheds, and Yakima indicated a 16% decrease in snow-dominated extent and a 25% increase in rain-dominated extent by the mid-21st century (2036-2050). Previous work predicted that snowmelt-dominated watersheds in the region will be reduced by 2080 [63]. Yakima River is one of three western US watersheds predicted to be almost entirely depleted in 1 April SWE [30] by the 2080s. The authors also reported that, by 2080, the Yakima watershed will experience a shift from a transient rain-snow basin to a rain-dominated basin [30].
In a larger part of the US, the major driver for severe spring flooding is snow meltwater [64,65]. ROS events are a large cause of peak flow events [64,66]. At higher elevations in snow-dominated areas, the projected future increase in ROS makes flood risk assessments complicated [67,68]. Our snow modeling simulations for the historical period of WY1980-WY2015 showed a similar outcome, where an increase of up to 4 days/decade was seen in some higher elevation areas.
Patterns in streamflow indicate a system in transition between the historical periods of 1970-1999, as well as into the future to 2099. Under a slightly drier future climate scenario, all streamflow at most sites is projected to experience declining summer minimum and maximum streamflow, with one high-elevation station (Methow) experiencing declines across all seasons and annually over the 1991-2099 period. Many systems are also projected to experience increasing spring minimum streamflow, due to higher winter streamflow influenced by increased rain in winter. On the other hand, winter streamflow is projected to increase strongly across the region. These finding have strong implications for increasing wintertime flooding events and drought scenarios in summer, highlighting the potential for increased extreme events of opposing direction that may occur throughout the water year.
When considering changes using the GEV approach, we observed a linear and strongly decreasing summer streamflow and increasing winter flows, again reinforcing the trends and the transition in the systems toward a rainfall-dominant winter streamflow pattern, with a shift in snowmelt streamflow peaks to an earlier and more drawn out melt, leading to lower streamflow availability during the summer. The lack of nonlinearity in responses highlights the combined impact of increasing temperature, reductions in snow, and increases in rainfall on the systems observed across the region.
In a study where the focus was to assess hydrology implications from climate change of Washington state, Elsner et al. [30] projected a similar shift to that reported in this study in future streamflow peak from summer to winter for Yakima River at Parker (see location in Figure 1) when compared to historical measurements. Our results are also comparable to a study by Vano et al. [69], who performed a hydrological seasonal sensitivity study to climate variations for several regions in the Pacific Northwest. They reported similar shifts in streamflow for Yakima River near Parker. In their study, out of the five basins studied, the Yakima River basin was the most sensitive to warming. The sensitivity is highlighted here in the observed trends in streamflow projected to occur, with similar responses noted for both the irrigated Yakima system and the adjacent, unimpacted river basins of Wenatchee, Stehekin, and Methow. These findings have implications for other systems adjacent to the sites we studied, including those systems draining into the Columbia River basin.
The changes observed in the multiple indicators have broad-reaching implications for flow. Water managers, agriculturalists, irrigation managers, and urban planners may require adjustments to planning and water retention scheduling to accommodate and retain winter rainfall storage for drier summer months in lieu of the reduced natural water storage in snowpack. Importantly, YTC will want to take into account the drier summer conditions, which may affect military missions. Hanford may want to consider the possibility of increased movement of solutes and sediments posed by the increase in winter runoff, as well as the overall increase in water moving through the surface, subsurface, and soils in the region. Additionally, the declines in summer streamflow will likely have an impact on other important components of the ecology and economy in the region, such as fisheries. Studies in the YRB indicate that a warming climate will decrease the salmonid habitat due to both increased streamflow during winter months and higher water temperatures [70].

Summary and Conclusions
Future climate predictions indicate that temperatures will continue to increase, and that this can negatively impact snow-dominated watersheds in regions that rely on the snowpack for water resources. Therefore, it is greatly important to investigate changes in climate during the design phase of new infrastructure. Moreover, military operations and missions can be affected by alterations in climate, making it necessary to account for these changes during future land management. This study focused on historical and future trends in climate, as well as their effect on streamflow for the Yakima River Basin, a snowdominated area located in the western US, in the same region as a large military installation (YTC). Our analyses included air temperature and SWE data from a total of 23 stations over a 30 year period. From the air temperature trends, we also projected by what year winter temperatures would be above freezing and would lead to a depleting snowpack. Additionally, we analyzed streamflow trends across five stations in the YRB and, because this watershed is heavily regulated, we also ran trend analyses for three close-by stations.
Our air temperature trend analyses agree with previous studies showing increasing air temperatures. The 30 year historical record of air temperature showed an increasing trend for almost each month starting in December through February, as well as for the 3 month average (December through February) air temperatures. Surprisingly, a few stations showed no increase in air temperature during the month of February, indicating this as the coldest of the months analyzed. More interesting, more than one-third of the stations (41%) are predicted to have above freezing air temperatures by 2040 in December through February. Most of these stations are located at lower elevations, but some stations are at higher elevations (i.e., Lost Horse and Green Lake). By 2060, a majority of the stations (82%) are projected to show air temperatures above freezing during this time period. The increase in air temperature will result in a shift in snow arrival and departure date, SWE peak flow, and maximum spring streamflow. The streamflow trend and GEV analyses indicated that the YRB basin is a system in transition, where streamflow is changing to earlier peak flows as a result of increasing temperatures. The findings from this study provide valuable input for planning infrastructure, water, and land management in the YRB region and locations with similar climate. Furthermore, this study highlights the local variability in climate, which proves the need for and importance of incorporating spatial differences during the design phase of hydraulic structures in order to avoid negative impacts. Table 1. Characteristics for the 22 automated SNOTEL stations from Natural Resources Conservation Service, including SNOTEL site name and number, latitude, longitude, elevation, and basin location. Table 2. Characteristics for the eight Washington state streamflow gages and sites analyzed in this study. Site name and abbreviated description, USGS ID, latitude, longitude, basin area, elevation of USGS gage, and basin location. Table S3. Seasonal and annual trends for minimum and maximum streamflow (m3/s) for 1991-2099 for GFDL-ES2M. Trends are values occuring over the 109 year time period. Bold values indicate significant trends (p-value <= 0.05). Figure S1. Precipitation comparison between 2006-2016 (current/historical) and 2070-2099 (future) for the CORDEX RCP 8.5 models that support the selection of our two future models. The CanESM2.CanRCM4 model shows the wettest winter and spring for the current climate compared to future projections and the GFDL-ESM2M.WRF model is the driest (winter and spring) in comparison to the future climate. Figure  S2. Left, GFDL-ESM2M for 1991-2099 maximum streamflow. Right, CanESM2 for 1991-