Concentration-Discharge Patterns Across the Gulf of Alaska Reveal Geomorphological and Glacierization Controls on Stream Water Solute Generation and Export

High latitude glacierized coastal catchments of the Gulf of Alaska (GoA) are undergoing rapid hydrologic changes in response to climate change and glacial recession. These catchments deliver important nutrients in the form of both inorganic and organic matter to the nearshore marine environment, yet are relatively understudied with respect to characterization of the solute generation processes and total yields. Using multiple linear regression informed by Bayesian Information Criterion analysis we empirically demonstrate how watershed characteristics aﬀect solute generation as represented by concentration-discharge relationships. We ﬁnd that watershed mean slope and relief control solute generation and that solute yields are inﬂuenced most by glacier coverage. We contribute a new ﬂux and concentration-discharge based conceptualization for understanding solute cycles across a hydroclimatic gradient of GoA watersheds that can be used to better understand future watershed responses to rapid hydrologic change.


Introduction
High latitude regions are warming at rates of two to three times the global average (IPCC, 2007).Precipitation regime changes associated with the increasing temperatures will result in increased precipitation, primarily in the form of autumn and winter rain (Beamer et al., 2017).Associated with these changes, global glacier volume will be 29-41% by 2100 compared to 2006 with a projected 20% decline of global glacier runoff (Radić et al., 2014).Ice fields and glaciers cover 18% of the 420,230 km 2 Gulf of Alaska (GoA) region and supply 47% of the freshwater water runoff (Neal et al., 2010).Looking to the future it is predicted that Alaska will experience a 30% decline in runoff by the end of the 21 st century (Bliss et al., 2014) accompanied by a forecasted decrease of glacier volume between 32±11 and 58±14% for RCP2.6 and RCP8.5 respectively (Huss & Hock, 2015).As such, the climate change predicted for the 21 st century will significantly alter the amount of freshwater discharging to the GoA along with changes in precipitation regimes.Currently, glacial fed streams have increased discharge and non-glacial fed streams have lower discharge.This paradigm will shift as coastal glacier coverage declines into the 21 st century.Today, glaciers act as a control on seasonal runoff variation within a catchment.Stream discharge within a glacierized basin varies little year to year and peak runoff is generally predictable.Precipitation fed streams, conversely, have higher interannual variations in discharge due to their susceptibility to interannual climate variability.For example, Fountain and Tangborn (1985) suggest that basins with glacial coverage around 36% have the lowest year-to-year variation in discharge.Over the seasonal to monthly timescales, variations are at a maximum in July and August for basins with less than 10% glacier cover (Fountain & Tangborn, 1985).
Glaciers contain vast stores of water as ice and the seasonal discharge from meltwater to the GoA delivers freshwater and vital nutrients (Milner et al., 2017;Neal et al., 2010;O'Neel et al., 2015).The flux of nutrients to the nearshore environment sustains many trophic levels (Arimitsu et al., 2018).Therefore, previous studies focused on river chemistry of the GoA primarily investigate major nutrients N and P, dissolved organic carbon (DOC) and particulate organic carbon (POC) (Fellman et al., 2014;Hood et al., 2015Hood et al., , 2020;;Hood & Berner, 2009).Past studies of micronutrients have focused on Al concentration in river plumes in the GoA (Brown et al., 2010) and Fe fluxes in glacierized and non-glacierized tributaries feeding the Copper River (Schroth et al., 2011).Brennan et al. (2014) examined grab samples of 61 streams across Alaska (13 from within the GoA region) for radiogenic Sr isotope ratios and concentrations of major cations, with the goal of broadly illustrating patterns in Sr isotope ratios and carbonate vs. silicate weathering.However, there remains a conspicuous lack of work examining the hydrogeochemical cycles of major elements, anions and micronutrients such as silica.Analysis of the major elements and anions in stream water can be used to reveal hydrologic properties within watersheds (Godsey et al., 2009).Additionally, investigating major element and anion concentrations may offer further insights into carbon cycling (Amiotte Suchet & Probst, 1993;A. F. White & Blum, 1995), chemical weathering-climate feedback (Eiriksdottir et al., 2013); and, specifically within the context of the GoA, how glacial coverage within a watershed affects physical and chemical weathering (Anderson, 2005;Sharp et al., 1995;Torres et al., 2017).
The chemical composition of glacier fed streams is distinct compared to non-glacier fed streams.Primarily, glacial melt runoff contains more K + and less Si compared to average global rivers (Anderson et al., 1997).The dominant ions within glacial fed streams are Ca 2+ , HCO 3 − , and SO 4 2− with high Ca 2+ :Si and low HCO 3 − : SO 4 2− (Tranter & Wadham, 2003).The major cation composition of glacier runoff is always Ca 2+ , controlled by the dissolution kinetics of carbonates, which are found within most bedrock (Raiswell, 1984).Further, specific discharge is a primary control on area-normalized weathering rates (yields) globally (Anderson et al., 1997;Maher & Chamberlain, 2014;Torres et al., 2017).
Rates of weathering, and both solute generation and flux, especially within glacierized catchments, are primarily controlled by stream discharge (Rose et al., 2018).Dilution curves, or concentration discharge (C-Q) relationships, are derived using a power law in the form of C=aQ b , where C is concentration, Q is discharge, and a and b are constants (Johnson et al., 1969).Different elements exhibit varying responses to increases in discharge (Godsey et al., 2009;Ibarra et al., 2016;Moon et al., 2014); this relationship is an important descriptor of the mechanisms controlling solute flux and weathering within a watershed.C-Q relationships can aid in elucidating future physical and chemical weathering regimes within the GoA watersheds which is particularly relevant under the current climate trajectory and projections of glacial recession because of more precipitation occurring as rain and increased air temperatures.
In this contribution we use physical watershed characteristics and USGS legacy water chemistry and discharge data from streams sites across the GoA to investigate how watershed characteristics affect solute yields and C-Q relationships.We explore the dissolved and suspended loads across stream sites to elucidate physical and chemical weathering patterns in the GoA.Additionally, we broadly show dominant weathering regimes and primary sources of solutes provided to stream water.This analysis is a unique approach to describe and explore controls on the generation and yields of solutes across the GoA.

Gulf of Alaska Region
The GoA watershed (Figure 1a.) and the seven sub-regions span a perhumid hydroclimate in the Southeast, to a subarctic hydroclimate in the northern regions.Precipitation varies considerably on a seasonal basis and across the region.Seasonal precipitation (Thornton et al., 2020) for the major regions of the GoA are shown in Figure 1c.The majority of precipitation occurs in the fall and winter with the Southeast, Central Coast and Prince William Sound regions receiving far greater precipitation than other regions.This is driven by orographic precipitation and primary storm track trajectories.Regions on the lee of the coastal ranges, Copper River, Knik Arm/Kenai Peninsula and Susitna River, receive far less precipitation compared to the windward.The glaciers of the GoA region occur primarily near the coast, however large glaciers also exist within the interior northern reaches of the Susitna and Copper River regions (Figure 1b).

Stream Sites and Watersheds
The USGS NWIS was queried to obtain all stream sites within Alaska that contain water chemistry analysis.We then spatially filtered the dataset to find stream sites that drain into the GoA.Each stream site was further spatially filtered and assigned a region based on the seven defined regions (Figure 1a).The stream sites contained within the GoA were filtered to those sites that have total dissolved solute (TDS) and total suspended sediment (TSS) paired with stream discharge.To explore concentration-discharge relationships we chose stream sites that have at least 12 paired concentration and discharge measurements for HCO 3 − , Ca 2+ , Mg 2+ , Na* (Cl-corrected; Ibarra et al., 2016;Moon et al., 2014), K + , SO 4 2− , SiO2 and TSS resulting in 34 stream sites total that fit our criteria.For watershed boundaries of the 34 stream sites with requisite C-Q data we use the watersheds delineated by Curran and Biles (2020).

Watershed characteristics
We calculated physical and climate watershed characteristics based on five separate datasets: landcover (CCRS, 2015), geology (Garrity & Soller, 2009), elevation statistics (Tadono et al., 2014), glacial coverage (RGI Consortium, 2017) and climate (Thornton et al., 2020).Each dataset was clipped using watershed boundaries and the specific parameters were grouped and calculated.A full description of the watershed characteristics and analysis is available in the supporting information (Text S1).

Weathering Regimes
We broadly investigate physical and chemical weathering regimes within watersheds of GoA by exploring the yield of TDS and TSS and the HCO 3 − :SiO2.Because we do not correct our data for rainfall contributions or carbonate weathering the TDS:TSS relationship illustrates a relative weathering index (RWI).For the calculation of TDS, we summed concentrations of Ca 2+ , Mg 2+ , Na + , K + , and SiO 2 .We compare mean TDS and TSS yields of global rivers (Gaillardet et al., 1999) along with the mean TDS and TSS yields of the GoA streams.Similarly, we show all  ] and [SiO2] with global rivers for comparison.For this study we do not correct for rainwater contribution because our main objective is to investigate solute generation and yields and to broadly appropriate sources of solutes.

Solute Yield and C-Q Relationships
We calculate solute yields by multiplying the solute concentrations by the instantaneous discharge measurements and dividing by the watershed area.Where we present mean yield values for each stream we first calculate the yields for each measurement and calculate the mean of each calculated yield value at each stream site.
To explore C-Q relationships we fit the chemistry and discharge data of the 34 sites to the power law function C=aQ b ('nls2' R package, Grothendieck, 2013; following code modified from Ibarra et al., 2016, 2017and Wymore et al., 2017) where C is the concentration, Q is discharge and a and b are constants.The value of b (the power law exponent or b-value), in logC-logQ space, is the slope of the resulting linear fit.The slope of the line (b-value) has a physical interpretation (Godsey et al., 2009).Slopes near -1 indicate simple dilution, and that concentration varies inversely with discharge.A slope around zero (typically -0.1 to +0.1) indicates chemostatic behavior within a watershed, meaning, as discharge increases the concentration remains relatively unchanged.Power law exponents greater than zero indicate enrichment of a solute as discharge increase.We calculated fits for HCO 3 − , Ca 2+ , Mg 2+ , Na*, K + , SO 4 2− , SiO2 and TSS.However, we do not use the b-values for K + and SO 4 2− in further analysis due to the poor quality of the modeled fits.

Statistical Analysis
To determine the relationships between physical watershed characteristics, solute yields, and C-Q relationships we use Bayesian Information Criterion (BIC) to select the minimum number of watershed characteristics that will best describe the variation in a solute yield or the power law exponent.Each of the physical and climate-based watershed characteristics within Table S2 were used in the BIC analysis.We then use the minimum numbers of parameters selected by the BIC in a multiple linear regression model to find how well the parameters describe the variation in solute yields and b-values.For each parameter used in the linear regression we show whether the solute yield or b-value is positively or negatively related and if the derived slope with respect to a given parameter is statistically significant (p<0.05).

Results and Discussion
3.1 Solute C-Q Relationships On a regional basis streams draining to the GoA exhibit variable C-Q relationships.The C-Q relationships of HCO 3 − , Ca 2+ , and SiO2 are used as representative solutes for rock-derived nutrients (Figure 2).The C-Q relationships of the solutes across the regions may be influenced by climate regimes (Godsey et al., 2019), lithology (Ibarra et al., 2016), geomorphology (Torres et al., 2015) and vegetation (Wymore et al., 2017).Median b-values at a regional level for HCO 3 − , Ca 2+ , and SiO2 range from -0.23 -0.01, -0.24 --0.07, and -0.25 --0.11 respectively (Table S2).The Southeast and Susitna River regions have the lowest median b-values, and the W. Cook Inlet/Kodiak region has the highest median b-values.Basins within the Susitna River region have 189 shallow slopes which may contribute to the lower overall b-values (e.g.Torres et al., 2015; discussed further below).In the Southeast region the basins receive a greater amount of annual precipitation, which may contribute to high mean specific discharge leading to lower b-values (e.g.Godsey et al., 2009Godsey et al., , 2019)).Within the W. Cook Inlet/Kodiak region the basins are relatively steep with high relief which may contribute to the higher observed b-values.Additionally, because historical USGS sampling efforts do not generally capture the highest discharge events, the available dataset is likely biased such that overall C-Q relationships may appear to be more chemostatic than if higher discharge events were more represented.

Weathering Regimes
Gulf of Alaska streams span a large range in dissolved solute and sediment yields.Figure 3a illustrates the TDS and TSS data for 34 of the GoA streams, with linear lines indicating the ratio of the load being exported in the dissolved phase vs. the suspended phase.The RWI describes the ratio between TDS and TSS yields that are not corrected for carbonate dissolution and represents the partitioning of landscape denudation into chemical vs. physical weathering, which on long timescales, at steady state, is set by the uplift rate (eg., (Dixon and von Blanckenburg, 2012;Gabet and Mudd, 2009;Hilley et al., 2010;Maher and Chamberlain, 2014;Waldbauer and Chamberlain, 2005).For example, in Figure 3a points plotting above the 0.75 line instantaneously export a majority (i.e., 75%) of the load as dissolved solids.Points plotting below the 0.1 line export a majority of the load as suspended sediment.Importantly, points plotting below 0.5 have lower inferred chemical weathering compared to physical weathering.Erosion by glaciers combined with steep topography produces enhanced physical weathering in the GoA watersheds.However, at a given erosion rate, chemical weathering is on average lower within catchments of the GoA compared to global rivers.
Mean TDS and TSS yield ratios of the selected streams of GoA are lower in overall fluxes (area-normalized) than many of the global rivers (Figure 3b).Compared to other northern high latitude rivers (Ob, Kuskokwim, Yukon, and Mackenzie) within the global dataset the GoA streams have similar RWI values.Globally, the suspended load of rivers dominates the material flux to the oceans (Walling & Webb, 1983).A positive trend between dissolved and suspended load exists, however a clear relationship is not systematic (Walling & Webb, 1983).Suspended sediment loads have been shown to have a negative relationship with basin area, a trend not shared by the dissolved load.Within the GoA dataset, and most notably within the Southeast region (Figure 3a-b), streams near the RWI = 0.1 line are generally the glacierized basins and streams near and above the RWI = 0.75 are non-glacierized (Figure 3b).Chemical weathering is limited within glacierized basins though there is an increase in physical erosion.We infer that increased generation of fresh mineral surfaces within glacierized basins is being outpaced by the time (i.e., kinetics) required for chemical weathering reactions to occur (Ferrier & Kirchner, 2008;Torres et al., 2017).Bicarbonate and SiO2 concentration ratios provide a mechanism to investigate the dominant chemical weathering sources within a basin at a broad scale.The GoA streams are dominated by high HCO 3 − :SiO 2 values (Figure 3c).Carbonate dissolution supplies the majority of ions to the dissolved load compared to silicate weathering regardless of primary bedrock lithology (Raiswell, 1984).Two non-glacial fed streams in the Knik Arm/Kenai Peninsula region and one non-glacial fed stream in the Susitna River region, are primarily influenced by silicate weathering (based on the dominant bedrock geology of these catchments; Table S2) with HCO 3 − :SiO2 values near 2:1, typical HCO 3 − :SiO2 of silicate net weathering reactions range from <1:1 to 4:1 (e.g., Bouchez and Gaillardet, 2014;Ibarra et al., 2016;Maher, 2011;Winnick and Maher, 2018).Furthermore, the ratio between HCO 3 − and Ca 2+ , illustrated in Figure 3d, is consistent with dissolution of carbonate by carbonation (Equation S1).Within glacierized basins the lithology is a secondary control on the chemical composition of stream water (Raiswell, 1984;White et al., 2001).Further, decreased temperatures within glacierized basins reduce silicate-weathering rates (Anderson, 2005).Additionally, silicate denudation rates in glacierized basins are lower compared to non-glacierized basins at similar specific discharge values (Anderson et al., 1997).Therefore, glacierized basins within the GoA do not export more area normalized SiO2 compared to non-glacierized basins.
Compared to mean global values, the streams of the GoA generally have lower mean HCO 3 − :SiO2 (Figure S5).Mean HCO 3 − :SiO2 values for high latitude northern rivers are also slightly elevated compared to GoA streams.We speculate that this could be driven by one of two factors.First, a lack of monolithologic carbonate catchments in the GoA.Second, lower overall primary productivity (due to colder temperatures and shorter growing season), leading to reduced soil CO2 values and thus less HCO 3 − generation thereby lowering the HCO 3 − :SiO2 values (all else being equal) (Raymond & Cole, 2003).

Geomorphology and Glacier Coverage
To explore controls on the above observations we rely on our BIC analysis.We find through BIC and multiple linear regression that solute yields are controlled by glacier coverage while C-Q relationships are primarily controlled by watershed geomorphology (Table S1).Earlier efforts have linked watershed geomorphology and C-Q relationships, suggesting the differences in transport vs. supply regimes are controlled by changes in fluid-residence times in the weathering zone (Maher, 2011;Maher & Chamberlain, 2014;McGuire et al., 2005;Torres et al., 2015).For example, we find that b-values for HCO 3 − , Ca 2+ , and SiO2 and watershed slope are positively correlated (Figure 4a).Fluid transit time distributions within watersheds with steeper topography are likely greater than watersheds with shallower slopes.Longer transit times may allow the concentration of the solutes to reach equilibrium with respect to the net weathering reactions (Maher, 2011;Maher & Chamberlain, 2014;Torres et al., 2015).The BIC analysis only identified glacier coverage as an important parameter with respect to the b-value for SiO2. Figure 4b illustrates the relationship between b-values and glacier coverage.We find this relationship to be undefined, however, there are some surprising observations.The median bvalues for HCO 3 − , Ca 2+ , and SiO2 of the glacier fed streams are -0.15,-0.19, -0.19 respectively, indicating solute concentrations are not significantly affected by changes in discharge (areanormalized).Further, there are several glacier fed streams in our dataset that exhibit at or near chemostatic behavior (b-values between -0.1 and 0.1) (e.g., Hindshaw et al., 2011Hindshaw et al., , 2014)).Broadly, the patterns in C-Q relationships within the glacierized basins generally conform to observations in non-glacier fed systems (e.g.Godsey et al., 2009).
Figure 4c illustrates the relationship between the yield of SiO2, Ca 2+ , and HCO 3 − and glacier coverage.The yield of these weathering derived solutes exhibits positive relationships with glacier coverage.Previous work examining weathering rates versus watershed glacier coverage showed a negative relationship across a global dataset (Torres et al., 2017).Within the GoA dataset, glacier coverages range from 0% to 40%, while the dataset presented by Torres et al. (2017) ranges from 0% to 95% (Figure 4d).Cation yields within glacierized basins of the GoA have a distinct positive relationship with glacier coverage (Figure 4d).Though we do not derive a relationship for the globally derived dataset, we propose that there may be two distinct trends when examining the cation yield vs. glacier coverage relationship.A threshold appears to exist at glacier coverages < 35-40% in which a positive relationship exists between cation yields and glacier coverage.Basins with >40% glacier coverage may have a slight negative relationship between cation yields and glacier coverage (and drives the relationship fit to the dataset by Torres et al., 2017;their Figure 1b).We suggest that the proportionally small area of proglacial environments in these catchments may be a limiting factor for solute generation in basins with high glacier coverage.As such, proglacial environments may be an important driver of continued weathering (Deuerling et al., 2017).Additionally, the lithology of the GoA watershed may contribute to the strong positive relationship observed between cation yields and glacier coverage.Metasedimentary rocks that are highly fractured are more easily weathered compared to granite or basalt bedrock (Bluth & Kump, 1994).Furthermore, future work is needed to partition solutes into specific lithologic source components (e.g. an inverse weathering model; Gaillardet et al., 1999).This will aid in further elucidating how modifications in the physical watershed characteristics due to glacial recession may alter fluxes of critical nutrients to the GoA.

Conclusions
We find that TDS and TSS yields of streams in the GoA cover a large range of values, implying variable rates of chemical and physical weathering, and variable partitioning of overall denudation.Streams fed by glaciers tend to have lower chemical weathering with higher amounts of physical weathering (lower RWI values).Conversely, streams with low to no glacier coverage have higher RWI values.In terms of dominant weathering regimes, carbonate dissolution provides the majority of solutes to streams in the GoA regardless of dominant lithology.
Based on our statistical analysis of solute yields, C-Q relationships, climate, and physical watershed characteristics we find that geomorphology and glacier coverage affect the yields and solute generation differently.Power law exponents are primarily explained by geomorphologyrelated landscape characteristics such as slope, elevation, and aspect.We infer that water transit times are the primary control on solute generation and influence the observed C-Q relationships (Torres et al., 2015).Solute yields on the other hand are controlled by the amount of glacier coverage within a catchment.This is unsurprising given that glaciers are a primary control on stream runoff within a given catchment (Fleming, 2005;Fountain & Tangborn, 1985) and that most solute yields increase with increases in runoff.However, a novel result of this work is that we demonstrate that physical watershed characteristics control solute generation (C-Q relationships) and glacier coverage control area normalized fluxes (yields).This new conceptualization for understanding solute delivery to the GoA indicates that glaciers ultimately control the absolute amount of solutes exported to the ocean due to runoff generation, however downstream characteristics control the generation of solutes displayed in C-Q patterns.Future work will need to aim at estimating annual fluxes of solutes and to forecast how solute yields and fluxes may change as glaciers continue to rapidly recede in the coming century.
characteristic.Text S2 explains the methods and results of the Bayesian Information Criterion (BIC) analysis and multiple linear regression.

Text S1.
The physical and climate-based watershed characteristics were calculated in QGIS 3.16.5 and with custom Python scripts.For watershed characteristics derived from raster datasets (elevation statistics, landcover, precipitation, and temperature) we used the Zonal statistics tool within the QGIS Raster analysis toolbox.To calculate watershed characteristics derived from vector data (geology and glacier coverage) a custom Python script was used to clip the input data with each watershed and calculate the percent area each parameter covers within the given watershed.Landcover, geology, and glacier coverage are reported in Table S2 as percent coverage.Climate parameters, mean precipitation and mean temperature, are averages within each watershed calculated from 1981 to 2010 using the DAYMET dataset (Thornton et al., 2020) accessed using Google Earth Engine.

Text S2.
Results of the BIC and multiple linear regression are shown in Table S2.The BIC analysis provides a list of the minimum number of parameters (watershed characteristics) that explain the variation of the given variable (b-values and solute yields).The resulting list of parameters are then given to a multiple linear regression (MLR) model.In Table S2 the column "Parameters Chosen" show the parameters used in the MLR along with the sign of the relationship with the given variable in parentheses.For each solute b-value and solute yield we provide the R 2 value for the MLR model.We also show what parameters have statistically significant (p>0.05)relationships in the final column.Dataset S1.Includes three shapefiles for the USGS watershed boundaries and USGS stream sites used in this study, and the regional boundaries of the GoA.

Figure 1 .
Figure 1.Regions of the Gulf of Alaska (a), regional glacier coverage (b), and associated monthly mean precipitation and air temperature (c).The coastal regions of Prince William Sound, Southeast and Central Coast receive the greatest amount of precipitation.Mean summer air temperatures across regions are similar while mean air temperature of the Copper River and the Susitna River regions are notably colder that the other regions.
] for the GoA streams and the mean [HCO 3 −

Figure 2 .
Figure 2. Concentration-discharge relationships for SiO2, Ca 2+ , and HCO 3 − for each of the seven regions of the GoA.An example plot (lower right) illustrates the physical interpretation of a given C-Q slope.Streams within the GoA region do not conform to simple dilution behavior based on intermediate slopes.This indicates that streams throughout the GoA behave similarly to streams within the conterminous United States (e.g.Godsey et al., 2009) why do we make this comparison here?All other comparisons are to global rivers.Regional difference of C-Q relationships may indicate different controls such as hydroclimate and geomorphology.

Figure 3 .
Figure3.TDS yield and TSS yield for all GoA sites with available data (a).Mean TDS and TSS yields for each stream site with data from global rivers(Gaillardet et al., 1999) (b).Symbols indicate the range of percent glacier coverage within each watershed.Watersheds with medium to high glacier coverage generally plot below the RWI = 0.1 implying that chemical weathering within moderately to highly glacierized catchments is relatively lower compared to physical weathering.Concentration of HCO 3 − vs. SiO2 for GoA stream sites (c).Three streams plot on or below the 2:1 line indicating that silicates are the primary weathering minerals.The rest of the streams plot above the 4:1 line suggesting carbonates as the primary weathering minerals.HCO 3 − vs. Ca 2+ concentration this is oppositely stated in the SI (d) for streams primarily plot near the 1:0.5 line indicating carbonate dissolution by carbonic acid.

Figure 4 .
Figure 4. Relationship between power law exponents and mean catchment slope (a) and glacier coverages (b) for SiO2, Ca 2+ , and HCO 3 − .Mean slope is used as an example geomorphic parameter.The b-values for SiO2, Ca 2+ , and HCO 3 − show a positive relationship with mean slope, while there is no relationship of these b-values and glacier coverage.However, calculated yields for SiO2, Ca 2+ , and HCO 3 − .vs. glacier coverage (c) show a positive relationship.Cation yields vs. glacier coverage (d) for the GoA streams with data compiled by Torres et al. (2017) for glacierized watersheds globally.These results indicate that to a certain extent, glaciers affect the solute yield, while watershed slope (and other geomorphic characteristics) control solute generation.Asterisks next to R 2 value denotes where p < 0.05.

Figure S1 .
Figure S1.Locations of USGS stream gage sites used for this study.Stream sites were chosen if concentrationdischarge measurements equal or exceeded 12 paired measurements.There is a notable lack of stream gages locations within the Central Coast region, likely due to the rugged terrain and abundance of tidewater glaciers.

Figure S2 .
Figure S2.Watershed boundaries used in this study to calculate watershed characteristics along with USGS stream gage sites.Several watersheds within the Susitna River region are nested within each other and do not appear within this map.

Figure S3 .
Figure S3.Concentration of HCO 3 − vs. SiO 2 for GoA stream sites plotted by region.Across the various regions of the GoA carbonate dissolution appears to be the primary weathering regime compared to silicate weathering.

Figure S4 .
Figure S4.Cross plots of TDS yield vs. TSS yield plotted separately by region.Streams within watersheds with medium to high glacier coverage generally have lower rates of chemical weathering compared to physical weathering.Streams with low to no glacier coverage have higher rates of chemical weathering compared to physical weathering.