Future streamflow simulation in a snow-dominated Rocky Mountain headwater catchment

Mountain systems are sensitive to climate change (Kohler & Maselli 2009) due to temperature and precipitation effects on snowpack development and, subsequently, snowmelt derived water supply (Barnett et al. 2005). The snow hydrology of mountains, primarily in the form of snowpack and glacier formation and loss, regulate freshwater that supplies the world's major river systems including the Indus, Ganges-Brahmaputra, and Colorado (Beniston et al. 1997; Water 2009; Ficklin et al. 2013). The impacts of climate variability on water resources can affect areas ranging from several hundred square kilometers (e.g. Milly et al. 2005; Lubini & Adamowski 2013) to as little as 8 km² (e.g. Wang et al. 2008). Climate change impacts on hydrologic processes in mountains has resulted in significant interest from scientists, and resource managers, and decision-making bodies due to concerns regarding changes in the reliability of water supplies (Li et al. 2013).

Projected changes in temperature and precipitation may influence the mean hydrologic processes of river basins that can affect the frequency and magnitude of extreme hydrologic events (Praskievicz & Chang 2009). Sustainable future water management of mountainous basins relies on accurate representation of climate variables where 50–70% of the total precipitation may fall in the form of snow (Serreze et al. 1999), and the seasonal snowmelt of the spring and early summer may account for 50–80% of the total annual runoff (Stewart et al. 2004).

The headwater basins in the Rocky Mountains of Canada are critical for freshwater resources. However, for the past century, river flows for these basins have declined by an average of 0.22%/yr (Rood et al. 2005). The reason for these declines is unclear, however snow and ice accumulation and melt timing are suspected as mechanisms affected by increased temperature and changes in precipitation (Hamlet & Lettenmaier 1999). Within the Canadian Rockies, glacier area has already been reduced by 25%, attributed to climate change beginning in the late 1800s (Luckman & Kavanagh 2000). The mean annual temperature has been increasing in the upper elevation of these mountain systems since the 1950s (Luckman 1990), with reduced snowpacks and shifts in seasonal release of meltwater from these basins (Barnett et al. 2005; Lapp et al. 2005). Changes in the low-order headwater basins of the Canadian Rocky Mountains may not be reflective of overall basin water budgets (Peters et al. 2013). However, precipitation changes in the region coupled with the ecologically sensitive nature of headwater catchments warrants investigation into possible future shifts in hydrology.

Simulated climate change effects on geographically low-elevation agriculture-dominated watersheds have shown changes in overall water supply rates, potentially affecting the reliability of flows during times of high human usage (e.g. Chien et al. 2013; Novotna et al. 2014; Neupane & Kumar 2015). In headwater basins of the North American Rocky Mountains, climate change impacts on hydrological processes have been generally shown to be associated with an earlier onset of melting (White et al. 1998; Cayan et al. 2001; Mote et al. 2005; Stewart 2009; Tinkham et al. 2015), and decreased mean annual streamflow (Zhang et al. 2001; Rood et al. 2005). In this study, we simulated the streamflow of the Rocky Mountain watershed with detailed incorporation of snow/glacier data using a process-based hydrologic model (i.e. Soil and Water Assessment Tool, SWAT). We then estimated the effects of potential climate variability on key hydrological processes including precipitation and snowmelt in the study watershed. This study explores the regional hydrologic response to climate change, in view of the impacts on ecosystem-services and the oil sands industry under a range of climate projections. These estimations may be important to assess the timing and source of future water availability, with the emphasis on expected changes related to sustained streamflow and potential ecosystem functioning in the watershed.

For this study, we focused on the headwater watershed of the Athabasca River, referred to as the Upper Athabasca watershed, located in south-western Canada. This watershed covers an area of about 3,500 km², with an elevation range from 922 to 3,736 m above mean sea level (Figure 1). The Athabasca River is the second largest river in Alberta, originating from the Athabasca glacier of the Rocky Mountains located in Jasper National Park that flows north-east emptying into Lake Athabasca downstream. The watershed is part of the Mackenzie River that eventually discharges into the Arctic Ocean. The watershed downstream is critical for diverse ecosystems including a staging area for a large number of waterfowl, primarily during spring and autumn seasons, and is recognized internationally as a RAMSAR wetland and UNESCO World Heritage site (Schindler et al. 2007; Pavelsky & Smith 2009).

Geologically, the watershed area is underlain almost entirely by sedimentary bedrock units ranging in age from Proterozoic to early Tertiary, with large variations in snow/glacier processes that cause landslides in the region (Jackson 2002; Selkowitz et al. 2002a). The climate of the western part of this region is mainly controlled by a stronger maritime influence originated from the North Pacific Ocean. The eastern part has a distinctly more continental climate (Selkowitz et al. 2002b). The annual temperature in the Lake Athabasca region for the period of 1971–2000 ranged between −3.5 and 7.6 °C with a mean value of 2.1 °C (http://climate.weather.gc.ca/climate_normals/index_e.html). Precipitation varies dramatically within the span of low and high elevation ranges. Mean annual precipitation is ∼504 mm, primarily in the form of snowfall, with higher precipitation in the north-eastern part of the region. These variations in climatic factors over relatively small distances affect microclimate, vegetation distribution, and ecosystem services of the region (Peterson et al. 1997).

Precipitation and stream discharge of the Athabasca River basin are shown in Figure 2. Analysis of 23 years (1979–2001) of precipitation data for eight different spatial locations of the basin (http://rda.ucar.edu/pub/cfsr.html) showed a higher amount of precipitation from May to August (Figure 2(a)). The minimum precipitation amount of 76 mm occurred in February and the maximum precipitation of 153 mm was observed in June. Stream discharge measured from May to September at Jasper in the Athabasca River basin (common outlet shown in Figure 1) showed a maximum value of 253 m³/s in July (http://wateroffice.ec.gc.ca/search/searchResult_e.html) potentially due to higher summer precipitation that increased surface runoff in these months. Seasonal analysis of precipitation and stream discharge data indicated that about 25% of total annual precipitation in the basin occurred during the spring season which corresponded to 11% of the mean annual stream discharge in the same season (Figure 2(b)). The minimum precipitation contribution of 20% occurred during the winter season, corresponding to a minimum discharge contribution of 4% during the same season. Water resources were most abundant during summer months with precipitation and snow/glacier-melt generally exceeding potential evapotranspiration and, therefore, the maximum discharge contribution of 67% occurred during this season.

Assessing climate change impacts on potential hydrologic processes of mountain basins is complicated due to climatic heterogeneity, lack of hydro-meteorological data, and uncertainty in snow/ice characteristics (Beniston 2003). There are various models developed for modeling watershed hydrologic processes subject to solid and liquid precipitation, including commercial software such as MIKE-SHE (http://mikebydhi.com) and public domain models such as HBV (Bergström 1992), Xinanjiang Model (Zhao et al. 1995), HEC-HMS (Azmat et al. 2016), SRM (Vafakhah et al. 2015), and SWAT (Neitsch et al. 2009). SWAT, a process-based hydrologic and water quality model developed for the USDA Agricultural Research Service (ARS), has gained international acceptance as a robust interdisciplinary model with a user-friendly interface set-up in a GIS framework that can integrate multiple environmental processes for development of better-informed policy decisions (Gassman et al. 2005). The model was originally developed to predict the impact of agricultural land management practices on water, sediment, and agricultural chemical yields in large-sized ungauged basins (Arnold et al. 1995; Neitsch et al. 2009). It has also been successfully calibrated and used to estimate the effects of potential climate variability on hydrologic processes for a number of diverse global basins (Stonefelt et al. 2000; Eckhardt & Ulbrich 2003; Song & Zhang 2012; Neupane et al. 2014, 2015; Awan et al. 2016).

A simple model was used to estimate stream water temperature from the SWAT simulations using weighted values derived from the contribution of different water sources to the stream water yield (Equation (1)). For each component, temperatures were assigned based on water sources. Rainfall was assumed to be the daily average temperature. Groundwater was assumed to be the annual mean temperature (≈4 °C) and snow/glacial meltwater was assumed to be 0.1 °C.

Spatial data required for SWAT simulations included topography, land use, and soil properties. Topographic data were acquired from the global digital elevation model (GDEM) sourced from the Advanced Spaceborne Thermal Emission and Reflection Radiometer (http://gdem.ersdac.jspacesystems.or.jp/search.jsp) with a 30 × 30 m resolution. The GDEM data were then used to delineate the watershed boundary, stream network, and topographic characteristics such as terrain length and slope of the stream channels. We used classified multispectral data derived from the Environmental Satellite (ENVISAT) Medium Resolution Imaging Spectrometer (MERIS) (GlobCover, http://due.esrin.esa.int/page_globcover.php) with a 300 × 300 m spatial resolution as land use input into the model (Arino et al. 2009). Based on the inclusion of this land use data, about 7% of the total basin area is covered by permanent snow/ice (Figure 3(a)).

The only soil data available for this basin were from the Food and Agriculture Organization of the United Nations (FAO 1995; Reynolds et al. 1999) with a broad spatial resolution of 10 km. To refine soil characteristics similar to the grain size of the topographic and landcover data, we first calculated the topographic saturation index (TSI) (Neupane et al. 2015) from the GDEM data using the flow accumulation and slope functions in ArcGIS. We then defined two soil layers (shallow soil with <10 cm from surface and deep soil >10 cm from surface; Figure 3(b)) from the TSI data using the maximum soil depth found from the original FAO soil data to help scale our TSI values (Neupane et al. 2014). Other important soil characteristics such as bulk density and nutrient characteristics were used ‘as is’ with values input from a user soil database based on the FAO/UNESCO Soil Map of the World (FAO/UNESCO 2003).

When soil data were combined with land use and slope values, we derived a total of 291 HRUs based on minimum area threshold values of 5, 10, and 10% for each land use, soil, and slope categories, respectively. From the land cover data, permanent glaciers were defined in nine HRUs from three different sub-basins with a total area of 236 km² and a total volume of 31 km³. To account for topographic gradient, we defined elevation bands in every 500 m elevation within each sub-basin (Fontaine et al. 2002). For this, the mean elevation of each elevation band and percentage of the sub-basin area within that band were entered as the SWAT model has the ability to include up to 10 elevation bands within each sub-basin. Details of the elevation gradient in each sub-basin with the number of elevation bands are shown in Figure 4.

For model calibration and confirmation, we used 23 years (1979–2001) of daily precipitation and daily maximum and minimum temperature data derived from grid-based observed meteorological data. The data were obtained from the National Centers for Environmental Prediction (NCEP) Climate Forecast System Reanalysis (CFSR) (http://rda.ucar.edu/pub/cfsr.html) for eight different spatial locations and the details are given in Table 1. These point data represent an area of ∼38 km² and have been successfully used for hydrologic simulations in other watershed systems (Fuka et al. 2013; Dile & Srinivasan 2014). Relative humidity, solar radiation, and wind speed data were calculated using a weather generator input (.wgn) file incorporated within SWAT using daily measured temperature and precipitation data. To account for glacier melt data, three spatial locations (Stations A, B, and C, as shown in Figure 1) were established at high elevations of glaciated regions to simulate mean daily air temperature applying the mountain climate simulator (MT-CLIM) (www.ntsg.umt.edu/project/mtclim) for the period 1982–2001 using daily meteorological inputs derived from the CFSR station grid-number 526-1178 and incorporated with the degree-day model. The derived glacier melt data were then incorporated with the SWAT model simulation outputs of stream discharge to compare with further simulations. Temperature and precipitation lapse rates required by SWAT for the distribution of meteorological variables based on elevation bands within each sub-basin were derived from previous studies conducted in US mountain systems (Fontaine et al. 2002; Pradhanang et al. 2011). Daily and monthly measured stream discharge data taken from the Canadian Hydrometric station number 07AA002 (http://wateroffice.ec.gc.ca/search/searchResult_e.html) located at Jasper (Figure 1) were compared with model simulation outputs.

Stream discharge simulated for the period 1979–1991 using SWAT default parameter values was first compared with measured stream discharge data taken from the gauging station at Jasper to analyze the un-calibrated model performance, referred to as the pre-calibration simulation. Subsequent simulations were organized as calibration (1979–1991) and confirmation (1992–2001) periods. The initial three years of calibrated simulation outputs were disregarded as a model spin-up period that allows the model to cycle multiple times to minimize the effect of user estimated parameter values (Zhang et al. 2007). For calibration and confirmation analysis of the model, we applied the SWAT-CUP (Abbaspour et al. 2007; Abbaspour 2012). The 95% prediction uncertainty (95PPU) (P factor) and the thickness of 95PPU (R-factor) were used to evaluate the accuracy of calibration and uncertainty analysis (Rostamian et al. 2008; Arnold et al. 2012). The Sequential Uncertainty Fitting (SUFI-2) algorithm, a semi-automatic inverse modeling approach, was used for this study due to its better capability of handling a large number of parameters with a lower number of model runs (Yang et al. 2008).

For this study, we chose the B1 (low), A1B (medium), and A2 (high) emission scenarios as the representative of all extreme conditions expected for the 21st century that are known as the Special Report on Emission Scenarios (SRES) developed by the Intergovernmental Panel on Climate Change (IPCC) (Maurer et al. 2010). These scenarios include both natural and anthropogenic drivers of climate change. We used daily Bias-Correction Constructed Analogue (BCCA) average temperature and precipitation data estimated for the SRES (http://gdo-dcp.ucllnl.org/downscaled_cmip_projections/) (Maurer et al. 2010; Brekke et al. 2013) for the end of the 21st century. These data are downscaled at 1/8 degree (∼12 × 12 km) spatial resolution that are suitable for hydrologic assessment studies. The bias-correction follows a basic approach of smoothening monthly mean values to avoid abrupt discontinuity to compensate for dry months. This generally helps to narrow the differences to obtain the best fit model. The general circulation model (GCM) structure is a major source of uncertainty for estimating the hydrologic impacts (Kay et al. 2009; Bennett et al. 2012); however, the Coupled Model Intercomparison Project Phase 3 (CMIP3) multimodel dataset (Meehl et al. 2007) were used in our study, based on their wider applicability with better performance, and more specifically over the western part of North America (Werner 2011). Also, these datasets were selected based on the hydrological impact studies recently conducted at the Pacific Climate Impacts Consortium (PCIC) for multiple watersheds of western Canada (Shrestha et al. 2012; Schnorbus et al. 2014). The average temperature and precipitation data derived from eight GCMs were incorporated into the SWAT for estimating their effects on annual and seasonal hydrologic processes including precipitation, surface runoff, stream discharge, water yield, evapotranspiration, soil water content, snowfall, and snowmelt of the Athabasca River basin.

To study the potential hydrologic changes, specifically for the end of the 21st century, the model simulations were run for the period 2081–2099 on a monthly basis. However, we ignored the initial three years (2081–2083) of simulation outputs for our analysis as a model spin-up period to compare with the baseline simulation results.

For this study, the key hydrologic parameters used for the Athabasca River basin SWAT simulations were selected based on the information derived from previous studies conducted in complex snow-dominated mountain basins (e.g. Pradhanang et al. 2011; Neupane et al. 2014, 2015) (Table 2). Model calibration and confirmation analyses are presented in Table 3. The correlation between measured and simulated discharge values during calibration and confirmation were improved when snow and glacier melt were specifically incorporated. The model showed the highest correlation with observed data for monthly simulations for the calibration of the band-glac scenario with NSE = 0.92, PBIAS = −14, RMSE = 27, and r² = 0.94. The statistical values were 0.85, −9, 34, and 0.87 for these indices, respectively for the band-glac scenario for the simulation confirmation period. These higher correlation indices indicate the importance of including glacier processes for predicting stream discharge simulation in the study watershed, similar to the findings of Abbaspour et al. (2010). Those researchers found that the highest correlation was obtained with the addition of glacier melt water, though their simulated area was larger than that used in our study. Rahman et al. (2013) also showed a similar high correlation between measured and model simulated stream discharge in the complex glacier-dominated upper Rhone River watershed specifically by the addition of elevation bands and glacier melt data in the basin; however, they were focused on a relatively larger glacier basin considering only nine hydrologic parameters for calibration and confirmation of the model.

Improved statistical values between measured and model simulated stream discharge data were also clearly represented by the hydrographs obtained after model calibration (Figure 5). However, the model showed overestimation during peak summer flows for daily simulations, potentially due to higher glacier melt estimation combined with measurement errors occurring during the same season (Rossi et al. 2009). The model also overestimated the winter baseflow components for both daily and monthly simulations that may be attributed to the difficulty in accounting sub-surface flow contribution to stream discharge, including a lack of proper representation for near stream saturation associated with excess runoff (Larose et al. 2007). Additional research is needed to address these discrepancies between measured and model simulated outputs, especially important for assessing seasonal changes in discharge.

In Rocky Mountain basins, the ablation of glaciers intensifies in the summer and this together with higher snowmelt at high elevations prolongs the higher flows during this season (Woo & Thorne 2003). There are only a limited number of studies conducted in parts of the North American Rocky Mountain system regarding changes in snowmelt hydrology (Ahl et al. 2008; Abbaspour et al. 2010; Watson & Putz 2013). However, ours is one of the first to identify the impacts of climate change at a smaller spatial scale.

The most sensitive input parameters regarding stream discharge simulation were identified on the basis of global sensitivity analysis, and are presented in Table 2. This analysis identified that the five most sensitive parameters for this study were SMFMX (maximum snowmelt rate), SMTMP (snowmelt base temperature), CN2 (surface runoff curve number), TLAPS (temperature lapse rate), and SFTMP (snowfall temperature), out of which four were observed as snow-related parameters. Betrie et al. (2015), in a study conducted for the entire Athabasca basin, also found the same parameters as the most important for SWAT simulations, however with different calibrated values reflective of the scale differences. The calibrated lower SMFMX value of 0.54 likely indicates minimal melting of snow at the beginning of the summer season in the basin. The SMTMP variable is the threshold temperature above which snowmelt occurs and therefore influences the simulated hydrographs' shape and peak flows. The negative calibrated SMTMP value indicates an early start of the melting process. We found the TLAPS as another influencing variable to influence the accuracy of simulated stream discharge. This is likely due to the extreme topography of the basin, similar to other Rocky Mountain watersheds affecting adiabatic lapse rates over short geographic distances influencing the amount and type of precipitation (Minder et al. 2010). The higher negative TLAPS value of −7.0 °C/km (Table 2) may indicate the influence of local climatic factors such as the presence of higher moisture-bearing winds in the basin coming from the Pacific Ocean. The CN2 variable is a runoff coefficient obtained by calculating the amount of surface runoff following a precipitation event and assigned to each HRU based on land use, soil type, and moisture content. As surface runoff is extremely sensitive to CN2, higher values increase surface runoff, reduce the infiltration rate, and decrease the groundwater recharge (Singh et al. 2005). The EPCO (plant uptake compensation factor), GWQMN (threshold depth of water in shallow aquifer required for return flow), ESCO (soil evaporation compensation factor), and SURLAG (surface runoff lag time) variables were found to be the least sensitive for simulating stream discharge. Uncertainty analysis of simulated stream discharge data also showed significant over- and underestimation values, mainly during the summer season (Figure 6), that may be attributed to measurement errors occurring during high flow seasons (Rossi et al. 2009).

Based on our simulation results, the effects of potential changes in climatic variables, specifically temperature and precipitation, are likely to bring substantial changes to the discharge regime of the Athabasca River. The most important results of this study indicate that streamflow for the Athabasca River in Jasper National Park may be 86% less (Figure 7(c)) than current flows by the end of this century, caused by climate change. The reason for this decline is attributed to decreased precipitation with minor effects associated with changes in snowmelt. The projected changes from the GCMs used for this study are not extremes, rather averages of 16 different models. These results support a continuing trend of lower water budgets for this catchment that appears to have begun in the 1950s at an approximate rate of −0.22% per year (Rood et al. 2005; Peters et al. 2013). Our results, based on the predicted 86% less discharge for the 2080s, indicated a decline of nearly 1.09% per year.

Monthly and annual changes in temperature and precipitation derived for B1, A1B, and A2 scenarios for the Athabasca River basin are presented in Table 4. Mean annual temperature was estimated to be higher than the current value with values ranging from 0.9 to 2.5 °C. The mean monthly temperature varied substantially, ranging from −2.7 to +5.2 °C. We found future precipitation decreased in all scenarios with a maximum change of approximately −86% in SRES-B1 (Figure 7(a)). In addition, maximum precipitation change was projected during the summer months and minimal change during the winter months. Lower precipitation corresponded with reduced surface runoff by 91% in the same scenario (Figure 7(b)). The SRES-B1 scenario had 35% lower ET than the baseline (Figure 7(e)) associated with a coincident reduction in the soil water (Figure 7(f)). The SRES-A2 scenario had the highest change in soil water of −70% compared to the baseline, likely as a result of increased temperature with resultant effects on latent heat exchange.

Higher temperatures combined with lower precipitation in the basin resulted in less snowfall and snowmelt for all scenarios with maximum changes of −89% (Figure 7(g)) and −79% (Figure 7(h)), respectively, for the SRES-B1 scenario. Future lower precipitation combined with reduced meltwater contribution will affect future water availability (Kundzewicz et al. 2007), exasperating current conflicts in human water use and ecosystem services in this drainage basin (Schindler & Donahue 2006; Squires et al. 2009).

Future seasonal changes in hydrologic processes of the Athabasca River catchment are shown in Figure 8. At present, the seasonal partitioning of contributed precipitation in the basin were estimated to be 19, 25, 33, and 23% for the winter, spring, summer, and autumn seasons, respectively (Figure 8(a)). Shifts in precipitation inputs based on the GCM data show, for example, that for the SRES-A2 scenario, more precipitation is expected in winter and autumn seasons of 29 and 27%, respectively, and less summer precipitation, 24%, as compared with the current climate. For the SRES-B1scenario, spring precipitation decreased from 25 to 19% for the basin.

Despite a lower precipitation contribution during the spring season, we estimated substantially higher surface runoff during this season in all scenarios, with a maximum contribution up to 47% in the SRES-A2 scenario (Figure 8(b)). This corresponded to a higher spring snowmelt contribution of 60% in the same scenario (Figure 8(h)). Lower summer surface runoff of 48% in SRES-A2 corresponded to minimal precipitation, leading to a reduced stream discharge contribution of 39% in the same scenario (Figure 8(c)). Due to potential higher spring surface runoff, total water yield of the basin also increased in all the scenarios up to 28% in SRES-A2 (Figure 8(d)). Similarly, reduced precipitation coupled with earlier spring snowmelt also reduced the water yield substantially during the summer season in all scenarios with a minimum contribution of 37% in SRES-A2. Due to changes in temperature and precipitation, projected snowfall in the basin also changed with a substantially lower contribution from melting during the summer season (Figure 8(g)). Therefore, the summer season was projected to be more dry with higher evapotranspiration (Figure 8(e)) and lower soil water content (Figure 8(f)). Climate change impact studies in recent years have begun for regional scale estimations. For example, Barnett et al. (2008) evaluated climate change impacts on water supply and regional hydrology in the western part of the United States using a high-resolution model and, similar to our estimations, they assessed increased runoff during spring season and a significant reduction during summer months at the peak of water demands.

For the baseline simulation, we estimated a maximum summer contribution of 58% for total mean annual stream discharge (Figure 8(c)). This flow includes the current higher summer precipitation combined with a higher glacier contribution during summer months in the basin. The large meltwater proportion is also due to abundant current snowpack levels in addition to seasonal glacier ablation from the high elevation of the Rocky Mountains, primarily during the summer season (Woo & Thorne 2003). However, this hydrologic peak during the summer is likely to diminish as glaciers in these basins continue to rapidly retreat in response to ongoing climate warming (Marshall & White 2010; Grover et al. 2014). Our simulations show decreased summer discharge with increased flows during winter and spring seasons in the future. This corresponded to the findings of Shepherd et al. (2010) who investigated similar seasonal changes in future water availability of Oldman River Basin in the North American Rocky Mountains.

Stream water temperature estimates for this study indicated that changes in water sources and potential increase in air temperature affected summer values most, with the largest difference between baseline estimated at 5.0 °C and the SRES-B1 scenario at 7.3 °C (Figure 9). For the SRES-A2 scenario, autumn stream water temperatures were elevated compared with all other simulations.

Lower streamflow and projected changes in seasonal flow caused by lower precipitation and a potential increase in temperature affects the snowpack of this Rocky Mountain basin and this will directly limit future water availability for any anticipated commercial appropriation, such as water for oil sands development or municipal water in this region (Pavelsky & Smith 2008; Squires et al. 2009). The shifts in quantity and timing of hydrologic discharge were estimated for this snow-dominated region over the coming century and which may influence freshwater ecosystem-services. Stream water quantity and flow timing influence critical water characteristics such as flow rate, dissolved chemical constituents, and water temperature affecting aquatic invertebrate and vertebrate life cycles (Fagre et al. 1997; Hauer et al. 2007). As water temperatures have been increasing at higher elevations in similar Rocky Mountain catchment reaches, sensitive aquatic invertebrate taxa may be running out of space with local extirpations, resulting in trophic shifts in mountain streams (Giersch et al. 2015).

Overall, warming is expected for the basin with higher temperature values significantly affecting both the water quality and quantity of the region (Chmura 2005). Water temperatures predicted here are crude estimates and do not include important factors such as energy balances associated with temperature effects on river ice or river hydrodynamics affecting ice formation and breakup. However, because changes in source and temperature increases during spring seasons are likely to occur, ice breakup in the river and subsequent flooding are likely to occur earlier in the season. Changes in seasonal hydrology can affect riparian ecosystems by changing the timing and frequency of flooding events (de Rham et al. 2008). These potential changes in temperature and precipitation combined with anthropogenic disturbance may increase the impacts on hydrologic components and surrounding ecosystem-services of the basin (Schindler & Donahue 2006).

Future hydrologic estimation in complex mountainous basins is complicated due to large climatic variations, and it is further elevated by the uncertainties of the General Circulation Models (GCMs) that are too coarse to adequately represent the orographic details of the mountains (Beniston 2003). Also, our research is based on a ‘static’ impact approach that does not include mesoclimatic interactions such as those from Pacific Decadal Oscillation, Arctic Oscillation, and the North Atlantic Oscillation which have been shown, in this basin, to greatly modify the timing and amount of streamflow (Burn 2008). While our simulations can be improved by using Representative Concentration Pathways (RCPs) laid under IPCC AR5 (IPCC 2013), these results should facilitate the consideration of implications for the possibility of future water limitations in this catchment and deliberation of climate change adaptation strategies.

The SWAT model was calibrated and tested for the snow-dominated upper Athabasca River watershed and the simulation results indicate that the model is a useful tool for assessing the effects of potential changes in temperature and precipitation on hydrological processes in the basin. Potential discharge of the Athabasca River was found to be dramatically lower for the 2080–2099 period, corresponding with lower precipitation, consequently leading to a reduced total annual water availability in the basin. Seasonally, winter discharge of the river was estimated to be higher, corresponding to similar changes in precipitation. Despite reduced future spring precipitation, discharge of the river was estimated to be higher, potentially indicating the importance of higher spring snowmelt in the basin. However, we estimated a substantial reduction in water availability in the basin during the summer season when there is a higher demand for water resources for agricultural, industrial, and domestic sectors. Finally, these assessments may improve our understanding of hydrological consequences of potential changes in temperature and precipitation of complex snow-dominated mountain basins, and provide better knowledge for future water resources management. Better representation of snowmelt processes through field-scale experiments should be investigated to refine model performance, primarily for simulating winter season base flow.

We would like to thank Drs P. Allen, J. Dunbar, S. Dworkin, and S. Alexander at Baylor University for their constructive suggestions during the development of this manuscript. Special thanks to the reviewers and editors for their constructive comments and suggestions, which improved the quality of the manuscript.