GLACIER CHANGE AND CONTRIBUTION TO STREAMFLOW IN THE CANOE BASIN, BRITISH COLUMBIA, 1948-2005 by Teresa Bessie Brewis B.Sc., University of Northern British Columbia, 2006 THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF SCIENCE IN NATURAL RESOURCES AND ENVIRONMENTAL STUDIES UNIVERSITY OF NORTHERN BRITISH COLUMBIA January 2013 © Teresa Bessie Brewis, 2013 1+1 Library and Archives Canada Bibliotheque et Archives Canada Published Heritage Branch Direction du Patrimoine de I'edition 395 Wellington Street Ottawa ON K1A0N4 Canada 395, rue Wellington Ottawa ON K1A 0N4 Canada Your file Votre reference ISBN: 978-0-494-94439-4 Our file Notre reference ISBN: 978-0-494-94439-4 NOTICE: AVIS: The author has granted a non­ exclusive license allowing Library and Archives Canada to reproduce, publish, archive, preserve, conserve, communicate to the public by telecommunication or on the Internet, loan, distrbute and sell theses worldwide, for commercial or non­ commercial purposes, in microform, paper, electronic and/or any other formats. L'auteur a accorde une licence non exclusive permettant a la Bibliotheque et Archives Canada de reproduire, publier, archiver, sauvegarder, conserver, transmettre au public par telecommunication ou par I'lnternet, preter, distribuer et vendre des theses partout dans le monde, a des fins commerciales ou autres, sur support microforme, papier, electronique et/ou autres formats. The author retains copyright ownership and moral rights in this thesis. Neither the thesis nor substantial extracts from it may be printed or otherwise reproduced without the author's permission. L'auteur conserve la propriete du droit d'auteur et des droits moraux qui protege cette these. Ni la these ni des extraits substantiels de celle-ci ne doivent etre imprimes ou autrement reproduits sans son autorisation. In compliance with the Canadian Privacy Act some supporting forms may have been removed from this thesis. Conform em ent a la loi canadienne sur la protection de la vie privee, quelques formulaires secondaires ont ete enleves de cette these. W hile these forms may be included in the document page count, their removal does not represent any loss of content from the thesis. Bien que ces formulaires aient inclus dans la pagination, il n'y aura aucun contenu manquant. Canada Abstract Glaciers provide late-summer runoff to the Columbia River Basin, but few studies quantify glacier runoff and its relation to recent climate variations. This study used digital photogrammetry to map glacier area and volume change in the Canoe Basin, headwaters of the Columbia River, from 1948 to 2005. Historic streamflow and glacier contribution to streamflow were simulated using the HBV-EC hydrologic model. From 1948 to 2005, glacier area and volume decreased by -3.8 ± 0.44 km2 (-6.8 ± 0.79%) and —0.23 + 0.14 km3, respectively. The majority of glacier change occurred from 1985 to 2005. Increasing streamflow from 1947 to 2006 is likely related to increased summer precipitation and not enhanced glacier wastage. Glacier contribution to streamflow decreased from 1947 to 2006, with a higher rate of decrease from 1985 to 2006 than from 1947 to 1985. From 1947 to 2006, glaciers contributed 14-16% (32-38%) of total annual (August) streamflow. Table of Contents Abstract................................................................................................................................................ ii Table of Contents............................................................................................................................... iii List of Tables....................................................................................................................................... v List of Figures................................................................................................................................... vii Acknowledgement.............................................................................................................................. x 1. 2. Introduction.................................................................................................................................1 1.1 Motivation..............................................................................................................................1 1.2 Research objectives.............................................................................................................. 3 1.3 Thesis outline........................................................................................................................4 Glacier change in the Canoe Basin, British Columbia, 1948-2005......................................5 2.1 Introduction........................................................................................................................... 5 2.2 Study A re a ............................................................................................................................ 7 2.3 Methods..................................................................................................................................9 2.3.1 Climatic Variability..................................................................................................... 9 2.3.2 Digital Photogrammetry............................................................................................ 13 2.4 Results..................................................................................................................................22 2.4.1 Climatic variability (1947-2006)............................................................................. 22 2.4.2 Stereo model verification and bias rem oval........................................................... 30 2.4.3 Glacier area and volume change...............................................................................32 2.5 Discussion........................................................................................................................... 39 2.5.1 Climatic variability.....................................................................................................39 2.5.2 Glacier terminus change............................................................................................41 2.5.3 Glacier area and volume change...............................................................................44 2.5.4 Sources of error and study lim itations..................................................................... 47 2.6 Conclusion.......................................................................................................................... 48 3. Quantifying glacier contribution to streamflow in the Canoe Basin, British Columbia, using the HBV-EC watershed model..............................................................................................50 3.1 Introduction......................................................................................................................... 50 3.2 Study A re a .......................................................................................................................... 52 3.3 Methods................................................................................................................................54 3.3.1 Data.............................................................................................................................. 54 3.3.2 The HBV-EC Hydrologic M odel..............................................................................57 iii 3.3.3 HBV-EC model set u p .............................................................................................. 64 3.3.4 HBV-EC calibration and validation........................................................................ 65 3.3.5 Data extraction and compilation..............................................................................68 3.3.6 Data analysis.............................................................................................................. 68 3.4 3.4.1 HBV-EC calibration and validation........................................................................ 71 3.4.2 Streamflow sensitivity to glacier area change in the HBV-EC m o d el............. 76 3.4.3 Simulated streamflow variability............................................................................ 78 3.4.4 Glacier contribution to streamflow......................................................................... 83 3.5 5. Discussion........................................................................................................................... 88 3.5.1 HBV-EC model perform ance.................................................................................. 88 3.5.2 Streamflow sensitivity to glacier area change in the HBV-EC m odel............. 90 3.5.3 S treamflow variability.............................................................................................. 90 3.5.4 Glacier contribution to streamflow......................................................................... 92 3.5.5 Sources of Error.........................................................................................................93 3.6 4. Results..................................................................................................................................71 Conclusion.......................................................................................................................... 94 Research summary and future studies................................................................................... 96 4.1 Summary of research.........................................................................................................96 4.2 Future studies...................................................................................................................... 97 Literature C ited.........................................................................................................................99 Appendix A: Comparison of spliced and rehomogenized Blue River climate d ata 108 Appendix B: Complete record of trends in the Blue River climate data................................ I l l Appendix C: Comparison of simulated streamflow generated by forcing the HBV-EC model with spliced and rehomogenized Blue River climate data......................................................... 113 List of Tables Table 2.1. Goodness-of-fit statistics for infilling monthly Blue River climate data for the entire record (n=279) and for the record divided into summer (May-September) and winter (October-April) seasons (n = 109)................................................................................................... 10 Table 2.2. Aerial photographs and aerial triangulation scans used in project.......................... 14 Table 2.3. Summary of trends in the Blue River mean, maximum and minimum temperature (°C yr_1) and total precipitation (mm y r-1) time series. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold...............................................................................23 Table 2.4. Summary of trends in Climate WNA (average Canoe Basin) mean, maximum and minimum temperature (°Cyr_1) and total precipitation (mm y r-1) time series. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold...................................................... 28 Table 2.5. Wilcoxon rank-sum test p-values for tests comparing Blue River climate variables during periods of calculated glacier change and known phase shifts of the PDO. Significant (p < 0.05) differences are shown in bold........................................................................................ 29 Table 2.6. Total RMSE for stereo models.....................................................................................31 Table 2.7. Checkpoint elevation statistics for 2005 and earlier stereo model comparisons.. 31 Table 2.8. Significance of linear biases present in stereo models. Bold biases were removed from the models................................................................................................................................. 31 Table 2.9. Morphological properties of select glaciers in the Canoe watershed, based on 2005 morphology...............................................................................................................................32 Table 2.10. Change in glacier terminus position and elevation in the ablation area of the glacier over the period of study. Negative values indicate retreat; positive values indicate advance............................................................................................................................................... 34 Table 2.11. Glacier area change from 1948 to 2005 in the Canoe Basin. Area is based on the fraction of glaciers that are covered in all years of photography................................................ 38 Table 2.12. Glacier volume change from 1948 to 2005 in the Canoe Basin............................ 38 Table 3.1. Area of land classes used in the HBV-EC simulations for each year..................... 55 Table 3.2. Linear regression equations for infilling daily climate data.....................................56 Table 3.3. Description of HBV-EC model parameters that require calibration (National Research Council, 2010).................................................................................................................. 59 Table 3.4. Previous studies that have used the HBV-EC model................................................ 60 Table 3.5. Parameter ranges used in each step of model calibration. Nash Sutcliffe statistic (E) and August error (AE) values are listed for the best-performing parameter set in each step...................................................................................................................................................... 72 Table 3.6. Range of calibrated parameters in the top eight parameter sets...............................73 Table 3.7. Nash Sutcliffe efficiency (E), absolute August Error (AE) and root mean square error (RMSE) values obtained during the calibration and validation periods for the top eight parameter sets. Simulated glacier volume change for the calibration period is also shown. v Simulated glacier volume change was not used in model validation; therefore, it is not shown.................................................................................................................................................. 73 Table 3.8. Summary of trends in the simulated streamflow record. Significant (p < 0.05) trends are shown in bold................................................................................................................... 78 Table 3.9. Mean, maximum and minimum glacier contribution to streamflow from 1947 to 2006.....................................................................................................................................................83 Table 3.10. Summary of model coefficients in the glacier contribution to streamflow record. Significant (p < 0.05) trends are shown in bold.............................................................................84 Table A .I. Summary of trends in the combined Blue River A and Blue River N stations for mean, maximum and minimum temperature (°C yr_1) and total precipitation (mm y r -1) time series. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold...............110 Table A.2. Summary of trends in the rehomogenized Blue River mean, maximum and minimum temperature (°Cyr_1) and total precipitation (mm y r-1) time series. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold.....................................................110 Table B .l. Summary of trends (°C yr_1) in the Blue River mean temperature data. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold.......................... 111 Table B.2. Summary o f trends (°Cyr_I) in the Blue River maximum temperature data. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold.......................... 111 Table B.3. Summary of trends (°Cyr_1) in the Blue River minimum temperature data. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold............................ I l l Table B.4. Summary of trends (mm y r-1) in the Blue River total precipitation data. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold.......................... 112 Table C .l. Performance statistics comparing streamflow simulated using the rehomogenized climate data to streamflow simulated using the spliced climate data. Statistics include root mean square error (RMSE), Nash Sutcliffe Efficiency (NSE) and relative volume error (%). Statistics were calculated for three periods in the simulated streamflow time series (19471006,1947-1968 and 1969-2006)................................................................................................ 116 Table C.2. Performance statistics comparing observed W ater Survey of Canada streamflow to streamflow simulated using the rehomogenized and spliced climate data sets. Statistics include root mean square error (RMSE), Nash Sutcliffe Efficiency (NSE) and relative volume error (%)..............................................................................................................................117 vi List of Figures Figure 2.1. Location of the Canoe Basin, British Columbia......................................................... 8 Figure 2.2. Glaciers in the Canoe Basin. The outlined glaciers have well-defined termini and were analyzed separately for morphology........................................................................................ 8 Figure 2.3. Scatter plots showing the overlapping period between the BRA and BRN climate stations (n=279). Regression lines are also shown (solid lines)................................................. 11 Figure 2.4. Glaciers covered by the 1948 photos overlain on a Landsat 5-4-3 composite (4 September, 1995).............................................................................................................................. 14 Figure 2.5. Glaciers covered by the 1955 photos overlain on a Landsat 5-4-3 composite (4 September, 1995).............................................................................................................................. 15 Figure 2.6. Glaciers covered by the 1970 photos overlain on a Landsat 5-4-3 composite (4 September, 1995).............................................................................................................................. 15 Figure 2.7. Glaciers covered by the 1985 photos overlain on a Landsat 5-4-3 composite (4 September, 1995)............................................................................................................................. 16 Figure 2.8. Glaciers covered by the 1996 photos overlain on a Landsat 5-4-3 composite (4 September, 1995). Additional area was covered in the 1996 photos; however, data collection was limited to the area shown due to high snow cover.................................................................16 Figure 2.9. Glaciers covered by the 2005 AT scansoverlain on a Landsat 5-4-3 composite (4 September, 1995).............................................................................................................................. 17 Figure 2.10. Mean annual, winter and summer temperature at Blue River from 1947 to 2006. The solid line is the 5-year running mean. Dashed lines show significant trends (p < 0.05; Kendall-Theil robust line)................................................................................................................ 24 Figure 2.11. Maximum annual, winter and summer temperature at Blue River from 1947 to 2006. The solid line is the 5-year running mean. Dashed line shows significant trend (p < 0.05; Kendall-Theil robust line)...................................................................................................... 25 Figure 2.12. Minimum, winter and summer temperature at Blue River from 1947 to 2006. The solid line is the 5-year running mean. Dashed lines show significant trends (p < 0.05; Kendall-Theil robust line)................................................................................................................ 26 Figure 2.13. Mean annual, winter and summer precipitation at Blue River from 1947 to 2006. The solid line is the 5-year running mean. Dashed line shows significant trend (p < 0.05; Kendall-Theil robust line)...................................................................................................... 27 Figure 2.14. Boxplots showing the difference in Blue River mean temperature and precipitation from 1947 to 1984 and from 1985 to 2005, periods o f calculated glacier volume change. Boxplots show the median and the 25th and 75th percentiles and the notches indicate the 95% confidence limits of the median. Whiskers extend to 1.5 times the interquartile range from the box and dots represent outliers..............................................................................29 Figure 2.15. Boxplots showing the difference in Blue River temperature and precipitation pre- and post-1977. Boxplots show the median and the 25th and 75th percentiles and the notches indicate the 95% confidence limits of the median. Whiskers extend to 1.5 times the interquartile range from the box and dots represent outliers....................................................... 30 vii Figure 2.16. An example of a systematic elevation bias present in the 1985 models before (left) and after (right) bias removal. Solid line shows the slope of the linear bias.................. 32 Figure 2.17. Glacier hypsometry mapped from (a) 1948 and (b) 2005 air photos................. 33 Figure 2.18. Elevation change of North Canoe glacier from 1970 to 1985............................ 35 Figure 2.19. North Canoe glacier terminus change from 1948 to 2005................................. 35 Figure 2.20. South Canoe and CC5/6 glacier terminus change from 1948 to 2005............... 36 Figure 2.21. Penny glacier terminus change from 1948 to 2005.............................................. 36 Figure 2.22. Zilmer W est and East glacier terminus change from 1948 to 2005................... 37 Figure 2.23. Difference in North Canoe glacier surface elevation from 1948 to 2005..........38 Figure 3.1. Location of the Canoe Basin, British Columbia.............................................. 53 Figure 3.2. Scatter plots showing the overlapping period between the BRA and BRN climate stations (n=5753). Regression lines are also shown (solid lines).............................................. 56 Figure 3.3. Land use classification used to set up the HBV-EC model............................ 64 Figure 3.4. Scatterplots of observed and simulated streamflow in the calibration (left) and validation (right) periods, for the top parameter set. One-to-one line is shown...................... 74 Figure 3.5. Scatterplot of observed and simulated annual flow volumes. One-to-one line is shown..................................................................................................................................................74 Figure 3.6. Comparison of observed and simulated total annual flow volumes based on the hydrologic year. Observed total flow volume is not available from 1995-1997......................75 Figure 3.7. Observed and simulated hydrographs for the years with (a) best and (b) worst model performance. Simulated discharge is displayed from the top performing parameter set. Day of year corresponds to the hydrologic year............................................................................76 Figure 3.8. Comparison of streamflow for the year 1998 simulated with 1955 and 2005 glacier extents, (a) Simulated streamflow using the best parameter set and difference between streamflow simulated with 1955 and 2005 glacier extents are shown, (b) Simulated streamflow for each glacier extent using the best eight parameter sets show uncertainty in the simulations. Day of year corresponds to the hydrologic year.................................................... 77 Figure 3.9. Simulated mean annual, summer and winter streamflow, based on the best parameter set. Solid black line shows the five-year running mean. Dashed black lines show significant trends (p < 0.05). Vertical lines indicate years when the glacier extent was updated in the model.........................................................................................................................79 Figure 3.10. Reconstructed hydrographs from end-points of the Kendall-Theil robust line, based on HBV-EC simulations from 1948 to 2006. Reconstructed hydrographs are shown for the top eight parameter sets....................................................................................................... 80 Figure 3.11. Reconstructed hydrographs from end-points o f the Kendall-Theil robust line, based on HBV-EC simulations for (left) 1948-1985 and (right) 1985-2005. Reconstructed hydrographs are shown for the top eight parameter sets.............................................................. 81 Figure 3.12. Time series of select 5-day mean periods. Time series are shown for the top eight parameter sets. The red lines show the 5-year running mean for the top parameter set. 82 Figure 3.13. Model output from the top eight parameter sets. Vertical lines indicate years when the glacier extent was updated in the HBV-EC model, (a) Simulated mean annual streamflow with (black) and without (grey) glacier cover, (b) Annual glacier contribution to streamflow, expressed as a fraction o f total streamflow............................................................. 84 Figure 3.14. Model output from the top eight parameters sets. Vertical lines indicate years when the glacier extent was updated in the HBV-EC model, (a) Simulated mean August streamflow with (black) and without (grey) glacier cover, (b) August glacier contribution to streamflow, expressed as a fraction of total streamflow............................................................. 85 Figure 3.15. Simulated hydrographs for the years with the lowest and highest glacier contribution to streamflow, which were 1999 (a) and 1958 (b). Plots show the simulated mean daily streamflow with (red) and without (blue) glacier cover. Day o f year corresponds to the hydrologic year. Simulated hydrographs are shown for the top eight parameter sets. 86 Figure 3.16. Annual (top) and August (bottom) glacier contribution to streamflow from 1947 to 2006. Solid black line shows the five-year running mean. The dashed black line represents statistically significant (p < 0.05) trends from 1948 to 2006. Vertical lines indicate years when the glacier extent was updated in the model..............................................................87 Figure 3.17. Day of onset of glacial melt relative to (left) mean annual temperature (R2 = 0.19, p = 0.0005) and (right) total annual precipitation (R2 = 0.17, p = 0.001). Solid black lines indicate significant linear trends. Days are based on the hydrologic year, 1 Oct - 30 Sep.......................................................................................................................................................88 Figure C .l. Simulated hydrographs for the years with (top) most similar and (bottom) least similar simulated streamflow using spliced and rehomogenized climate data. Day of year corresponds to the hydrologic year................................................................................................115 ix Acknowledgement I would like to thank my supervisor Dr. Brian Menounos for his advice, support and encouragement throughout this project. I would also like to thank my committee members Dr. R. Dan Moore and Dr. Peter Jackson for their guidance and insightful contributions to the project. I would also like to thank my external examiner Dr. Paul Whitfield for participating in my defence and for providing valuable and constructive comments on my thesis. I wish to thank Christina Tennant for her assistance with digital photogrammetry, and Dr. Roger Wheate for his assistance with air photo formatting, GIS and remote sensing. I am grateful to Georg Jost for providing HBV-EC model optimization codes adapted specifically for this study, and to Scott Emmons and Natalie Saindon for their assistance with model calibration. Finally, I would like to thank my husband and my parents for their love and support throughout the thesis process. This research was supported by the National Sciences and Engineering Research Council of Canada and the Columbia Basin Trust. x 1. Introduction 1.1 Motivation In British Columbia (B.C.) and the northwest United States (U.S.), water from the Columbia River is used for hydroelectric power generation, irrigation, agriculture, recreation and tourism. Unprecedented population growth in western North America pressures this already heavily managed freshwater resource (Bamett et al., 2005). The Columbia River Treaty between Canada and the United States sets the framework for the cooperative management of freshwater resources in the Columbia River Basin. Although the treaty does not expire, 2024 is the earliest date that either the U.S. or Canada can terminate the Treaty provided that at least 10 years advance notice are given (Columbia Basin Trust, 2012). Successful planning and management of the Columbia River for future decades requires an understanding of the complex relation between climate, glaciers and surface runoff. In glacierized watersheds, snowmelt produces high flows in the spring and glacier melt contributes flow throughout the summer. Glacierized watersheds typically have high flows during late summer and their annual peak flows may occur later than their nival counterparts (Young, 1990). Glaciers can also have a moderating effect on annual streamflow by storing water during cool and wet periods and releasing it in warm and dry periods (Fountain and Tangbom, 1985; Chen and Ohmura, 1990; Collins and Taylor, 1990). This moderating effect is particularly important to low-elevation areas that are prone to water shortages and drought conditions, such as the southern portion of the Columbia Basin. Glacier wastage in the last century was extensive worldwide (Dyurgerov and Meier, 2000) and significant throughout B.C. and western North America (e.g. Granshaw and Fountain, 2006; Moore et al., 2009; Bolch et al., 2010). In particular, glacier area decreased by 15% in both the southern interior of B.C. (Bolch et al., 2010) and the Columbia Basin 1 from 1985 to 2005 (Jost et al., 2012). As a glacier retreats in response to climate warming, streamflow initially increases in flow. However, as glaciers retreat below a critical area, the streamflow generated begins to decrease (Moore et al., 2009). Glaciers in many watersheds in south and central B.C. have already retreated below this critical area and, as a result, streamflow in these basins shows a decreasing trend (Stahl and Moore, 2006; Moore et al., 2009). In 2007, the Inter-governmental Panel on Climate Change (IPCC) summarized projected global surface temperature increases for the 21st century (Solomon, et al., 2007). Using emissions held constant to AD 2000 levels, global surface temperatures are likely to increase by 0.3-0.9°C by 2100. More realistic emissions scenarios project a surface temperature increase of 1.1-6.4°C (Solomon et al., 2007). In the Pacific Northwest, average annual temperature is predicted to increase by up to 2.1 °C by 2100 (Payne et al. 2004). It is expected that future warming will result in further glacier wastage which will negatively affect freshwater resources. More research is required to quantify glacier contribution to streamflow and examine how glacier retreat has and will continue to impact water resource availability. In British Columbia, glaciers are typically located at high elevation in remote locations and field work is resource-intensive, time-consuming and expensive. However, glacier change can often be studied remotely using aerial photographs and satellite imagery. Aerial photographs can often provide a historic record o f glacier change spanning several decades. Digital photogrammetry is an efficient tool to calculate glacier area and volume change (e.g. Luckman et al., 1987; DeBeer and Sharp, 2007; Bolch et al., 2010) and it can be successfully used with limited training (Gao and Liu, 2001). 2 Hydrologic models are a valuable tool for resource managers, and can be used to understand historic streamflow processes and predict streamflow response to climate variations under various warming scenarios. Model calibration is improved by using multiple criteria in parameter optimization, and glacier mass balance is particularly beneficial for calibration of models in glaciated watersheds (e.g. Braun and Aellen, 1990; Konz and Seibert, 2010). Taken together, digital photogrammetry and hydrologic modeling provide a reasonable method to study historic and future impacts of climatic variability and glacier change on surface flows in remote locations. The goal of this thesis is to improve our understanding o f the role that glaciers play in the surface hydrology of mountain environments. Geodetic data for glaciers in the Canoe River Basin is used in conjunction with a hydrologic model to quantify the importance of glacier runoff in the headwaters of the Columbia Basin. 1.2 Research objectives The goal of this study is to quantify glacier area and volume change in the Canoe Basin and examine trends in streamflow relative to glacier change and climatic variability. The specific objectives of this thesis are to: a) use digital photogrammetry to determine glacier area and volume change from 1948 to 2005; b) simulate streamflow and quantify glacier contribution to streamflow from 1947 to 2006 using the HBV-EC model; and c) examine streamflow and glacier contribution to streamflow relative to glacier change and climate variability. 3 1.3 Thesis outline The objectives of this study are addressed in the following chapters of this thesis. Chapters 2 and 3 are written as individual studies, each with separate methods, results and discussion sections. Chapter 2 reports on the use of digital photogrammetry to derive glacier area and volume change in the Canoe Basin from 1948 to 2005. The results are compared to observed climatic change from 1948 to 2005. Chapter 3 discusses calibration and use of the HBV-EC model to simulate streamflow and glacier contribution to streamflow in the Canoe Basin from 1947 to 2006. Trends in streamflow and glacier contribution to streamflow are assessed relative to glacier change and climate variability. Chapter 4 synthesizes major findings of the study and provides suggestions for future research opportunities. 4 2. Glacier change in the Canoe Basin, British Columbia, 1948-2005. 2.1 Introduction Glaciers have a moderating effect on streamflow. Glaciers store water during cool, wet years and release it during warm, dry years (Krimmel and Tangbom, 1974; Collins and Taylor, 1990). This buffering effect is particularly important in low elevation and drought susceptible areas, including the southern portion of the Columbia Basin (Hamlet and Lettenmaier, 1999). Glaciers worldwide substantially retreated during the 20th century (Dyurgerov and Meier, 2000; Moore et al., 2009), and this recession impacted the timing and quantity of glacier runoff in many mountain basins (Kaser et al., 2010). Glaciers respond to long-term changes in temperature and precipitation and largescale ocean-atmosphere climate regimes (such as the El Nino Southern Oscillation and the Pacific Decadal Oscillation) by adjusting their physical properties. Mass balance records are essential to understand how climate affects glacier mass balance; these records can also improve the calibration of hydrologic models (Braun and Aellen, 1990; Konz and Seibert, 2010; Jost et al., 2012). Hydrologic models can then be used to study the impacts of glacier and climate change on streamflow and water resource availability. Glacier change or mass balance can be measured using field techniques (surveying the glacier terminus and ablation stakes and digging snow pits to measure snow depth and density) or geodetic analyses (aerial photographs and satellite imagery). Glaciers are typically remote and difficult to access and thus field surveys are time-consuming, expensive and often dangerous. Traditional methods of measuring glacier mass balance in the field can yield the mass balance gradient or show changes in density throughout the snow/ice pack. Long term glacier mass balance records using traditional methods are only available at a few sites in western North America 5 (including Place, Peyto and Helm Glaciers; W orld Glacier Monitoring Service, 2012). Alternatively, aerial photogrammetry is an inexpensive and efficient way to study glacier area and volume change (Fox and Nuttall, 1997) even though geodetic methods are unable to recover the mass balance gradient. Geodetic data, including satellite imagery and aerial photographs, provide repeated coverage of much o f British Columbia. Several studies examine glacier change in the Canadian portion of the Columbia Basin; however, few studies examine these changes repeatedly and over multi-decadal time scales. Glacier volume change in the Columbia Mountains from 1985 to 1999 was quantified using aerial photography and data from the Shuttle Radar Topography Mission (Schiefer et al., 2007). Glacier area and volume change (derived from volume-area scaling) in the Columbia Mountains from 1951/2 to 2001 was estimated using aerial photographs and Landsat imagery (DeBeer and Sharp, 2007). Glacier area change in western Canada was calculated between 1985 and 2005 using Landsat imagery (Bolch et al., 2010). Glacier change studies with both short temporal resolution and long-term records are available at few sites in the vicinity of the Columbia Basin. Terminus behavior at Illecilliwaet Glacier, in the Columbia Basin, was recorded since 1887 on a near-decadal basis, using photographs, field surveys and remotely sensed imagery (Champoux and Ommaney, 1986). Luckman et al. (1987) examined terminus behavior of glaciers in the Premier Range, which includes the Canoe Basin, from 1948 to 1985 using aerial photography and field investigations. Terminus behavior from 1959 to 2007 at nearby Castle Creek Glacier was examined by dating annual push moraines (Beedle et al., 2009). The Intergovernmental Panel on Climate Change (IPCC) recognizes the importance of long-term studies of glacier change in developing our understanding of hydrologic response to climate change (Bates et al., 2008). The objective o f this study is to use digital 6 photogrammetry to measure decadal-scale changes in glacier area and volume from 1948 to 2005 in the Canoe Basin, British Columbia. Patterns of glacier change are compared among glaciers in the basin and overall glacier behavior is examined in relation to climatic variability. Glacier volume change and mapped glacier extents are used in Chapter 3 of this thesis, which uses the HBV-EC watershed model to simulate glacier contribution to streamflow. 2.2 Study Area The Canoe Basin is located in the Premier Range of mountains, west of Valemount, British Columbia, in the headwaters of the Columbia River (Figure 2.1). For the purpose of this study, the Canoe Basin includes the watershed area above the Water Survey of Canada (WSC) streamflow gauge 08NC004 “Canoe River below Kimmel Creek” (52.728° N, -119.408° W). The watershed encompasses 307 km2 o f mountainous terrain, has a median elevation of 1957 m above sea level (a.s.l.) and glaciers occupied 18% of the catchment in 2005 (Figure 2.2). At Blue River, approximately 65 km south of the outlet of the Canoe Basin, mean annual precipitation is 1100 mm of which 38% falls as snow. Mean annual temperature is 4.5°C with a monthly minimum and maximum of -9°C (January) and 16.4°C (July), respectively. As further explained in section 2.3.1 of this thesis, climate data from Blue River were used in this thesis because it provided the longest complete record of temperature and precipitation within proximity to the Canoe Basin. Although several meteorological stations are closer to the Canoe Basin, the length of those station records were shorter than those from Blue River (e.g. Cariboo Lodge, McBride 4SE) or they had extensive periods with missing data (e.g. Valemount East, Cariboo Lodge). 7 120* W 119* W 118* W BC/Alberta Bordei Glaciers Major W atersheds Rivers/Lakes himbla •Basin North Tl Basin 7-%C . ,--------------------------a—-lfj-------- ——— 120* W *"■* 119* W *t 118" W Figure 2.1. Location o f the Canoe Basin, British Columbia. n9"4o,o"w 119‘20'CTW ) Penny North Canoe South Canoe / CC6 ^J^Zilm er west Zilmer East r ' -52*40’0"N 52‘40'0'N2.5 l__ 75 10 Km i_____ | 119"40'crw 119"20'0T/V Figure 2.2. Glaciers in the Canoe Basin. The outlined glaciers have well-defined termini and were analyzed separately for morphology. 8 2.3 Methods 2.3.1 Climatic Variability Data Monthly mean climate data (mean temperature, maximum temperature, minimum temperature and total precipitation) from the Blue River North (BRN) and Blue River A (BRA) meteorological stations were downloaded from the National Climate Data and Information Archive (Environment Canada, 2010). These stations are located approximately 65 km south of the watershed outlet. Although several meteorological stations are closer to the Canoe Basin, the length of those station records was shorter than BRN and BRA (e.g. Cariboo Lodge, McBride 4SE) or they had extensive periods with missing data (e.g. Valemount East, Cariboo Lodge). A complete climate record was obtained by combining the records from BRN (1947-1985) and BRA (1969-2006) stations. Least squares regression based on monthly values was used to determine the relationship between climate variables (n = 279; Table 2.1, Figure 2.3) at the two stations during the period of overlapping records (1969-1985). The resulting equations were used to infill missing monthly values (0.5% of the climate record) in the BRA data, and were then applied to BRN data to extend the climate record back to 1947. The combined record of BRN and BRA stations is referred to as Blue River climate data in this study. For comparison, regression coefficients (R2) and standard error (SE) were also calculated using the monthly BRN and BRA data divided into seasons (summer: May-Sep; winter: Oct-Apr; n = 109; Table 2.1). While the regression coefficients are similar in the winter and summer, all standard error values were lower for the summer than the winter data. To ensure splicing the BRN and BRA records did not introduce spurious trends to the climate record, the monthly climate record created by combining the BRN and BRA stations was compared to Environment Canada’s rehomogenized monthly 9 climate data (Mekis and Vincent, 2011; Vincent et al., 2012) for Blue River (see Appendix A for comparison). Monthly mean temperature, maximum temperature, minimum temperature and total precipitation were also obtained from the Climate WNA v4.62 (http://www.genetics.forestry.ubc.ca/cfcg/ClimateWNA/ClimateWNA.html). program The program extracts and downscales PRISM data and historical and future climate data simulated by various global circulation models (Daly et al, 2002; Mitchell and Jones, 2005), and methodologies are described in Wang et al. (2012). Climate WNA model outputs were analyzed for comparison to the Blue River climate data. Climate WNA model outputs were extracted for the geographic location corresponding to an approximate average elevation (2000 m) and geographic position (52.758° N, -119.681° W) of the Canoe Basin. Table 2.1. Goodness-of-fit statistics for infilling monthly Blue River climate data for the entire record (n=279) and for the record divided into summer (May-September) and winter (October-April) seasons (n = 109). Climate Variable Mean Temperature Max Temperature Min Temperature Precipitation Entire Record R2 Standard Error 0.99 0.36 °C 0.99 0.49 °C 0.46 °C 0.99 0.88 14.4 mm R2 0.99 0.99 0.99 0.91 10 W inter Standard Error 0.34 °C 0.45 °C 0.42 °C 16.8 mm R2 0.99 0.99 0.98 0.92 Summer Standard Error 0.27 °C 0.45 °C 0.37 °C 12.3 mm o (M O o O o o -15 -10 ■5 0 5 10 15 20 0 -10 BRN Mean Monthly Temperature (*C) 10 20 BRN Maximum MotfWy Tem perature (*C) o ur> fx o o o (N o o o o o E 5 2 2 o m o 20 -10 0 0 10 50 100 150 200 250 BRN Total Monthly Precipitation (mm) BRN Minimum Monthly Temperature f C ) Figure 2.3. Scatter plots showing the overlapping period between the BRA and BRN climate stations (n=279). Regression lines are also shown (solid lines). Data Analysis The Mann-Kendall rank-based test (Mann, 1945 and Kendall, 1975) and the KendallTheil robust line (Theil, 1950) were used to detect trends in the climate data (mean/max/min temperature and total precipitation). The Mann-Kendall test is commonly used in trend detection of hydrometeorological parameters (e.g. Bum et al., 2004; Rood et al., 2005) 11 because it is non-parametric and so an assumption of normality is not required. A rank-based test considers only the relative rank of each data point and is therefore more resistant to outliers in the dataset than least squares regression. The Mann-Kendall test assesses all possible pairs of data points to determine whether successive values are lower or higher. Kendall’s tau (x) correlation coefficient is calculated as follows: S = P —M (2.1) T=^75 < 2 ' 2 ) where S is the test statistic, P is the number of concordant data pairs (number of times the y values increase as the x values increase), M is the number of discordant pairs (number of times the y values decrease as the x values increase), and n is the number o f data points in the time series. The Mann-Kendall test is not recommended for use on datasets where seasonality or autocorrelation exists (Bum et al., 2004). Seasonality was avoided in the datasets by using one mean value per year, for each parameter of interest. To avoid the potential impacts of autocorrelation on trend detection, this study used a modified version of the Mann Kendall test that automatically corrects the power of the trend test by adjusting the number of independent samples in the time series (Yue and Wang, 2004). The modified Mann-Kendall used the 95% significance level to test for trends. Trend magnitude was calculated as the slope of the Kendall-Theil robust line which is more resistant to outliers than linear regression. Linear trends were examined over the entire record and over shorter time periods within the record aligning with photographic records, glacier volume change calculations and known periods of climatic fluctuation. For each period, data analysis was completed on an 12 annual (hydrological year from 1 Oct to 30 Sep), winter and summer (Oct-Apr and May-Sep) and seasonal (Jan-Mar, Apr-Jun, Jul-Sep, Oct-Dec) basis. Blue River climate data were tested using the non-parametric Wilcoxon rank-sum test (Wilcoxon, 1945) to determine whether temperature or precipitation were statistically different in the years 1948-1984 and 1985-2005, and the years 1948-1977 and 1978-2005. 2.3.2 Digital Photogrammetry Data Aerial triangulation (AT) scans from 2005 and aerial photographs from 1948, 1955, 1970, 1985 and 1996 (Table 2.2) were obtained from the British Columbia Government and the Canadian National Air Photo Library. Photos from 1948 provided almost complete glacier coverage, lacking only Mt. Kimmel and a slice o f the North Canoe glacier accumulation area (Figure 2.4), and had good contrast with no fresh snow. Photos from 1955 were missing the eastern part of the watershed (Figure 2.5), had fresh snow in many areas and poor contrast particularly in the accumulation areas. Photos from 1970 provided complete coverage (Figure 2.6); however, although contrast was good at lower elevations, contrast was poor at higher elevations and the glacier surface could not be distinguished. Photos from 1985 covered all glaciers (Figure 2.7) with good contrast and no fresh snow. Photos from 1996 were missing the eastern part of the watershed (Figure 2.8) and had extensive fresh snow cover. AT scans from 2005 provided almost complete coverage of the glaciers in the watershed, with the exception of Mt. Kimmel in the southeast (Figure 2.9) and had good contrast with no fresh snow. 13 Table 2.2. Aerial photographs and aerial triangulation scans used in project. Acquisition date(s) Roll Number of photos used Scale Percent of glaciers covered 1948 1948-09-10 1948-08-17 1948-08-17 1955 1955-08-30 1970 1970-08-21 1985 1985-07-09 1996 1996-07-13 2005 2005-08-05 A11668 A l 1610 A11565 A 14930 15BC5394 15BC85071 15BCB96011 30BCC5108 30BCC5107 48 12 15 13 16 33 1:35,000 1:70,000 1:80,000 1:60,000 1:40,000 1:20,000 96 % 93 % 100% 100 % 24% 99 % Figure 2.4. Glaciers covered by the 1948 photos overlain on a Landsat 5-4-3 composite (4 September, 1995). 14 Figure 2.5. Glaciers covered by the 1955 photos overlain on a Landsat 5-4-3 composite (4 September, 1995). ▲ Watershed Outlet • Checkpoint Patches x Ground Control Poin Figure 2.6. Glaciers covered by the 1970 photos overlain on a Landsat 5-4-3 composite (4 September, 1995). Figure 2.7. Glaciers covered by the 1985 photos overlain on a Landsat 5-4-3 composite (4 September, 1995). I Figure 2.8. Glaciers covered by the 1996 photos overlain on a Landsat 5-4-3 composite (4 September, 1995). Additional area was covered in the 1996 photos; however, data collection was limited to the area shown due to high snow cover. 16 Figure 2.9. Glaciers covered by the 2005 AT scans overlain on a Landsat 5-4-3 composite (4 September, 1995). Stereo Model Development The VR Mapping software suite was used for digital photogrammetry. Aerial Triangulation (AT) scans from 2005 were imported into VR AirTrig, and oriented into photo strips. The 2005 AT scans were the basis of external orientation for the air photos from all previous years. Aerial photographs from 1948, 1955, 1970, 1985 and 1996 were imported into a VR AirTrig file for each year, and oriented into photo strips according to the flight path. A coordinate file was created for use with all the models. Interior orientation was completed in VR AirTrig using six tie points between adjacent photographs and several pass points between photograph strips, to align the photo pairs and strips relative to each other. External orientation was completed in relation to the 2005 AT scans. Ground control points (GCPs) were collected in the 2005 AT scans and matched to points in earlier photos, at locations where the elevation should not change from year to year (Kaab and Vollmer, 2000; Schiefer and Gilbert, 2007). These GCPs represented stable features on bedrock and in areas free of snow/ice, vegetation and signs of erosion. The GCPs were collected to provide distribution across space and elevation within and close to the watershed (Figures 2.4-2.8). Due to the small portion of glaciers covered in the 1996 photos, only 11 GCPs were collected; the other years had between 19 and 24 GCPs. Air photo footprints, snow/ice cover and contrast/shadow in the photos were different for each year; however, an attempt was made to use the same GCPs wherever possible. After GCP collection, model bundle adjustment was run in VR AirTrig to determine the root mean square error (RMSE) in the tie points, pass points and GCPs. Points with RMSE higher than 1.0 m were re-sampled to ensure correct placement. Upon successful internal and external orientation of air photos for each year, the stereo models were built using VR Orientation. VR One was used to view photo pairs in stereo and collect glacier extents and elevation mass points. The elevation bias between 2005 AT scans and older photographs was calculated using a series of checkpoint patches. At each checkpoint patch, surface elevations were sampled on a 5 x 5 m grid. Checkpoint patches were placed in areas that should not change in elevation over time (Schiefer and Gilbert, 2007); therefore, any difference in elevation between years would indicate a bias between the models. Ten checkpoint patches were used in the 1996 photography because of the small portion of glacier area covered in the photos. The other years used 17 or 18 checkpoint patches. For each year, the differences in checkpoint elevations were plotted against elevation, easting and northing coordinates, and linear models were developed to determine the presence and statistical significance of any systematic biases in the data. Using the coefficients from the linear models, significant trends were removed from each year of 18 photography (Schiefer et al., 2007). In cases where more than one significant trend was identified, the trend with the largest correlation coefficient (R2) value was removed. Only one trend was removed from each dataset to avoid over-manipulating the data. Data extraction Three-dimensional glacier extents were mapped from stereo photo pairs in the VR One software. Where the glacier termini were partially obscured by debris cover, models from previous and subsequent years were compared; if a change in elevation was apparent, that portion of the glacier was mapped as debris-covered. It is also common for some glacier areas to be obscured by shadows and steep cliffs (Fox and Nuttall, 1997). In these cases, adjacent photo pairs were compared to ensure accurate glacier mapping. Elevation mass points were collected over the glacier surface for each year on a 100 m grid. Full mass point coverage was collected from the 2005 and 1985 photos, and almost complete coverage was obtained for the 1948 photos. In the northwest area of the 1948 photos, 110 mass points in an elevation band from 2500-2800 m could not be collected because the photography was unavailable (Figure 2.4). An attempt was made to develop a linear model between elevation change from 1948 to 1985 and the 1985 elevation to estimate the elevation change of the missing data; however, a strong correlation did not exist. Therefore, the missing values were replaced using the average difference in elevation (19851948) for the elevation band, which was 1.03 m. Because 1948, 1985 and 2005 had complete coverage, these years were used in glacier volume change calculations. The 1955 and 1970 photos have incomplete coverage due to poor contrast in the accumulation area; as a result, these photos were not used in volume change calculations. Mass point coverage in the 1955 photos was limited to the glacier snouts, with few mass points collected in the accumulation 19 area. In the 1970 photos, approximately half of the accumulation area has coverage with mass points. Abundant fresh snow limited mass point collection in the 1996 photos. Linear biases identified in the checkpoint analyses were removed from each mass point dataset. Glacier morphology and area/volume change Glacier extents were imported into ArcMap for analysis. Flowsheds of the six glaciers with well-defined snouts were delineated and the mass point data were used to derive glacier hypsometry for each flowshed in 1948 and 2005. The change in the terminus position was mapped and morphological properties were analyzed for each glacier. Glacier terminus change in this study is defined as the greatest distance between glacier termini in successive years of photography and is measured perpendicular to the glacier terminus and parallel to the glacier flowline (Andreassen et al., 2005). For all years except 1996, total glacier area was calculated by summing the areas of all glacier polygons. The error term associated with glacier area for each year was calculated by creating a buffer around the glacier polygon, similar to methods used in Granshaw and Fountain (2006) and Bolch et al. (2010), with a width equal to the RMSE from the X and Y positional errors calculated in VR for stereo model batch adjustment. The area of the buffer is the error term for the glacier area. The glacier areas were standardized to include only the glacier areas common in all years to ensure a change in glacier area was not simply a result of missing data (Table 2.2). To determine the change in glacier area between subsequent years, the older extent was subtracted from the newer extent. The RMSE associated with the change in glacier area between years was calculated as: RMSE = V(bufferA)2 + (bufferB )2 20 (2.3) where bufferA and bufferB represent the buffer areas (m2) in the two years of photography. No area was calculated for 1996; photography in 1996 was only used to map the North Canoe terminus. Bias-corrected mass points were used to analyze elevation change and calculate volume change. Mass point surface elevations from older years were subtracted from more recent years to calculate change in glacier surface elevation through time. As complete mass point coverage was limited to 1948, 1985 and 2005, volume change calculations were based on these years o f data. Glacier volume change was calculated by multiplying the average elevation difference over the glacierized area by the largest glacier cover of the two years. The glacier area for volume change calculations was different than that for area change calculations, as each calculation required the glacier area common to all years of interest (i.e. glacier area change was calculated between 1948, 1955, 1970, 1985 and 2005, whereas volume change was calculated between 1948, 1985 and 2005). The multiplicative error term for glacier volume change was calculated as a combination of the errors associated with measuring the glacier area change (buffer areas) and measuring the difference in glacier surface elevation (mean elevation difference after bias removal). Volume error was calculated as: Volume E rror = Volume Change x / VAC/ + ( —-—) V Z m ean-' (2.4) where AC is the area change between the two years (m ), AE is the error term for area change (m2), zmean is the mean elevation difference between the two years (m) and z is the bias-corrected checkpoint standard deviation (m). Glacier volume change was converted to snow water equivalent with an average density of 690 kg m '3, based on densities of 900 kg m '3 and 550 kg m '3 in the ablation and accumulation areas, respectively, and an Accumulation Area Ratio of 0.6 (Schiefer et al., 21 2007). Although calculation o f volume change was limited to 1948-1985-2005, mass points collected from the 1955 and 1970 photos provided an indication of thickening or thinning of glacier snouts. 2.4 Results 2.4.1 Climatic variability (1947-2006) A complete record of trends in the climate data is provided in Appendix B. For Blue River, mean annual and winter (Oct-Apr) temperatures show similar increasing trends, while mean summer (May-Sep) temperature did not change significantly through time (Table 2.3, Figures 2.10-2.13). Mean annual temperature significantly increased at a rate of 0.024 °C yr-1 over the period of record (Figure 2.10), which equates to a total increase of 1.44 °C over 60 years. The increase in mean annual temperature is primarily due to increases in mean monthly temperatures in January, February and March (Appendix B, Table B .l). Much yearto-year variability is apparent throughout the record, with 1958 and 1998 the warmest years on record and 1950 and 1956 the coldest years on record. Minimum annual temperature increased at a significant rate of 0.034 °C yr-1. Maximum (minimum) winter temperature increased at a significant rate of 0.034 °C yr-1 (0.050 °C yr-1), which equates to a total increase of 2.0 °C (3.0 °C) over 60 years. No significant trends were observed in maximum and minimum summer temperatures. Total summer precipitation at Blue River for the period 1947-2006 increased at a rate of 2.664 mm yr-1 (Table 2.3, Figure 2.13). The highest annual precipitation was recorded in 1999. 22 Table 2.3. Summary o f trends in the Blue River mean, maximum and minimum temperature (°C yr_1) and total precipitation (mm y r -1) time series. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold. 1947-2006 1947-1985 1985-2006 Tau________ KT Slope________Tau______ KT Slope________ Tau_______ KT Slope Mean Temperature Annual Winter Summer 0.285 0.294 0.101 0.024 0.036 0.005 0.148 0.126 0.001 0.020 0.027 0.000 0.160 0.212 -0.052 0.016 0.053 -0.013 -0.016 0.045 -0.016 -0.003 0.013 -0.003 0.026 -0.013 0.026 0.017 -0.009 0.017 0.269 0.207 0.358 0.051 0.057 0.048 0.195 0.229 -0.035 0.039 0.060 -0.007 0.134 -0.123 0.312 2.162 -1.578 0.229 0.169 0.195 6.146 2.800 5.542 Maximum Temperature Annual Winter Summer -0.047 0.250 -0.047 -0.010 0.034 -0.010 Minimum Temperature Annual Winter Summer 0.332 0.322 0.169 0.034 0.050 0.011 Total Precipitation Annual Winter Summer 0.156 -0.060 0.296 1.793 -0.522 2.664 23 4.169 Mean Annual Temperature CC) <75 o Minimum Summer Temperature (*C) o p in p 1950 Figure 2.12. Minimum, winter and summer temperature at Blue River from 1947 to 2006. The solid line is the 5-year running mean. Dashed lines show significant trends (p < 0.05; Kendall-Theil robust line). 26 Total Annual Precipitation (mm) 8 8 O 8 a> 8 CO o o © <£> 1960 1970 1960 1960 1980 1990 2000 1970 1990 2000 1970 1990 2000 Total Winter Precipitation (mm) 8 s<0 8 o o Total Summer Precipitation (mm) <0 o Q V g 1960 Figure 2.13. Mean annual, winter and summer precipitation at Blue River from 1947 to 2006. The solid line is the 5-year running mean. Dashed line shows significant trend (p < 0.05; Kendall-Theil robust line). 27 A comparison of trends between the Blue River data (Table 2.3) and Climate WNA (average Canoe Basin) model output (Table 2.4) indicates that the recorded mean, maximum and minimum temperature trends are similar; however, not all trends are significant in both climate records. Table 2.4. Summary o f trends in Climate WNA (average Canoe Basin) mean, maximum and minimum temperature (°Cyr~ ) and total precipitation (mm y r - ') time series. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold. 1947-2006 1947-1985 1985-2006 Tau________ KT slope________Tau______ KT slope_________ Tau________ KT slope Mean Temperature Annual W inter Summer 0.327 0.288 0.189 0.026 0.036 0.011 0.130 0.109 0.032 0.015 0.023 0.002 0.177 0.199 0.004 0.029 0.044 0.010 -0.013 0.007 -0.013 -0.024 0.038 -0.024 -0.011 0.013 -0.011 0.026 0.061 0.026 0.008 -0.013 0.008 0.070 0.070 0.025 0.058 0.058 0.250 0.040 0.040 0.216 0.216 -0.090 0.133 0.133 -0.050 2.511 1.500 0.875 0.224 0.151 0.235 -0.273 -0.260 -0.095 -12.500 -7.714 -1.857 Maximum Temperature Annual W inter Summer -0.080 0.065 -0.080 Minimum Temperature Annual Winter Summer 0.224 0.224 0.259 0.030 Total Precipitation Annual W inter Summer 0.121 0.106 0.110 6.273 3.050 2.750 Blue River climate variables were compared during periods of calculated glacier change (1947-1984 and 1985-2005) and during periods before/after known phase changes of the Pacific Decadal Oscillation (1947-1977 and 1977-2005; Mantua, et al., 1997) to assess whether significant differences exist in temperature and precipitation between the two periods. Boxplots graphically show the distributions between the two periods (Figures 2.14 and 2.15) and the results of the Wilcoxon rank-sum test identify statistically significant differences between periods (Table 2.5). Based on results of the Wilcoxon rank-sum test, 28 winter temperatures were significantly higher from 1985 to 2005 than from 1947 to 1984. Summer precipitation and winter temperature were also significantly higher post-1977. Table 2.5. W ilcoxon rank-sum test p-values for tests comparing Blue River climate variables during periods of calculated glacier change and known phase shifts o f the PDO. Significant (p < 0.05) differences are shown in bold. Periods of glacier change PDO phase shifts _______________ (1948-1985 compared to 1985-2005) (1947-1977 compared to 1978-2005) ________________Temperature________ Precipitation_________ Temperature________ Precipitation Summer 0.13 0.070 0.16 0.00063 W inter 0.0071 0.59 0.011 0.19 © Io 8«■> o «o o o © a a to 8 1947-1964 1986-2005 1965-2005 Figure 2.14. Boxplots showing the difference in Blue River mean temperature and precipitation from 1947 to 1984 and from 1985 to 2005, periods of calculated glacier volume change. Boxplots show the median and the 25* and 75* percentiles and the notches indicate the 95% confidence limits o f the median. W hiskers extend to 1.5 times the interquartile range from the box and dots represent outliers. 29 tf) m 8 <0 O o v> tn o oo £ tD <30 O 8 rt U) 1947-1977 1978-2005 1947-1977 1978-2005 1947-1977 1978-2005 o o o 16 IE o CM North Canoe South Canoe Penny CC5/6 Zilmer West Zilmer East Ui o o ID O O O 06 0.2 00 08 1.0 Fraction area o o o CO (b) 2005 oo o CM North Canoe South Canoe Penny CC5/6 Zilmer West Zilmer East Ui o o ID O o o 00 0.2 0.4 06 0.8 1.0 Fraction area Figure 2.17. Glacier hypsometry mapped from (a) 1948 and (b) 2005 air photos. 33 The termini of glaciers in the Canoe Basin responded differently over the period of record: four glaciers underwent a minor advance from 1955 to 1985 and two glaciers retreated continuously from 1948 to 2005 (Table 2.10). All glaciers experienced net retreat from 1948 to 2005. North Canoe and South Canoe glaciers retreated and thinned in the ablation area from 1948 to 1955, advanced and thickened from 1955 to 1970, advanced and thinned from 1970 to 1985 (Figure 2.18) and retreated and thinned from 1985 to 2005 (Figures 2.19 and 2.20). In contrast, Penny and CC5/CC6 glaciers advanced and thinned from 1948 to 1955, advanced and thickened from 1955 to 1970, advanced with no thickness change from 1970 to 1985 and retreated and thinned from 1985 to 2005 (Figures 2.21 and 2.20). Zilmer West and East glaciers retreated and thinned continuously over the period of record, with no evidence of thickening or terminal advance (Figure 2.22). Table 2.10. Change in glacier terminus position and elevation in the ablation area of the glacier over the period o f study. Negative values indicate retreat; positive values indicate advance. Change in glacier terminus position (m) and ablation area elevation Glacier North Canoe South Canoe CCS and CC6 Zilmer West Zilmer East Penny 1948-1955 1955-1970 1970-1985 1985-2005 - 85 m thinned + 115 m thickened + 120 m thickened + 20 m thickened - 830 m thinned + 210 m thinned + 210 m thinned + 90 m no change - 2 10m thinned - 125 m thinned + 200 m thickened + 90 m no data -4 5 0 m thinned -7 0 0 m thinned -2 3 0 m thinned -3 0 0 m thinned -2 5 0 m thinned - 350 m thinned - 580 m thinned -2 1 0 m thinned -2 1 0 m thinned - 200 m thinned -7 5 m insufficient data 34 - 60 m thinned Eltvatlon Chwige (m) 005-1970 ■ -45 t o -25 H -25 t o -10 -10 to *3 *3 to 3 m 3*010 I 10 to 25 ■ 25 to 45 Figure 2.18. Elevation change of North Canoe glacier from 1970 to 1985. CD 1 0 4 6 4 - I j 1QSS I 1 1070 C^J-toes I Iiooe I 1 2005 1.000 Meters Figure 2.19. North Canoe glacier terminus change from 1948 to 2005. 35 Figure 2.20. South Canoe and CC5/6 glacier terminus change from 1948 to 2005. j | 1046 r ~ i ioss CZ3,e70 CZI 1085 1 ,0 0 0 M e te rs l f 2005 Figure 2.21. Penny glacier terminus change from 1948 to 2005. 36 250 500 750 Figure 2.22. Zilmer W est and East glacier terminus change from 1948 to 2005. As a result of the small area changes from 1948 to 1955 and from 1955 to 1985, the total area change from 1948 to 1985 is within the bounds of the error term (Table 2.11). Glacier area decreased by -3 .8 ± 0.44 km2 (-6.8 ± 0.79%) from 1948 to 2005, with the majority of retreat between 1985 and 2005. Glacier volume change from 1948 to 1985 was small and within the bounds of the error term (Table 2.12). The error term is larger for the 1948-1985 calculation than the 19852005 calculation because the 1948 and 1985 photos are all referenced to the 2005 AT scans, which have relatively small error terms. Glacier volume change from 1985 to 2005 was -0.32 + 0.06 km3 (Table 2.12), or -0.22 ± 0.04 km3 water equivalent. Glacier volume change cannot be expressed as a percentage because the total volume o f glacier ice is unknown. From 1948 to 2005, there was minimal overall elevation change in the accumulation area of any glaciers (e.g. North Canoe glacier, Figure 2.23), although some areas showed changes of + 10 m. 37 Table 2.11. Glacier area change from 1948 to 2005 in the Canoe Basin. Area is based on the fraction of glaciers that are covered in all years o f photography. Year Glacier Area (km2) 2005 5 2 .1 + 0 .0 2 1985 1970 55.6 ± 0 .1 6 55.7 ± 0 .2 6 1955 1948 56.5 ± 0.25 55.9 ± 0 .1 9 Glacier Area Change (%) (km2) Cumulative Glacier Area Change (km2) (%) -3.8 ± 0 .4 4 -6.8 ± 0.79 -3.5 ± 0.16 -0.07 ± 0.30 -6.3 ± 0.29 -0.12 ± 0 .5 4 -0.32 ± 0.44 -0.57 ± 0.79 -0.78 ± 0.36 -1.4 ± 0.65 -0.25 ±0.41 -0.45 ± 0.733 0.52 ±0.31 0.94 ± 0.55 0.52 ±0.31 0.94 ± 0.55 - - - Table 2.12. Glacier volume change from 1948 to 2005 in the Canoe Basin. Years Glacier Volume Change (km3) 1948-2005 1985-2005 -0.23 ±0.14 -0.32 ± 0.06 1948-1985 0.089 ±0.13 E lev atio n C hang* (m) 2009-1948 ■ -50 t o -25 M -25 t o -10 -10 to -3 -3 to 3 £ 3 to 10 ■ 10 to 25 ■ 25 to 40 Figure 2.23. Difference in North Canoe glacier surface elevation from 1948 to 2005. 38 2.5 Discussion 2.5.1 Climatic variability Mean annual temperature at Blue River increased by 1.44°C from 1947 to 2006, at a rate of 0.24°C per decade (equivalent to 2.4°C per century). This is a higher rate of warming than that observed in the central Canadian Rocky Mountains (1.43°C per century; Luckman, 1998), the Canadian portion of the Columbia River Basin (0.7°C - 1.7°C from 1901 to 2004; Murdock and Werner, 2011), the Pacific Northwest (0.7°C - 0.9°C per century; Mote, 2003) and the global warming trend (0.6°C + 0.2°C per century; Houghton et al., 2001). The higher rate of warming at the Blue River stations could be attributed to the period of record encompassing only the recent six decades, as it has been observed that the rate of warming increased in the second half of the twentieth century. The rates of change at Blue River are based only on the combined Blue River stations, and thus may not be representative of a regional rate of change. Globally, the rate of warming from 1956 to 2005 was almost twice the rate of warming from 1906 to 2005, and the period from 1995 to 2006 had eleven of the twelve warmest years on record globally since 1850 (IPCC, 2007). The largest warming trends at Blue River were observed during January, February and March. Similarly, winter was the season with the largest, positive trend in the Canadian Rocky Mountains (Luckman, 1998), in the Canadian portion o f the Columbia River Basin (Murdock and Werner, 2011) and throughout the Pacific Northwest (Mote, 2003). The years 1958 and 1998 were the warmest on record at Blue River and 1998 and 2005 were the warmest years globally (Solomon et al., 2007). Extreme temperatures of 1958 and 1998 occurred during large El Nino events. While the trend in annual precipitation was not significant, total summer precipitation at Blue River increased by 160 mm over the 60 year period which is equivalent to a 3% 39 increase per decade. Precipitation is typically more variable than temperature across the region (e.g. Luckman, 1998) and this variability is also evident in the different precipitation trends identified in the Blue River and Climate WNA model outputs. In the Canadian portion of the Columbia Basin, total annual precipitation increased by as much as 4% per decade (Murdock and Werner, 2011) and some stations in southeast British Columbia recorded precipitation increases greater than 60% over the last century (Mote, 2003). The greatest increase in precipitation at Blue River was during the months of April, May and June. Similarly, precipitation increased most during the spring (3-7%) and summer (1-4%) in the Canadian portion of the Columbia Basin (Murdock and Werner, 2011) and from April to July in the Pacific Northwest (Mote, 2003). On average, there were cooler winters and drier summers from 1948 to 1977 compared to 1978-2005, with no significant difference in summer temperature or winter precipitation. These results coincide with a shift of the Pacific Decadal Oscillation (PDO) from negative to positive in 1977 (Mantua, et al., 1997). The positive (negative) PDO phase generally results in above (below) normal temperatures and below (above) normal precipitation in British Columbia (Fleming and Whitfield, 2010) and the Pacific Northwest (Mantua et al., 1997). In B.C. specifically, the positive PDO phase is generally associated with warm, dry winters (Stahl et al., 2006). Inter-annual to inter-decadal climate variability that is linked to ENSO or the PDO complicate trend detection in hydrometeorology, as the upward leg of a cycle may be misinterpreted as a trend, particularly over time scales too short to discern the cycle (Woo et al., 2006). It is thus possible that the climatic trends observed at the Blue River stations are due to natural cyclical patterns, rather than longer term (century or millennium) changes in climate. It is also possible that longer term climatic trends could be masked by the quasi40 periodic patterns attributed to ENSO or the PDO. Longer climate records are necessary to differentiate between long-term trends and natural cycles in climate data. 2.5.2 Glacier terminus change Glaciers in the Canoe Basin retreated from 1948 to 1955, slowed their recession with some termini advancing from 1955-1985, and retreated again from 1985 to 2005. Luckman et al. (1987) noted that 80% of glaciers in the Premier Range retreated from 1948 to 1955 and over 50% of glaciers advanced from 1955 to 1985. According to Luckman et al. (1987), morphological evidence suggests that the maximum downvalley terminus position from 1955 to 1985 occurred a few years prior to 1985, followed by glacial retreat. A dendrochronological age for this advance at Penny Glacier, for example, is AD 1980-1981 (Luckman et al., 1987). Glacier behavior in the Canoe Basin mirrors records of glacier change in western North America over the 20th century, with glacier retreat from 1900 to 1950, slowed recession or advance from 1950 to 1980 and rapid retreat since the 1980’s (e.g. Moore et al., 2009). The retreat-advance-retreat pattern of the 20th century is partly attributed to largescale changes in climate arising from a shift in the PDO (e.g. Koch et al., 2009). In the Pacific Northwest, the negative PDO phase typically results in below normal temperatures and above normal precipitation (Mantua et al., 1997). In B.C., the warm, dry winters during the positive phase of the PDO (Stahl et al., 2006) favor negative glacier mass balance and vice versa. Luckman et al. (1987) speculate that the glacier advance culminating in the midto late-1970’s was a direct response to colder than normal summers from 1954 to 1968 and higher than normal winter precipitation from 1951 to 1976. The Blue River climate data supports a period of cool, wet conditions from 1947 to the late 1970’s. 41 Glaciers in many mountain ranges in western North America advanced from the 1950’s to the 1980’s. Many glaciers in Garibaldi Park in the southern B.C. Coast Mountains expanded in the 1960’s and 1970’s (Koch et al., 2009). Boundary Glacier in the Canadian Rocky Mountains advanced 100 m from 1966 to 1970 (Gardner and Jones, 1985). Within the Columbia River Basin, Illecillewaet Glacier, for example, advanced more than 90 m from 1972 to 1984 (Champoux and Ommanney, 1986). In the northwest United States, several glaciers advanced at Mount Hood, Oregon, from 1947 to the late 1970’s (Lillquist and Walker, 2006), at Mount Baker, Washington, from the late 1970’s to the early 1980’s (Harper, 1993) and in Glacier National Park, Montana, from 1966 to 1979 (Carrara, 1987). Results of glacier terminus mapping in the current study and those from Luckman et al. (1987) are similar, but some discrepancies exist for North Canoe and Zilmer glaciers. For example, both studies note that the North Canoe glacier advanced from 1955 to 1985, but Luckman et al. (1987) indicate that the terminus advanced +25 m from 1948 to 1955. In the Canoe Basin study, however, the North Canoe glacier retreated -8 5 m over the same period of time. The Canoe Basin study also found that Zilmer West and East glaciers continuously retreated over the period of record. In contrast, Luckman et al. (1987) report that Zilmer glacier did not change from 1973 to 1976 and subsequently advanced +70 m from 1976 to 1985. Luckman et al. (1987) did not identify any recent terminal moraines at Zilmer glacier; therefore, they did not determine a date for the maximum extent of recent terminal advance. The difference in terminal change mapped in this thesis and in Luckman et al. (1987) is likely related to the methods used in each study. Luckman et al. (1987) used aerial photography (1948-1976), and field observations and oblique photography (1984-1985) to map glacier change. In contrast, the current study mapped glacier extents using digital photogrammetry, which can increase mapping accuracy. 42 A glacier’s sensitivity to climatic change is influenced by glacier physical properties including slope, aspect, elevation and hypsometry (Arendt et al., 2002; Granshaw and Fountain, 2006; DeBeer and Sharp, 2007; Jiskoot et al., 2009). As a result, some variability in glacier response to climate change should be expected in a region (e.g. DeBeer and Sharp, 2007). Even if the large-scale patterns in climate variables are consistent across a region, the micro-climate over each individual glacier could differ because of glacier elevation and hypsometry (e.g. Oerlemans et al., 1998; Arendt et al., 2009). While some glaciers respond rapidly to changes in climate, others take much longer (Raper and Brathwaite, 2009). Therefore, glaciers may still be responding to climatic influences that occurred before a given study period. Complicating the issue, the response time of a glacier to climate change varies among glaciers. Although factors influencing glacier response times are not completely understood, Pel to and Hedlund (2001) found that glacier physical properties (slope, degree of crevassing and glacier terminus velocity) are important variables. Glaciers in different parts of the Canoe Basin responded differently to climate change. Zilmer West and East glaciers retreated continuously over the period o f record, while the other glaciers in the Canoe Basin advanced from 1955 to 1985. The varied response of these glaciers may relate to variability in regional precipitation, as well as differences in the sensitivities and response times of the glaciers. Zilmer W est and East glaciers have the lowest median elevations and elevation ranges of glaciers in the Canoe Basin (Table 2.9). The small elevation range could indicate that a greater proportion of the glacier area is near the ELA. Thus, a small increase in the ELA due to a rise in temperature could result in a large increase in the area of the glacier below the ELA (Lillquist and Walker, 2006) and enhanced ablation. 43 Studies have shown that sensitivity to temperature is greater when a glacier has more area at low elevations (Larsen et al., 2007). Zilmer W est and East glaciers have the highest toe elevations of glaciers in the Canoe Basin (Table 2.9). If a large proportion of a glacier is at low elevation, even a small increase in the ELA could result in substantially more negative mass balance and terminal retreat. Typically, high elevations are colder than low elevations, which means that a greater proportion of precipitation falls as snow. Thus, during warm or wet years, glaciers at a higher elevation would have a more positive mass balance than those at a lower elevation (Harper, 1993). The high toe elevations of Zilmer W est and East glaciers could help to explain the sustained terminus retreat observed for these glaciers. The shape of a glacier and the relative size of the accumulation and ablation areas can also result in different responses to climatic variability. For example, in a glacier that is shaped like a triangle, with a large accumulation area and a long narrow ablation area, a small mass input in the accumulation area will transfer into a large mass input in the ablation area due to a funneling effect. This means that the snout might gain a large amount of mass relative to a glacier that is shaped like a rectangle (Furbish and Andrews, 1984). To remain in equilibrium with a large mass input at higher elevations, a glacier of this shape might also have a terminus with a lower elevation. 2.5.3 Glacier area and volume change There was a small increase in total glacier area in the Canoe Basin from 1948 to 1955 despite terminal retreat and thinning in the ablation zones of the glaciers. Although most glacier termini advanced and thickened from 1955 to 1970, the significant retreat of Zilmer W est and East glaciers resulted in a net decrease in glacier area in the watershed. 44 Glacier area change from 1970 to 1985 was within the bounds of the error term. There was a substantial decrease in glacier area from 1985 to 2005. From 1970 to 1985, the snouts of North Canoe and South Canoe glaciers advanced but thinned. Arendt et al. (2002) studied 67 glaciers in Alaska, o f which they found that 10% of the glaciers advanced during times when these glaciers thinned. They also found that there was a low degree of correlation between a glacier’s length change and the rate of thickness change. Glacier mass readjustment indicates that flow dynamics play a role in glacier change over time scales around 10-40 years (Arendt et al., 2002). Because of this, it is problematic to infer glacier volume change or mass balance based solely on glacier area or length change (Arendt et al., 2002). From 1948 to 1985, there was a small increase in glacier volume (within the bounds of the error term) despite an overall decrease in glacier area. Glacier volume (and area) in the basin substantially decreased from 1985 to 2005. The significant reduction in glacier volume from 1985 to 2005 accords with the results of the following glacier change studies in western North America. Schiefer et al. (2007), DeBeer and Shaip (2007), Beedle et al. (2009) and Bolch et al. (2010) discuss area and volume change for glaciers in the vicinity of the Canoe Basin. Schiefer et al. (2007) calculated glacier volume change in the Columbia Mountains (including the Canoe Basin) from 1985 to 1999 using digital terrain models and Shuttle Radar Topography Mission (SRTM) data. They calculated that glacier volume decreased at a rate of -1.23 + 0.31 km 3yr-1 over a total glacier area of 2,326 km2, which is a thinning rate of -0.53 ± 0.13 myr-1. In comparison, glaciers in the Canoe Basin thinned at a rate of -0.29 ± 0.07 myr~! from 1985 to 2005. DeBeer and Sharp (2007) assessed glacier area and volume change (derived using volume-area scaling) in the Columbia Mountains from 1951/2 to 2001 using aerial 45 photographs and Landsat ETM+ (Enhanced Thematic Mapper +) imagery. They found that glacier area (volume) in the Columbia Mountains study area decreased by -5 % (—1.1 ± 0.2 km3). In comparison, Canoe Basin glacier area (volume) decreased by - 6.8 ± 0.79% (-0.23 ± 0 . 1 km3) from 1948 to 2005. Although the volume calculations cannot directly be compared between studies, the percent change in glacier area in the Canoe Basin was similar to that in the Columbia Mountains over similar time periods. Beedle et al. (2009) reconstructed the terminus position of Castle Creek Glacier, located in the Premier Range, by mapping a series of annual push moraines. The study found that the glacier retreated continuously from 1959 to 2007. The study also found that lower rates of retreat were recorded in the mid- to late-1970’s, in comparison to the rest of the record. Bolch et al. (2010) calculated glacier area change in western Canada from 1985 to 2005 using Landsat Thematic Mapper imagery and aerial photography. The Canoe Basin lies within the Southern Interior region, as defined in their study. In the Southern Interior, glacier area change from 1985 to 2005 was -1 5 .2 ± 4.4% (Bolch et al., 2010). Similarly, glacier area decreased by approximately -15% in the Columbia Basin from 1985 to 2005 (Jost et al., 2012). In comparison, glacier area decreased by -6.3 ± 0.29% in the Canoe Basin over the same time period. Observations of glacier volume change in B.C. and western North America indicate an accelerated rate of glacier decline in the late 1990’s and early 2000’s, similar to the accelerated glacier decline in the Canoe Basin from 1985 to 2005. Bolch et al. (2010) found that the percent rate of decline in glacier area was higher from 2000 to 2005 than from 1985 to 2000 in all regions where data allowed the comparison. In the Southern Interior, specifically, the rate of glacier area change increased from -0.27% yr-1 in the period 198546 2000 to -0.44% yr~‘ in the period 2000-2005 (Bolch et al., 2010). Arendt et al. (2009) found that the retreat rates for 76% of glaciers studied in Alaska and northwestern Canada increased from the m id-1990’s to the mid-2000’s. Similarly, the rate of glacier area change in the Canoe Basin increased from -0.005% yr-1 in the period 1970-1985 to -0.18% yr-1 in the period 1985-2005. 2.5.4 Sources of error and study limitations Climate record Climate varies substantially with small changes in geographic position and elevation, particularly in mountain environments. Precipitation, in particular, is dramatically different on the windward and leeward sides of mountain ranges, and elevation influences the phase and amount of precipitation. Topography in mountainous terrain makes it difficult to extrapolate a point measurement of temperature and especially precipitation across an entire watershed, even one as small as the Canoe Basin. Climate data used in this study are from the Blue River meteorological stations, located approximately 65 km south of the Water Survey of Canada streamflow gauge at the outlet of the Canoe Basin. This is a large distance over which to extrapolate climate data, particularly in mountain environments. Also, the climate station is at 689 m a.s.l., well below the median elevation of the Canoe Basin of 1957 m a.s.l. However, the choice of climate station is justified in this study because it provides the longest complete record of temperature and precipitation. Closer meteorological stations at Valemount (“Valemount” and “Valemount East”) might provide a better approximation of climate in the Canoe Basin, as they are at 797 m a.s.l. and only 13 km from the watershed outlet. However, the climate records are incomplete and have long periods of missing data. 47 D igital photogram m etry The accuracy of digital photogrammetry depends on the quality o f air photos and the subjective ability of the researcher to interpret the photos. Poor contrast is a frequent problem in snow/ice covered regions if light conditions are not ideal (Keutterling and Thomas, 2006) and early season snowfall can mask glacier termini (Granshaw and Fountain, 2006). Relief and shadow are problematic and can result in incomplete or inaccurate data collection (Surazakov and Aizen, 2006). In this study, challenges in air photo interpretation resulted from photo coverage and scale, poor contrast and seasonal snow cover. The larger scales of photos from 1948, 1996 and 2005 (1:20,000-1:40,000) likely resulted in more accurate data collection in comparison to the smaller scales of the 1955, 1970 and 1985 photos (1:60,000-1:80,000). Fresh snow and poor contrast in the accumulation area of the 1955 photos and poor contrast in the accumulation area of the 1970 photos limited mass point collection; therefore 1955 and 1970 air photos were not used in volume change calculations despite near complete (1955) and complete (1970) watershed coverage in the air photos. Photos from 1996 were missing much of the watershed and fresh snow completely obscured glacier termini with the exception of the North Canoe glacier. 2.6 Conclusion In this study, digital photogrammetry was used to examine glacier area and volume change in the Canoe Basin from 1948 to 2005 in relation to climate variability. Although every glacier terminus retreated from 1948 to 2005, several glaciers advanced from 1955 to 1985. Other studies have likewise noted advances of glaciers in western North America during this time, and this advance coincides with cool, wet conditions during the negative 48 phase of the PDO. The different responses of glaciers in the watershed to climate are likely related to differences in glacier morphology. No simple relationship between glacier morphology and climatic response was apparent in this study, however. Negligible area and volume changes occurred from 1948 to 1985, but glaciers dramatically retreated and lost volume from 1985 to 2005. These results accord with the general consensus that the rate of glacier retreat increased during the latter decades of the 20th century. Mean annual and winter temperature, and total summer precipitation increased from 1947 to 2006. The decrease in glacier volume over the period of record is largely attributed to the increase in temperature. 49 3. Quantifying glacier contribution to streamflow in the Canoe Basin, British Columbia, using the HBV-EC watershed model. 3.1 Introduction Glaciers are an important source of late-summer streamflow, and the impact of glaciers on streamflow is evident in watersheds with as little as 2-5% glacier cover (Stahl and Moore, 2006). Glaciers moderate seasonal to interannual variability in streamflow by storing water in the winter and during cool years and releasing it in the summer and during warm years (Krimmel and Tangbom, 1974). This moderating effect is particularly important for drought susceptible regions, including the southern Columbia Basin, since glacier melt produces streamflow during the summer months which are typically hot and dry. Glaciers throughout British Columbia significantly retreated during the past century (Bolch et al., 2010) and streamflow in many glacierized watersheds has changed as a result (e.g. Stahl and Moore, 2006). In particular, glacier area in the Mica Basin, a large watershed in the Canadian portion of the Columbia Basin, decreased by 14% from 1985 to 2005 (Jost et al., 2012). Continued warming is expected to further reduce glacier area in the Columbia Basin and will thus directly affect the timing and magnitude of streamflow. Glacier contributions to streamflow may be defined as "glacier melt" or "glacier wastage" (Comeau et al., 2009). "Glacier melt" refers to the ice melt from the ablation area during the summer that is offset by the water equivalent of snow falling in the accumulation area in the winter. "Glacier wastage" is the volume of ice melt in excess of the annual input of snow in the accumulation area, and results in a negative glacier mass balance (Comeau et al., 2009). Glacier wastage represents the additional volume of streamflow that is associated with the loss of mass from glacier storage. 50 Glacier contribution to streamflow is quantified for areas in western North America using various methods. Granshaw and Fountain (2006), for example, estimated that glaciers contribute up to 6 % of August-September streamflow in the North Cascades National Park Complex, Washington. This is based on an estimate of glacier volume loss and streamflow records from 1958 to 1988. Comeau et al. (2009) used a hydrologic model to demonstrate that glacial wastage (melt) could contribute as much as 10% (27%) of streamflow from July to September in sub-basins of the North and South Saskatchewan river basins, Alberta. Nolin et al. (2010) analyzed stream chemistry to conclude that up to 41-73% of late summer streamflow (melt and wastage) originates from glaciers in a watershed on Mount Hood, Oregon. Jost et al. (2012) quantified glacier contribution to streamflow in the Mica Basin, which is located in the Canadian portion of the Columbia Basin and includes the Canoe Basin. The study concluded that 6 % of mean annual streamflow and 25% (35%) of August (September) streamflow is derived from glacier wastage in the Mica Basin, which currently contains less than 5% glacier cover. Hydrologic modeling provides one way to quantify glacier contribution to streamflow. Hydrologic model calibration typically optimizes model parameters to best reproduce historic streamflow during the calibration period. However, it has been shown that the same streamflow can be reproduced using different combinations of values for the model parameters; this phenomenon is known as equifinality. The implication of equifinality is that one set of parameters might reproduce streamflow while poorly representing the hydrologic processes that are responsible for the streamflow response (Braun and Aellen, 1990). Model calibration is therefore greatly improved using a multi-criteria approach. Calibrating a model to reproduce streamflow and glacier mass balance helps to ensure that the model accurately 51 represents glacier accumulation and ablation while still reproducing historic streamflow (e.g. Braun and Aellen, 1990; Konz and Seibert, 2010). In this study the HBV-EC watershed model was used to simulate streamflow and quantify glacier contribution to streamflow from 1947 to 2006 in the Canoe Basin, in the headwaters of the Columbia River Basin. Simulated streamflow and glacier contribution to streamflow were compared to observed glacier change and climate variability. Glacier contribution to streamflow is defined as the difference in streamflow between a watershed with glacier cover and that same watershed when all glaciers are converted to open alpine land class. Glacier contribution to streamflow is therefore a measure o f glacier wastage (as defined by Comeau et al., 2009), because streamflow resulting from the melting seasonal snowpack is included in model simulations both with and without glacier cover. The HBVEC model was calibrated with a multi-criteria approach using observed streamflow and records of glacier volume change. This study also includes an evaluation of whether updating glacier extents has a significant effect on the reproduction of historic hydrographs. 3.2 Study Area The Canoe Basin is located in the Premier Range of mountains, west of Valemount, British Columbia, in the headwaters of the Columbia River (Figure 3.1). For the purpose of this study, the Canoe Basin includes the watershed area above the Water Survey of Canada (WSC) streamflow gauge 08NC004 “Canoe River below Kimmel Creek” (52.728° N, -119.408° W). The watershed encompasses 307 km 2 of mountain terrain, has a median elevation of 1957 m a.s.l., and glaciers occupied 18% of the catchment in 2005. At Blue River, approximately 65 km south of the outlet of the Canoe Basin, mean annual precipitation is 1100 mm of which 38% falls as snow. Mean annual temperature is 4.5°C 52 with a monthly minimum and maximum of -9 °C (January) and 16.4°C (July), respectively. As further explained in section 3.3.1 of this thesis, climate data from Blue River were used in this thesis because it provided the longest complete record of temperature and precipitation within proximity to the Canoe Basin. Although several meteorological stations are closer to the Canoe Basin, the length of those station records was shorter than those from Blue River (e.g. Cariboo Lodge, McBride 4SE) or they had extensive periods with missing data (e.g. Valemount East, Cariboo Lodge). 120* W 119* W 118* W BC/Alberta Border Major W atersheds Rivers/Lakes fumbia Basin Basin 50 Km j 120* W 5 V' -% South ;fhf>mp^on Basin 119’ W 118° W Figure 3.1. Location o f the Canoe Basin, British Columbia. 53 3.3 Methods 3.3.1 Data A 25 m digital elevation model (DEM) produced under the Terrain Resource Inventory Management (TRIM) program was the basis o f watershed and stream delineation in the HBV-EC model using the Green Kenue interface version 3.1.55 (National Research Council, 2010). Forest cover polygons were extracted from the Vegetation Resources Inventory (VRI) dataset (B.C. Ministry o f Forests, Lands and Natural Resource Operations (BCMFLNRO), 2012), with ‘forest’ defined as all areas with crown closure greater than 20%. No lakes are present in the watershed. Glacier polygons for 1948, 1955, 1970, 1985 and 2005 were obtained from digital photogrammetry (please refer to Chapter 2 of this thesis). All areas not included in ‘glacier’ or ‘forest’ polygons were considered to be ‘open’ (Table 3.1). ‘Open’ areas have less than 20% crown closure, and also include rocky areas and plains. Glacier areas and calculated volume changes differ in this Chapter from those expressed in Chapter 2 because of the areas used to derive those estimates. In Chapter 2, glacier areas were based on the areas common to all years o f photography. To ensure complete glacier cover in the HBV-EC model for Chapter 3, missing glacier extents were replaced by polygons from the next closest year of photography. The glacier volume change reported in Chapter 2 was -0.22 ± 0.04 km 3 water equivalent; this was calculated using glacier areas common to 1948, 1985 and 2005. For HBV-EC calibration, glacier areas common only to 1985 and 2005 were used, and the resulting glacier volume change is -0.24 + 0.05 km 3 water equivalent over the twenty year period (1985-2005). 54 Table 3.1. Area of land classes used in the HBV-EC simulations for each year. Land Use Class Glacier (km2) Open (km2) Forest (km2) 1948 60.2 167.4 79.5 1955 61.0 166.6 79.5 1970 60.0 167.6 79.5 1985 59.9 167.7 79.5 2005 55.9 171.7 79.5 Streamflow and climate data are available at monitoring stations within proximity to the watershed. Streamflow data used in model calibration and validation were from the Water Survey of Canada’s streamflow gauge 08NC004 “Canoe River below Kimmel Creek,” from 1972 to 2006. The location of the streamflow gauge was used as the watershed outlet in the hydrologic model. Climate data (daily precipitation and temperature) were combined from Environment Canada’s meteorological stations Blue River North (BRN) and Blue River A (BRA) to form a complete record from 1947 to 2006 (Environment Canada, 2010). Although several meteorological stations are closer to the Canoe Basin, the length of those station records was shorter than BRN and BRA (e.g. Cariboo Lodge, McBride 4SE) or they had extensive periods with missing data (e.g. Valemount East, Cariboo Lodge). A complete climate record was obtained by combining the records from BRN (1947-1985) and BRA (1969-2006) stations. Least squares regression based on daily values was used to determine the relationship between climate variables at the two stations during the period of overlapping records (1969-1985; Figure 3.2). The resulting equations (Table 3.2) were used to infill missing values (0.5% of the climate record) in the BRA data and were then applied to BRN data to extend the climate record back to 1947. The combined record o f BRN and BRA stations is referred to as Blue River climate data in this study. To ensure splicing the BRN and BRA stations did not introduce spurious trends to the climate record, the monthly climate record created by combining the BRN and BRA records was compared to Environment 55 Canada’s rehomogenized monthly climate data (Mekis and Vincent, 2011; Vincent et al., 2012) for Blue River (see Appendix A for comparison). Table 3.2. Linear regression equations for infilling daily climate data. Regression Equation BRA = (1.0283 x BRN) - 0.085 BRA = (0.6615 x BRN) + 0.4754 Climate Variable Mean Temperature Precipitation Sample Size 5753 5753 R 0.99 0.65 Standard Error 1.1 °C 2.8 mm For comparison, regression coefficients (R2) and standard error (SE) were also calculated using the BRN and BRA data divided into seasons (summer: May-September, sample size = 2448; winter: October-April, sample size = 3335). The correlation between BRN and BRA mean temperature was similar in the winter (R2 = 0.98, SE = 1.3°C) and in the summer (R 2 = 0.96, SE = 0.84°C). Precipitation was more highly correlated between BRN and BRA in the summer (R 2 = 0.70, SE = 2.7 mm) than in the winter (R 2 = 0.67, SE = 2.9 mm). O • ro • o _ o (H t I O_ Ia . o __ o £ a > . ~ » * jjf lp * • / ♦ %* * «- o „ - • • • • ** . • J • • ## * •• • • « jT / / • 4* 4 / * y / • ! A 3 2 < QC o o C j» - < C L C D o t ^ 'jU iS S r .y ■ o c~> — I ^ K o o - . i ? * •- * jB E & U 'i; \ -30 20 -10 0 20 10 30 40 BRN Total Daily Precipitation (mm) BRN Mean Daily Temperature (*C) Figure 3.2. Scatter plots showing the overlapping period between the BRA and BRN climate stations (n=5753). Regression lines are also shown (solid lines). 56 Monthly potential evapotranspiration (PET) for the Canoe Basin was estimated using the Hamon (1961) temperature-based method: PET = 0.55 D2Pt (3.1) where D is daylight hours in 12-hour units, and Pt is the saturated water vapour density, calculated using: A QC *0.062 Ta Pt = (3.2) where Ta is the air temperature. Several other methods are available to estimate potential evapotranspiration, including mass-transfer methods, radiation based methods and combination methods. However, these methods require measures of humidity, radiation and wind speed (Xu and Singh, 2002), which are not available for the Blue River meteorological stations. The Hamon method was chosen in this study because the minimal input requirements (air temperature and daylight hours) are readily available. 3.3.2 The HBV-EC Hydrologic Model The HBV-EC hydrologic model is based on the 1997 version of the Hydrologiska Byrans Vattenbalansavdelning model (HBV-97) developed at the Swedish Meteorological and Hydrological Institute (Lindstrom et al., 1997). The HBV-EC model is a semi­ distributed conceptual model that simulates snow/glacier melt and runoff hydrographs using a day as a maximum time step. The Environment Canada (EC) version of the HBV model uses the user-friendly Green Kenue graphic user interface (National Research Council, 2010). HBV-EC requires a digital elevation model (DEM), land use classification (glacier, open, forest, lake) and climate data (daily temperature and precipitation, and monthly potential 57 evapotranspiration). There are also 34 model parameters (Table 3.3) which require calibration. Although the model simulates glacier runoff, it does not incorporate changes in glacier area or volume through time. The HBV-EC model has been applied in several studies (Table 3.4), and has proven useful for simulating runoff in mountain watersheds. 58 Table 3.3. Description of HBV-EC model parameters that require calibration (National Research Council, 2010). Parameter Description RFCF SFCF TLAPSE TT TTI PLAPSE (low) PLASPE (high) EMID EPGRAD ETF TFRAIN TFSNOW Snow AM Snow TM Snow CMin Snow DC Snow M RF Snow CRFR Snow WHC Snow LWR Soil FC Soil BETA Soil LP Glacier MRG Glacier AG Glacier DKG Rainfall correction factor Snowfall correction factor Temperature lapse rate (°C rri1) Threshold temperature (°C) Threshold temperature interval(°C) Fractional increase in precipitation with elevation below EMID (mm m '1) Fractional increase in precipitation with elevation above EMID Threshold elevation(m) Fractional rate o f decrease of potential evaporation per metefrm'1) Potential evaporation temperature anomaly Fraction of rain reaching the ground Fraction of snow reaching the ground Slope and aspect influence on melt factor Snowmelt threshold temperature (°C) Melt factor for open areas at the winter solstice (mm °C 1 day'1) Melt factor increase from winter to summer solstice (mm °C 1 day ') Ratio of melt factors in forest and open areas Rate o f water refreezing in the snowpack (mm °C 1 d a y 1) Snowpack water holding capacity Amount of water that can be retained by a snowpack (mm) Soil field capacity (mm) Relationship between soil water release and infiltration Soil moisture content below which evaporation is supply-limited (proportion of FC) Ratio between glacier ice melt and seasonal snow melt Relationship between runoff coefficient and glacier snowpack water equivalent (m m 1) Difference between the minimum and maximum glacier water storage outflow coefficients (d a y 1) Minimum glacier water outflow coefficient Recession coefficient Proportion o f fast reservoir released per day Fast reservoir release rate adjustment (exponent) Proportion o f slow reservoir released per day Runoff fraction sent to fast reservoir Initial rate o f discharge from the fast reservoir (1000 L s '1) Glacier KGMin Glacier KGRC Runoff KF Runoff Alpha Runoff KS Runoff FRAC Initial Fast Reservoir Discharge Initial Slow Reservoir Discharge Initial rate o f discharge from the slow reservoir (1000 L s'1) 59 Table 3.4. Previous studies that have used the HBV-EC model. Basin Area (km2) Glacier Cover (%) Lillooet River Basin, B.C. 2160 17 Simulated streamflow to apply a new version of the model. Hamilton et al. (2000) M’Clintock Basin, Yukon 1700 0 Examined discharge under ice-covered conditions. Stahl et al. (2008) Bridge River Basin, B.C. 152.4 61.8 Coupled a glacier response model with the HBV-EC model to simulate the impacts of warming trends on glacier mass balance and streamflow. Hirshfield (2010) Goathom Creek and Moffat Creek Watersheds, B.C. 126; 539 0.7; 0 Simulated the impacts of climate and land use change on streamflow. Jost et al. (2012) Mica Basin, in the Canadian portion of the Columbia River Basin 20,742 5 Used the model to quantify glacier contribution to streamflow. Study Watershed Moore (1993) 60 Use of HBV-EC model HBV-EC Structure and Routing The HBV-EC model divides the watershed into climate zones based on elevation bands and horizontal climate gradients. Each climate zone is sub-divided into Grouped Response Units (GRUs), which are areas of similar land cover, elevation, slope and aspect (Stahl et al., 2008). Water is routed through the fast and slow reservoirs, and water release from those reservoirs is controlled by several parameters. Runoff FRAC is the fraction of runoff that the model directs to the fast reservoir; if a watershed response to precipitation is quick (delayed), Runoff FRAC will be high (low). The rate of release of water from the fast reservoir is determined by Runoff KF and Runoff Alpha, and the release rate from the slow reservoir is determined by Runoff KS (CHC and WSC, 2010). The total streamflow modelled at the basin outlet is the combination o f outflows from the glacier reservoirs, the fast reservoir and the slow reservoir (Stahl et al., 2008). Water Balance The HBV-EC model is based on the following water balance equation: „ d [S P + S M + U Z + L Z + lakes] F — t — y — -----------------^ ----------------- ( j -j ) where P is precipitation, E is evaporation, Q is runoff, SP is snowpack, SM is soil moisture, UZ is the upper groundwater zone, LZ is the lower groundwater zone and lakes is the volume of the lakes in the watershed (Bergstrom, 1995). Precipitation The difference between the amount of precipitation recorded at the climate station and the amount of precipitation in the watershed is accounted for by changing the Snowfall 61 Correction Factor (SFCF) and Rainfall Correction Factor (RFCF) and by changing the precipitation gradient with elevation (PGRADH and PGRADL). The phase of precipitation (rain or snow) falling in different parts of the watershed is determined by the threshold temperature (TT; °C). Mixed precipitation occurs within the threshold temperature interval (TTI; °C) above or below the threshold temperature. The temperature lapse rate (TLAPSE; in °C m_1) and precipitation gradients above (PGRADH) and below (PGRADL) the threshold elevation (EMID) are used to adjust the input climate data according to basin elevation. The effects of canopy interception on the fraction of precipitation reaching the ground are accounted for by the unitless rain (TFRAIN) and snow (TFSNOW) throughfall parameters (Stahl et al., 2008). Snowmelt The HBV-EC model calculates snowmelt using a temperature index approach (Bergstrom, 1995) according to a melt factor (MF) which varies between minimum and maximum values (CMIN and CMIN + DC) on the winter and summer solstices, respectively (Hamilton et al., 2000). MF is adjusted for slope (s) and aspect (b) according to the equation MF = MFflat x [1 —AM x sin(s) x cos(b)] (3.4) where MFflat is the melt factor for flat terrain (mm d~'), and AM is a dimensionless model parameter (CHC and WSC, 2010). Forest canopy shading and sheltering, which reduce snowmelt in forested areas, are accounted for by the melt reduction factor (MRF). When the air temperature is below the melt threshold, liquid water can refreeze in the snowpack at a rate controlled by CRFR. Water release from the snowpack to soil moisture storage or the glacial ice surface occurs when the water holding capacity (WHC) of the snowpack is 62 exceeded (Moore, 1993). The maximum amount of liquid water that a snowpack can retain is defined by LWR (CHC and WSC, 2010). Soil Routing and Evaporation The HBV-EC model calculates soil moisture differently for open and forested areas in each elevation band. The threshold soil moisture below which actual evapotranspiration falls below potential evapotranspiration (LP) and the soil field capacity (FC) are model parameters. The amount of water in excess of the soil field capacity drains out of the soil and into the fast and slow reservoirs. The accumulated runoff from forest and open areas in all elevation bands is routed to the fast and slow reservoirs. Evapotranspiration is calculated according to the monthly average potential evaporation values with actual evapotranspiration being less than the potential rate, and occurs only when the entire snowpack has depleted (Hamilton et al., 2000). Glacier Routing The model treats each glacier GRU separately. The outflow from each glacier GRU is calculated based on the outflow coefficient, which varies from KGMIN when there is deep snow and poorly developed glacier drainage channels, to KGMIN + DKG in the late-summer conditions of bare ice and efficient glacier drainage channels (Stahl et al., 2008). The parameter MRG accounts for an enhanced melt factor due to a reduction in surface albedo after the seasonal snowpack has melted (Stahl et al., 2008). 63 3.3.3 HBV-EC model set up The HBV-EC model was set up for the Canoe Basin using the Green Kenue interface (National Research Council, 2010). The watershed extent was extracted from the DEM, and the location of the WSC streamflow gauge defined the watershed outlet node. Land use polygons (Figure 3.3) were overlain on the watershed using the updated glacier extent derived for each decade from the photogrammetric analysis. The watershed was divided into grouped response units using ten elevation bands, two slope bands, and four aspect bands. One climate zone was used for the entire watershed. The daily Blue River climate data are distributed over the GRUs according to the model parameters discussed previously. Glacier Forest Open Figure 3.3. Land use classification used to set up the HBV-EC model. 64 3.3.4 HBV-EC calibration and validation The HBV-EC model was calibrated using observed streamflow from 1986 to 2005 because this period aligns with a period of calculated glacier change in the Canoe Basin (Chapter 2). The model was validated against the remaining observed streamflow record from 1972 to 1985. Observed streamflow data were missing from 13 April 1995 to 25 August 1997; therefore, these years were not included in the assessment of model performance during the calibration period. Model calibration was also constrained with the glacier volume change for the period 1985-2005 of between —0.19 km 3 and —0.29 km 3, derived from digital photogrammetry (see Chapter 2). Although the HBV-EC model does not include provisions to change the glacier extents in subsequent years, the spatial change in glacier extent was incorporated into the model by dividing the calibration period into two time-steps. The first time-step, 1986-1995, was simulated using the 1985 glacier extents. The second time step, 1996-2005, was simulated using the 2005 glacier extents. A five-year model spin-up period from 1965 to 1971 was used to ensure model reservoirs were in equilibrium with climate. The methods used in model calibration were based on the optimization routine provided by Georg Jost (B.C. Hydro). Model parameters were optimized using a Latin hypercube sampling of 100,000 iterations. Latin hypercube sampling is an effective method to systematically sample across a parameter space (McKay et al., 1979). Computing requirements of Latin hypercube sampling for 13 parameters is a considerable limitation, and one computer would take approximately 16 days to complete 10,000 iterations for the Canoe Basin. To minimize time requirements, 13 computers were used, each running the HBV-EC optimization sampling code on two CPUs. Each computer ran 385 iterations on each CPU. This method of high performance computing was justified because, while 10,000 iterations 65 had to be run, each iteration was independent of all others. Using this approach, 10,000 iterations were completed in less than 24 hours. To conduct the Latin hypercube sampling, an initial range of values was specified for each of the 13 parameters identified by Jost et al. (2012) as correlated with glacier mass balance and Nash Sutcliffe efficiency. Snow LWR was adjusted to 250 mm, as recommended by Georg Jost (B.C. Hydro). Default values were used for the remaining 20 parameters that were not included in the optimization routine. Initial ranges for the 13 calibrated parameters were based on values used in previous studies that used HBV-EC (Moore, 1993; Hamilton, 2000; Stahl et al., 2008; Hirshfield, 2010; Jost et al., 2012). Goodness-of-fit or model performance criteria adopted for this study include the Nash Sutcliffe efficiency (E; Nash and Sutcliffe, 1970), percent error in mean August streamflow (August error; AE) and root mean square error (RMSE). The Nash Sutcliffe efficiency, August error, and root mean square error are calculated as follows: E = (3.5) 1 — E t= i(Q o -Q o )i _ £ t = i Qs ~ £ t = i Qo (3.6) (3.7) RM SE where Q0 is observed discharge, Qs is simulated discharge, t is time and the overbar denotes the average value. August error was used as a performance measure to ensure the model accurately reproduced late summer streamflow, which is influenced by the presence of glaciers in a watershed. Model performance was constrained to reproduce glacier volume change over the calibration period (1986-2005) to -0.24 ± 0.05 km 3, as calculated via digital photogrammetry (Chapter 2 of this thesis). 66 Model calibration was run in three steps, with each step based on a smaller range of parameter values than the previous. In the first calibration step, 40,000 iterations of the optimization routine were run using the initial parameter ranges. Each iteration was divided in two time periods (1986-1995 and 1996-2005), and the goodness-of-fit statistics were averaged over the two time periods to get the goodness-of-fit statistic for the entire calibration period. Hydrographs were generated for each parameter set that met criteria to visually examine the hydrograph in comparison to the observed streamflow record. All parameter sets that met criteria (glacier volume change -0.24 + 0.05 km3, August error (AE) < 10%, and Nash Sutcliffe efficiency (E) > 0.85) were extracted and new parameter ranges were compiled using the maximum and minimum values. In the second calibration step, 50,000 iterations of the optimization routine were completed using the refined parameter ranges. To further improve model calibration, the parameter ranges were then reduced according to the maximum and minimum values in the top 10 parameter sets from the second step. The calibration routine was run for an additional 10,000 iterations using the refined parameter ranges. After the third optimization step, the performance criteria were refined to E > 0.90 and AE < 5%. All parameter sets that met the revised criteria were used to define uncertainty bounds for simulated model outputs. The model was validated by comparing simulated streamflow to observed WSC streamflow from 1972 to 1985 using each of the top parameter sets. Glacier change was not used as a test of model validation because the poor quality of the 1970 air photos did not allow for a glacier volume change calculation. Therefore, model validation was only based on E and AE from 1972 to 1985. 67 3.3.5 Data extraction and compilation Following model calibration and validation, the HBV-EC model was forced with Blue River climate data from 1947 to 2006 to reconstruct the streamflow record. For each of the top parameter sets, the HBV-EC model was run with each year of glacier extents, using a five-year spin up period (1942-1946) with each simulation. Time series (total discharge in m3s_l and ice melt in mm d_1) were compiled using the 1948 glacier area for years 19471951, the 1955 glacier area for years 1952-1962, the 1970 glacier area for years 1963-1977, the 1985 glacier area for years 1978-1995 and the 2005 glacier area for years 1996-2005. Thus, the compiled time series accounted for changing glacier area throughout the simulations. Discharge time series were also extracted using only 2005 glacier extents (i.e. without changing glacier extent through time) to examine the effect of using static and dynamic glacier cover in the HBV-EC model. A simulation was also run using a land classification that assumed zero glacier cover (by converting ‘glacier’ land classification to ‘open’) in the watershed, for use in calculating glacier contribution to streamflow. The HBV-EC model was also forced with Environment Canada’s rehomogenized daily climate data (Mekis and Vincent, 2011; Vincent et al., 2012) for Blue River, using the best performing model parameter set. Simulated streamflow was compared to streamflow simulated using the combined climate records from Environment Canada’s BRA and BRN stations (see Appendix C for model output and comparison). 3.3.6 Data analysis Climate data from Blue River were analyzed in Chapter 2 of this thesis (2.3.1 and 2.4.1 discuss data/methods and results, respectively). 68 Stream flow sensitivity to gla cier area change in the H B V -E C m odel The effect of updating glacier extents in historic streamflow simulations was examined by comparing streamflow simulated with static and dynamic glacier cover. The model was used to simulate the year with the maximum August glacier contribution to streamflow from 1996 to 2005 using 1955 and 2005 glacier extents, which were the largest and smallest glacier extents on record, respectively. As streamflow in August and September would likely be impacted by changing glacier area, monthly relative volume errors were calculated for these months. Uncertainty bounds were assigned using the top eight parameter sets. Streamflow variability The Mann-Kendall test (Mann, 1945, and Kendall, 1975) provides a means to assess significance of monotonic trends (95% significance level). Because autocorrelation reduces the number of independent observations, a modified Mann-Kendall test was used to correct the power of the test by adjusting the degrees of freedom (Yue and Wang, 2004). The Kendall-Theil robust line (Theil, 1950) was used to determine the magnitude of observed trends. Time series were analyzed over several time periods aligning with the photographic record and with recorded climatic fluctuations. Simulated mean streamflow was tested on an annual (hydrological year from 1 Oct to 30 Sep), winter and summer (Oct-Apr and May-Sep) and monthly (June, July, August, September) basis. Simulated streamflow was also compared over periods of calculated glacier volume change (1948-1985 and 1985-2005). Researchers commonly examine changes in the timing o f streamflow through time using the date of center of volume of flow (e.g. Stewart et al., 2004; Stewart et al., 2005; Bum, 2008) and the date of peak flow (e.g. Bum and Elnur, 2002; Bum et al., 2004). 69 However, trends in the timing of the center of volume of flow can be misleading, particularly when glacier melt or intense late-summer precipitation are significant sources of runoff (Dery et al., 2009). The date of peak flow, which is the timing of maximum runoff in snowmelt dominated basins, can be masked or overshadowed by precipitation events that occur in autumn or winter, particularly in small basins (Dery et al., 2009). As an alternative method, Dery et al. (2009) suggest dividing annual daily records of flow into pentads (5-day averages) and examining trends in these records. The five day period is long enough to reduce the impact of short intense precipitation or warming events, yet short enough to reproduce small changes in the timing of flow (Dery et al., 2009). Endpoints of the Kendall-Theil robust line for each 5-day period were used to reconstruct initial and final hydrographs over a period of record, which can then be visually examined for changes in the timing and volume of flow (Dery et al., 2009). This method was implemented on the simulated streamflow record from 1948 to 2006 and the shorter periods of 1948-1985 and 1985-2006. Time series of select 5day mean periods were also analyzed for the presence of trends. Glacier contribution to streamflow Percent glacier contribution to streamflow was calculated by subtracting discharge without glacier cover from discharge with glacier cover and dividing by the total discharge. Glacier contribution to streamflow was calculated on an annual, seasonal and monthly basis. Mean, maximum and minimum glacier contributions to streamflow from 1947 to 2006 were calculated and expressed as uncertainty bounds within the top eight parameter sets. Time series produced from the HBV-EC model using the best performing parameter set were used for trend analyses. All time series of glacier contribution to streamflow were tested for the presence of autocorrelation prior to trend analysis. 70 Significant (p < 0.05) trends were identified using the modified Mann-Kendall equation, and the magnitude of trends was calculated as the slope of the Kendall-Theil robust line. Trends were examined on an annual, seasonal and monthly basis, and also in relation to periods o f calculated glacier volume change. The date of onset of glacial melt was also identified from the time series of the top parameter set, defined as the first day when glacial melt occurred each year. Days were assigned according to the hydrologic year (1 Oct - 30 Sep). 3.4 Results 3.4.1 HBV-EC calibration and validation In the first set of 40,000 iterations, 21 parameter sets met the glacier volume change criteria with less than 10% August error (AE) and Nash Sutcliffe statistic (E) greater than 0.85. The best parameter set resulted in E = 0.89, AE = 2.2% and glacier volume change of -0.24 km3. In the second step of model calibration, 225 of the 50,000 parameter sets met criteria, with the best parameter set achieving E = 0.90, AE = 0.53% and glacier volume change -0.21 km3. After the third step of model calibration, 8 parameter sets out of all 100,000 iterations met the refined criteria (Table 3.5); they were subsequently used to define uncertainty bounds for the simulated data (Table 3.6). In the calibration period, the best parameter set out of the 100,000 iterations yielded a glacier change of -0 .2 8 km3, with E = 0.91, AE = 2.1% and RMSE = 5.21 m3s_1 (Table 3.7, Figure 3.4). During the model validation period, the top parameter set reproduced streamflow with E = 0.91, AE = 4.9% and RMSE = 5.18 m3s-1 (Table 3.7, Figures 3.4). In the validation period, the top eight parameter sets had E values ranging from 0.89 to 0.91 and AE values ranging from 4.9% to 10.8% (Table 3.7). Total annual flow volumes were well reproduced by the model (Figures 3.5 and 3.6). 71 Table 3.5. Parameter ranges used in each step of model calibration. Nash Sutcliffe statistic (E) and August error (AE) values are listed for the bestperforming parameter set in each step. Glacier DKG Glacier KGMin Runoff KF Runoff Alpha Runoff KS Runoff FRAC E (best) Absolute AE (best) Glacier Volume Change (km3) 0 .0 - 0 .7 0 .0 1 -0 .5 5 0.07 - 0.49 0.31430 3.64869 2.65 - 4.0 0 .0 - 3 .0 1.73919 0 .4 -0 .8 5 0 .4 - 0 .8 0 .4 -0 .9 5 0.60630 1.78445 q oi I o Glacier MRG 0.00050 O cn 1 00 o Snow MRF 0.0004 - 0.0006 o rS I Snow DC 0.00025 - 0.00065 b Snow CMin 0.00002 - 0.0022 © Snow AM Best Parameter Set (of 100,000 iterations) 0.00744 0 1 PLAPSE Step 3 (10,000 iterations) 0.007 - 0.008 O Temperature lapse rate (°C m' ) Fractional increase in precipitation with elevation (mm m '1) Slope and aspect influence on melt factor Melt factor for open areas at the winter solstice (mm “C 1d a y 1) Melt factor increase from winter to summer solstice (mm °C'1 d a y 1) Ratio of melt factors in forest and open areas Ratio between glacier ice melt and seasonal snow melt Difference between the minimum and maximum glacier water storage outflow coefficients (d ay 1) Minimum glacier water outflow coefficient Proportion of fast reservoir released per day Fast reservoir release rate adjustment (exponent) Proportion of slow reservoir released per day Runoff fraction sent to fast reservoir Step 2 (50,000 iterations) 0.006-0.008 © 04 1 O TLAPSE Step 1 (40,000 iterations) 0.004 - 0.0095 1 Description © Parameter 1.1 - 2 .0 0 .0 1 -0 .2 5 0.05 - 0.25 0 .0 5 -0 .2 5 0.22932 0.004 - 0.06 0.004 - 0.06 0.005 - 0.06 0.04598 0 .0 1 -0 .3 5 0 .0 2 -0 .3 5 0 .0 8 -0 .3 2 0.21942 0 .0 -0 .5 0.03 - 0.45 0 .0 3 -0 .4 0.19083 0.0005 - 0.08 0 .0 0 6 -0 .0 8 0.007 - 0.035 0.00722 0 .2 5 -0 .9 5 0 .2 5 -0 .8 0 0 .4 - 0 .7 0.60348 0.89 0.90 0.95 % 0.91 0.91 2.1 % 2.1 % 0.21 -0.28 -0.28 2 .2 % -0.24 72 - Table 3.6. Range o f calibrated parameters in the top eight parameter sets. Parameter Min Max Range (Max-Min) TLAPSE (°C m '1) 0.00741 0.00782 0.00040 PLAPSE (mm m ') 0.00046 0.00055 0.00009 Snow AM 0.12437 0.42031 0.29353 Snow CMIN (mm ° C ' day ') 2.89884 3.94393 1.04509 Snow DC (mm °C '‘ day'1) 1.71657 2.91963 1.20306 Snow M RF 0.45663 0.70029 0.24367 Glacier MRG 1.60732 1.94752 0.34020 Glacier DKG (d a y 1) 0.10459 0.24591 0.14133 Glacier KGMin 0.03558 0.05846 0.02288 Runoff KF 0.13244 0.30410 0.17166 Runoff Alpha 0.03917 0.27002 0.23085 Runoff KS 0.00722 0.01976 0.01254 Runoff FRAC 0.54965 0.67105 0.12140 Table 3.7. Nash Sutcliffe efficiency (E), absolute August Error (AE) and root mean square error (RMSE) values obtained during the calibration and validation periods for the top eight parameter sets. Simulated glacier volume change for the calibration period is also shown. Simulated glacier volume change was not used in model validation; therefore, it is not shown. Calibration Parameter Set 1 2 3 4 5 6 7 8 MAX MIN E AE RMSE 0.91 0.90 0.90 0.90 0.90 0.90 0.90 0.90 0.91 0.90 2.1 % 1.2 % 2.7 % 0.5 % 2.3 % 1.9% 0.2% 0.5 % 2.7 % 0.2% 5.21 5.24 5.26 5.28 5.31 5.34 5.38 5.38 5.38 5.21 Validation Glacier Volume Change (km3) -0.28 -0.19 -0.24 -0.21 -0.20 -0.24 -0.26 -0.20 -0.19 -0.28 73 E AE RMSE 0.91 0.90 0.90 0.90 0.90 0.89 0.90 0.90 0.91 0.89 4.9 % 8.7 % 10.4 % 8.6% 10.8 % 9.7 % 9.0 % 8.9 % 10.8 % 4.9 % 5.18 5.31 5.31 5.28 5.32 5.47 5.42 5.36 5.47 5.18 (m Simulated Streamflow Figure 3.4. Scatterplots of observed and simulated streamflow in the calibration (left) and validation (right) periods, for the top parameter set. One-to-one line is shown. CM CM CO CO 14 16 18 20 22 Observed Annual Flow Volume (1 011 m3) Figure 3.5. Scatterplot of observed and simulated annual flow volumes. One-to-one line is shown. 74 op 1975 1960 1990 1995 2000 2005 Figure 3.6. Comparison of observed and simulated total annual flow volumes based on the hydrologic year. Observed total flow volume is not available from 1995-1997. Visual examination of the simulated hydrographs indicates that the HBV-EC model accurately simulates the streamflow response to snow and ice melt and the seasonal peak flows. Winter low flows are also well simulated. The HBV-EC model, however, lacks the ability to accurately simulate the streamflow response to intense precipitation events, particularly those that occur between August and October. M odel performance was highest in 1990 (E = 0.94, AE = 0.90%; Figure 3.7) and lowest in 1987 (E = 0.76, AE = 4.40%; Figure 3.7). Percent relative volume errors for 1990 and 1987 were 0.2% and 8.0%, respectively. The years with highest and lowest model performance were both in the calibration period. 75 (a) Best year = 1990 o o 0 200 100 300 Day of Year o o o GO <£> O V o 0 200 100 300 Day of Year Figure 3.7. Observed and simulated hydrographs for the years with (a) best and (b) worst model performance. Simulated discharge is displayed from the top performing parameter set. Day o f year corresponds to the hydrologic year. 3.4.2 Streamflow sensitivity to glacier area change in the HBV-EC model Simulated streamflow time series using static and dynamic glacier extents were compared to determine the effect of changing glacier cover in historic streamflow simulations. This analysis was based on simulations for the year 1998, which had the highest August glacier contribution to streamflow between 1996 and 2005. From the best parameter set, simulated late summer streamflow was lower using the 2005 (smaller) glacier extent 76 (Figure 3.8, top) with monthly relative volume errors of 7.4% and 6.8% for August and September, respectively. However, the difference between the simulations with 1955 and 2005 glacier cover is within the uncertainty bounds of the top eight parameter sets (Figure 3.8, bottom). s oIG 2005giaaers 1955gteciers 1955-2005 glaciers o o o 0 200 100 300 Day of Year o at o 1960 1960 1960 1990 2000 1970 1960 1990 2000 1970 1980 1970 o to m 1950 1960 o s o 1960 1960 1980 1970 1990 2000 o CD Rsconttructed Di•charge ( rnJ « 1 ) © 5 day mean penod 6 ? 2 Aug - 6 Aug o accessed 31 October 2012. 99 Bum, D.H. 2008. Climatic influences on streamflow timing in the headwaters of the Mackenzie River Basin. Journal o f Hydrology, 352: 255-238. Bum, D.H. and Elnur, M.A.H. 2002. Detection of hydrologic trends and variability. Journal o f Hydrology, 255:107-122. Bum, D.H., Cunderlik, J.M. and Pietroniro, A. 2004. Hydrological trends and variability in the Liard River Basin. Hydrological Sciences, 49(1): 53-67. Carrara, P.E. 1987. Holocene and latest Pleistocene glacial chronology, Glacier National Park, Montana. Canadian Journal o f Earth Sciences, 24: 387-395. Champoux, A and Ommanney, C.S.L. 1986. Evolution of the Dlecillewaet Glacier, Glacier National Park, B.C., using historical data, aerial photography and satellite image analysis. Annals o f Glaciology, 8: 31-33. Chen, J. and Ohmura, A. 1990. On the influence of Alpine glaciers on runoff. International Association o f Hydrological Sciences, 193: 117-125. Collins, D.N. and Taylor, D.P. 1990. Variability of runoff from partially-glacierised Alpine basins. International Association o f Hydrological Sciences, 193: 365-372. Columbia Basin Trust. 2012. Columbia River Treaty. Website < http://www.cbt.org/crt> accessed 22 March 2012. Comeau, L.E.L., Pietroniro, A., and Demnuth, M.N. 2009. Glacier contribution to the North and South Saskatchewan Rivers. Hydrological Processes, 23: 2640-2653, doi: 10.1002/hyp.7409. Daly, C., Gibson, W.P., Taylor, G.H., Johnson, G.L. and Pasteris, P. 2002. A knowledgebased approach to the statistical mapping of climate. Climate Research, 22: 99-113. DeBeer, C.M. and Sharp, M J. 2007. Recent changes in glacier area and volume within the southern Canadian Cordillera. Annals o f Glaciology, 46: 215-221. Dery, A.J., Stahl, K., Moore, R.D., Whitfield, P.H., Menounos, B. and Burford, J.E. 2009. Detection of runoff timing changes in pluvial, nival, and glacial rivers of western Canada. Water Resources Research, 45, doi: 10.1029/2008WR006975. Dyurgerov, M.B. and Meier, M.F. 2000. Twentieth century climate change: Evidence from small glaciers. Proceedings o f the National Academy o f Sciences, 97(4): 1406-1411. Environment Canada. 2010. National Climate Data and Information Archive. Website c http://www.climate.weatheroffice.gc.ca/Welcome e.html> accessed 07 December 2012 . 100 Fleming, S.W. and Clarke, G.K.C. 2003. Glacial control of water resource and related environmental responses to climatic warming: empirical analysis using historical streamflow data from northwestern Canada. Canadian Water Resources Journal, 28(1): 69-86. Fleming, S.W. and Whitfield, P.H. 2010. Spatiotemporal mapping o f ENSO and PDO surface meteorological signals in British Columbia, Yukon, and southeast Alaska. Atmosphere-Ocean, 48 (2): 122-131. Fountain, A.G. and Tangbom, W.V. 1985. The effect of glaciers on streamflow variations. Water Resources Research, 21(4): 579-586. Fox, A.J. and Nuttall, A.M., 1997. Photogrammetry as a research tool for glaciology. Photogrammetric Record, 15, pp. 725-737. Furbish, D.J. and Andrews, J.T. 1984. The use of hypsometry to indicate long-term stability and response of valley glaciers to changes in mass transfer. Journal o f Hydrology, 30(105): 199-211. Gao, J. and Liu, Y. 2001. Applications of remote sensing, GIS and GPS in glaciology: a review. Progress in Physical Geography, 25(4): 520-540. Gardner, J.S. and Jones, N.K. 1985. Evidence for a Neoglacial advance o f the Boundary Glacier, Banff National Park, Alberta. Canadian Journal o f Earth Science, 22: 1751 1755. Giovannini, G. and Lucchesi, S. 1983. Effect of fire on hydrophobic and cementing substances of soil aggregates. Soil Science, 136(4): 231-236. Giovannini, G., Lucchesi, S. and Giachetti, M. 1988. Effect of heating on some physical and chemical parameters related to soil aggregation and erodibility. Soil Science, 146(4): 255-261. Granshaw, F.D. and Fountain, A.G. 2006. Glacier change (1958-1998) in the North Cascades National Park Complex, Washington, USA. Journal o f Glaciology, 52(177): 251-256. Greene, A. M. 2005. A time constant for hemispheric glacier mass balance. Journal o f Glaciology, 51(174): 353-362. Hamilton, A.S., Hutchinson, D.G. and Moore, R.D. 2000. Estimating winter streamflow using conceptual streamflow model. Journal o f Cold Regions Engineering, 14: 158175. Hamlet, A.F. and Lettenmaier, D.P. 1999. Columbia River streamflow forecasting based on ENSO and PDO climate signals. Journal o f Water Resources Planning and Management, 125(6): 333-341. 101 Hamon, W. 1961. Estimating potential evapotranspiration. Journal o f Hydraulics Division, Proceedings o f the American Society o f Civil Engineers, 871: 107-120. Harper, J.T. 1993. Glacier terminus fluctuations on Mount Baker, Washington, U.S.A., 1940-1990, and climatic variations. Arctic and Alpine Research, 25(4): 332-340. Hirshfield, F. 2010. The impact of climate change and harvest of mountain pine beedle stands on streamflow in northern British Columbia. Master of Science Thesis. University of Northern British Columbia. Houghton, J.T., Y. Ding, D.J. Griggs, M. Noguer, P.J. van der Linden, X. Dai, K. Maskell, C.A. Johnson, editors. 2001. Climate change 2001: the scientific basis. Contribution of Working Group I to the Third Assessment Report of the Intergovernmental Panel on Climate Change. Cambridge University Press, Cambridge and New York. 881 pp. Intergovernmental Panel on Climate Change. 2007. Climate Change 2007: Synthesis Report. Contribution of Working Groups I, II, and II to the Fourth Assessment Report of the IPCC. Geneva, Switzerland. 104pp. Jiskoot, H., Curran, C.J., Tessler, D.L. & Shenton, L.R., 2009. Changes in Clemenceau Icefield and Chaba Group glaciers, Canada, related to hypsometry, tributary detachment, length-slope and area-aspect relations. Annals o f Glaciology, 50, pp. 133143. Jost, G., Moore, R.D., Menounos, B., and Wheate, R. 2012. Quantifying the contribution of glacier runoff to streamflow in the upper Columbia River basin, Canada. Hydrological and Earth Systems Science, 16: 849-860. Doi: 10.5194/hess-16-8492012 . Kaab, A. and Vollmer, M., 2000. Surface geometry, thickness changes and flow fields on creeping mountain permafrost: automatic extraction by digital image analysis. Permafrost Periglacial Processes, 11, pp. 315-326. Kaser, G., Groshauser, M. and Merzeion, B. 2010. Contribution potential of glaciers to water availability in different climate regimes. Proceedings o f the National Academy o f Sciences, 107(47): 20223-20227. Kendall, M.G. 1975. Rank Correlation Measures. Charles Griffin, London. Keutterling A. and Thomas, A. 2006. Monitoring glacier elevation and volume change with digital photogrammetry and GIS and Gepatschfemer glacier, Austria. International Journal o f Remote Sensing, 27(19): 4371-4380. Koch, J., Menounos, B., and Clague, J. 2009. Glacier change in Garibaldi Provincial Park, southern Coast Mountains, British Columbia, since the Little Ice Age. Global and Planetary Change, 66: 161-178. 102 Konz, M. and Seibert, J. 2010. On the value of glacier mass balances for hydrological model calibration. Journal o f Hydrology, 385: 238-246. Krimmel, R.M. and Tangbom, W.V. 1974. South Cascade Glacier: The moderating effect of glaciers on runoff: Western Snow Conference, 42d Annual Meeting, Anchorage, Alaska, 1974, Proceedings, p. 9-13. Larsen, C.F., Motyka, R.J., Arendt, A.A., Echelmeyer, K.A. and Geissler, P.E. 2007. Glacier changes in southeast Alaska and northwest British Columbia and contribution to sea level rise. Journal o f Geophysical Research, 112, F01007, doi: 10.1029/2006JF000586. Leung, L.R. and Wigmosta, M.S. 1999. Potential climate change impacts on mountain watershed in the Pacific Northwest. Journal o f the American Water Resources Association, 35(6): 1463-1471. Lillquist, K. and Walker, K. 2006. Historical glacier and climate fluctuations at Mount Hood, Oregon. Arctic, Antarctic and Alpine Research, 38(3): 399-412. Lindstrom, G., Johansson, B., Persson, M., Gardelin, M., Bergstrom, S. 1997. Development and test of the distributed HBV-96 hydrological model. Journal o f Hydrology, 201: 272-288. Luckman, B.H., Harding, K.A., and Hamilton, J.P. 1987. Recent glacier advances in the Premier Range, British Columbia. Canadian Journal o f Earth Sciences, 24: 11491161. Luckman, B.H., 1998. Landscape and climate change in the central Canadian Rockies during the 20th century. The Canadian Geographer, 42 (4): 319-36. Mann, H.B. 1945. Non-parametric tests against trend. Econometrica, 13: 245-259. Mantua, N., Hare, S., Zhang, Y., Wallace, J.M., and Francis, R. 1997. A Pacific interdecadal climate oscillation with impacts on salmon production. Bulletin o f the American Meteorological Society, 78: 1069-1079. McKay, M.D., Beckman, R.J. and Conover, W.J. 1979. A comparison of three methods for selecting values of input variables in the analysis of output from a computer code. Technometrics, 21(2): 239-245. Mekis, E. and Vincent, L.A. 2011. An overview of the second generation adjusted daily precipitation dataset for trend analysis in Canada. Atmosphere-Ocean, 49(2): 163177. Miles, E.L., Snover, A.K., Hamlet, A.F., Callahan, B., and Fluharty, D. 2000. Pacific Northwest Regional Assessment: The impacts of climate variability and climate 103 change on the water resources of the Columbia River Basin. Journal o f the American Water Resources Association, 36 (2): 399-420. Mitchell, T.D. and Jones, P.D. 2005. An improved method o f constructing a database of monthly climate observations and associated high-resolution grids. International Journal o f Climatology, 25: 693-712. Moore, R.D. 1993. Application of a conceptual streamflow model in a glacierized drainage basin. Journal o f Hydrology, 150: 151-168. Moore, R.D., Fleming, S.W., Menounos, B., Wheate, R., Fountain, A., Stahl, K., Holm, K. and Jakob, M. 2009. Glacier change in western North America: influences on hydrology, geomorphic hazards and water quality. Hydrological Processes, 23: 4261. Mote, P.W. 2003. Trends in temperature and precipitation in the Pacific Northwest during the twentieth century. Northwest Science, 77(4): 271-282. Murdock, T.Q. and A.T. Werner. 2011. Canadian Columbia Basin Climate Trends and Projections: 2007-2010 Update. Pacific Climate Impacts Consortium, University of Victoria, Victoria, BC, 43 pp. Nandakumar, N. and Mein, R.G. 1997. Uncertainty in rainfall-runoff model simulations and the implications for predicting the hydrologic effects of land-use change. Journal o f Hydrology, 192:211-232. Nash, J.E. and Sutcliffe, J.V. 1970. River flow forecasting through conceptual models. Part I: a discussion of principles. Journal o f Hydrology, 10: 282-290. National Research Council. 2010. Green Kenue Version 3.1.55. Website < http://www.nrccnrc.gc.ca/eng/solutions/advisorv/green kenue index.html> accessed 7 January 2013. Nolin, A.W., Phillippe., J., Jefferson, A., and Lewis, S. L. 2010. Present-day and future contributions of glacier runoff to summertime flows in a Pacific Northwest watershed: Implications for water resources. Water Resources Research, 46, W12509, doi: 10.1029/2009WR008968. Oerlemans, J., Naderson, B., Hubbard, A., Huybrechts, Ph., Johannesson, T., Knap, W.H., Scheits, M., Stroeven, A.P., van de Wal, R.S.W, Wallinga, J., and Zuo, Z. 1998. Modelling the response of glaciers to climate warming. Climate Dynamics, 14: 267274. Payne, J.T., Wood, A.W., Hamlet, A.F., Palmer, R.N., and Lettenmaier, D.P. 2004. Mitigating the effects of climate change on the water resources of the Columbia River Basin. Climatic Change, 62: 233-256. 104 Pelto, M.S. and Hedlund, C. 2001. Terminus behavior and response time of North Cascade Glaciers, Washington, USE. Journal o f Glaciology, 47(158): 497-506. Pelto, M.S. and Reidel, J. 2001. Spatial and temporal variations in annual balance of North Cascade glaciers, Washington 1984-2000. Hydrological Processes, 15: 3461-3472. Raper, S.C.B. and Brathwaite, R J. 2009. Glacier volume response time and its links to climate and topography based on a conceptual model of glacier hypsometry. The Cryosphere,3: 183-194. Rood, S.B., Samuelson, G.M., Weber, J.K. and Wywrot, K.A. 2005. Twentieth-century decline in streamflows from the hydrographic apex of North America. Journal o f Hydrology, 306: 215-233. Rood., S.B., Pan, J., Gill, K.M., Franks, C.G., Samuelson, G.M., and Shepherd, A. 2008. Declining summer flows of Rocky Mountain rivers: Changing seasonal hydrology and probable impacts on floodplain forests. Journal o f Hydrology, 349: 397-410. Schiefer, E. and Gilbert, R., 2007. Reconstructing morphometric change in a proglacial landscape using historical aerial photography and automated DEM generation. Geomorphology, 88, pp. 167-178. Schiefer, E., Menounos, B. and Wheate, R. 2007. Recent volume loss of British Columbia glaciers, Canada. Geophysical Research Letters, LI 6503, doi: 10.1029/2007GL030780. Scott, D.F. 1993. The hydrological effects of fire in South African mountain catchments. Journal o f Hydrology, 150: 409-432. Seibert, J. McDonnel, J.J. and Woodsmith, R.D. 2010. Effects of wildfire on catchment runoff response: a modeling approach o detect changes in snow-dominated forested catchments. Hydrological Research, 41: 378-390. Shen, Z.Y., Chen, L. and Chen, T. 2012. Analysis of parameter uncertainty in hydrological and sediment modeling using GLUE method: a case study of SWAT model applied to Three Gorges Reservoir Region, China. Hydrology and Earth System Sciences, 16: 121-132. Solomon, S., Qin, D., Manning, M., Chen, Z., Marquis, M., Averyt, K., Tignor, M. and Miller, H. (eds.) 2007. Contribution of Working Group I to the Fourth Assessment Report of the Intergovernmental Panel on Climate Change. Cambridge University Press, Cambridge, United Kingdom and New York, NY, USA. 996. Stahl, K. and Moore, R.D. 2006. Influence of watershed glacier coverage on summer streamflow in British Columbia. Water Resources Research, 42. W06201, doi: 10.1029/2006WR005022. 105 Stahl, K.., Moore, R.D. and McKendry, I.G. 2006. The role o f synoptic-scale circulation in the linkage between large-scale ocean-atmosphere indices and winter surface climate in British Columbia, Canada. International Journal o f Climatology, 26: 541-560. Stahl, K., Moore, R.D., Shea, J.M., Hutchinson, D. and Cannon, A.J. 2008. Coupled modeling of glacier and streamflow response to future climate scenarios. Water Resources Research, 44, W02422, doi:10.1029/2007WR005956. Stewart, I.T., Cayan, D.R and Dettinger, M.D. 2004. Changes in snowmelt runoff timing in western North America under a ‘business as usual’ climate change scenario. Climatic Change, 62: 217-232. Stewart, I.T., Cayan, D.R. and Dettinger, M.D. 2005. Changes toward earlier streamflow timing across western North America. Journal o f Climate, 18: 1136-1155. Surazakov, A.B. and Aizen, V.B. 2006. Estimating volume change of mountain glaciers using SRTM and map-based topographic data. IEEE Transactions on Geoscience and Remote Sensing, 44(10): 2991-2995. Theil, H. 1950. A rank-invariant method of linear and polynomial regression analysis. Indagationes Mathematicae, 12: 85-91. Vincent, L.A., Wang, X.L., Milewska, E.J., Wan, H., Yang, F. and Swail, V. 2012. A second generation of homogenized Canadian monthly surface air temperature for climate trend analysis. Journal o f Geophysical Research, 111, D18110, doi: 10.1029/2012JDOI7859. Wang, T., Hamann, A., Spittlehouse, D., and Murdock, T.N. 2012. Climate WNA - HighResolution Spatial Climate Data for Western North America. Journal o f Applied Meteorology and Climatology, 61: 16-29. Wilcoxon, F. 1945. Individual comparisons by ranking methods. Biometrics, 1: 80-83. Woo, M., Thome, R. and Szeto, K. 2006. Reinterpretation of streamflow trends based on shifts in large-scale atmospheric circulation. Hydrological Processes, 20: 3995-4003. Doi: 10.1002/hyp.6590 World Glacier Monitoring Service. 2012. Website < www.geo.uzh.ch/microsite/wgms/ > accessed 31 October 2012. Xu, C.-Y. and Singh, V.P. 2002. Cross comparison of empirical equations for calculating potential evapotranspiration with data from Switzerland. Water Resources Management, 16: 197-219. Young, G.J. 1990. Glacier Hydrology. In Northern Hydrology Canadian Perspective. 1990. NHRI Science Report No. 1. Eds. Prowse, T.D. and Ommaney, C.S.L. Ministry of Supply and Services Canada, 308 p. 106 Yue, S. and Wang, C. 2004. The Mann-Kendall test modified by effective sample size to detect trend in serially correlation hydrological series. Water Resources Management, 18: 201-218. 107 Appendix A: Comparison of spliced and rehomogenized Blue River climate data Environment Canada’s monthly rehomogenized climate data (Mekis and Vincent, 2011; Vincent et al., 2012) for Blue River were analyzed for trends. The objectives of this additional analysis were to: • Compare the magnitude and significance of trends in the rehomogenized data to trends in the climate record from the combined BRA and BRN stations (spliced data). • Assess whether splicing data from the BRA and BRN stations introduced spurious trends in the time series. Trend identification and comparison Trend analysis was completed using monthly data from each climate record. To avoid complications with autocorrelation in trend detection, trends were identified using the modified Mann-Kendall test (95% significance level) that corrects the power of the trend test by adjusting the number of independent samples in the time series (Yue and Wang, 2004). Trend magnitude was calculated as the slope of the Kendall-Theil robust line. Significant trends in temperature are generally present in the rehomogenized and spliced datasets in the periods 1947-2006 and 1947-1985 (Tables A .l and A.2). However, additional significant trends are present in the rehomogenized temperature record. Trends in the rehomogenized temperature record have higher Kendall-Theil slopes than those present in the spliced temperature record (Tables A .l and A.2). Specifically, the trend in mean annual (summer) temperature is 28% (35%) greater in the rehomogenized data than the spliced data. 108 Significant trends in summer precipitation are present in both datasets (1947-2006 and 1947-1985; Tables A .l and A.2). The magnitude of significant trends is higher in the spliced record than in the rehomogenized record. Specifically, the trend in summer precipitation is 17% greater in the spliced data than the rehomogenized data. Summary Where significant trends are present in the spliced data, they are also present in the rehomogenized data. However, there are additional trends present in the rehomogenized data that are not present in the spliced data. It is unlikely that the method used to splice BRA and BRN station records together introduced false trends into the climate record. In summary, the magnitude of significant trends present in the rehomogenized data is similar or higher (lower) to the magnitude of significant trends present in the spliced temperature (precipitation) data. The results presented in Chapter 2 of this thesis may thus underestimate the magnitude of significant trends in the temperature time series. 109 Table A .l. Summary o f trends in the combined Blue River A and Blue River N stations for mean, maximum and minimum temperature (°C yr_1) and total precipitation (mm y r _1) time series. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold. 1947-2006 1947-1985 1985-2006 Tau________ KT Slope________Tau______ KT Slope_________ Tau________ KT Slope Mean Temperature Annual W inter Summer 0.285 0.294 0.101 0.024 0.036 0.005 0.148 0.126 0.001 0.020 0.027 0.000 0.160 0.212 -0.052 0.016 0.053 -0.013 -0.016 0.045 -0.016 -0.003 0.013 -0.003 0.026 -0.013 0.026 0.017 -0.009 0.017 0.269 0.207 0.358 0.051 0.057 0.048 0.195 0.229 -0.035 0.039 0.060 -0.007 0.134 -0.123 0.312 2.162 -1.578 0.229 0.169 0.195 6.146 2.800 5.542 Maximum Temperature Annual W inter Summer -0.047 0.250 -0.047 -0.010 0.034 -0.010 Minimum Temperature Annual W inter Summer 0.332 0.322 0.169 0.034 0.050 0.011 Total Precipitation Annual W inter Summer 0.156 -0.060 0.296 1.793 -0.522 2.664 4.169 Table A.2. Summary of trends in the rehomogenized Blue River mean, maximum and minimum temperature (°C yr_1) and total precipitation (mm y r"1) time series. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold. 1947-2006 1947-1985 1985-2006 Tau________ KT Slope________Tau______ KT Slope_________Tau________ KT Slope Mean Temperature Annual W inter Summer 0.406 0.380 0.266 0.037 0.050 0.017 0.306 0.233 0.069 0.030 0.042 0.014 0.148 0.150 -0.073 0.041 0.061 0.016 1.362 -0.434 0.039 0.050 0.160 0.212 -0.048 0.016 0.053 -0.013 -0.014 0.108 0.247 -0.039 0.009 0.033 -0.033 0.327 0.250 0.418 0.067 0.073 0.052 0.173 0.229 -0.039 0.025 0.060 -0.011 0.119 -0.045 0.263 3.112 -0.614 0.221 0.108 0.177 5.450 2.857 5.708 0.019 Maximum Temperature Annual W inter Summer -0.047 0.250 -0.047 0.017 0.033 Minimum Temperature Annual W inter Summer 0.373 0.365 0.233 Total Precipitation Annual Winter Summer 0.101 -0.035 0.228 2.200 110 3.389 Appendix B: Complete record of trends in the Blue River climate data Table B. 1. Summary of trends (°Cyr "') in the Blue River mean temperature data. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold. Annual W inter Summer JFM AMJ JAS OND 1947-2006 1947-1985 1985-2006 Tau________KT Slope________Tau______ KT Slope_________ Tau_______ KT Slope 0.285 0.148 0.020 0.024 0.160 0.016 0.294 0.036 0.126 0.027 0.212 0.053 0.005 0.001 0.101 0.000 -0.052 -0.013 0.200 0.318 0.073 0.178 0.060 0.050 0.190 0.021 0.077 0.013 -0.048 -0.033 0.046 0.000 -0.007 0.000 0.048 0.000 0.038 0.007 -0.146 -0.031 0.009 0.000 Table B.2. Summary of trends (°Cyr ) in the Blue River maximum temperature data. Significant (p < Oi trends (Kendall-Theil slope; KT) are shown in bold. Annual Winter Summer JFM AMJ JAS OND Tau -0.047 0.250 -0.047 0.190 -0.091 -0.005 0.049 1947-2006 KT Slope -0.010 0.034 -0.010 0.030 -0.024 -0.002 0.005 1947-1985 Tau KT Slope -0.016 -0.003 0.045 0.013 -0.016 -0.003 0.045 0.014 -0.178 -0.057 0.031 0.005 0.009 0.011 1985-2006 Tau KT Slope 0.026 0.017 -0.013 -0.009 0.026 0.017 -0.082 -0.042 -0.130 -0.083 0.048 0.030 -0.078 -0.067 Table B.3. Summary of trends (°Cyr ') in the Blue River minimum temperature data. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold. Annual Winter Summer JFM AMJ JAS OND 1947-2006 1947-1985 1985-2006 Tau________ KT slope________ Tau______ KT slope________ Tau_______KT slope 0.332 0.034 0.269 0.051 0.195 0.039 0.322 0.050 0.207 0.057 0.229 0.060 0.358 0.169 0.011 -0.035 -0.007 0.048 0.321 0.244 0.079 0.160 0.110 0.043 0.366 0.037 0.329 -0.171 0.048 -0.040 0.000 0.004 0.035 0.313 0.039 -0.006 0.098 0.020 -0.062 -0.014 0.056 0.015 111 Table B.4. Summary of trends (mm y r-1) in the Blue River total precipitation data. Significant (p < 0.05) trends (Kendall-Theil slope; KT) are shown in bold. Annual W inter Summer JFM AMJ JAS OND 1947-2006 1947-1985 1985-2006 Tau________ KT slope________Tau________ KT slope_________Tau_____ KT slope 0.134 0.156 1.793 2.162 0.229 6.146 -0.060 -0.522 -0.123 -1.578 2.800 0.169 0.296 0.312 4.169 2.664 0.195 5.542 -0.162 -0.696 -0.120 -0.911 0.108 1.262 1.645 0.359 0.213 0.165 1.392 1.597 0.171 1.221 0.212 2.400 0.078 2.100 -0.100 0.017 0.149 -0.929 0.238 2.350 112 Appendix C: Comparison of simulated streamflow generated by forcing the HBV-EC model with spliced and rehomogenized Blue River climate data The HBV-EC model was forced with daily rehomogenized climate data (Mekis and Vincent, 2011; Vincent et al., 2012) using the best performing model parameter set. The objectives of this analysis were to: • Compare streamflow simulated using the rehomogenized climate data to streamflow simulated using the combined Blue River A and Blue River North climate data (spliced data), and • Compare the streamflow simulated using the rehomogenized and spliced data to the observed Water Survey of Canada streamflow. Comparison o f time series generated from rehomogenized and spliced climate data Streamflow simulated using rehomogenized climate data was compared to streamflow simulated using spliced climate data. The simulated hydrographs were visually observed and simulated annual peak flows were compared. Goodness-of-fit statistics were also calculated (Table C .l) to compare the streamflow simulated with rehomogenized data to the streamflow simulated using spliced data. Goodness-of-fit statistics were calculated from 1947 to 2006, 1947 to 1968 and 1969 to 2006. Over the period of record, rehomogenized and spliced climate data resulted in similar simulated streamflow (RMSE = 1.43 mV; NSE = 0.99; AE = 2.1%). The two climate data sets simulated the most similar streamflow in 1996 (Figure C .l; RMSE = 0.19 mV1; NSE = 0.99; AE = 0.2%) and the least similar streamflow in 1968 (Figure C .l; RMSE = 2.45 mV1; 113 NSE = 0.98; AE = 8.9%). Upon visual observation of the hydrographs, there was a greater difference between streamflow simulated using the rehomogenized and spliced climate data from 1947 to 1968 than from 1969 to 2006. From 1947 to 1968, the use of rehomogenized climate data resulted in an average of 8% higher simulated peak flows in the summer and higher baseflow in the winter. From 1969 to 2006, there was less difference between both the winter baseflows and the annual peak flows (average of 0.4%) simulated using the two climate datasets. There is less difference between the streamflow simulated using the rehomogenized and spliced data sets during the summer months (May to September; Table C .l) than during the winter months (October to April; Table C .l). 114 8 Discharge ( m 8 o O c« o 0 100 200 300 Day of Year o Discharge 1n r V ) CO o o 0 100 200 300 Day of Year Figure C .l. Simulated hydrographs for the years with (top) most similar and (bottom) least similar simulated streamflow using spliced and rehomogenized climate data. Day o f year corresponds to the hydrologic year. 115 Table C .l. Performance statistics comparing streamflow simulated using the rehomogenized climate data to streamflow simulated using the spliced climate data. Statistics include root mean square error (RMSE), Nash Sutcliffe Efficiency (NSE) and relative volume error (%). Statistics were calculated for three periods in the simulated streamflow time series (1947-1006,1947-1968 and 1969-2006). Performance Statistic RMSE NSE RVE (%) Monthly RVE (%) Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec Period of Comparison 1947-2006 1.43 0.991 2.4 1947-1968 1.94 0.985 7.1 1969-2006 1.02 0.996 0.2 7.0 7.0 7.0 0.8 0.9 1.1 2.1 2.1 3.3 6.2 6.7 6.9 19.0 19.0 19.0 3.0 0.6 4.5 6.7 5.7 7.6 16.1 19.0 19.5 0.3 0.3 0.6 0.2 1.0 0.8 0.5 0.2 1.0 0.8 0.2 0.2 Model performance using the rehomogenized climate data Hydrographs generated using the rehomogenized and spliced climate data sets were compared to the observed Water Survey o f Canada streamflow record. Model performance using each climate data set was analyzed using various goodness-of-fit measures (Table C.2). Overall model performance in both the calibration and validation periods was similar between the two climate data sets. Both climate data sets result in higher model performance in the summer months than in the winter months. 116 Table C.2. Performance statistics comparing observed W ater Survey o f Canada streamflow to streamflow simulated using the rehomogenized and spliced climate data sets. Statistics include root mean square error (RMSE), Nash Sutcliffe Efficiency (NSE) and relative volume error (%). Performance Statistic______________________________ Period of comparison RMSE NSE RVE (%) Monthly RVE (%) Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec Calibration (1986-2005) Spliced Rehomogenized 5.21 5.35 0.901 0.906 1.1 1.6 52.1 37.6 5.8 4.1 5.2 0.7 2.9 2.4 1.3 3 16.8 40.1 53.9 39.2 6.4 4.1 6.6 0.5 2.4 2.1 0.6 2.8 17.7 41.6 Validation (1971-1985) Rehomogenized Spliced 5.185 5.180 0.905 0.905 5.1 4.9 50.0 48.5 25.4 6.8 0.6 6.1 1.3 5.6 4.3 0.4 27.6 57.3 49.4 47.9 24.9 6.9 1.3 6.4 1.1 4.9 3.2 2.1 27.6 57.0 Summary The spliced and rehomogenized climate data sets result in similar streamflow simulated in the summer months (RVE < 7.6% from 1947-1968 and RVE < 1 % from 19692006). As the main objective of this thesis to simulate and analyze summer streamflow as it relates to glacier change in the watershed, the major findings of this thesis do not seem to arise from the use of the spliced climate data. 117