INTRODUCTION
The glacier surge cycle has been attributed to dynamics that are internal to the glacier system, with contrasts in the surge and quiescent velocities controlled by sliding at the glacier bed (e.g. Meier and Post, Reference Meier and Post1969). The surge phase is characterized by rapid sliding, which is thought to arise from instabilities in the subglacial drainage system (Kamb, Reference Kamb1987), changes in thermal conditions at the bed (Clarke and others, Reference Clarke, Collins and Thompson1984), changes in the hydraulic conductivity or mechanics of the subglacial till (Truffer and others, Reference Truffer, Harrison and Echelmeyer2000), or various combinations thereof (e.g. Fowler, Reference Fowler1987; Murray and others, Reference Murray, Strozzi, Luckman, Jiskoot and Christakos2003). Globally, only a small fraction of glaciers surge, with some mountain ranges having a high density of surge-type glaciers (e.g. Jiskoot and others, Reference Jiskoot, Boyle and Murray1998). Though these mountain ranges are found to span various climatic and geological regimes (Cuffey and Paterson, Reference Cuffey and Paterson2010), surge-type glaciers seem to cluster within a given temperature – precipitation envelope (Sevestre and Benn, Reference Sevestre and Benn2015). The observations that glacier surges are rare on a global scale and are triggered by instabilities at the glacier bed have led to the suggestion of an underlying causative relationship with bedrock geology (e.g. Post, Reference Post1969). Large-scale univariate statistical analyses aimed at testing this hypothesis have been successful in correlating glacier surges with limestone bedrock in Spitsbergen (Hamilton and Dowdeswell, Reference Hamilton and Dowdeswell1996), and to basalts and sedimentary rocks in East Greenland (Jiskoot and others, Reference Jiskoot, Murray and Luckman2003), while a multivariate analysis suggests a correlation between shale/mudstone bedrock and surging for all of Svalbard (Jiskoot and others, Reference Jiskoot, Boyle and Murray1998). Where correlations exist between bedrock lithology and glacier surges, explanations tend to centre around the ways, in which bedrock erodibility influences till production (e.g. Hamilton and Dowdeswell, Reference Hamilton and Dowdeswell1996), or how bedrock can control the grain size distribution of till and thus its mechanics (Jiskoot and others, Reference Jiskoot, Boyle and Murray1998). Other geological variables, like proximity to major faults or the number of subglacial geological contacts, have been unsuccessful at predicting the geographical occurrence of surge-type glaciers (e.g. Clarke and others, Reference Clarke, Schmok, Ommanney and Collins1986; Jiskoot and others, Reference Jiskoot, Boyle and Murray1998).
Building on the hypothesis that glacier surges are related to the underlying bedrock geology, we analyze the suspended sediment size distributions (SSSDs) from the proglacial streams of both surge-type and non-surge-type glaciers that are underlain by different fractions of two distinct lithologies. In glacial systems, the SSSD likely depends on the size distribution of the sediment source and the transport capacity of the stream (e.g. Karlsen, Reference Karlsen1991). Seasonal glacier dynamics may also play a role as changes in the subglacial drainage morphology lead to changes in the SSSD (e.g. Humphrey and Raymond, Reference Humphrey and Raymond1994; Thayyen and others, Reference Thayyen, Gergan and Dobhal1999; Pandey and others, Reference Pandey, Singh and Hasnain2002). We therefore design our study to control for both bedrock lithology and glacier-type, where the term ‘glacier-type’ denotes surge-type versus non-surge-type. We search for differences in SSSDs using tests of significance on the distribution statistics and principal component analysis (PCA).
FIELD SITE
We sampled a total of 20 glaciers in the Donjek Range and Maxwell Group of the St. Elias Mountains, Yukon, Canada (Fig. 1). The study area comprises rocks from the Kaskawulsh Group, which belongs to the Alexander terrane (Wheeler, Reference Wheeler1963; Dodds and Campbell, Reference Dodds and Campbell1988; Israel and Cobbett, Reference Israel, Cobbett, Emond, Bradshaw, Lewis and Weston2008). The study glaciers are underlain by Paleozoic metasedimentary (MS) rocks intruded by various fractions of Jurassic aged felsic plutonic rock. In mapping outcrops at the glacier margins of the MS basins, we observed rocks ranging from slate to schist, with hornfels near the contacts of the plutonic rocks. An outcrop of amphibolite appeared near the termini of two glaciers at the northern end of the study area (Glaciers 12 and 13; Fig. 1), but Dodds and Campbell (Reference Dodds and Campbell1988) show that this rock type is limited to the terminus area. Wheeler (Reference Wheeler1963) and Dodds and Campbell (Reference Dodds and Campbell1988) have identified sandstones and volcanics in the MS units. We, however, did not encounter these rock types in any of the sampled basins, where we visited an average of ~3 outcrops per basin. Thin section analysis indicates that the MS rocks are almost all rich in quartz, plagioclase, calcite and biotite. Previous work indicates that the plutonic rocks range in composition from granodiorite to tonalite (Wheeler, Reference Wheeler1963), which was also confirmed by our thin section analysis. Although we find a range of compositions for the MS rocks, the thin sections show that the MS rocks have a narrow distribution of grain sizes with a smaller mean grain size than the more coarse grained plutonic rocks.
Glaciers sampled in the Donjek Range and Maxwell Group, St. Elias Mountains, Yukon, Canada. Surge index given by Clarke and others (Reference Clarke, Schmok, Ommanney and Collins1986) for all glaciers, with large symbols representing the classification groups used in this study. We use ‘glacier-type’ to indicate whether the glacier is surge-type or non-surge-type. Bedrock geology compiled by Gordey and Makepeace (Reference Gordey and Makepeace1999). Note that the ice extents from geological map are not current. Landsat 7 imagery used for grayscale under-shading.
Most of the study glaciers are unnamed and henceforth referred to as Glaciers 1–20. Maxwell Glacier (Glacier 17) has been previously studied (Kruszynski, Reference Kruszynski1989), while Glaciers 1 and 2 have been previously documented in the literature as South- and North Glaciers, respectively (e.g. MacDougall and Flowers, Reference MacDougall and Flowers2011). Glaciers 1 and 2 differ in glacier-type, length and orientation, yet ice temperature measurements indicate that these glaciers share a similar thermal structure, characterized by temperate ice in the accumulation area and cold ice in the lower ablation area (Wilson and others, Reference Wilson, Flowers and Mingo2013). Trapridge Glacier also displays a similar thermal structure (Clarke and others, Reference Clarke, Schmok, Ommanney and Collins1986), despite being ~70 km away. We therefore assume a relatively common thermal structure for all glaciers within the sample area, as further supported by the modelling work of Wilson and Flowers (Reference Wilson and Flowers2013).
We assign glacier-type based on the surge index (S.I.) of Clarke and others (Reference Clarke, Schmok, Ommanney and Collins1986), where the likelihood that a glacier is surge-type is ranked from 0–5. The Clarke and others (Reference Clarke, Schmok, Ommanney and Collins1986) scheme is as follows: ‘i = 0, no special features suggesting surges (very probably normal); i = 1, uncertain surge characteristics; i = 2, possible surge characteristics (one or two features); i = 3, probable surge characteristics (several features); i = 4, very probable surge characteristics (strong surface evidence); i = 5, definitive surge characteristics…’ For our study, we use a binary classification of glacier-type, where glaciers with an index value <3 are taken as non-surge-type, and >3 as surge-type (see Table 1 and Fig. 1). We classify the 8 glaciers with an index of 3 based on time-lapse satellite imagery, aerial photography and field observations, and consider the implications of our classification in the discussion. The time-lapse imagery was compiled from georeferenced USGS Landsat preview images dating from 1975 to present. A binary scheme is also used to group glaciers by basin lithology. Glacier 1 (South Glacier) is underlain solely by granodiorite (felsic plutonic), while 10 of the glaciers are underlain solely by undifferentiated MS rocks. The remaining glaciers have a mixed lithology (MX), being underlain partially by MS rocks and partially by felsic plutonic rocks. The plutonic rocks are confined to the upper basins of these glaciers. Based on the designation of glacier-type and lithology described above, we classify each glacier as a member of one of four groups: (1) metasedimentary surge-type (MS-S), (2) metasedimentary non-surge-type (MS-NS), (3) mixed lithology surge-type (MX-S) and (4) mixed lithology non-surge-type (MX-NS), where Glacier 1 has been included in the MX group. Although we discuss the geology in terms of the bedrock lithology, we are not implying that these glaciers rest on a hard bed. Rather, we observe widespread till in the forefield of all glaciers, regardless of glacier-type. We are not able to differentiate between hard and soft bedded glaciers, but our observations indicate the likely presence of a soft bed proximal to the terminus for all glaciers in this study.
Characteristics of sampled glaciers. Glaciers 1–20 are classified in this study by MS-S, MS-NS, MX-S and MX-NS based on the S.I. of Clarke and others (Reference Clarke, Schmok, Ommanney and Collins1986) and bedrock lithology. MS % indicates the estimated fraction of MS rocks underlying the glacier, where values <100% indicate that the remaining fraction is plutonic. Discharge (Q) was measured at time and location of suspended sediment sample collection, with no measurements collected at Glaciers 3 and 6.
* Denotes that a small portion of the terminus overrides an amphibolite unit.
METHODS
Sample collection and analysis
Glaciers 2–20 were all sampled over 5 d between 1 and 7 August 2013. One sample per glacier was obtained by submerging a 250 mL Nalgene bottle in the proglacial stream within 50 m of the glacier terminus. The effect of stream sampling location is discussed below. Dilution gauging was carried out at each site with a slug injection of salt following the methods of Moore (Reference Moore2004). The conductivity response of the salt injection was measured 50–100 m downstream, and the discharge calculated from the dispersal characteristics of the conductivity plume, as more thoroughly described in Crompton and others (Reference Crompton, Flowers, Kirste, Hagedorn and Sharp2015). Suspended sediment concentrations were measured separately on the filtrate of vacuum pumped samples.
All samples were analyzed by laser diffraction through a Malvern Mastersizer 2000 particle size analyzer, yielding measurements of sediment volume in each size class. Grain sizes greater than fine sand were excluded from the analysis by wet-sieving the samples through a 90 µm mesh. Composite grains that may have flocculated in the glacial environment (Woodward and others, Reference Woodward, Porter, Lowe, Walling and Evans2002; Chanudet and Fillela, Reference Chanudet and Fillela2006) were dispersed using a 0.05% sodium hexametaphosphate solution, and after 24 h of soaking were placed in a sonic bath leading into the loading cell. Error in the size distribution arises from the assumption that all particles are spherical. Furthermore, the distribution is likely finer than measured given the presence of inorganic colloids that may not have completely deflocculated (Chanudet and Fillela, Reference Chanudet and Fillela2006) and the presence of platy particles that align perpendicular to the Mastersizer laser (e.g. Konert and Vandenberghe, Reference Konert and Vandenberghe1997).
To quantify the temporal variation in SSSD in relation to discharge and suspended sediment concentration, 28 samples were collected throughout the 2013 melt season at Glacier 1. Multiple dilution gauging experiments were performed to construct a rating curve, from which we estimated discharge at the time of sample collection (methods described in Crompton and others (Reference Crompton, Flowers, Kirste, Hagedorn and Sharp2015)). In addition to sampling by submerging a bottle in the proglacial stream, samples were also collected with an automatic water sampler. For a 5 min sampling period on 23 July 2013, during which the discharge can be assumed constant, we tested a variety of sampling procedures to characterize the uncertainty associated with sampling methodology. Results indicate that there is negligible difference between the samples collected by submerging the bottle in the middle of the stream (as done for Glaciers 2–20), submerging the bottle in the middle of the stream but 50 m upstream from the usual sampling site, and using the automatic water sampler. An integrated depth sampling technique resulted in a higher mean grain size, while submerging the bottle in an eddy resulted in a lower mean grain size. We therefore avoided using these techniques at Glaciers 2–20.
Distribution statistics and hypothesis testing
Mastersizer measurements are output as a cumulative sediment size distribution for each sample, from which we construct a probability density function (the SSSD) by numerical differentiation using a two-point centred difference stencil. We then interpolate the SSSD using the Matlab® piecewise cubic interpolation pchip function, which preserves the location and magnitude of all maxima and minima within the data. For each distribution, we compute a mean, standard deviation, skew and kurtosis, hereafter referred to as the ‘distribution statistics’. The skew and kurtosis are computed as normalized moments, and the kurtosis is normalized by a value of 3 (e.g. Joanes and Gill, Reference Joanes and Gill1998). We compute the median grain size, but do not include it in the analysis as it furnishes information similar to that captured by the mean. We also analyze the distribution quantiles, but find these less effective than the distribution statistics in characterizing the samples. We perform all calculations using arithmetic rather than geometric quantities, as the latter lead to a poor representation of the distribution at the fine end of the spectrum where we observe significant structure.
For each group (MS-S, MS-NS, MX-S, MX-NS), we compute the mean and std dev. of each distribution statistic (i.e. the mean of the means, the std dev. of the kurtoses, etc.). In an effort to address the question of whether the groups are statistically distinguishable, we test the null hypothesis that the distribution statistics have means that are equivalent across groups by carrying out an analysis of variance (ANOVA). Should the ANOVA F-test obtain a decidedly small p-value (e.g. p < 0.05), we conclude that at least one group is unlike the others in terms of the statistic in question. In order to further determine which groups are different from one another, we employ a Tukey honest significant difference (HSD) test (with α = 0.1). This multiple-comparison procedure determines the pairwise subsets that are significantly different in terms of the chosen statistic. The Tukey HSD test is a more conservative form of the t-test, having the advantage of accounting for the error rate, or probability that two groups are different by chance variability (Type I error) (e.g Dowdy and others, Reference Dowdy, Wearden and Chilko2011). Before significance testing is carried out, we first ensure that the sample statistics are normally distributed within groups using the χ 2 test at the 95% significance level following the methods of Davis and Sampson (Reference Davis and Sampson2002).
Principal component analysis
We compare the shapes of the SSSDs using PCA, a commonly used multivariate technique in the natural sciences with applications ranging from meteorology (e.g. Hardy and Walton, Reference Hardy and Walton1978) to geochemistry (e.g. Mitchell and others, Reference Mitchell, Brown and Fuge2006). PCA identifies linear combinations of variables that explain the greatest variance within a dataset (e.g. Jolliffe, Reference Jolliffe2002). In this application, we take the variables to be the grain sizes as binned according to the standard
$\mathop {\log} \nolimits_2 $
or phi (ϕ) scale, with the value of the kth bin taking on the average grain size within that bin (e.g. Davis and Sampson, Reference Davis and Sampson2002). Binning the grain sizes in this way allows for a more dense sampling of the finer grain sizes, where we observe significant structure in the measured distributions. To find the principal components (PCs), we construct a matrix X of dimensions m × n, where there are n discrete grain sizes ranging from ϕ = 3 to 12, and m samples. We then compute the eigenvectors (ν) and eigenvalues (λ) of the variance-covariance matrix of X. The normalized eigenvalues give the percent variance explained by each eigenvector (latency), while the eigenvectors (loads) provide the directions along which the greatest variance occurs. Projecting the data matrix X onto the eigenvectors yields a matrix of scores (s) as
${\bf X}\nu = {\bf s}$
, creating a transformation that allows the data to be visualized in PC space.
RESULTS
SSSDs from all glaciers show a peak at ϕ = 10.5 (Fig. 2). Although each sediment sample has a large fraction of grains around ϕ = 10.5, the precise location of the peak in this region is fixed by the interpolation method. Some of the distributions are unimodal, with the most pronounced examples in the MX-S group (Fig. 2b). A second peak at ϕ ≈ 8 creates a bimodal structure for several distributions, and is most apparent in the MS-NS group (Fig. 2c). Given our analytical equipment and the fact that we are only analyzing grains smaller than medium sand, our sampling and analysis methodologies are not entirely consistent with previous work (e.g. Pandey and others, Reference Pandey, Singh and Hasnain2002; Haritashya and others, Reference Haritashya, Kumar and Singh2010). As a result, we avoid direct comparison with other studies, as suggested by Karlsen (Reference Karlsen1991).
SSSDs collected in proglacial streams. Thick black lines represent the average distribution within groups for (a) MS-S, (b) MX-S, (c) MS-NS and (d) MX-S glaciers. The shaded region in panel (b) shows ± 1σ computed from the mean of 28 samples collected at Glacier 1 during the 2013 melt season. The distributions are plotted using the
$\phi {\kern 1pt} {\kern 1pt} (log_2 )$
scale. Note that distribution statistics are computed in absolute space, while PCA uses data equally spaced on the ϕ scale.
A comparison of the distribution statistics (Fig. 3) proves useful in differentiating the four groups (MS-S, MS-NS, MX-S, MX-NS), where distinctions are not necessarily obvious when examining the distribution curves themselves (Fig. 2). Small p-values from the ANOVA (Table 2) indicate that at least one group is significantly different than the rest, which is the case for each distribution statistic. At the α = 0.1 significance level, the Tukey HSD test shows that, for all distribution statistics, the MX-NS (4) and MS-NS (2) groups are significantly different, the MS-NS (2) and MX-S (3) groups are significantly different and the MS-S (1) and MS-NS (2) groups are significantly different with the exception of the skew (Fig. 3; Table 2). We cannot reject the null hypothesis for equivalence of means at the α = 0.1 significance level for the MS-S (1) and MX-S (3) groups for any of the distribution statistics (Fig. 3; Table 2, columns 1 and 3), indicating that there are notable similarities among all surge-type glaciers, regardless of bedrock lithology. Lastly, we find that the MS basins are significantly different than MX basins for all distribution statistics, regardless of glacier type. The same cannot be said for surge-type versus non-surge-type glaciers for a test that does not account for lithological class (not shown).
Statistics of SSSDs for each group: MS-S, MS-NS, MX-S and MX-NS. Horizontal bars indicate mean values of the groups for the given distribution statistic. (a) Mean. (b) Std dev. (c) Skew. (d) Kurtosis.
ANOVA using the distribution statistics, with group means shown in first four columns, where the groups are MS-S and MS-NS, and MX-S and MX-NS. Small p-values (from ANOVA F-test) suggest significant difference across groups for each distribution statistic. The last column shows pairwise groups considered significantly different from one another (at α = 0.10, i.e. p-value < 0.010) as computed by the Tukey HSD test.
The results above indicate that SSSDs (1) differ as a function of lithology, especially for the non-surge-type glaciers, (2) differ as a function of glacier-type in the MS basins and (3) are statistically indistinguishable for surge-type glaciers with different lithologies. Although the MX-S group has a lower mean and std dev. than the MX-NS group (Fig. 3; Table 2, columns 3 and 4), these differences are not statistically significant. We speculate on reasons for this lack of differentiation in the discussion.
PCA performed on the entire dataset (Fig. 4) indicates that 67.5% of the variance can be explained by the first PC (PC 1), 22.4% by the second (PC 2) and 5.7% by the third (PC 3). With 95% of the variance explained by the first three components, the remaining components represent noise within the data (Jolliffe, Reference Jolliffe2002). Similar to the distribution statistics, the PCA shows that the MS-NS group is clearly distinct from the others (Fig. 4a), corroborating the results of the distribution statistics in discriminating glacier-type in the MS basins, and lithology for the non-surge-type glaciers. This discrimination is most clear in the scores of PC 1 (x-axis in Fig. 4a). Based on the loads (Fig. 4b), PC 1 is characterized by the contrasting contributions of finer (ϕ = 6 − 12) and less fine (ϕ = 3 − 6) sediment, while PC 2 is characterized by contrasting contributions of silt (ϕ = 5 − 9) versus clay and sand (ϕ < 5 and 9 < ϕ). The fact that no individual load stands out in PC 1 indicates that no single grain size can account for most of the variance within the data (Davis and Sampson, Reference Davis and Sampson2002). The PCA was also performed using subsets of the data that comprised two groups at a time (e.g. MS-S and MX-S) but did not yield additional insight.
Results of principal component analysis on full dataset. (a) Data projected onto first and second PCs (scores), with groups distinguished as shaded polygons by solid fill. Symbol types are consistent with Figures 1 and 3 as described in text. (b) Eigenvectors (loads) for PC 1 (filled bars) and PC 2 (open bars). Grain sizes shown by ϕ scale and binned by fine clay (fc), medium clay (mc), coarse clay (cc), very fine silt (vfsi) to coarse silt (csi) and very fine sand to fine sand (fs).
DISCUSSION
Glacier classification and experimental robustness
The distinctive characteristics of the SSSDs from MS-NS glaciers are important to the outcomes of this study. To ensure that the influence of the MS-NS group is not dominated by processes unaccounted for, we examined the influence of proglacial discharge, suspended sediment concentration and glacier length, area and orientation. At Glacier 1, we find that there is notable variation in the mean SSSD as a function of discharge (Fig. 5 inset). This relationship is best fit (r = 0.64) with a quadratic function. The dependence of distribution mean on discharge at Glacier 1 raises some concern about the discharge dependence of samples collected at Glaciers 2–20, though this concern is partially allayed by the lack of observed correlation between SSSD statistics and discharge across all glaciers (Fig. 5), or within any of the groups. The mean discharge is roughly equal across three of the groups, while MS-NS has a higher mean discharge (Table 1) but lower mean SSSD. This further supports the limited role of discharge in controlling the SSSD as high discharge does not correlate with high mean SSSD in this case. We also calculate that, for all samples, the maximum measured grain size of 90 µm would be in suspension. We do this by conservatively estimating a stream slope of 0.01 and a water depth of 0.1 m based on the geometry of the stream at Glacier 5, which had the lowest measured discharge. With an estimated particle density of 2.7 g cm−3, we calculate a particle Reynolds number of ~5 and a Shields stress of ~6 (e.g. Garcia, Reference Garcia2008). When plotting these values on the Shields-Parker diagram, we see that the maximum grain size is well above the threshold for suspension. These calculations account for shear velocity and turbulence at the grain boundary. To investigate the effect of stream turbulence or sampling location, recall that we conducted experiments on sampling methodology at Glacier 1 that indicated no significant difference in the SSSD statistics for different sampling locations, provided the sample was not collected in an eddy. The same observation was made by Humphrey and Raymond (Reference Humphrey and Raymond1994), but with sediment size truncated at 50 µm, instead of 90 µm. Lastly, we find that the SSSD statistics show no significant correlations with the suspended sediment concentration, glacier length or glacier orientation.
Suspended sediment size distribution means. Symbol types consistent with Figures 1, 3 and 4. Average of 28 samples at Glacier 1 shown in the middle of cross bars, with bar lengths of ± 1σ for discharge and distribution mean over all 28 samples. The inset shows the mean of SSSDs from Glacier 1 samples with a quadratic fit to the discharge (solid black line) as μ = 0.91Q2 − 0.39Q + 17.79 (r = 0.64). Dashed lines show 95% confidence bounds. Axes on inset share units with outer plot. The box encloses the distribution means of five samples collected at Q = 3.0 m3 s−1, using different sampling methodologies as described in text.
We classified basins as MX if any plutonic rocks were outcropping at the glacier margin or found within the proglacial till, leaving little ambiguity in the lithological classification. The assignment of glacier-type is perhaps more ambiguous. Glaciers 13–16 have a S.I. of 3, and their glacier-type classification is not entirely obvious. Field evidence and time-lapse imagery do not show any signs of surge activity for glaciers 14–16, but we cannot rule out the possibility that these glaciers surged before 1975, from which we acquired the first satellite image. However, since these glaciers are situated in MX basins, where there is little difference in SSSD between surge-type and non-surge-type glaciers, the glacier-type does not change the outcome of the study. Glaciers 1, 3, 8, 17 and 18 also have a S.I. of 3, but field evidence and satellite imagery strongly indicate that they are surge-type. If focusing solely on the surge-index from Clarke and others (Reference Clarke, Schmok, Ommanney and Collins1986), we see that there is a high concentration of glaciers with a surge index of 3 in the MX group, as opposed to the more clearly separated surge indices in the MS group (Table 1). The latter may indicate that surges are more evident on MS rock, while the former may be an indication of the ambiguity in classifying glacier-type within the MX-NS group.
Although we verify the S.I. with satellite imagery spanning from 1975 to present, surge indices from Clarke and others (Reference Clarke, Schmok, Ommanney and Collins1986) were published in 1980 using aerial photography that may have dated back decades before publication. Surge indices assigned based on modern imagery could potentially be different, in part due to the possibility of changes in glacier behaviour over time. For example, Budd (Reference Budd1975), Dowdeswell and others (Reference Dowdeswell, Hodgkins, Nuttall, Hagen and Hamilton1995) and Hoinkes (Reference Hoinkes1969) hypothesize that a significant change in mass balance can result in a transition in glacier-type. Evidence for this phenomenon is thought to exist in the St. Elias Mountains (Frappé and Clarke, Reference Frappé and Clarke2007; Flowers and others, Reference Flowers, Roux, Pimentel and Schoof2011). Lastly, another potential source of uncertainty in our study arises from not knowing the precise phase of each glacier surge cycle, with Glacier 6 being the only glacier to show visual signs of surging during the time of sampling. Glacier 1 is likely undergoing a slow surge (De Paoli and Flowers, Reference De Paoli and Flowers2009), but this interpretation was only possible through more extensive field work than was carried out on Glaciers 2–20.
Interpretation of correlations
In the following discussion we consider possible geological and glaciological conditions that may give rise to similarities or differences in SSSDs among groups. We then speculate on correlations between SSSDs and glacier-type. We take caution not to generalize the relationship between glacier-type and SSSD, given the similar SSSDs across MX basins and the resulting indistinguishability of SSSDs by glacier-type for the dataset as a whole.
Bedrock geology
In explaining how bedrock geology controls the SSSD, we might hypothesize that the SSSD is solely a function of the size distribution of the underlying bedrock mineralogy (e.g. Karlsen, Reference Karlsen1991). For example, Dreimanis and Vagners (Reference Dreimanis, Vagners, Yatsu and Falconer1972) argue that the finest grain size in till stems from the terminal or mineral grain size within the bedrock. Alternatively, bedrock geology could influence the SSSD through chemical processes like mineral dissolution or precipitation (e.g. Pandey and others, Reference Pandey, Singh and Hasnain2002; Crompton and others, Reference Crompton, Flowers, Kirste, Hagedorn and Sharp2015). The degree of comminution will also affect the sediment size distribution within the till (e.g. Hooke and Iverson, Reference Hooke and Iverson1995). Bedrock geology may have an ultimate control on the shape, hardness and initial size of clasts within the till, and thus the mechanics of till deformation and the resulting comminution.
Geological controls can be invoked to explain the difference in SSSDs between MX and MS basins, but might also be used to resolve differences between the MS-NS and MS-S groups, given that lithology alone does not represent the fracture characteristics or geochemistry of the bedrock. An alternative geological classification may help to differentiate SSSDs between MX-S and MX-NS groups, which are indistinguishable with the current classification scheme. The lack of differentiation between SSSDs from MX-S versus MX-NS basins may also be explained by noting that the classification scheme fails to represent the fraction of MS rock in MX basins. However, in defense of the lithological classification, thin section analysis shows that all the MS rocks have smaller grain sizes than the plutonic rocks, indicating that lithology exerts at least some control on grain size of the source rock.
Glacier processes
As mentioned above, the SSSD could be a function of till deformation processes; till deformation in turn depends on glaciological variables such as basal shear stress and effective pressure. Alternatively, the SSSD may simply reflect the natural sampling of sediment by the evolving subglacial drainage system (Karlsen, Reference Karlsen1991; Humphrey and Raymond, Reference Humphrey and Raymond1994; Thayyen and others, Reference Thayyen, Gergan and Dobhal1999). During the 1982/83 surge of Variegated Glacier, a fining of grain sizes was observed with an increase in sliding speed (Humphrey and Raymond, Reference Humphrey and Raymond1994). If fines are purged during surge events, we might expect surge-type glaciers in early-to-mid quiescence to have grain size distributions characterized by higher means. Such a situation would be consistent with the higher mean SSSDs that we observe in MS-S than MS-NS basins.
If the SSSD is predominately controlled by natural sampling mechanisms of the subglacial stream, then the SSSD may not reflect the grain size distribution of subglacial sediment. In this case, bedrock lithology may still play an indirect role in determining SSSD. For example, lithology may broadly control valley shape or bed roughness, thus influencing the resistance to ice flow at a large scale and basal water flow at a smaller scale, both contributing to the morphology and evolution of the drainage system and thus sediment evacuation. Logical steps forward in differentiating between geological and glaciological controls on SSSD would be to analyze the grain size distribution by mineral, to compare the mineral saturation states with the chemical composition of the water, and to perform laboratory comminution experiments on the source rocks.
Relationships between SSSD and glacier-type
With the understanding that the SSSD may be a product of geological or glaciological influences, we now consider correlations between glacier-type and grain size, under the assumption that the SSSD reflects the subglacial sediment size distribution. There are three possibilities: (i) SSSD and glacier-type are correlated but are not causative in either direction, (ii) glacier-type affects subglacial sediment size or (iii) subglacial sediment size affects glacier-type.
As an example of case (i), SSSD could be a proxy for bedrock erodibility or till production rate, and thus a control on till thickness or spatial extent, which may govern glacier-type. In this case, SSSD and glacier-type are both functions of bedrock characteristics, but do not depend directly on each other. Although the extent and thickness of subglacial till has not been quantified for the study glaciers, we observe that both surge-type and non-surge-type glaciers in this study reveal till in the forefield as they recede. As an example of case (ii), we might expect that sediment near the terminus of a surge-type glacier will have been cyclically exposed to subaerial processes at a much higher frequency than sediments underlying non-surge-type glaciers. Grain size sorting and comminution in the subaerial environment depend on the amount of freeze-thaw, fluvial erosion from rain and snowmelt, and subaerial chemical weathering. Therefore, we might expect that sediments sourced from the lower elevations of surge-type glaciers will have different size distributions than sediments sourced from subglacial environments alone. As an example of case (iii), sediment size could influence glacier-type through the hydro-mechanical properties of the till, such as hydraulic conductivity (Freeze and Cherry, Reference Freeze and Cherry1979). The hydraulic conductivity impacts basal water pressure, which is important not only in controlling the morphology of the basal drainage system, but the mechanical behaviour of the till itself (e.g. Tulaczyk and others, Reference Tulaczyk, Kamb and Engelhardt2000; Iverson, Reference Iverson2010). The fine fraction of till may also control rates of regelation and freezing front migration into the till (e.g. Rempel, Reference Rempel2008). A complex interplay of these processes could lead to a decoupling of the ice from its bed (e.g. Iverson, Reference Iverson1999), failure within the till (e.g. Truffer and others, Reference Truffer, Harrison and Echelmeyer2000) or an inhibition of sliding altogether, all of which impact temporal variations in basal velocity. At present, we do not have sufficient information to favour any one of cases (i)–(iii) over the others.
CONCLUSION
We have attempted to probe the relationship between bedrock geology and glacier surges in and around the Donjek Range of the St. Elias Mountains, Yukon. Our results show that the SSSDs from the proglacial streams of glaciers underlain solely by MS rocks can be distinguished from those of glaciers underlain by a mix of MS and felsic plutonic rocks, with the differences being more pronounced when considering only non-surge-type glaciers. The most intriguing results of our study are (1) the significant differences in the SSSDs between surge-type and non-surge-type glaciers underlain solely by MS rocks, and (2) the remarkably similar distributions among all surge-type glaciers, regardless of bedrock lithology. The results of the SSSDs can be explained with equal plausibility by considering either geological or glaciological processes. Where relationships exist between SSSD and glacier type, the existence and direction of causality remains unknown. Analysis of the hydrochemistry and mineralogy of stream sediments, along with crushing experiments on the source rocks, should help resolve these ambiguities and enable us to better identify the environments that predispose glaciers to surging.
ACKNOWLEDGEMENTS
We thank the Kluane First Nation (KFN), Parks Canada and the Yukon Territorial Government for granting us permission to work in traditional KFN territory and Kluane National Park and Reserve. We are grateful for financial support provided by the Natural Sciences and Engineering Research Council of Canada, the Garfield Weston Foundation of the Association of Canadian Universities for Northern Studies, the Yukon Geological Survey, Simon Fraser University, the Northern Scientific Training Program and the Polar Continental Shelf Project. We kindly acknowledge Trans North Helicopter pilot Dion Parker, and the Arctic Institute of North America's Kluane Lake Research Station for facilitating field logistics. We are grateful to Flavien Beaud and Laurent Mingo for all aspects of field assistance, and to Steve Israel for collaboration on bedrock mapping. We are grateful to Harlan Campbell for assistance with the statistics including the ANOVA and Tukey HSD test. We thank Peijun (Perry) Sang and Jiguo Cao for useful discussions about PCA. John Clague provided access to the Mastersizer, and Rita Dubman, Rachel Hay and Hazel Wong assisted with sample analysis. We thank two anonymous reviewers and Scientific Editor Ian Hewitt for their insightful comments.
Open access
