HYDROLOGIC CONTROLS OF A WETLAND COMPLEX IN BRITISH
COLUMBIA’S
INLAND TEMPERATE RAINFOREST
by
Jeremy Edric Morris
B.Sc., Simon Fraser University, 2018

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
August 2021
© Jeremy Morris, 2021

Committee Members:
Supervisor: Stephen Déry, Ph.D.
Department of Geography, Earth, and Environmental Sciences
University of Northern British Columbia, Prince George, Canada
Committee Member: Darwyn Coxson, Ph.D.
Ecosystem Science and Management, Faculty of Environment
University of Northern British Columbia, Prince George, Canada
Committee Member: Jun Yin, Ph.D.
BC Ministry of Forests, Lands, Natural Resource Operations and Rural Development
Government of British Columbia, Prince George, Canada

ii

Abstract
Ancient Forest/Chun T’Oh Whudujut Provincial Park (AFP) is part of British Columbia’s
(BC) inland temperate rainforest that receives high total annual precipitation amounts
(>1000 mm) from orographic enhancement. AFP hosts a ~4 km2 complex of valley bottom
wetlands which is supported by this precipitation, a good portion of which falls as snow.
This study examines the hydrology of the wetland complex to determine the primary
water sources and their influence on flow paths and recharge. Water levels and
meteorological conditions were monitored for 2019 and 2020, and water samples were
collected for isotopic data during the 2020 snow free period. Kriging of water level data
revealed a northeastward nearly flat hydraulic gradient that shifted orientation during
wet and dry periods through the summertime. Rainfall amounts were above average at
700.2 mm and 676.4 mm while snow water equivalent varied at 762.8 mm (below average)
and 1082.8 mm (above average) for 2019 and 2020, respectively. The reduced snowpack
of 2019 yielded lower water levels through the summertime when compared to those of
2020. End Member Mixing Analysis (EMMA) of stable water isotope data indicates that
rainfall is not a significant enough recharge source to alter the wetland groundwater
composition, while snowmelt is likely the dominant source of input for the groundwater.
Cross correlation between rainfall and water level data indicates however that water
levels do respond to rainfall events with lag times ranging from 13 – 30 hrs. This
observation leads to the conclusion that rainfall serves to flush stored snowmelt
iii

generated water from the soils into the wetland complex, though not in a significant
enough volume to replace snowmelt as the dominant soil water source. Since previous
research has concluded that by 2050 as much as 50 % of the current snow contribution to
total annual precipitation will be replaced by rainfall, this study indicates that this shift
in precipitation regime may result in lower water levels in the wetland complex.
Numerical modelling of the relationship between wetland water levels and precipitation
phases and amounts would improve the understanding of the climate resiliency of valley
bottom wetlands in the Robson Valley.

iv

Table of Contents
Chapter 1 - Introduction ............................................................................................................ 1
1.1

Overview ........................................................................................................................ 1

1.2

Research Questions ....................................................................................................... 3

1.3

Thesis Structure ............................................................................................................. 4

Chapter 2 - Background ............................................................................................................. 5
2.1

What are wetlands? ....................................................................................................... 5

2.2

Function of freshwater wetlands in mountain regions............................................ 6

2.3

Existing Research........................................................................................................... 9

2.3.1

Previous studies ..................................................................................................... 9

2.3.2

Snowmelt Hydrology in Wetlands .................................................................... 10

2.4

Climate change and Hydrology in Northern BC .................................................... 11

Chapter 3 - Data and Methods ................................................................................................ 14
3.1

Study Site - Ancient Forest/Chun T’Oh Whudujut Provincial Park .................... 14

3.2

Hydrometric Data ....................................................................................................... 18

3.2.1

Piezometers and Loggers .................................................................................... 18

3.2.2

Fraser River Data.................................................................................................. 23
v

3.2.3
3.3

Climate Data ......................................................................................................... 23

Isotope Data ................................................................................................................. 24

3.3.1

Sampling process and sites ................................................................................. 24

3.3.2

Lab Analysis ......................................................................................................... 25

3.4

Analysis ........................................................................................................................ 26

3.4.1

Mapping the Groundwater Table ...................................................................... 26

3.4.2

Slim Creek and Fraser River Interpolation....................................................... 27

3.4.3

Hydraulic Gradient Calculation ........................................................................ 28

3.4.4

Cross Correlation ................................................................................................. 28

3.4.5

Isotope Interpretation .......................................................................................... 31

3.4.6

End Member Mixing Analysis ........................................................................... 31

Chapter 4 - Results .................................................................................................................... 33
4.1

Study Context .............................................................................................................. 33

4.1.1

Geomorphology and Wetland occurrence ....................................................... 33

4.1.2

Climate Conditions .............................................................................................. 36

4.2

Water level data ........................................................................................................... 41

4.3

Kriging .......................................................................................................................... 45
vi

4.4

Hydraulic Gradient ..................................................................................................... 47

4.5

Cross Correlation of Rainfall and Water Level ....................................................... 50

4.6

Stable Water Isotopes.................................................................................................. 54

4.7

End Member Mixing Analysis (EMMA) .................................................................. 56

Chapter 5 - Discussion ............................................................................................................. 59
5.1

Controls on Hydraulic Head and Gradient ............................................................. 59

5.1.1

Wetlands and flow paths .................................................................................... 59

5.1.2

Hydraulic Gradients in Wet and Dry Conditions ........................................... 61

5.2

Recharge Sources and Volumes ................................................................................ 64

5.2.1

Snowmelt and rainfall volumes ......................................................................... 64

5.2.2

Recharge sources .................................................................................................. 66

5.3

Implication Under Changing Climate ...................................................................... 71

5.4

Study Limitations ........................................................................................................ 73

Chapter 6 - Conclusions and Future Research..................................................................... 75
6.1

Summary ...................................................................................................................... 75

6.2

Future Work ................................................................................................................. 77

6.2.1

Hydrological Modelling ...................................................................................... 77
vii

6.2.2

Wetland Ecohydrology ....................................................................................... 77

6.2.3

Remote Sensing Dataset ...................................................................................... 78

6.2.4

Summary ............................................................................................................... 78

References................................................................................................................................... 79

viii

List of Tables
Table 4.1 Monthly precipitation data from Dome Mountain including the 14-year mean,
2018/2019 and 2019/2020 water year data and their respective ratios to the mean. ......... 40
Table 4.2 Peak correlation coefficients and their corresponding lags. The lags represent
water table response times to rainfall inputs where the peak correlation coefficient
exceeded the 95 % confidence interval ................................................................................... 53
Table 4.3 Two component EMMA of groundwater samples from piezometer B2 as a
mixture of rainfall and snowmelt ............................................................................................ 57
Table 4.4 Two component EMMA of creek water samples from West Creek as a mixture
of rainfall and snowmelt ........................................................................................................... 58

ix

List of Figures
Figure 1.1 The Robson Valley has an array of climate, hydrometric and snow pillow
stations in the area surrounding AFP. The valley bottom varies from 5 to 15 km in
width. ............................................................................................................................................. 3
Figure 2.1 Idealized wetland schematic. Source: http://patti-isaacs.com ............................ 6
Figure 2.2 Schematic of the major components of the carbon cycle (from Kayranli et al.,
2010) ............................................................................................................................................... 8
Figure 3.1 Hillshade DEM from LIDAR survey over the study area with ground surface
contours plotted at 2 m intervals. Included on this map are the locations of the
piezometer, creek stage and precipitation monitoring sites. ............................................... 15
Figure 3.2 Conceptual stratigraphic model of AFP (not to scale). Unit thickness
estimates are derived from regional borehole logs, while the wetland sediment types
included in the table were collected through field observations. The Area of Interest
(AOI) sits within the lowland area. The three dashed red lines refer to the relative
positioning of the three piezometer transects, which the included soils table refers to. . 16
Figure 3.3 Site map outlining AFP and the study area within including water level
(piezometers) and creek stage monitoring sites, the precipitation sampler site and the
Ancient Forest climate station. WCRK and ECRK refer to the West and East Creeks, the
two primary inflow creeks for the wetland complex. The teal dashed line indicates the
area over which hydraulic gradient calculations were conducted. .................................... 19
x

Figure 3.4 Water level data continuity used for this study collected with Odyssey
Capacitance Loggers during the 2019 and 2020 growing seasons. See Figure 3.3 for site
locations. ...................................................................................................................................... 20
Figure 3.5 Photo of PVC piezometers (bottom left) and Odyssey Capacitance Level
Loggers (top left) used in the study. Typical piezometer site with logger installation
(right). .......................................................................................................................................... 21
Figure 3.6 Schematic of terms used to calculate hydraulic head and water level relative
to ground surface. ...................................................................................................................... 22
Figure 3.7 Water sampling methods used. Left two images show equipment used for
pumping water from piezometers. Right image shows rainfall sampler. ......................... 24
Figure 3.8 Hypothetical cross correlation example (from Lee et al., 2006) displaying a
rainfall event with following rise in water table (top) and corresponding cross
correlation plot (bottom). .......................................................................................................... 30
Figure 4.1 LIDAR DEM with locations of various wetland types (excluding swamps
which are distributed amongst the marshes), and 2 m ground elevation contours. The
dashed blue circle indicates the position of a fan depositional feature. The yellow
diamonds indicate piezometers, with the labels A, B, and C referring to the primary
transects. ...................................................................................................................................... 35
Figure 4.2 Snow water equivalent, total precipitation, and air temperature from Dome
Mountain snow pillow for the study period. Each plot displays the 2018/2019 (red) and
xi

2019/2020 water years, and the 14-year mean of the datasets. The SWE data are plotted
at the daily timescale, while precipitation and air temperature data are plotted at
monthly timescales. ................................................................................................................... 38
Figure 4.3 Daily average water levels relative to ground surface (dashed black lines) for
the growing seasons from wetland piezometers for 2019 (top) and 2020 (bottom). Each
plot is split into transect within which each colour reflects a different piezometer of the
transect. Note different y-axis scales. ...................................................................................... 42
Figure 4.4 Differences in average monthly water level in piezometers between 2019 and
2020. Positive values indicate that water levels were higher in 2020 than 2019 (negative
indicates lower). Note that the y-axis varies by transect...................................................... 44
Figure 4.5 Map of kriged water level data from 1 July 2020. Light green values are
higher elevations, dark blue colours are lower elevations. Black dots represent the data
points used in the kriging process. Black lines are ground elevation contours in a 2 m
contour interval. White dashed lines are manually interpreted groundwater flow paths
based on the kriged groundwater table contours. ................................................................ 46
Figure 4.6 Bimonthly plots of kriged water table for the growing season of each year.
Data from 2019 are displayed on the top row, while 2020 data are along the bottom row.
Note that the number of data points used in the interpolation vary by plot due to data
loss during the study. Data points along the Fraser River and Slim Creek were
estimated by extrapolating single water level values to the length of the two rivers. .... 47
xii

Figure 4.7 Weekly hydraulic gradient (dh/dx) and Fraser River stage fluctuations
through the growing seasons of 2019 and 2020. .................................................................... 49
Figure 4.8 Cross correlation of rainfall (at AF station) and wetland water level data
including 95 % confidence levels (dashed blue lines) for 1 June to 1 November 2019. A
12 hour rolling mean was applied prior to correlation. Lags are in 15-minute
increments. Note the y-axis scale varies by site..................................................................... 52
Figure 4.9 Meteoric Water Line plot displaying results from stable water isotopic
analysis of samples from this study. Colours represent month of sampling in 2020.
Shapes represent different sample types. The Global Meteoric Water Line (GMWL) is
constructed here with the equation δ2H = 8.79× δ18O + 13.416. The Local Meteoric Water
Line (LMWL) here is defined by the line of best fit through the 2020 rainfall samples.
Grey samples are from Knudsvig (2015). ............................................................................... 56
Figure 5.1 Map displaying the spatial occurrence of wetland in light of ground
elevations, estimated surficial flow paths and the hydraulic gradient from 26 August
2019. Forested swamps constitute much of the unmapped areas between marshes. ...... 61
Figure 5.2 Water table maps from kriging representing the water table after dry periods
(left column) and the following wet periods (right column). .............................................. 63
Figure 5.3 Temporal plot of stable water isotope data during 2020. Note that the
snowmelt sample in each plot is the same sample collected in March, repeated for
comparison purposes. ............................................................................................................... 69
xiii

Figure 5.4 Stable water isotopic data (first two rows) with daily meteorological data
from the AF weather station, including air temperature, rainfall and snow depth in
2020............................................................................................................................................... 71

xiv

Acknowledgements
This project would have not been possible were it not for the expertise and support from
my supervisor Dr. Stephen Déry, and committee members Dr. Darwyn Coxson and Dr.
Jun Yin. Were it not for Dr. Coxson’s biological interest in the area, Dr. Yin’s expertise in
groundwater investigations, and Dr. Déry’s tireless support for my scatter brain, this
project would not have happened. Lastly, thank you to Dr. Jason Leach of NRCAN for
providing external examination of this thesis.
Thanks to the funders of this project: BC Parks through the Living Labs program, Natural
Sciences and Engineering Research Council of Canada through the Discovery Grant
program and UNBC for the Research Project Award.
Next those that provided logistical support: BC Ministry of Forests, Lands, Natural
Resource Operations and Rural Development (FLNRORD) and Water Survey of Canada
(WSC) for lending equipment and time, and for creating high-quality, publicly available
datasets. An enormous thank you to the Hakai Institute and Dr. Brian Menounos for
fitting this project into a busy LIDAR flight plan and doing so at a discounted rate.
Additional thanks to the InnoTech Alberta laboratory for affordable and rapid stable
water isotope analysis.
Thanks to my friends, family and colleagues who volunteered along the way: Hunter,
Selina, Hadleigh, Aseem, Kelly, Ivy, Hilary, Robin, Harlan, Margot, Danny, Natalya,
xv

Pierrick. I’m not sure they knew what they were getting into with this but now they’ve
all experienced the joys of bushwhacking through swamps.
Finally, I would not have been able to enter this Master’s program were it not for the
unyielding support from my partner Kirsti. Life continued throughout this project but
her support kept me going through thick and thin.

xvi

“How do these catchments store water for weeks or months, but then release it in minutes or
hours in response to rainfall inputs?” – Dr. James W. Kirchner

xvii

Chapter 1 - Introduction
1.1 Overview
Wetlands comprise 14% of Canada’s land area (National Wetlands Working Group, 1988)
and provide essential ecosystem services for both human and natural functions
(Richardson, 1994). For example, wetlands function as biodiversity hotspots by providing
unique habitats for both plant and animal species (Tews et al., 2004). Wetlands are also
effective terrestrial sediment (Phillips, 1989), carbon (Chmura et al., 2003)

and

contaminant (O’Geen et al., 2007) sinks. These storage abilities are increasingly being
utilized by engineers as constructed wetlands, and are becoming common solutions to
treat wastewaters and remediate damaged ecosystems (Ketcheson et al., 2016; Vymazal,
2010). With increasing global surface temperatures expected to surpass 1.5°C by 2040
compared to preindustrial times (IPCC, 2018), climate change poses a risk of significantly
altering the distribution and ecosystem services of wetlands within Canada, especially
those in northern latitudes where climate change is amplified (Dawson et al., 2003).
The Ancient Forest/Chun T’oh Whudujut Provincial Park (AFP) of the Robson Valley
hosts part of British Columbia’s (BC’s) inland temperate rainforest including a complex
of wetlands on the banks of the Fraser River (Figure 1.1). Wetlands in the region can
become inundated during the spring as the melting of seasonal snowpacks occurs, and
fall rainstorms once again provide water inputs prior to snowfall and snow cover during
1

winter. This seasonal water influx provides sufficient water across the complex to
support an array of rare and at-risk aquatic species (D. Coxson, personal communication,
20 October 2019). It has yet to be determined whether this water influx occurs
predominantly by delivery via stream inflows, lateral groundwater influx from Slim
Creek or the Fraser River, or by direct rainfall input or in situ snowmelt. With no prior
research on the hydrology of wetlands in the Robson Valley, and the fact that wetlands
have been found to either enhance or inhibit processes such as groundwater recharge, or
downstream discharge (Bullock & Acreman, 2003) there is no existing knowledge in the
area to draw upon. Moreover, the seasonal snowpacks of the Cariboo Mountains
surrounding the Robson Valley are anticipated to decline with climate change. Therefore,
it is imperative to understand the hydrologic function of spring snowmelt events in the
wetlands to anticipate downstream impacts of reduced meltwater input to the wetlands
and subsequently to the Fraser River (Déry et al., 2014; Islam & Déry, 2017). Thus,
developing an understanding of how influential snowmelt runoff is for the wetland, the
resultant flow patterns and water residence time within the wetland is important to
foresee ecosystem change as a result of climate change.

2

Figure 1.1 The Robson Valley has an array of climate, hydrometric and snow pillow stations in the area surrounding
AFP. The valley bottom varies from 5 to 15 km in width.

1.2 Research Questions
This research will address the following questions:
1. What are the primary drivers of water level fluctuations in the wetland complex
during the snow free period (April-November)? Is this system more dependent on
precipitation, streamflow, or groundwater?
2. How reactive is the wetland hydrological system to varying inputs? Do water
levels respond to individual rain events?
3

3. How do fluctuations in inputs affect the water levels or the hydraulic gradient?

1.3 Thesis Structure
This thesis is organized as follows: Chapter 2 provides background through a review of
previous research on wetlands, their hydrology, and climate change impacts on these
systems. Chapter 3 includes a description of the study site and explores the methods used
to collect hydrometric and geochemical data, and the methods behind the various
analyses. Chapter 4 outlines the results from the collection and analysis of the
aforementioned datasets. Chapter 5 provides a holistic discussion of the results in light
of one another to address the defined objectives. Lastly, Chapter 6 concludes the thesis
with the major facts learned and how this information can guide future research in AFP.
The main period of interest is between the months of April and November, referred to
throughout the thesis as the growing season, for the years of 2019 and 2020.

4

Chapter 2 - Background
2.1 What are wetlands?
Wetlands are areas that remain wet for long durations due to high water tables resulting
in the development of hydric soils, hydrophyte plants and biological processes that are
tolerant of wet and anaerobic conditions (Figure 2; National Wetlands Working Group,
1997). It is well established that wetlands perform various hydrologic functions including
the modification of flood volumes, dry season runoff, and groundwater recharge rates
(Bullock & Acreman, 2003). In Canada, wetlands comprise 14% of land area (National
Wetlands Working Group, 1997) and provide ecosystem services for both human and
natural functions (Richardson, 1994).
Wetland types vary based on their substrate and their hydrogeomorphologic setting.
National Wetlands Working Group (1997) identifies two broad categories of wetlands as
being either organic (wherein organic soils have developed on peat) or mineral (mineral
soil dominated with limited or no peat development). In general terms, bogs, fens and
peat rich swamps are organic wetlands whereas mineral wetlands consist of marshes and
swamps dominated by trees (National Wetlands Working Group, 1997). These different
wetland classes are typically associated with disparate hydrological systems. In
freshwater (non-tidal) systems, marshes and swamps can receive water from
precipitation, groundwater or surface waters, whereas fens are more frequently
5

associated with groundwater seepage and bogs can be solely reliant on precipitation
(Mitsch & Gosselink, 2007; National Wetlands Working Group, 1997). This variety in
water sources results in different soil, nutrient, and vegetative conditions for each
wetland class.

Figure 2.1 Idealized wetland schematic. Source: http://patti-isaacs.com

2.2 Function of freshwater wetlands in mountain regions
Mountainous areas can host freshwater wetlands that develop due to landscape position
among glacial and fluvial landforms such as kettles or terraces near valley bottoms. Aside
6

from being biodiversity hotspots, wetlands provide a plethora of ecosystem services;
acting as sinks for organic material, plant detritus and sediments due to slow rates of
decomposition and water flow (Mitsch & Gosselink, 2007). The variety of redox processes
that can occur in wetlands (Figure 2.2) means many of their services cannot be
generalized. Wetlands store a significant portion of global carbon with wetlands of the
United States (US) alone storing 11.52 gigatonne of carbon (four years’ worth of total US
emissions); over 90% of this carbon being stored in freshwater wetlands specifically
(Nahlik & Fennessy, 2016). However, due to the variety of hydrogeomorphic and
biogeochemical conditions in which they occur, freshwater wetlands can be sources or
sinks of carbon dioxide and methane (Kayranli et al., 2010; Mitsch et al., 2013). A study
of wetlands located on glacial sediments in Michigan, US found wetlands with
groundwater input resulted in higher rates of methanogenesis and sulphate reduction
compared to wetlands reliant on precipitation (Whitmire & Hamilton, 2008).

7

Figure 2.2 Schematic of the major components of the carbon cycle (from Kayranli et al., 2010)

Wetlands have also been found to modify discharge of the water to downstream areas as
well as regional climates. Glenn & Woo (1997) found that a valley bottom fen on
Ellesmere Island, NWT acted as a conduit for snowmelt waters upon reaching its
saturation but produced negligible outflows during dry spells. Investigations of an
upland swamp with predominantly organic sapric soils overlying sands in the Budderoo
Plateau of SW Australia revealed the swamp stored rainwater during dry periods but
during wet periods could act as a rapid conduit for through flow and overland flow
(Fryirs et al., 2014). By retaining water through otherwise dry periods, wetlands can
modify local climates. In the alpine of the Zoige Plateau of China, wetlands were found

8

to cool and humidify local climates, and that their removal led to reduced precipitation
and surface runoff (Bai et al., 2013).

2.3 Existing Research
2.3.1 Previous studies
A previous study in the Robson Valley sought to determine the influence of precipitation
in the area on shallow groundwater levels and soil moisture. Knudsvig (2015) determined
that soil moisture data from upland sites around AFP had significant correlations with
rainfall following 3-8 days of lag time but only weak correlations between rainfall and
groundwater levels. Using stable water isotope values from snow and groundwater, the
study determined that 40% of the shallow groundwater in the area could be attributed to
snow waters (Knudsvig, 2015). Another similar project that took place in northern BC
sought to establish a water budget for the Coles Lake watershed to facilitate numerical
modelling. Abadzadesahraei (2018) deployed an array of tipping bucket rain gauges,
piezometers with water level loggers, and measured river discharge at sites above and
below Coles Lake over two growing seasons in 2013-2014. The study showed that water
budgets for these types of system can be reasonably modelled numerically
(Abadzadesahraei, 2018). Mercer (2018) conducted a water budget and isotope study on
wetlands in alpine environments over two years in Banff National Park, Alberta. The
research found that that regardless of a three-week difference in snowmelt timing
9

between the two years, the water table positions remained stable. Isotopic end member
mixing analysis during the study found that up to 53% of the wetland water was
attributed to groundwater contributions in one of the peatlands, suggesting that
groundwater influx provided some measure of climate resilience to the alpine wetlands
(Mercer, 2018).
2.3.2 Snowmelt Hydrology and Wetlands
Research on wetlands in the Boreal Plain region of Alberta has shown that the
downstream runoff from wetlands in summer months depends on the antecedent
saturation caused by spring snowmelt (Wells et al., 2017). Other studies from the
Canadian Prairies have illustrated the snowmelt storage capacity of wetlands and
significant water level responses to snowmelt input allowing for “fill and spill” processes
to allow summertime output of waters from wetland complexes (Hayashi et al., 2016; Liu
& Schwartz, 2012; Shaw et al., 2012). Snowmelt contributions to wetlands are therefore
well established to be significant drivers to fluctuations in water levels. The importance
of snowmelt will be even more significant within the AFP wetland complex as the site
receives substantially more precipitation than within the prairies. Of the 1200 mm of
average annual precipitation that occurs in the Robson Valley, 800 mm of it can occur as
snow in the region of the AFP (Déry et al., 2014). Annual snowfall within the prairie

10

research sites, where much of the existing research has taken place, can be as low as 100
mm (Ketcheson & Price, 2016).

2.4 Climate change and Hydrology in Northern BC
Climate change studies in the broader region have examined elevation-dependent air
temperature trends (Sharma & Déry, 2016) and the effects of reduced snowpacks on
streamflow timing (Kang et al., 2016) and snowmelt contribution to streamflow (Kang et
al., 2014) over the last half century. Further research has examined the future impacts of
climate change on the precipitation regimes and runoff volumes and timing within the
Fraser River Basin (FRB) (Islam et al., 2019).
Air temperature trend analysis by Sharma & Déry (2016) has shown that the Cariboo
Mountains Region has warmed by minimum and maximum air temperatures of 1.9°C
and 1.2°C, respectively over the period of 1950-2010, a faster rate than regionally and
globally, and that trends are amplified at higher altitudes and during wintertime, most
likely due to snow-albedo feedbacks. Kang et al. (2016) conducted Variable Infiltration
Capacity (VIC) modelling within the FRB of the period 1949 to 2006 showing that
although there is a lack of trend in total precipitation amounts, negative trends in
snowmelt contribution to runoff are a result of changes in precipitation phase as the nival
regime shifts towards a pluvio-nival regime. The largest decline in snow water equivalent
(SWE) (-200 mm) over the period was found in the Upper Fraser Basin (UF), driving
11

a -19% change in a snowmelt-runoff ratio (Kang et al., 2014). For the same time period,
further VIC modelling has shown that peak discharge of the Fraser River has advanced
by ten days, and that the UF and the Thompson-Nicola Basins are the primary drivers of
the Fraser River’s spring freshet (Kang et al., 2016). Additional modelling has applied
climate projections to assess precipitation and runoff regimes going into the future. VIC
modelling using inputs from Representative Concentration Pathway 8.5 climate data to
midway through the 21st century suggests that the negative trend in snowmelt
contribution to the hydrology of the Fraser River will continue as a 50% reduction in the
fraction of precipitation falling as snow is expected by 2050 (Islam et al., 2019).
Climate impacts on wetlands have been investigated in modelling and observational
studies. Modelling of wetland extent in the Mid-Atlantic Highlands of the US has shown
that under current climate change trajectories the spatial extent of wetlands in the region
will decrease, potentially leading to higher flooding hazards (Pitchford et al., 2012).
Further modelling using the Soil Water Assessment Tool conducted for the prairiepothole region of Canada suggests that when all wetlands were removed from the
landscape, streamflow would increase by 55% in winter and 16% in summer for the
region (Muhammad et al., 2018). Mercer (2018) found that in three alpine wetlands near
Banff, Alberta, water table elevation did not vary significantly between a dry year and a
wet year indicating that these systems could be more resilient to climate change than

12

expected. However, given the variability of geologic and climatic conditions of alpine
wetlands, these results are difficult to generalize. The proposed research within the
Robson Valley would therefore expand the existing research on climate resilience of
wetlands in mountain regions.

13

Chapter 3 - Data and Methods
The following sections provide a description of the study site, and list all the datasets
used in this study. The datasets are split into physical and chemical sections, within each
are descriptions of methods used in both data collection and analysis.

3.1 Study Site - Ancient Forest/Chun T’Oh Whudujut Provincial Park
The AFP is located approximately 115 km east of Prince George, BC in the Cariboo
Mountains (Figure 1.1) within the Rocky Mountain Trench (RMT) that defines the
boundary between the Columbia Mountain and Northern Rocky Mountain ranges.
Elevations in this area range from roughly 600 m (asl) at the base of the RMT to 3250 m
(asl) at Mt Sir Alexander roughly 50 km to the northeast. The RMT itself is a tectonically
generated valley caused by regional scale geologic faulting, which has resulted in a wide
valley bottom that can reach up to 15 km in width.
A fixed wing airborne laser scanning survey was conducted on 29 August 2019 by the
Hakai Institute to collect Light Detection and Ranging (LIDAR) point cloud data of the
study site within AFP. The survey was conducted using a Riegl VQ-580 infrared (1024
micron) laser scanner with dedicated Applanix POS AV Global Navigation Satellite
System inertial measurement unit. The aircraft flew at an average altitude of 1470 m with
a flight line separation of 500 m with the LIDAR unit operating at 400 kHz, resulting in

14

an average ground point density of 2.71 pts m-2. LIDAR point cloud data were processed
to create a high resolution, bare earth digital elevation model (DEM) (Figure 3.1).

Figure 3.1 Hillshade DEM from LIDAR survey over the study area with ground surface contours plotted at 2 m
intervals. Included on this map are the locations of the piezometer, creek stage and precipitation monitoring sites.

To interpret the geology of the study site, several well logs from within 20 km of the site
in the valley bottom were pulled from the BC well database (Figure 1.1). The well logs,
though not very detailed, provide the thickness of several surficial sediment units in the

15

bottom of the valley. In addition, the LIDAR DEM was used to interpret the relative
depositional sequence of units, in tandem with regional soil maps from Maxwell (1987).
Furthermore, soil sample descriptions were recorded during piezometer deployment
(discussed below). These data were compiled into a conceptual stratigraphic model of the
valley near the study zone (Figure 3.2).

Figure 3.2 Conceptual stratigraphic model of AFP (not to scale). Unit thickness estimates are derived from regional
borehole logs, while the wetland sediment types included in the table were collected through field observations. The
Area of Interest (AOI) sits within the lowland area. The three dashed red lines refer to the relative positioning of the
three piezometer transects, which the included soils table refers to.

16

The AFP is bordered by the ridgetop of Driscoll Ridge to the southwest and the Fraser
River to the northeast. At the base of Driscoll Ridge is a glaciolacustrine terrace up to 150
m thick of silt and sand likely deposited in Glacial Lake Fraser following the last glacial
maximum (Clague et al., 2021). Following deglaciation, the modern Fraser River incised
into the glaciolacustrine sediment forming the terraces and leaving lowland floodplains
where the river likely formed and abandoned different meanders and backwater
floodplains. During this process a mixture of fluvial and backwater or oxbow lake
sediment deposited in a layer up to 50 m thick. In present day, these lowlands are where
wetland complexes occur and drain into the Fraser River. A variety of wetland types
occur along the lowlands, including black spruce bogs, forested swamps, fens, and sedge
marshes.
The area of interest for this study is bordered by Slim Creek to the southeast, the Fraser
River to the northeast, and a slope break between the terrace and the lowlands (Figure
3.1). The lowland area comprises a wetland complex including forested swamps, fens,
marshes, and bogs. Two creeks flow into the complex from an overlying terrace, which
receive runoff from several tributaries that originate on Driscoll Ridge. As water flows
into the wetlands, topography flattens substantially leading to a loss in flow velocity.
Eventually this causes the creeks to lose their channel definition and the water begins to
flow in a distributed manner (thus formation of forested swamps). Several beaver ponds

17

lie within the flow path, and around their shores lie the marshlands. Surface water
eventually flows out to the Fraser River or Slim Creek through several small channels at
the base of the wetland complex. This area was chosen for study due to the existence of
several long-term ecological monitoring plots located within the wetland complex and
out of interest of determining the nature of water supplies to the system.

3.2 Hydrometric Data
3.2.1 Piezometers and Loggers
A network of 18 water level monitoring sites was deployed using Odyssey Capacitance
water level loggers (Figure 3.3). The data loggers were configured to record water levels
at a 15-minute interval. Fifteen 1’’ PVC piezometers with measurement points ~1.2 m
below natural ground surface were installed. Twelve were deployed in a grid of three
transects, with two additional piezometers installed at the bases of East and West Creeks,
and one installed between Slim Creek and the main grid. An additional three stream stage
monitors were established in East Creek, West Creek and Slim Creek, though only the
data from the Slim Creek logger were used in analysis. The three piezometer transects
that make up the grid are labelled “A”, ”B”, and “C” from ordered from the bottom end
to the upper end of the wetland complex (Figure 3.3).

18

Figure 3.3 Site map outlining AFP and the study area within including water level (piezometers) and creek stage
monitoring sites, the precipitation sampler site and the Ancient Forest climate station. WCRK and ECRK refer to the
West and East Creeks, the two primary inflow creeks for the wetland complex. The teal dashed line indicates the area
over which hydraulic gradient calculations were conducted.

These deployments were staggered over the two-year field campaign and the respective
time series are summarized in Figure 3.4. During site visits, manual water level

19

measurements were taken using a Heron Instruments Dipper-T water level tape with a
precision of 0.5 mm to quality check the level logger data.

Figure 3.4 Water level data continuity used for this study collected with Odyssey Capacitance Loggers during the
2019 and 2020 growing seasons. See Figure 3.3 for site locations.

The Odyssey water level loggers consist of a conductor wire surrounded by a Teflon
insulator (Figure 3.5). As the wire and insulator are submerged in water (e.g.

20

groundwater in a piezometer), capacitance builds across the insulator and the logger
records the capacitance. A linear relationship can therefore be established between the
capacitance and the water level in a piezometer. Prior to field deployment, calibration
data were collected for each level logger by submerging them to known depths and
recording the resultant capacitance value to apply to a linear equation. The resulting
instrument precision is 0.8 mm water depth.

Figure 3.5 Photo of PVC piezometers (bottom left) and Odyssey Capacitance Level Loggers (top left) used in the study.
Typical piezometer site with logger installation (right).

To tie this relative water level value to a known datum, the height between the top of the
level logger casing to the ground surface was recorded with a measuring tape once
loggers were installed in piezometers. The hydraulic head (ℎ) is subsequently obtained
via:
21

ℎ = 𝐶𝐸 − 𝐶𝑇𝑊

(1)

where 𝐶𝐸 represents the elevation of the top of casing above sea level and 𝐶𝑇𝑊 is the
measured depth of water from the top of casing as measured by the level loggers (Figure
3.6). Top of casing elevation was determined by adding the height of the top of casing at
a given piezometer to the elevation of the ground surface at that piezometer as
determined by the LIDAR survey.

Figure 3.6 Schematic of terms used to calculate hydraulic head and water level relative to ground surface.

22

Ground elevation data from the LIDAR survey allow for a common sea level datum to be
applied across the piezometer network. Where the term water level is used in this thesis,
it refers to depth of water relative to ground surface, calculated by subtracting the casing
to ground height from the CTW value measured by the level loggers. This term is used
instead of groundwater level as water levels in the piezometers were often above ground
surface.
3.2.2 Fraser River Data
River stage data from the Water Survey of Canada’s Real-Time Hydrometric Database
were used to estimate the water level of the Fraser River during the study. Daily average
data were acquired for the Fraser River at Hansard site (08KA004) for the duration of the
study period (Figure 1.1).
3.2.3 Climate Data
A climate station with a record from 2009 to present is located roughly 2.2 km southwest
of the wetland complex (Figure 3.3) in a 20 m diameter clearing amongst a grove of old
growth cedar trees. The climate station has a near continuous 15-minute record for
several atmospheric parameters including rainfall, snow depth, and air temperature (see
Hernández-Henríquez et al., 2018 for instrument models and precision).
Additionally, 23 km southeast of the wetlands lies a snow pillow station operated by the
BC Ministry of Environment & Climate Change Strategy (MOE) since 20 October 2006
23

located at Dome Mountain (Figure 1.1). In addition to SWE data generated by the snow
pillow, the site also measures cumulative precipitation, air temperature and snow depth
data at an hourly interval.

3.3 Isotope Data
3.3.1 Sampling process and sites
For stable water isotopic analysis, a water sampling program was adopted to collect
samples from Slim Creek, West Creek, and the Fraser River, along with groundwater and
surface ponding samples from piezometer site B2.

Figure 3.7 Water sampling methods used. Left two images show equipment used for pumping water from piezometers.
Right image shows rainfall sampler.

Groundwater samples were collected by hand pump or Solinst peristaltic pump. A single
well volume would be purged before collecting the sample. Samples were collected in 20
ml glass vials that were subsequently taped shut with electrical tape. Sample bottles were
24

triple rinsed with the sampling water prior to use and samples were filtered with 0.45 µm
glass fibre filters within 24 hours of collection. An approximation of snowmelt water was
sampled by collecting a bulk shallow groundwater sample at the end stage of snowmelt.
Rainfall was collected via an evaporation-proof precipitation collector (tube-dip-in-water
style) following the design of Gröning et al. (2012). Rainfall samples were retrieved
monthly.
3.3.2 Lab Analysis
Water samples were submitted to the InnoTech Alberta lab in Victoria, BC, for analysis
of stable water isotopes δ2H, and δ18O. Sample analysis was conducted on a Delta V
Advantage mass spectrometer with HDevice and GasBench II peripherals for δ2H and
δ18O analyses, respectively. For δ2H measurements, the samples were vapourized into H2
gas, which was then sampled for H2 isotopic composition. For δ18O, CO2 gas is first
injected into the sample to force an exchange of oxygen between the gas and the sample.
The gas is then sampled and analyzed for its oxygen isotope composition. Both isotope
measurements are then reported vs a Vienna Standard Mean Ocean Water (VSMOW) lab
sample analysis and reported as δ18O or δ2H vs. VSMOW to an accuracy of ± 1 permil and
± 0.2 permil for δ2H and δ18O, respectively as per methods described in Nelson (2000).
Additional data were incorporated into this study from the previous effort in the area by
Knudsvig (2015). This study collected samples from shallow piezometers at the AFP and
25

Lunate Creek weather stations (Figure 1.1), from several groundwater seeps near Lunate
Creek, and from the snowpack. Snowpack samples were collected using a federal snow
sampler to collect snow cores that were subsequently deposited into plastic bags. As the
samples melted, any excess air in the bags was removed to prevent sample
contamination.

3.4 Analysis
3.4.1 Mapping the Groundwater Table
To map the groundwater table across the wetland, ordinary kriging was conducted across
the piezometer sites using the R-programming language. First, the hydraulic head data
were averaged to a daily timescale and combined to a master dataset. This dataset was
then split into daily shapefiles, which included the site geometry data and the mean
hydraulic head for a given day. Kriging was conducted by following the methods
outlined in Hengl (2007) and requires weighting for neighbouring values. The weights in
this case were determined through a fitted semivariogram model. The data were first
detrended via a first order polynomial, and then used to compute a semivariogram cloud.
A sample semivariogram plot provided estimates for the Sill, Range and Nugget to be
used when fitting the model semivariogram. Selecting the appropriate theoretical
function to use while fitting the model semivariogram was done by running the model
with several different functions (e.g. spherical, exponential, gaussian) until the model
26

best fit the sample semivariogram. Once the fitted semivariogram was calculated for the
weights in the interpolation, ordinary kriging was conducted using the “krige” command
from the “gstat package” (Pebesma, 2004).
3.4.2 Slim Creek and Fraser River Extrapolation
Two geographic point datasets were created to facilitate extrapolation of the two
waterways bordering the wetlands. Both datasets were generated by tracing a line feature
along the lengths of the two rivers, converting the lines to points spaced at 100 m
intervals. The elevations of the points were extracted from the DEM. Upon generating the
point datasets, water level time series were joined to each river point feature and the
water level adjusted relative to the elevation of each point. At Slim Creek, the time series
data were generated from an in-situ water level logger (Figure 3.3). For the Fraser River,
the river stage time series was extracted from Water Survey of Canada data at Hansard
(Figure 1.1). The Hansard time series was then normalized to water levels observed on
the morning of 29 August 2019 when the Lidar flight took place, and the zeroed time
series was added to the elevations of each simulated point. These water level (head)
values were used in the kriging process.
Since both methods used elevation data from the DEM at each control point, the
resulting time series at each point fluctuates relative to the observed ground elevation at
that site, therefore uncertainty is minimized through the process. In using Fraser River
27

stage data from Hansard, 50 km downstream of the study area, uncertainty is introduced
in the magnitude and timing of the river stage fluctuation. The Bowron River forms one
substantial inflow between the study site and Hansard; therefore a future improvement
on this process could be achieved by removing the Bowron River influence from the
Fraser River time series prior to use in the extrapolation process.
3.4.3 Hydraulic Gradient Calculation
To estimate the hydraulic gradient (dh/dx) through the snow free periods of the study,
slope calculations were conducted on the kriged groundwater table maps. To
accomplish this the “terrain” command from the “raster package” (Hijmans, 2020) was
used to calculate mean slopes from the raster output from the kriging process at a
weekly time step. Slope rasters were generated from the water table rasters, mean
slopes were calculated from all values in the slope rasters. For this process, the area of
interest (AOI) was reduced to an area directly surrounding the piezometer grid to
minimize errors caused by topography (Figure 3.3).
3.4.4 Cross Correlation
To understand the reactivity of the wetland water levels to rainfall during the snow free
period, cross correlation was conducted between the rainfall data from the weather
station and the raw water level time series. Figure 3.8 displays a single hypothetical
rainfall event resulting in a rise in the corresponding water table. Cross correlation
28

provides a measure of the statistical significance of the correlation between two time
series and, in this case, the lag of the water level response to rainfall input. Correlations
at the 95 % confidence interval are deemed to be those exceeding the standard error (as
~2/√N, where N is the number of observations in the dataset; Lee et al., 2006). Peak
positive correlations that exceed the 95 % confidence interval reflect the response time
between two time series that are in phase. In Figure 3.8, this occurs at day 7 wherein the
water level is rising fastest following the rainfall input. Negative correlations are ignored
as this indicates the time series are out of phase.

29

Figure 3.8 Hypothetical cross correlation example (from Lee et al., 2006) displaying a rainfall event with following
rise in water table (top) and corresponding cross correlation plot (bottom).

Prior to cross correlating a 12 hour (96 timestamp) rolling mean was applied to both
datasets, and both time series were trimmed down to the period of 1 June to 1 November
2019 as this period removes any effect of potential snowmelt or snowfall on the wetland
water levels. The 12 hour rolling mean was conducted to filter instrument noise that
interfered with the analysis. The 12 hour period was selected following exploration of the
effects of different rolling mean periods on the outputs; less than or greater than 12 hours

30

would result in low correlations. The correlation itself was conducted using the “ccf”
command from the R-programming language.
3.4.5 Isotope Interpretation
Isotope data were plotted on a meteoric water line plot (δ18O vs δ2H) and compared with
results from Knudsvig (2015). A global meteoric water line (GMWL) was added to the
plot using the equation
𝛿 2 𝐻 = 8.0𝛿 18 𝑂 + 10 ‰

(2)

defined by Craig (1961). Additionally, a local meteoric water line (LMWL) was estimated
by plotting a line of best fit through the rainfall data points.
3.4.6 End Member Mixing Analysis
To determine proportions of input waters in the wetland groundwater, a two end
member mixing model analysis (EMMA) was undertaken. The analysis involves selecting
samples with end member compositions of some conservative tracer (in this case δ18O
and δ2H) and a mixed sample with some mid-range concentration of that tracer. For this
study, the groundwater samples taken from the wetland were considered to be the mixed
samples, while snow waters and rainfall were chosen as the end members. Samples from
Slim Creek and the Fraser River were excluded for this analysis as the hydraulic gradient
maps from kriging did not support any input of these water sources into the wetlands.
An equation was then set up as follows:
31

𝑓𝐺𝑊 = 𝑓𝑆𝑊 + 𝑓𝑅𝑁 = 1,

(3𝑎)

𝛿 18 𝑂𝐺𝑊 = 𝛿 18 𝑂𝑆𝑊 𝑓𝑆𝑊 + 𝛿 18 𝑂𝑅𝑁 𝑓𝑅𝑁 ,

(3𝑏)

𝑓𝑆𝑊 =

(𝛿 18 𝑂𝐺𝑊 − 𝛿 18 𝑂𝑅𝑁 )
(𝛿 18 𝑂𝑆𝑊 − 𝛿 18 𝑂𝑅𝑁 )

(3𝑐)

where 𝑓𝐺𝑊 (assumed to equal 1), 𝑓𝑆𝑊 , and 𝑓𝑅𝑁 are the fractions of groundwater, snow
water, and rainfall in the mixed sample, and 𝛿 18 𝑂𝐺𝑊 , 𝛿 18 𝑂𝑆𝑊 , and 𝛿 18 𝑂𝑅𝑁 are their 𝛿 18 𝑂
compositions, respectively. This process is carried out for several sets of samples where
the groundwater and rainfall samples share similar time periods. A single snowmelt
composition was used as it is assumed that upon reaching the groundwater system that
the snowmelt waters from that year would not significantly change in isotopic
composition.

32

Chapter 4 - Results
This chapter discusses the results from data collection and analysis. This discussion starts
with descriptions of the geomorphology of the area of interest as revealed by LIDAR and
local climate conditions during the study. Next is a description of water level time series
collected and the results of analysis involving this dataset starting with results from
kriging interpolation of the water level data and hydraulic gradient calculations,
followed by results from cross correlation analysis between water level and rainfall time
series. Lastly, the stable water isotope dataset is described prior to discussion of EMMA
using that dataset.

4.1 Study Context
4.1.1 Geomorphology and Wetland occurrence
Water flows by gravity from high elevations to low elevations with the rate of flow being
determined by hydraulic gradients (often controlled by local topography) and geologic
materials that water interacts with along its flow path. Thus, to understand how water
flows through an area it is necessary to investigate the local topographic and geologic
conditions in the study area. Below is a description of the geomorphology within and
surrounding the wetland complex to understand likely flow paths for water and the
sediments the water flows through. Also included is a brief description of locations of

33

different wetland types as these landscape features can also exert control over the flow of
water.
The LIDAR DEM shows that the wetland complex is flat to gently sloping with a subtle
depositional fan at the mouths of the two inflow creeks (Figure 4.1). In this upper part of
the wetland complex water from the inflow creeks retains some velocity, allowing the
creek water to maintain channel-forms that are visible in the DEM. In situ, these channels
have been observed to be typically less than 1 m in width and up to 1 m depth. These
channels persist until the slope reduces between transects C and B, where diffuse
overland flows and occasional open ponding occur. From this area towards the Fraser
River that the wetland complex broadly comprises forested swamp and marsh wetland
types.

34

A
B
C

Figure 4.1 LIDAR DEM with locations of various wetland types (excluding swamps which are distributed amongst
the marshes), and 2 m ground elevation contours. The dashed blue circle indicates the position of a fan depositional
feature. The yellow diamonds indicate piezometers, with the labels A, B, and C referring to the primary transects.

Soil compositions transition with the changing slope and wetland occurrence (Figure
3.2 includes a table of soil types). There is a general fining of sediments from the inflow
creeks to the base of the wetland complex near transect A. At the creek mouths,
sediments are generally a sandy loam with low fines (silt and clay) and dominantly fine
sands associated with low to moderate energy creek deposition. Laterally this likely
35

varies due to past channel migration along the subtle fan with some finer sediment
zones where overbank flooding has occurred. Near the base of the fan at transect C the
finer fraction increases to roughly 30-40% and the dominant soil is silt clay loam with
minor organic soil development. Transitioning to the main wetland area of transects B
and A there is a clear shift to organic soil dominance with peat occurring where
persistent overland water does not occur. This organic soil appears to overlie a dense
silty-clay unit that has evidence of reduction (gleying) and oxidation (rust colour)
occurring.
4.1.2 Climate Conditions
To understand wetland hydrology during the study period it is crucial to identify the
meteorological conditions that affect the supply of water to the system. The following
description is based on datasets collected by MOE at the Dome Mountain snow pillow
site as this station represents the regional precipitation (rain and snow), air temperature
and snowpack conditions. These climatic conditions are the drivers of the hydrologic
cycle within the study site, and therefore hold significant impact on the results of this
study.
Snow water equivalent (SWE) and precipitation during winter of 2018/2019 was drier
than average (Figure 4.2). A peak SWE of 762.8 mm (89.2 % of the average 864.4 mm)
occurred on 10 May 2019 (average peak occurs on 14 May) after which the snow pillow
36

(at 1772 m elevation) was snow free by 27 June 2019 (average snow off occurs 22 July).
Air temperatures during this winter were above average aside from a cold snap during
February caused by dominance of continental arctic air that resulted in an average
temperature of -17.1 °C for the month. Conversely, winter 2019/2020 was substantially
wetter with peak SWE of 1082.7 mm (127.5 % of average) occurring on 13 May 2020,
melting off the snow pillow by 17 July 2020. Air temperatures for the full water year of
2019/2020 were generally below average aside from November, December and
September.

37

Figure 4.2 Snow water equivalent, total precipitation, and air temperature from Dome Mountain snow pillow for the
study period. Each plot displays the 2018/2019 (red) and 2019/2020 water years, and the 14-year mean of the datasets.
The SWE data are plotted at the daily timescale, while precipitation and air temperature data are plotted at monthly
timescales.

38

Summer conditions for both 2019 and 2020 were relatively cool and wet. From June
through September precipitation was above average for both water years (Figure 4.2).
Summer 2020 was exceptionally wet with rainfall exceeding the 14-year mean (20 October
2006 – 31 December 2020) by at least 23.8% each month with June and July yielding
precipitation amounts of 159.6 % and 184.3 % of normal, respectively. June, July and
August 2019 also exceeded the monthly mean with amounts at 127.6 %, 129.9 % and 149.7
% of normal, respectively. These high precipitation amounts tapered off into the fall
months of 2019; however, in 2020 precipitation continued to exceed the 14-year mean into
the winter. Overall, the 2018/2019 water year yielded an average total precipitation of
1493 mm (102.9 % of the mean) whereas the 2019/2020 water year proved to be
substantially wetter than normal with a cumulative precipitation amount of 1759 mm
(123.7 % of the mean).

39

Table 4.1 Monthly precipitation data from Dome Mountain including the 14-year mean, 2018/2019 and 2019/2020
water year data and their respective ratios to the mean.

Totals

14 Year Mean

2018/2019

Month

(mm)

(mm)

% Mean

(mm)

% Mean

Oct

112.5

49.9

44.3

143.1

127.2

Nov

170.8

164.1

96.1

223.6

130.9

Dec

111.8

134.9

120.7

123.5

110.5

Jan

168.9

141.4

83.7

159.5

94.5

Feb

145.6

179.6

123.3

137.1

94.2

Mar

124.9

112.0

89.6

150.6

120.6

Apr

79.2

88.4

111.6

66.4

83.8

May

74.8

39.0

52.2

98.5

131.7

Jun

97.2

124.0

127.6

155.1

159.6

Jul

139.1

180.7

129.9

256.3

184.3

Aug

92.8

139.0

149.7

114.9

123.8

Sep

104.8

110.2

105.1

130.5

124.5

Total

1422.4

1463.0

102.9

1759.2

123.7

Total Rain

558.0

700.2

125.5

676.4

121.2

Total Snow

864.4

762.8

88.2

1082.8

125.3

40

2019/2020

4.2 Water level data
Climatic factors are main drivers of water level fluctuations. This dataset reveals the
interplay between water inputs and topography as controls on water flow, and what the
resultant effects to water levels are in various locations throughout the growing season.
Water levels across the wetland complex are near or above ground surface for the
duration of the study. Water levels generally rise during spring snowmelt before falling
to a relative base level for the remainder of the growing season (Figure 4.3). Variability
in the water level responses to snowmelt, precipitation, and drying can be attributed to
landscape position and sediment composition. Transect C is the highest in the wetland
complex and contains the coarsest sediments (silty clay loam) resulting in the greatest
variation in water levels. The lower transects (A and B) show very stable water level time
series owing to fine organic soil compositions and landscape positioning among marsh
type wetlands. Both transects lie within swamps dominated by either black spruce trees
(transect B) or brush (transect A) and are therefore commonly submerged in shallow
overland flow. Aside from spring snowmelt, fluctuations in water levels are a result of
rainfall or dry periods. As expected, water levels climb during wet periods, and fall
during dry spells.

41

Figure 4.3 Daily average water levels relative to ground surface (dashed black lines) for the growing seasons from
wetland piezometers for 2019 (top) and 2020 (bottom). Each plot is split into transect within which each colour reflects
a different piezometer of the transect. Note different y-axis scales.

42

Water levels were different between 2019 and 2020 (Figure 4.4). Generally, 2020 saw
higher water levels than 2019 with sites A2 and B4 being exceptions for part of the
growing season and A4 having lower levels for all of the 2020 growing season. The
magnitude of difference varies by site with most differences falling within 10 cm of the
previous year. Site C4 has particularly large differences with a maximum of 50 cm higher
water levels on average for the month of June in 2020 likely due to the coarser sediments
in the area. The higher water levels for 2020 are likely due to the higher than average
snowpack given that rainfall amounts and summertime air temperatures were similar
across the two years (Table 4.1).

43

Figure 4.4 Differences in average monthly water level in piezometers between 2019 and 2020. Positive values indicate
that water levels were higher in 2020 than 2019 (negative indicates lower). Note that the y-axis varies by transect

Water levels at sites A2, A4, and B4 are lower in 2020 than the year prior. This may be
due to variations in evapotranspiration between the two years, a greater reliance on
rainfall than initially thought, or unique microtopography allowing draining to occur.

44

4.3 Kriging
After determining the hydraulic head in various locations of the wetland complex it is
then important to map the groundwater table across the area. This is undertaken by using
the hydraulic head at each observation point (piezometer) and interpolating the water
table across the entire area of interest. The resulting maps of the water table reveal flow
directions for water across the wetland complex for various times of the year providing
insights into how topography and water inputs impact flow patterns in the system.
Interpolation by kriging of hydraulic head data across the wetland complex displays
hydraulic gradients across the area of interest (Figure 4.5 and Figure 4.6). In general,
groundwater flows northeastward toward the Fraser River from the base of the cutbanks.
The flow directions are partly radial from the top of the wetland complex as the two
primary outflows are to the north (direct to the Fraser) and west (to Slim Creek).

45

Figure 4.5 Map of kriged water level data from 1 July 2020. Light green values are higher elevations, dark blue colours
are lower elevations. Black dots represent the data points used in the kriging process. Black lines are ground elevation
contours in a 2 m contour interval. White dashed lines are manually interpreted groundwater flow paths based on the
kriged groundwater table contours.

Examining the groundwater table at several times of the growing season, the hydraulic
gradient is generally stable through the year (Figure 4.6). Note that the distribution of
data points varies through time. Topography clearly exerts a strong control on flow
46

gradients throughout growing seasons here. Some subtle fluctuations in water table
gradient are visible by comparing the water table at different times of the year likely due
to variations of precipitation conditions through the growing season though these are
most apparent at the top of the zone.

Figure 4.6 Bimonthly plots of kriged water table for the growing season of each year. Data from 2019 are displayed on
the top row, while 2020 data are along the bottom row. Note that the number of data points used in the interpolation
vary by plot due to data loss during the study. Data points along the Fraser River and Slim Creek were estimated by
extrapolating single water level values to the length of the two rivers.

4.4 Hydraulic Gradient
As the delivery of sediments, nutrients and water to the wetland complex and beyond
are controlled by the changes in water flow, an examination of how the hydraulic
gradient evolves through time is done.
47

The hydraulic gradient, or mean groundwater table slope, in the area of the piezometer
grid shows subtle changes through time (Figure 4.7). During the 2019 growing season,
the gradient ranged between 0.0102 (1.02 cm m-1) on 8 April and 0.0093 (0.93 cm m-1) on
8 July. The 2020 growing season saw a greater range in gradient from 0.0184
(1.84 cm m-1) on 14 April to 0.0083 (0.83 cm m-1) on 20 October. The hydraulic gradients
for both growing seasons are shallower than the 0.0694 (6.94 cm m-1) average slope of the
ground surface over the same area. Both years have steeper gradients in spring with
additional peaks in June 2019 and August 2020. Comparing the shifts in gradient with
the water level in the Fraser River shows a partial inverse relationship between the two.

48

Figure 4.7 Weekly hydraulic gradient (dh/dx) and Fraser River stage fluctuations through the growing seasons of
2019 and 2020.

In 2020 during the rising limb of the Fraser River, the hydraulic gradient diminishes as
the hydraulic head at the base of the study area climbs. The gradient remains low until
mid-July when the Fraser reaches the falling limb, causing the hydraulic gradient to
climb to another peak. July 2020 saw above average precipitation (Figure 4.2), in
additional to marking the end of the melt period on Driscoll ridge. Thus, the sharp
49

increase in hydraulic gradient is likely a result of compounding the increased surface
runoff to the wetland zone with the drop in Fraser River stage at the bottom of the
system. This relationship is less clear in 2019, however following a subdued peak river
stage in the Fraser, the hydraulic gradient does steepen during the initial sharp drops in
river level during the falling limb. Further analysis would be needed to determine other
events which impact the hydraulic gradient such as rainfall events or increased stream
inputs to the wetland.

4.5 Cross Correlation of Rainfall and Water Level
To understand how responsive the wetland water levels are to precipitation inputs, cross
correlation between rainfall and water levels was conducted. Cross correlation
determines whether two time series are correlated with one another and, if they are
correlated, whether the relationship has any time dependency. In this context, the cross
correlations between rainfall and water levels reveal how reactive the water levels in the
wetland complex are to rainfall inputs, thus determining whether rainfall plays a
significant role in subsurface recharge.
Cross correlation analysis of rainfall and water level data yielded mixed results. Figure
4.8 displays the cross correlograms between the AF station rainfall data and each water
level logger dataset with the 95 % confidence intervals. Correlograms for A1 and B4
appear to have asymptotic behaviour, and the correlation coefficients for A4 never cross
50

zero, therefore these sites are omitted from analysis here. This is likely caused by the fine
sediment compositions at these sites reducing infiltration capacity and thus preventing a
significant water level response to single rainfall events. Conversely, correlograms for
A2, B1, B2, C1, C2, C3 and C4 show asymmetry with significant negative correlations
before lag zero and clear peaks of significant positive correlations at positive lags.

51

Figure 4.8 Cross correlation of rainfall (at AF station) and wetland water level data including 95 % confidence levels
(dashed blue lines) for 1 June to 1 November 2019. A 12 hour rolling mean was applied prior to correlation. Lags are
in 15-minute increments. Note the y-axis scale varies by site.

Table 2 summarizes results from the correlograms, displaying the peak correlation
coefficients and the lag at which they occurred. Peak correlations for sites A2, B2, C1, C2,
C3 and C4 occur between 13 and 30 hours after rainfall events. The correlogram for B1
peaks at a lower lag of 3.75 hours. This relatively short lag time at B1 may indicate that
52

some annular spaces may exist around this piezometer and the annular spaces create
preferential paths to allow water to reach the groundwater table. Otherwise, the
variability of response times at this site likely reflects differing local sediment
compositions and the typical depth to water prior to rainfall as this would influence
infiltration timing.
Table 4.2 Peak correlation coefficients and their corresponding lags. The lags represent water table response times to
rainfall inputs where the peak correlation coefficient exceeded the 95 % confidence interval for 1 June to 1 November
2019.

Site

Lag (15 min)

Lag (hours)

Correlation Coeff.

A2

73

18.25

0.29

B1

15

3.75

0.18

B2

52

13.00

0.31

C1

88

22.00

0.25

C2

107

26.75

0.06

C3

105

26.25

0.12

C4

120

30.00

0.09

A limitation to this analysis is the correlation was carried out on the entire time series
between 1 June to 1 November 2019. This approach assumes that the response of water
levels to rainfall occurs on a uniform timespan. The antecedent soil conditions and
groundwater table elevation play a large role in the water level response rate to rainfall
53

events (Beiter et al., 2020). Lower water tables and dry vadose zone conditions would
likely result in higher lag times than is presented here, and the opposite if water tables
were high and vadose zones wet. Further study of the seasonal variation in peak lag
could determine how the water table response times would vary under differing
climate conditions.

4.6 Stable Water Isotopes
Although the above analysis provides information on how water flows through the
wetland complex and how it relates to precipitation inputs, it does not provide insight
into the provenance of the water entering the system. For this, stable water isotope data
were collected from 40 water samples collected around the study area. Stable water
isotopes, δ2H and δ18O, can be used as conservative tracers to determine sources of water
reaching a given point. For this process, EMMA was conducted. Additionally, comparing
stable isotopic ratios from samples collected at different sites and times of year provides
insights into how dominant water sources may change through time.
Figure 4.9 summarizes results from stable water isotope analysis from this study. Rainfall
samples range from -18.17 δ18O ‰, -146.73 δ2H ‰ to -9.22 δ18O ‰, -67.69 δ2H ‰. A local
meteoric water line (LMWL) defined by the line of best fit through the six rainfall samples
yields an equation of δ2H = 8.79× δ18O + 13.416. The dominant cluster of samples (Figure
4.9 inset) sits within -16.33 to -19.84 δ18O ‰ and -126.04 to -148.89 δ2H ‰ consisting of
54

samples collected from river sites (Fraser River, Slim Creek, West Creek), groundwater
(piezometer B2), the pond, and snowmelt. Within this cluster, the groundwater samples
are generally more enriched in heavy isotopes relative to the other sample types apart
from three surface water samples (West Creek, overland, pond) collected in October. The
three October samples are all from the local surface waters and are likely enriched in
heavy isotopes due to evaporation. The isotope sample results match reasonably those of
Knudsvig (2015; grey samples in Figure 4.9) with the groundwater samples of that study
plotting among the dense cluster of 2020 samples. The snow samples from that study are
more depleted in heavy isotopes than that of the 2020 snowmelt sample as these samples
were not exposed to the same degree of evaporative processes.

55

Figure 4.9 Meteoric Water Line plot displaying results from stable water isotopic analysis of samples from this study.
Colours represent month of sampling in 2020. Shapes represent different sample types. The Global Meteoric Water
Line (GMWL) is constructed here with the equation δ2H = 8.79× δ18O + 13.416. The Local Meteoric Water Line
(LMWL) here is defined by the line of best fit through the 2020 rainfall samples. Grey samples are from Knudsvig
(2015).

4.7 End Member Mixing Analysis (EMMA)
Results from EMMA vary slightly between δ2H and δ18O (Table 4.3 and 4.4). Table 4.3
indicates a waning influence of rainfall on groundwater composition as the growing
season progressed in 2020. With δ2H as the tracer, Rain:Snow ratios in groundwater
decrease from 50 % (45 % for δ18O) in June to 12 % (14 %) in October. This result does not
reflect the ratio of input volumes given that the snowmelt influx occurs essentially as

56

plug flow in May, more so it reflects both the stability of groundwater composition
through the growing season in comparison to the continual enrichment of heavy isotopes
by evaporation in rainwater as the growing season progresses to the fall.
Table 4.3 Two component EMMA of groundwater samples from piezometer B2 as a mixture of rainfall and snowmelt

Groundwater

Rain:Snow

Snow:Rain

Month (2020)

δ2H

δ18O

δ 2H

δ18O

Jun
Aug
Oct

50.93%
35.60%
12.25%

44.49%
31.31%
13.97%

49.07%
64.40%
87.75%

55.51%
68.69%
86.03%

Table 4.4 shows that unlike the wetland groundwater, the creek water flowing into the
wetland through West Creek retains either a stable, snowmelt dominated composition,
or a slightly increasing rainfall influence on composition. Here the results vary dependent
on whether δ2H or δ18O is the tracer. δ2H EMMA results show a stable rain to snow ratio
between 15.20 % in June, 16.97 % in August, and 15.40 % in October. In contrast, δ 18O
EMMA shows an increasing trend with the rain to snow ratio increasing from 8.69 % in
June to 14.40 % in October. This result indicates an increasing influence of evaporation
on the creek water, either through direct evaporation from the creek, or by the addition
of evaporation affected rainwater via surface runoff or interflow.

57

Table 4.4 Two component EMMA of creek water samples from West Creek as a mixture of rainfall and snowmelt

Creek Water

Rain:Snow

Snow:Rain

Month (2020)

δ 2H

δ18O

δ 2H

δ18O

Jun
Aug
Oct

15.20%
16.97%
15.40%

8.69%
10.32%
14.40%

84.80%
83.03%
84.60%

91.31%
89.68%
85.60%

58

Chapter 5 - Discussion
This study aimed to describe the flow patterns of groundwater through the AFP wetland
complex, determine the reactivity of water levels to precipitation inputs, and reveal the
predominant sources of water for the wetland complex. Results presented in the previous
sections show a stable flow system exists through the growing season wherein snowmelt
is a dominant source of water, and that in general water levels have a 13 - 30 hr lagged
response to rainfall inputs. Additionally, the winter of 2019 held lower snowpack than
that of 2020, with similar summer rainfall amounts.

5.1 Controls on Hydraulic Head and Gradient
5.1.1 Wetlands and flow paths
Given that water flow and wetland presence are closely related, a discussion of their
relationship is required. Marshes and ponds typically occur perpendicular to the
groundwater flow direction (Figure 5.1), aside from where they appear along drainage
routes (such as at the eastern end of the study area). Ponds in the site develop due to
beaver activity and the marshes subsequently develop around their up-gradient edges.
This orientation matches the general shape of the groundwater table contours where the
marshes and ponds occur. This suggests that the ponds may alter the natural hydraulic
head and control the orientation of the flow field across the complex.

59

Most of the ponds identified have what appear to be historical beaver dam berm
features that can range to 1 m tall, and up to 100 m long with at least one active beaver
den identified in the area. No beaver dam activity was observed to alter the flow
patterns during this study; however, it is likely that continued beaver presence will
result in further hydrological changes.
The two bogs in the site emerge in areas outside of where the primary flow paths are.
This has likely played a role in their formation as peat development requires very stable
water tables and they thus often develop in flat areas with poor drainage. Swamps,
although not specifically mapped in Figure 5.1, occur throughout much of the wetland
complex surrounding the marsh areas. This is the most spatially dominant type of
wetland in the complex and they likely serve to transmit overland flow to down
gradient areas.

60

Figure 5.1 Map displaying the spatial occurrence of wetland in light of ground elevations, estimated surficial flow
paths and the hydraulic gradient from 26 August 2019. Forested swamps constitute much of the unmapped areas
between marshes.

5.1.2 Hydraulic Gradients in Wet and Dry Conditions
Water level plots and kriged water table maps show that water levels and flow patterns
are stable throughout respective growing seasons. 2020 experienced marginally higher
61

water levels along the C transect (higher in the gradient) than in 2019, whereas sites
within B and A transects (lower) had very stable water levels within the same elevation
range throughout the year. Comparing the water table surface to the ground contours on
Figure 5.1, topography exerts a strong control on flow patterns through the study zone.
Direct comparisons of the two study years are challenging due to the varying distribution
of instrumentation; however, hydraulic gradients across the wetland complex are
generally stable throughout growing seasons. Fluctuations in the groundwater table
correspond to differing precipitation inputs through the season (Figure 5.2) with some
modification by Fraser River levels. During dry periods, it appears that more eastward
gradients emerge possibly causing water to follow a shorter flow path to Slim Creek.
Conversely, wet periods result in longer flow paths in the northeastward direction to the
Fraser River.

62

Figure 5.2 Water table maps from kriging representing the water table after dry periods (left column) and the
following wet periods (right column).

In addition to the shifts in the hydraulic gradient in the wetland complex, cross
correlation analysis suggests that water levels are correlated with fluctuations in
precipitation but with time lags potentially lasting up to ~30 hours. This may indicate that

63

substantial rainfall events are needed to impact water levels as the lag time allows for
evapotranspiration to remove water prior to recharge.
The above discussion on wet vs dry groundwater table configurations and cross
correlation describes the celerity of wetland hydrologic systems (McDonnell & Beven,
2014). As the isotope EMMA results indicate, rainwater decreases in its influence on
groundwater composition shifting later into the year. This led to the conclusion that rain
events serve to flush stored snowmelt water through the watershed, as most of the water
in the wetland complex has flowed down to that location from higher in the watershed.
Thus, what the cross correlation results, and the wet/dry period comparisons describe,
are short term perturbations of a hydrologic system that covers a larger area and flows
over a longer temporal scale than that of the perturbations themselves. Thus, larger scale
factors such as the height of the Fraser River, or the volume of higher elevation snowpack,
likely exert higher degrees of control on the wetland watershed velocity while the rainfall
events cause perturbations in the wetlands that can be described as the wetland
watershed celerity.

5.2 Recharge Sources and Volumes
5.2.1 Snowmelt and rainfall volumes
Examining Figure 4.4 suggests that seasonal snowpacks influence wetland water levels
late into the growing season. Since 2019 and 2020 yielded similar amounts of rainfall
64

during the growing seasons (700 and 676 mm, respectively; Figure 4.2), the below average
snowpack of 2019 would have resulted in reduced snowmelt input to the wetland during
the springtime. In addition, the reduced snowpack, and higher air temperatures of 2019
during early spring resulted in a shorter melt period. In 2019 the snowmelt period started
on 10 May at Dome Mountain (Figure 4.2) and ended by 27 June. Conversely, the 2020
melt period was 15 days longer starting 13 May and ending on 17 July. The comparison
of water levels between the two years suggests that this reduction in snowmelt may have
resulted in lower water tables into the late summertime in 2019. Since a significant portion
of the annual precipitation, and therefore wetland source water, is composed of snowfall
(52 % and 62 % in 2019 and 2020, respectively) these data indicate that reduced winter
snowpack can result in lower wetland water levels through to the following winter.
During the snowmelt process, climate and site conditions can influence the amounts of
snowmelt water that become runoff or recharge (Corps of Engineers, 1998). Prior to
production of meltwater, snowpacks enter a “ripening” phase wherein loose and dry
snow densifies and becomes moist. The snowpack becomes ripe for snowmelt once it
reaches an isothermal state at 0 °C. At this stage, liquid water can flow through the
ripened snowpack, although the presence of ice lenses within the snowpack leads to nonuniform flow paths. The non-uniformity of flow paths through the snowpacks, and the
presence of frozen soil or ground ice result in variability in how much snowmelt water

65

flows vertically into the soil to become recharge or flows horizontally out of the snowpack
as runoff.
On several wintertime visits to the wetland complex, ice was encountered at the base of
the snowpack and ice lenses were present within the snowpack as a result of mid-winter
rainfall events. Ice at ground level is likely a result of overland flow waters freezing in
place before the onset of snow and would prevent in situ groundwater recharge where
present. The extent of the ground ice would be difficult to quantify, and its formation is
likely dependent on the average air temperature at the onset of snow. Thus, determining
the amount of in situ snowmelt recharge to the groundwater remains a challenge for this
study, however snowmelt generated runoff and groundwater from above the wetland
complex likely provide a significant flux to the complex. It has been previously observed
that soils do not typically freeze over winter due to thick and rapid onsets of snow in
early wintertime in non-wetland parts of the area. Therefore, it is very likely that during
snowmelt, ample runoff generated from upslope would reach the wetland complex to
provide recharge and soils along upper slopes would be saturated with snowmelt water
to be released to the wetlands as stormflow runoff throughout the snow free period.
5.2.2 Recharge sources
Isotope data suggest that rainfall serves to flush out stored soil water generated from
snowmelt. Creek water constitutes a partial mix of rain and snowmelt; however,
66

snowmelt remains the dominant source late into the growing season. Results from mixing
model analysis presented conflicting results with respect to the changing groundwater
sources (Table 4.3). An examination of the temporal shifts in isotope composition of
varying water types (Figure 5.3) reveals that groundwater isotopic ratios remain fairly
stable throughout the season, yet rainfall isotope ratios shift throughout the growing
season to become more enriched in heavy isotopes by evapotranspiration as the summer
progresses. The mixing model results suggest that groundwater becomes less influenced
by rainwater as the season progressed, which would seem counter-intuitive as
cumulative rainwater input would increase later into the season. Though rainfall
(enriched in heavy isotopes by evapotranspiration) is being added to the system, if the
soils between the ridge and the wetland are saturated with snowmelt (heavy isotope
depleted) water, then the addition of rainfall will flush the stored snowmelt water down
gradient to the wetlands. In theory, with enough rainfall this flushing would result in
rainwater replacing snowmelt water as the dominant water type in the soil. However, it
is possible that a combination of low infiltration rates due to fine sediment compositions,
and elevated evapotranspiration rates prevent rainfall from reaching groundwater
sampling depths before the rainwater is either evaporated away or flows away via diffuse
overland flow. It is also worth noting that plant transpiration activity likely increases into
the summer season resulting in more water uptake by vegetation later in the summer
than in the spring. Additionally, in late summer, vegetation is fully exposed (snow free)
67

and thick with foliage, likely intercepting a sizeable component of precipitation. Though
cross correlation did find a water level response time of 13 - 30 hr to rainfall events, the
amount of water level response in much of the wetland (transects A and B) is subdued,
suggesting that the rainfall may not provide as significant a recharge mechanism as the
snowmelt process.

68

Figure 5.3 Temporal plot of stable water isotope data during 2020. Note that the snowmelt sample in each plot is the
same sample collected in March, repeated for comparison purposes.

Examining the isotope data in the context of the meteorological conditions explains some
of the patterns in that dataset (Figure 5.4). Proceeding into the summertime, rainfall
samples clearly become enriched in the heavy isotopes. This correlates well with the
increasing air temperature, suggesting higher evaporation rates have removed lighter
69

isotopes first in transport. Note also that the river samples from West Creek (blue “+” in
Figure 5.3) become slightly more enriched in heavy isotopes later in the season. This
reflects the results of EMMA analysis of creek water (Table 4.4), which showed either an
increasing presence of rainwater in the creek as the season progressed or increasing
evaporation of water from the creek. Since surface waters can receive rainfall waters
through surface runoff and interflow, it follows that this EMMA result accurately
conveys an increasing influence of rain on the surface waters.

70

Figure 5.4 Stable water isotopic data (first two rows) with daily meteorological data from the AF weather station,
including air temperature, rainfall and snow depth in 2020.

5.3 Implication Under Changing Climate
These findings suggest that the AFP wetland complex has some reliance on the slow
release of snowmelt. Within the Fraser River Basin total annual precipitation is expected
to slightly increase, however the snow contribution to total annual precipitation is
expected to decline by 50 % by 2050 under RCP 8.5 (Islam et al., 2017). Due to the
71

severity of this shift in precipitation phase, it is likely that water levels within the
wetland complex will be affected. Islam et al. (2017) note that this reduction in SWE will
likely be greater at low to mid-elevations (under 1200 m), which could result in a
significant reduction in snow accumulations at the AFP since Driscoll Ridge is only 1803
m at its highest point with the wetland complex itself at ~640 m.
Though this study showed that reduced winter snowpack affected summertime water
levels in a year with below average winter precipitation, it does not directly address
reduced snowpacks where the loss of SWE is compensated by increased rainfall.
However, given that the snowmelt process occurs over a period of 1-2 months it is
likely that this slow release of water to the wetland complex allows for water to reach
the wetland over an extended period without being dried out by evapotranspiration.
Conversely, if the spring snowmelt were replaced by rainfall inputs, much of the rain
input would be reduced by evapotranspiration between rainfall events. Analysis
conducted by Berghuijs et al. (2014) indicates that catchments with high fractions of
precipitation falling as snow are associated with higher long term stream flows than
catchments that are not, and that a precipitation phase change from snow to rain results
in a significant decrease in streamflow. Though more research would be needed to
confirm these effects in the context of wetlands, it is likely that the shift from nival to
pluvial precipitation regimes will alter the hydraulic head fluctuations of the wetland

72

complex. Numerical modelling with downscaled climate projection data could confirm
whether the changing phase of precipitation in the region would significantly impact
valley bottom wetlands in the Robson Valley.

5.4 Study Limitations
This study had several limitations. A number of challenges were encountered with the
Odyssey Capacitance logger deployments. First a piezometer diameter of 1” potentially
caused instrument errors as the sensor cabling could contact the insides of the
piezometers resulting in erroneous measurements. Secondly, several battery failures
occurred during the study resulting in data loss and ultimate omission of several sites
from the analysis. This included a logger placed in West Creek with the intention of
building a rating curve for creek discharge data that was made impossible with the data
loss.
For the cross correlation analysis limitation lies in the choice to carry out the analysis on
the entire time series between 1 June to 1 November 2019. This approach assumes that
the response of water levels to rainfall occurs on a uniform timespan. In reality,
antecedent soil conditions and groundwater table elevation play a large role in the water
level response rate to rainfall events (Beiter et al., 2020). Lower water tables and dry
vadose zone conditions would likely result in higher lag times than is presented here,
and the opposite if water tables were high and vadose zones wet.
73

The limited number of samples used in EMMA introduced some uncertainty to the
results. The use of a single snowmelt sample in the analysis creates uncertainty due to a
lack of knowledge on sample representativeness to the area’s snowmelt waters more
broadly. The composition of snowmelt waters can vary spatially and temporally as a
result of many different processes such as melt and refreeze within the snowpack,
differential accumulation of snow over an area, sublimation and other processes related
to site properties like slope and aspect (Penna et al., 2019). If the snowmelt sample used
for EMMA was not representative of the broader snowmelt signal, then it could directly
skew the results leading to artificially high or low snowmelt (or rainfall) contribution
estimates in the groundwater.
The study would have benefitted from information on vertical hydraulic gradients within
the wetland to understand whether any groundwater upwelling occurred. It would have
also helped to have estimates on rates of overland flow across the swamps of the wetland;
however, this is notoriously challenging to accurately measure. Lastly, this study
occurred during two relatively wet growing seasons, limiting the ability to examine the
wetlands in wet and dry conditions. Hydrological investigations can always use longer
study durations for this reason.

74

Chapter 6 - Conclusions and Future Research
6.1 Summary
Using several different approaches this study represents a broad scale attempt at
characterizing the hydrology of the wetland complex in AFP. Beginning with the
description of the geomorphology and wetland occurrence it indicates that the wetland
complex lies on gently sloping to flat terrain comprising dominantly organic soils
overlying most likely flood plain deposits left by past meanders of the Fraser River. The
wetland complex sits below glaciolacustrine terraces, through which two channels
provide streamflow inputs sourced from Driscoll Ridge above. Outflow from the wetland
complex is facilitated by several small channels that drain to the Fraser River, however
much of the output likely also occurs diffusely via overland flow or through the
subsurface to Slim Creek.
Different wetland types are distributed across the study site; with a series of fens
occurring along the toe slopes of the surrounding terraces, swamps where slowing creek
inflows lose their channel form and flow diffusely through forested areas, marshes where
beaver activity has promoted ponding, and bogs where sufficient peat has developed to
prevent free standing water.
Climate conditions for the two study years were variable. 2019 shows a below average
winter snowpack with above average rainfall for the following summer and fall. 2020
75

shows above average precipitation for the entire year. This variation in precipitation
during the study length allowed for an examination of water levels. Section 5.2.1 shows
that the relatively deeper snowpack in 2020 likely resulted in persistently higher water
levels, that were not seen in water level loggers in 2019. This observation draws the
conclusion that growing season water levels depend on the winter snowpack, with
summer rainfall serving as a flushing mechanism for the stored snowmelt waters.
Kriging interpolation of water level data showed that flow paths and the hydraulic
gradient are strongly controlled by the flat topography resulting in a dominantly
northeastward flow direction with a low gradient. This is modified during wet vs. dry
periods for the site where drying events steepened the gradient higher in the system
resulting in flow directions shifting to rainwater more eastward directly to Slim Creek.
Rainfall was found to influence water levels and the hydraulic gradient during the
growing season. Cross correlation of water levels and rainfall found that in some sites
rainfall events would result in a water level response with significant lag times. Lag times
generally ranged from 13 to 30 hours (omitting the spuriously short 3.75-hour response
in B1). Results from isotope analysis suggest that rainfall may serve to flush stored
snowmelt derived soil water into the wetland complex, though further analysis should
be done to confirm this finding.

76

6.2 Future Work
This project provided a first attempt at characterizing the transport of water in the
wetland complex of AFP, but many aspects of this hydrological system could be better
understood with further investigation.
6.2.1 Hydrological Modelling
Hydrology research within AFP would greatly benefit from a modelling exercise by
improving our understanding of the dominant processes in the site and providing
insights into the local hydrology under future climate scenarios. Given that the
wetlands are reliant on snowmelt, and that snowpacks are declining in the region (Islam
et al., 2019), modelling could provide insight into how the shifting precipitation regime
impacts water delivery to the wetlands and growing season water levels. Long term
monitoring would likely be required as calibration data for a modelling exercise to be
possible. Further monitoring would ideally capture wet and dry years allowing for
comparisons
6.2.2 Wetland Ecohydrology
Since this project was initiated to better understand the wetlands due to the unique
ecosystems they support, an ecohydrological study would close the gap between this
study and the ecological research in AFP. Such a study could include researching how
the flow of water affects nutrient transport or how flow patterns relate to occurrence of
77

various vegetation types. More detailed tracer research could better reveal flow paths
and transit times for snowmelt waters generated on Driscoll Ridge to reach the wetland
complex as this would affect nutrient transport. Lastly, further research opportunities
exist here on the influence of beaver activity on the wetland hydrology, including the
impacts of dams on the hydraulic gradients.
6.2.3 Remote Sensing Dataset
In addition to expanding the hydrological investigations, the airborne dataset from Hakai
Institute presents several opportunities for further research. The LIDAR dataset could be
used for forest biomass surveys, while the 5-band orthomosaic imagery could allow for
high resolution ecosystem mapping across the complex to potentially delineate wetland
types across the complex with a higher degree of accuracy than presented here.
6.2.4 Summary
Expanding the hydrology research in AFP would yield a greater understanding of this
complex system which supports a diverse array of plant life and megafauna. Better
understanding the hydrological and ecosystem services of the wetland complex in AFP
can provide insights into similar systems within the Robson Valley and the broader
inland temperate rainforest including how resilient or at risk these ecosystems are in the
face of changing climate conditions.

78

References
Abadzadesahraei, S. (2018). Quantifying the water balance of two northeastern boreal
watersheds, British Columbia. [Doctoral dissertation, University of Northern British
Columbia]. https://doi.org/10.24124/2018/58834
Bai, J., Lu, Q., Zhao, Q., Wang, J., & Ouyang, H. (2013). Effects of alpine wetland
landscapes on regional climate on the Zoige Plateau of China. Advances in
Meteorology, 2013(1). https://doi.org/10.1155/2013/972430
Beiter, D., Weiler, M., & Blume, T. (2020). Characterising hillslope-stream connectivity
with a joint event analysis of stream and groundwater levels. Hydrology and Earth
System Sciences, 24(12), 5713-5744. https://doi.org/10.5194/hess-24-5713-2020
Berghuijs, W. R., Woods, R. A., & Hrachowitz, M. (2014). A precipitation shift from
snow towards rain leads to a decrease in streamflow. Nature Climate Change, 4(7),
583–586. https://doi.org/10.1038/nclimate2246
Bullock, A., & Acreman, M. (2003). The role of wetlands in the hydrological cycle.
Hydrology and Earth System Sciences, 7(3), 358–389. https://doi.org/10.5194/hess-7358-2003
Chmura, G. L., Anisfeld, S. C., Cahoon, D. R., & Lynch, J. C. (2003). Global carbon
sequestration in tidal, saline wetland soils. Global Biogeochemical Cycles, 17(4), 1111.
79

https://doi.org/10.1029/2002GB001917
Clague, J. J., Roberts, N. J., Miller, B., Menounos, B., & Goehring, B. (2021). A huge flood
in the Fraser River valley, British Columbia, near the Pleistocene Termination.
Geomorphology, 374, 107473. https://doi.org/10.1016/j.geomorph.2020.107473
Dawson, T. P., Berry, P. M., & Kampa, E. (2003). Climate change impacts on freshwater
wetland habitats. Journal for Nature Conservation, 11(1), 25–30.
https://doi.org/10.1078/1617-1381-00031
Déry, S. J., Knudsvig, H. K., Hernández-Henríquez, M. A., & Coxson, D. S. (2014). Net
snowpack accumulation and ablation characteristics in the Inland Temperate
Rainforest of the Upper Fraser River Basin, Canada. Hydrology, 1(1), 1–19.
https://doi.org/10.3390/hydrology1010001
Fryirs, K., Gough, J., & Hose, G. C. (2014). The geomorphic character and hydrological
function of an upland swamp, Budderoo plateau, southern highlands, NSW,
Australia. Physical Geography, 35(4), 313–334.
https://doi.org/10.1080/02723646.2014.890092
Glenn, M. S., & Woo, M. (1997). Spring and summer hydrology of a valley-bottom
wetland, Ellesmere Island, Northwest Territories, Canada. Wetlands, 17(2), 321–329.
https://doi.org/10.1007/BF03161420

80

Gröning, M., Lutz, H. O., Roller-Lutz, Z., Kralik, M., Gourcy, L., & Pöltenstein, L. (2012).
A simple rain collector preventing water re-evaporation dedicated for δ 18O and δ
2H analysis of cumulative precipitation samples. Journal of Hydrology, 448–449, 195–
200. https://doi.org/10.1016/j.jhydrol.2012.04.041
Hayashi, M., van der Kamp, G., & Rosenberry, D. O. (2016). Hydrology of prairie
wetlands: Understanding the integrated surface-water and groundwater processes.
Wetlands, 36, 237–254. https://doi.org/10.1007/s13157-016-0797-9
Hengl, T. (2007). A Practical Guide to Geostatistical Mapping of Environmental
Variables. European Commission - Joint Research Centre, EUR 22904, 143.
Hernández-Henríquez, M. A., Sharma, A. R., Taylor, M., Thompson, H. D., & Déry, S. J.
(2018). The Cariboo Alpine Mesonet: Sub-hourly hydrometeorological observations
of British Columbia’s Cariboo Mountains and surrounding area since 2006. Earth
System Science Data, 10(3), 1655–1672. https://doi.org/10.5194/essd-10-1655-2018
Hijmans, Robert J. (2020). raster: Geographic Data Analysis and Modeling. R package
version 3.4-5. https://CRAN.R-project.org/package=raster
IPCC. (2018). Global warming of 1.5°C. An IPCC Special Report on the impacts of global
warming of 1.5°C above pre-industrial levels and related global greenhouse gas emission
pathways, in the context of strengthening the global response to the threat of climate

81

change,.
Islam, S. U., Curry, C. L., Déry, S. J., & Zwiers, F. W. (2019). Quantifying projected
changes in runoff variability and flow regimes of the Fraser River Basin, British
Columbia. Hydrology and Earth System Sciences, 23(2), 811–828.
https://doi.org/10.5194/hess-23-811-2019
Islam, S. U., & Déry, S. J. (2017). Evaluating uncertainties in modelling the snow
hydrology of the Fraser River Basin, British Columbia, Canada. Hydrology and Earth
System Sciences, 21(3), 1827–1847. https://doi.org/10.5194/hess-21-1827-2017
Islam, S. U., Déry, S. J., & Werner, A. T. (2017). Future climate change impacts on snow
and water resources of the fraser river basin, British Columbia. Journal of
Hydrometeorology, 18(2), 473–496. https://doi.org/10.1175/jhm-d-16-0012.1
Kang, D. H., Gao, H., Shi, X., Islam, S. U., & Déry, S. J. (2016). Impacts of a Rapidly
Declining Mountain Snowpack on Streamflow Timing in Canada’s Fraser River
Basin. Scientific Reports, 6, 19299. https://doi.org/10.1038/srep19299
Kang, D. H., Shi, X., Gao, H., & Déry, S. J. (2014). On the Changing Contribution of
Snow to the Hydrology of the Fraser River Basin, Canada. Journal of
Hydrometeorology, 15(4), 1344–1365. https://doi.org/10.1175/jhm-d-13-0120.1
Kayranli, B., Scholz, M., Mustafa, A., & Hedmark, Å. (2010). Carbon storage and fluxes
82

within freshwater wetlands: A critical review. Wetlands, 30(1), 111–124.
https://doi.org/10.1007/s13157-009-0003-4
Ketcheson, S. J., & Price, J. S. (2016). Snow hydrology of a constructed watershed in the
Athabasca oil sands region, Alberta, Canada. Hydrological Processes, 30(14), 2546–
2561. https://doi.org/10.1002/hyp.10813
Ketcheson, S. J., Price, J. S., Carey, S. K., Petrone, R. M., Mendoza, C. A., & Devito, K. J.
(2016). Constructing fen peatlands in post-mining oil sands landscapes: Challenges
and opportunities from a hydrological perspective. Earth-Science Reviews, 161, 130–
139. https://doi.org/10.1016/j.earscirev.2016.08.007
Knudsvig, H. (2015). Characterization of Regional Air Temperature and Precipitation and its
Relationship to Shallow Groundwater Resources in the Ancient Forest, British Columbia
[Master’s Thesis, University of Northern British Columbia].
https://doi.org/10.1016/b978-0-12-347255-7.50009-8
Lee, L. J. E., Lawrence, D. S. L., & Price, M. (2006). Analysis of water-level response to
rainfall and implications for recharge pathways in the Chalk aquifer, SE England.
Journal of Hydrology, 330(3–4), 604–620. https://doi.org/10.1016/j.jhydrol.2006.04.025
Liu, G., & Schwartz, F. W. (2012). Climate-driven variability in lake and wetland
distribution across the Prairie Pothole Region: From modern observations to long-

83

term reconstructions with space-for-time substitution. Water Resources Research,
48(8), W08526. https://doi.org/10.1029/2011WR011539
Maxwell, R. (1987). Biophysical Soil Resources and Land Evaluation of the Northeast Coal
Study Area, 1977-1978 Jarvis Creek-Morkill River Area.
McDonnell, J. J., & Beven, K. (2014). Debates-The future of hydrological sciences: A
(common) path forward? A call to action aimed at understanding velocities,
celerities and residence time distributions of the headwater hydrograph. Water
Resources Research, 50(6), 5342-5350. https://doi.org/10.1002/2013WR015096
Mercer, J. J. (2018). Insights into mountain wetland resilience to climate change: An evaluation
of the hydrological processes contributing to the hydrodynamics of alpine wetlands in the
Canadian Rocky Mountains [Master’s Thesis, University of Saskatchewan].
https://doi.org/10.13140/RG.2.2.12825.47205
Mitsch, W. J., Bernal, B., Nahlik, A. M., Mander, Ü., Zhang, L., Anderson, C. J.,
Jørgensen, S. E., & Brix, H. (2013). Wetlands, carbon, and climate change. Landscape
Ecology, 28(4), 583–597. https://doi.org/10.1007/s10980-012-9758-8
Mitsch, W. J., & Gosselink, J. G. (2007). Wetlands (5th ed.). John Wiley & Sons, Inc.
Muhammad, A., Evenson, G. R., Stadnyk, T. A., Boluwade, A., Jha, S. K., & Coulibaly, P.
(2018). Assessing the importance of potholes in the Canadian Prairie Region under
84

future climate change scenarios. Water (Switzerland), 10(11), 1657.
https://doi.org/10.3390/w10111657
Nahlik, A. M., & Fennessy, M. S. (2016). Carbon storage in US wetlands. Nature
Communications, 7, 13835. https://doi.org/10.1038/ncomms13835
National Wetlands Working Group. (1997). The Canadian Wetland Classification
System. In The Wetland Book: I: Structure and Function, Management, and Methods.
https://doi.org/10.1007/978-90-481-9659-3_340
Nelson, S. T. (2000). A simple, practical methodology for routine VSMOW/SLAP
normalization of water samples analyzed by continuous flow methods. Rapid
Communications in Mass Spectrometry, 14(12), 1044–1046.
https://doi.org/10.1002/1097-0231(20000630)14:12<1044::AID-RCM987>3.0.CO;2-3
O’Geen, A. T., Maynard, J. J., & Dahlgren, R. A. (2007). Efficacy of constructed wetlands
to mitigate non-point source pollution from irrigation tailwaters in the San Joaquin
Valley, California, USA. Water Science and Technology, 55(3), 55–61.
https://doi.org/10.2166/wst.2007.072
Pebesma, E. J. (2004). Multivariable geostatistics in S: the gstat package. Computers &
Geosciences, 30(7), 683–691. https://doi.org/10.1016/j.cageo.2004.03.012
Penna, D., & van Meerveld, H. J. (2019). Spatial variability in the isotopic composition of
85

water in small catchments and its effect on hydrograph separation. Wiley
Interdisciplinary Reviews: Water, 6(5), e1367. https://doi.org/10.1002/wat2.1367
Phillips, J. D. (1989). Fluvial Sediment Storage in Wetlands. Journal of the American Water
Resources Association, 25(4), 867–873. https://doi.org/10.1111/j.17521688.1989.tb05402.x
Pitchford, J. L., Wu, C., Lin, L., Petty, A. T., Thomas, R., Veselka, W. E., Welsch, D.,
Zegre, N., & Anderson, J. T. (2012). Climate change effects on hydrology and
ecology of wetlands in the mid-atlantic highlands. Wetlands, 32(1), 21–33.
https://doi.org/10.1007/s13157-011-0259-3
Richardson, C. J. (1994). Ecological functions and human values in wetlands: A
framework for assessing forestry impacts. Wetlands, 14(1), 1–9.
https://doi.org/10.1007/BF03160616
Sharma, A. R., & Déry, S. J. (2016). Elevational dependence of air temperature
variability and trends in British Columbia’s Cariboo Mountains, 1950-2010.
Atmosphere - Ocean, 54(2), 153–170. https://doi.org/10.1080/07055900.2016.1146571
Shaw, D. A., Vanderkamp, G., Conly, F. M., Pietroniro, A., & Martz, L. (2012). The FillSpill Hydrology of Prairie Wetland Complexes during Drought and Deluge.
Hydrological Processes, 26(20), 3147–3156. https://doi.org/10.1002/hyp.8390

86

Tews, J., Brose, U., Grimm, V., Tielbörger, K., Wichmann, M. C., Schwager, M., &
Jeltsch, F. (2004). Animal species diversity driven by habitat
heterogeneity/diversity: the importance of keystone structures. Journal of
Biogeography, 31(1), 79–92. https://doi.org/10.1046/j.0305-0270.2003.00994.x
U.S. Army Corps of Engineers (1998). Engineering and Design: Runoff from Snowmelt.
Engineer Manual 1110-2-1406.
Vymazal, J. (2010). Constructed wetlands for wastewater treatment. Water, 2(3), 530–
549. https://doi.org/10.3390/w2030530
Wells, C., Ketcheson, S., & Price, J. (2017). Hydrology of a wetland-dominated
headwater basin in the Boreal Plain, Alberta, Canada. Journal of Hydrology, 547, 168–
183. https://doi.org/10.1016/j.jhydrol.2017.01.052
Whitmire, S. L., & Hamilton, S. K. (2008). Rates of anaerobic microbial metabolism in
wetlands of divergent hydrology on a glacial landscape. Wetlands, 28(3), 703–714.
https://doi.org/10.1672/06-126.1

87