in riparian groundwater linked hydrological pathways in the boreal forest?

. The riparian zone, or near-stream area, plays a fundamental role for the biogeochemistry of headwaters. Here wet, carbon-rich 10 soils can change water chemistry before it enters the stream. In the boreal forest, the riparian zone plays an especially important role in the export of dissolved organic carbon (DOC) to streams. However, the riparian zone is not uniform and spatial variability of riparian groundwater hydrology and chemistry can be large. Terrestrial topographic depressions create hydrological pathways towards focal points in the riparian zone, which we refer to as Discrete Riparian Inflow Points (DRIPs). Combining the chemical function of the riparian zone and the convergence of hydrological pathways, we hypothesize that 15 DRIPs play a disproportionally large role in conveying DOC to small streams. Earlier work has demonstrated that runoff from DRIPs can make up the majority of riparian flow contributions to streams, but so far it is unknown how their groundwater chemistry differs from the rest of the riparian zone. We therefore ask the question: are DOC concentrations in riparian groundwater linked to hydrological pathways in the boreal forest? To answer this question we sampled riparian groundwater during six campaigns across three boreal headwater streams in Sweden. The groundwater wells were distributed in ten DRIP 20 and non-DRIP pairs (60 wells), following transects from upland (20 meters lateral distance from the stream bank) to near stream area (<5 meters lateral distance). The variability in dissolved organic carbon (DOC), pH and electrical conductivity (EC) was analyzed using linear mixed effect models (LMM). We explained the variability using three factors: distance from the stream, seasonality and DRIP/non-DRIP. Our results showed that DRIPs provided DOC rich water (34


Introduction
Headwater streams can be seen as the capillaries of the landscape: although small in appearance, collectively they make up the majority of a stream network. The rich variety in hydrology, biology and chemistry of headwaters is tightly connected to processes in their catchments (Bishop et al., 2008;Hunsaker and Levine, 1995). Lateral groundwater inputs account for a large part of the streamflow of small streams, magnifying groundwater controls on stream CO 2 emissions (Hotchkiss et al., 2015). 5 These controls are governed by groundwater-surface water exchange in the last interface between the landscape and stream ecosystems (Hayashi and Rosenberry, 2002). This near-stream area, so called riparian zone (RZ), holds important functions such as chemical transformation of hillslope water (Cirmo and McDonnell, 1997), thermal regulation (Davies-Colley and Rutherford, 2005) and erosion control (Smith, 1976). A few characteristics of the boreal RZ that leads to its unique ecosystem functions are high groundwater levels, dynamic redox potential, build-up of soil organic matter, and diverse vegetation (Grabs 10 et al., 2012;Kuglerová et al., 2014b;Lidman et al., 2017). In terms of the hydrological role of the RZ, it has been demonstrated that riparian water dominates streamflow generation, instead of event-based water contributions from hillslopes (McGlynn and McDonnell, 2003). Combined with the chemical transformation of water in the riparian zone, stream biogeochemistry is therefore largely controlled by riparian zones (Ledesma et al., 2018;Lidman et al., 2017). However, RZ's are not homogenous strips surrounding surface waters, but contain an array of heterogeneities in hydrogeology, soil development and vegetation 15 across small spatial scales (Buttle, 2002;Kuglerová et al., 2014b). Moreover, wetness state changes the chemical function of the RZ (Vidon, 2017). It is therefore important to further investigate which parts of the riparian zone matter most for element transport, stream flow generation and associated biogeochemical processes.
In hydrological models streamflow generation has often been conceptualized as a diffuse process, which limits the ability to 20 express points of focused groundwater discharges (Briggs and Hare, 2018). Some models, such as the RIM model and DSL concept, have considered the vertical heterogeneity in riparian groundwater fluxes to boreal streams (Ledesma et al., 2015;Seibert et al., 2009). But also longitudinally along streams reaches it is necessary to account for hydrological and biogeochemical heterogeneity within the RZ. For example, permanently saturated riparian areas have been identified as main stream flow generators (Penna et al., 2016), and have been associated with denitrification, as well as retention and 25 transformation of (labile) OM, compared to drier, oxic, riparian soils (Blackburn et al., 2017;Burgin and Groffman, 2012;Ledesma et al., 2018). In terms of vegetation, groundwater discharge zones are hotspots for diversity (Kuglerová et al., 2014a).
Although these studies show that heterogeneity in the saturation or wetness conditions could be good predictor for heterogeneity in soil chemistry, the connection between spatial variability in groundwater chemistry and hydrological pathways within the riparian-upland continuum has not been demonstrated. The hydrological connection between the upslope 30 catchment, riparian zones and consequently the stream network are highly variable: where some parts of the riparian zone only drain small individual hillslopes, others function as main hydrological flow paths funneling subsurface water through riparian input zones (Leach et al., 2017). Combining their chemical signature and hydrological upslope connectivity, contributions of such focused riparian inputs could therefore function as important control points in the landscape (Bernhardt et al., 2017). The difficulty is that incorporating these control points into models or practical applications means that they have to be 35 characterized in order to explain stream dynamics. Especially for informing distributed models that overpass catchment scale, determination and characterization of these control points remains one of the challenges for the scientific community (Briggs and Hare, 2018).
For the hydrological characterization of riparian inputs, various approaches can be used across scales. Although subsurface 40 pathways do not entirely follow surface topography (Devito et al., 2005), it has been demonstrated that topographic depressions are a good indicator for accumulation areas of water, ponding, shallow groundwater tables and concentrated flow paths in the near-stream area (Ågren et al., 2014;Jencso et al., 2009;Wallace et al., 2018). As such, topographic models can predict where along a stream network disproportionally large amounts of groundwater connect with the stream. Mixing models using water temperature and chemistry can further depict whether the topography-based predictions of focused riparian inputs to streams are in line with reality (Leach et al., 2017). However, these so called Discrete Riparian Inflow Points (DRIPs, Fig. S2), provide continuous flows of subsurface water during low flow periods, but have also been observed to be highly dynamic in their activation during hydrological events (Ploum et al., 2018). Contrary to the incorporation of water from ephemeral streams in 5 perennial stream networks, or the connection of intermittent sections of a stream network (Ågren et al., 2015), DRIPs are dominated by subsurface flowing water and the discharge to the stream is the first exposure to an open channel. A recent study demonstrated temporal dynamics in increasing greenhouse gas evasion from the stream reach in close downstream proximity of DRIPs (Lupon et al., 2019). Also in Arctic systems has the presence of riparian wet areas partially explained stream CO2 evasion (Rocher-Ros et al., 2019). The latter suggests that both the hydrological fluxes as well as biogeochemical reactions in 10 the stream are associated with the hydrological activity of DRIPs. However, in order to determine whether DRIPs matter for stream biogeochemistry, chemical characterization of the discharging groundwater is needed.
Characterizing groundwater chemistry is an especially challenging task. Previously this challenge has been by-passed by assuming that groundwater is a well-mixed source of water (Kirchner, 2003), or by inferring groundwater chemistry from base 15 flow chemistry of streams (Peralta-Tapia et al., 2015). Also the RIM model has provided a framework to infer groundwater chemistry profiles from stream chemistry (Seibert et al., 2009). However, even at the local scale spatial variability in groundwater chemistry overrules temporal variation and requires regular sampling of extensive well networks (Kiewiet et al., 2019). Within meters of each other, groundwater signatures can vary greatly (Penna et al., 2016). Three key parameters for chemical characterization of groundwater in boreal forests are dissolved organic carbon (DOC), pH and ionic strength. DOC 20 concentrations in groundwater is the result of interaction between water and carbon rich materials in the shallow subsurface environment that are associated with paludification (Lavoie et al., 2005). More specifically for near stream areas, the width of the riparian zone is associated with the size of the potential carbon pool and the subsequent DOC concentrations (Ledesma et al., 2015). Apart from its role in food-web structures and carbon transport, DOC also increases the acidity (decrease pH) of soils and surface waters . Electrical conductivity (EC) can be used as a proxy for the ionic strength, or 25 total amount of dissolved ions in water (Corwin and Lesch, 2005). Water contact time with minerals and weathering processes are important factors determining EC (Saarenketo, 1998), with increasing EC indicating longer interactions (Hayashi, 2004;Peralta-Tapia et al., 2015).
In the context of spatial variability of riparian groundwater chemistry, it can be expected that DOC, pH and EC differ between 30 DRIP and non-DRIP riparian areas (Fig. 1). DRIPs are associated with high groundwater levels and wet, organic rich soils with vegetation that favors wet conditions, while non-DRIPs have drier top soils and deeper groundwater levels (Kuglerová et al., 2014a). Inherent to their topographic setting, DRIPs drain a large upland area, while non-DRIPs typically drain only a small surrounding area of the riparian zone or they are recharge zones for adjacent DRIPs. Moreover, the water in DRIPs travels a longer distance horizontally; in presumably wet, highly permeable, organic rich soil. Non-DRIP water, on the other 35 hand, is likely to infiltrate vertically through an oxic, organic rich top soil, before being transported a relative short distance horizontally through supposedly more mineral substrate. This implies that the contact time of the water with wet, organic soil and drier, mineral soil is different for both cases, which should lead to contrasting water chemistry. In this study we characterize groundwater in a paired well network that is specifically designed to incorporate (saturated) riparian areas with large contributing areas (DRIPs) and drier parts of the riparian zone with small contributing areas (non-DRIPs). We hypothesize 40 that groundwater in DRIPs has higher DOC concentrations and lower pH compared to non-DRIPs. The deeper groundwater levels in non-DRIP areas, and longer contact times with mineral soil relative to organic soil, leads us to expect that EC will be higher in non-DRIP water compared to DRIPs.

Material and methods
To test our hypothesis we collected DRIP and non-DRIP groundwater across a riparian gradient during different seasons. 10 Using linear mixed effect models (LMM's) we analyzed the role of DRIPs on biogeochemical composition of riparian groundwater, in relation to spatial and temporal variability. We performed our study in Krycklan, a boreal forested catchment in northern Sweden.

Study area
The Krycklan catchment is situated near the town of Vindeln, Sweden (64°14´N, 19°46´E, Fig. 2). The bedrock is 15 predominately Svecofennian metasediments and metagreywacke. Quaternary deposits consist mostly of till (51%) and sorted sediments (30%). Land cover is dominated by forest (87%), and there is 9% mire cover. Furthermore there are sporadically thin soils and bedrock, and a small fraction of arable land (2%). The climate is characterized as cold humid temperate type, with almost 6 months of snow cover. The yearly average temperature is 1.8 °C, and annual precipitation is 614 mm, and the annual mean runoff approximates 311 mm (Laudon et al., 2013). The well network is situated along streams referred to as C4, 20 C6, and C8 (Laudon et al., 2013), with a drainage area of respectively 18, 110, and 230 ha. Catchments C4 and C6 have been widely studied in regard of lateral flow and groundwater and surface water interaction and can be referred to in other studies as Kallkälsmyrsbäcken and Stortjärnsbäcken (Laudon et al., 2004b. Flows vary from a few liters per second baseflow to 200 l/s peak flows (Ploum et al., 2018). The yearly hydrograph is characterized by sustained baseflow throughout the winter months, followed by spring snowmelt floods in April and May (Fig. S1). In summer and autumn low flow conditions are 25 common with occasional rain-induced flow events. From November onwards, flow reduces as temperatures fall below 0 °C and baseflow conditions set in.

Site selection and sampling well infrastructure
Discrete Riparian Inflow Points (DRIPs) were selected by considering wet areas, based on a topographic wetness index, and selecting large step changes in catchment area along stream networks using flow accumulation algorithms (Ågren et al., 2014;Beven and Kirkby, 1979;O'Callaghan and Mark, 1984). The DRIPs (n = 10) were selected with contributing upslope area varying from 0.6 to 7.7 ha, with a mean contributing area of 2.7 ha. Non-DRIPs had an upslope contributing area between 4 10 and 80 m 2 (on average 17 m 2 ). The DRIPs have been field-validated and surveyed on species richness (Kuglerová et al., 2014a).
For some sites chemical and thermal signatures further corroborated the location were riparian water discharged into the stream (Leach et al., 2017).
The setup of this study consists of a well network with a total of 60 fully screened PVC wells (30 mm diameter) arranged in 15 10 paired transects. Each transect consisted of a riparian well, situated typically between 1 and 5 meter from the stream, a transition well at approximately 10 meters from the stream, and an upland well 20 meters from the stream. Transects followed the local topography, to approximate local hydraulic gradients and flow paths. The non-DRIP transects were installed close (<50 m) to each DRIP transect to ensure similarity in local conditions. All wells were drilled until resistance, or an aquitard layer. Riparian wells had a mean depth of 95 cm, transition wells 99 cm, and upland wells 121 cm. We assumed that the water 20 sampled from the well is a weighted average of the phreatic aquifer, down to the depth of the well. Given the exponentially decaying hydraulic conductivity with depth, this assumption would imply that, under fully saturated conditions, the majority of the water is therefore from the upper soil layers, referred to as the dominant source layer (Ledesma et al., 2015). Given that context, lateral flow below the well bottom was considered negligible compared to the flow in the vertical domain of our well installations. For a small subset of riparian sampling wells, water levels were available from directly neighboring wells (<2 m apart). Figure S1 shows an exemplar time series of these wells and a hydrograph for 2018. The mean depth to water table for those time series was 9.6 cm for DRIP wells, and 54.5 cm for non-DRIP wells during the year of 2018.

Groundwater sampling and chemical analysis
The well network was sampled using suction cup lysimeters and vacuumed glass bottles (Blackburn et al., 2017). The wells 5 were pumped before installing the suction cups to ensure water from the aquifer was sampled and without any stagnant well water. The bottles were collected after approximately 24 hours and subsampled, filtered and analyzed within 48 hours. In addition, a more intensive sampling campaign was conducted for a series of riparian wells only. These were sampled following a similar protocol, but instead of suction cup lysimeters, a peristaltic pump was used for the collection of water samples.

10
Water samples were collected during spring, summer and autumn of the 2016 and 2017 hydrological years. In total 359 samples were analyzed from six sampling campaigns, of which 200 from DRIP wells and 159 from non-DRIP wells. Non-DRIP wells occasionally had too low water level to collect a representative water sample. For analysis of dissolved organic carbon (DOC), a subsample was filtered (0.45 µm) into acid-washed high-density polyethylene bottles (rinsed three times) and kept at 4 °C before laboratory analysis. DOC was measured by acidifying the sample and combustion using a Shimadzu TOC-VPCH. The 15 pH and EC were subsampled without headspace into acid-washed high-density polyethylene bottles (rinsed three times) and kept at 4 °C before laboratory analysis. Samples were analyzed using a Mettler Toledo DGi117-water probe for pH and Mettler Toledo InLab741 probe for electrical conductivity.

Statistical analysis
We used linear mixed-effect models (LMM) to analyze patterns in DOC, pH and EC. The analysis was performed in R using 20 lmer models from the R-package lme4 (Bates and Maechler, 2009;Bates et al., 2014). The LMM's provided a non-parametric approach to explain variability in the response variables by fixed effects (factors that were included in the study design) and random effects. Random effects account for factors which were not part of the study design, but possibly affected variability in DOC, pH and EC. The fixed effects considered in this study were the hydrological pathways (HP -DRIP, non-DRIP), position in the landscape relative to the stream (POS -riparian, transition, upland), season when the samples were taken summer,autumn), and the two-way interaction between HP and POS and TIME respectively. The included random effects were the stream identity and the transect identity along which the wells were situated. In this way we accounted for specific catchment and hillslope properties. The model structure selection was based on the lowest AIC (Akaike's Information Criterion). We evaluated the model performance using Type II Wald F tests with Kenward-Roger degrees of freedom (since all 30 explanatory variables are factors). F statistics indicate the explained variance as a ratio of unexplained variance. An effect was considered significant if p-values <0.05. We evaluated the assumption of Gaussian distribution of errors by inspecting residuals and quantile distributions. For DOC five outliers, and for pH two outliers were removed from the upper quantile. For EC one in the lowest tail and two in the highest tail of the distribution. For comparing contrasts of levels within explanatory factors (for example DRIP vs. non-DRIP comparisons), we investigated least square means using R-package lsmeans, including Tukey 35 adjustment to account for potential differences in sample size (Lenth and others, 2016). Furthermore, the marginal and conditional coefficient of determination (R 2 mar and R 2 con ) was presented to compare explained variance by the fixed effects, and the variance explained by the fixed and random effects together (R-package MuMln) (Barton, 2014).

DOC
The water collected in wells situated in the DRIPs had a higher mean DOC concentration (33.9 mg/l) compared to non-DRIP wells (19.9 mg/l, Fig. 3). DOC concentrations in DRIPs increased upland wells towards the riparian wells from 29.2 mg/l to 36.3 mg/l, while in non-DRIP riparian wells DOC concentration increased from 16.4 to 20.1 mg/l (DF=19, p=0.03). When we 5 only accounted for gradients between upland and riparian wells (without distinction between DRIP or non-DRIP sites) differences were not as large, but still significant (from 22.8 mg/l to 28.2 mg/l, DF=327, p=0.0001). Although DOC concentrations in the upland groundwater of DRIPs was already high, the overall gain in DOC concentrations from the upland to the riparian wells was most accountable to DRIPs (DF=326, p=0.0003). Average DOC concentrations were contrasting in the upland wells (29.2 mg/l and 16.4 mg/l for DRIP and non-DRIPs, respectively), but were statistically not significant 10 (p=0.1844, Fig. 4 upper left panel). In summer and autumn, DOC concentrations in DRIP groundwater (36.4 and 33.3 mg/l) were twice as high as non-DRIP groundwater (18.0 and 17.7 mg/l, Fig. 5, upper panels). However, during snowmelt in spring, this difference reduced. This change was a result of an average 20% decrease in DOC concentrations in DRIPs (28.5 mg/l) compared to the summer average. In non-DRIP areas there were no significant contrasts, although there was a small, increase in spring (21.6 mg/l) compared to summer and autumn (18.0 and 17.7 mg/l, p=0.4986 and p=0.3019). Overall, the fixed effects 15 alone explained 22% of the variance in DOC found in the groundwater well network. With the random effects included the explained variance was 68%.

pH
Although typically associated with DOC, the pH was not as distinctly different between DRIP and non-DRIP water as DOC (Fig. 3). Overall the fixed effects accounted for 13% of the variance, and 55% including random effects (Table 1). Mean pH levels were 5.38 for DRIPs and 5.66 for non-DRIPs (DF=16, p=0.2). Instead, position in the landscape had more effect on the variability in pH: the upland pH was similar at DRIPs and non-DRIPs and decreased towards the riparian area (5.66 to 5.40, 25 P<0.0001). Although no significant effect was found for interaction between the landscape position and hydrological conditions (Table 1), the least square means analysis showed a pronounced decrease in pH from upland to riparian wells in the DRIP areas (5.57 to 5.19, P<0.0001, Fig. 4 middle panels). The second important explanatory variable was seasonality (TIME in Table 1). The most notable was the increasing pH from the summer to autumn (5.37 to 5.70, P<0.0001), both in DRIP and non-DRIP areas (Fig. 5, center panel and center-right panel). In the transition to spring, pH decreased again (pHspring=5.48), 8 mostly due to a shift in the DRIPs (p=0.04). Furthermore the variability in pH in non-DRIP water was high compared to DRIP areas, especially during summer (Fig. 5, center plot).

EC
Mean electrical conductivity from DRIP water was 36.2 µS/cm, which was lower (p=0.08) compared to the mean of non-DRIP water (51.6 µS/cm, Fig. 3). The variance in EC was mostly explained by POS and TIME, and the interaction between HP and 5 POS (Table 1). Overall the conductivity increased from the upland to the riparian wells (39.3 to 48.0 µS/cm) and increased as well from spring to autumn (39.7 to 48.7 µS/cm, Fig. 4 lower panels). The interactions between groundwater conditions and the position relative to the stream were mostly related to two specific contrasts. The variability in EC in non-DRIP groundwater increased from the upland to riparian wells, while in DRIP areas the EC remained stable (Fig. 4, bottom row). Moreover large differences were found between DRIP and non-DRIP in the riparian wells, where the EC in non-DRIP riparian areas was twice 10 as high as the EC in DRIPs (63.6 µS/cm compared to 32.4 µS/cm). In the upland areas, the DRIP and non-DRIP water was similar. Non-DRIP water increased from 40.5 µS/cm to 62.4 µS/cm from the upland wells towards the riparian wells, while DRIPs even decreased in conductivity (38.2 and 32.4 µS/cm for upland and riparian wells). Over the different seasons (TIME), the contrasts between DRIP and non-DRIP chemistry were consistent (Fig. 5, lower panels). The interaction between groundwater and seasonality (TIME) was not found to have an effect on EC. The only specific contrasts for both DRIP and 15 non-DRIP was a 5µS/cm decrease from autumn to spring (PDRIP=0.05, P non-DRIP =0.0007). Overall, the explained variance of our LMM was 70% for EC, compared to 22% when only accounted for fixed effects (Table 1).

Discussion
Our riparian groundwater sampling campaigns demonstrated that on average DOC concentrations in DRIPs were almost twice 5 as high compared to the less hydrologically active riparian areas (non-DRIPs). Groundwater chemistry of DRIPs was more constant from the upland to the riparian zone, and moreover concentrations remained stable across seasons. The groundwater chemistry of non-DRIPs was characterized by 40% higher electrical conductivity than DRIPs, and increasing variability towards the stream and across seasons. Differences in pH were less distinct, and mostly accountable to seasonal changes. These results confirm our hypothesis that DRIPs have a more organic DOC-rich groundwater chemistry, while non-DRIP water can 10 be associated with more mineral chemistry. However, apart from the commonly tested factors, we found that site specific properties remain to play a major role in explaining spatiotemporal variability in groundwater chemistry.
We found that during spring flood conditions, the high DOC concentrations in DRIP groundwater decreased 20%, and became less spatially variable. We believe that snowmelt dilution is a likely cause for the decreased DOC concentrations during spring. 15 Furthermore, ice sheet formation in the DRIP areas has been reported previously, which can route water over the ice surface instead of the organic rich subsurface flow paths, such as the DSL (Ploum et al., 2018). These overland-flow findings are similar to dilution effects and soil frost effects reported for wetland dominated streams during spring floods (Laudon et al., 2004a(Laudon et al., , 2011. In contrast to DRIPs, riparian groundwater in non-DRIP areas increased in DOC and in variability during the snowmelt season (from 17.7 to 21.6 mg/l). This is likely associated with the increase in groundwater level (Fig. S1), and the 20 activation of the dominant source layer in the upper section of the soil (Ledesma et al., 2015). The increased variability could be related to different timing of rising groundwater levels, for example due to local conditions that affect snow melt rates on hillslopes such as shading or sun exposure. As such our sampling campaigns provided a snapshot of the elapse of the snowmelt flood.
With comparison of riparian groundwater chemistry through the DRIP/non-DRIP concept, we have studied two extreme riparian hydrological connectivity types: DRIPs had hydrological connection with large upslope contributing areas (on average 2.7 ha), and mostly saturated soil conditions (Fig. S1), while non-DRIPs were characterized by draining individual hillslopes (on average 17 m 2 ) and having lower groundwater levels in the riparian zone (Fig. S1). Earlier work in the study area has demonstrated that the extend of the riparian zone and contributing area play an important role in the available soil carbon pool 5 and the related DOC export from riparian zones to streams (Ledesma et al., 2015). However the latter covers riparian zones with contributing areas that range between 2.5 and 1500 m 2 . Between such riparian hillslope contributing areas and initiation of streams (e.g. 10-20 ha), there is a wide range of features that focus water towards the perennial network. Where ephemeral streams are often clear extension of the stream channel, which activate mostly during hydrological events (Ågren et al., 2015), DRIPs have no such stream-like features and should be more associated as a part of the terrestrial landscape than the stream 10 network. Such features have been represented in different landscapes across the world and highlight specific processes such as: groundwater discharge zone, groundwater hotspots, cryptic wetlands, swales, focused seepage, discrete seepage, springs, upwelling zones, preferential discharge, and zero-order basins (Creed et al., 2003;Hayashi and Rosenberry, 2002;Tsuboyama et al., 2000). With the term DRIPs we aimed to fill the gap between riparian hillslopes and (fractal) stream networks as riparian landscape features that have hydrological connection to a large upland contributing areas, but lack stream channel formation. 15 Given the stable DOC concentrations and the large role in stream flow generation (Leach et al., 2017;Ploum et al., 2018), the DRIP concept could potentially be used to scale riparian contributions to headwaters on catchment level . A preliminary analysis showed that 57% of the Krycklan catchment is draining into the stream network through DRIPs, spatially covering only 12% of the riparian zone (supplementary material). However, the topography driven 20 approach behind our DRIP concept might miss certain contributions that are not necessarily related to surface topography, especially in areas where phreatic aquifers are not underlain by till deposits, or on scales that surpass the headwater basins (Devito et al., 2005). Previous work has demonstrated that in boreal catchments the input of deeper/older groundwater (with high EC) increases with drainage area, up to a threshold where old and new groundwater input reach a balance (Peralta-Tapia et al., 2015). Future work could be directed towards further chemical analysis of DRIPs and non-DRIPs and their role in 25 groundwater-surface water interactions throughout the catchment. Although the visible effect of DRIPs on streams likely decreases in higher order streams, links between hydrological pathways on groundwater chemistry dynamics have been found to significantly affect the chemistry of a fifth order river (Carlyle and Hill, 2001). Further, flow paths known as watertracks have been shown as important biogeochemical controls on higher order Arctic rivers (Harms and Ludwig, 2016;McNamara et al., 1999). 30 Spatial characterization of groundwater chemistry has been studied as an integrated signal of the phreatic aquifer (Kiewiet et al., 2019), but also using piezometers or lysimeters at specific depths, to depict vertical chemical profiles (Grabs et al., 2012;Lidman et al., 2017). Our approach was considered to represent a mixture of riparian groundwater that is likely to be incorporated in the stream during various hydrological conditions. Where the aforementioned studies relate vertical water 35 chemistry profiles to water level fluctuations to obtain process-based understanding, our study focused on finding patterns in generalizable factors such as spatial distributions (upland to riparian), different seasons (spring, summer and autumn) and hydrological connectivity. In that way, our study can be contextualized as an approach that potentially allows characterization of control points in the landscape with use minimal information. The relative contributions and biogeochemical characteristics of DRIPs and non-DRIP riparian zones in the longitudinal dimension, can potentially be combined with models that specify 40 vertical profiles of groundwater chemistry, such as the RIM model (Seibert et al., 2009). As such we can identify within the riparian zone which parts exert a large control on stream water quality and quantity.
Along the stream networks, the delineation of DRIP/non-DRIP areas in the riparian zone can help to implement hydrologically adapted buffers in forest management, ensuring that waterbodies maintain a good water quality (Kuglerová et al., 2014b;Tiwari et al., 2016;Wallace et al., 2018). Traditionally, forest practices considered fixed width buffers even though the riparian function is not homogeneous around all water bodies (Buttle, 2002). Besides the extend of the riparian soils (Ledesma et al., 5 2018), species richness within the RZ (Kuglerová et al., 2014a), and the extend of (ephemeral) stream networks (Ågren et al., 2015), our results support that variable widths should be considered in riparian buffer management. We found that DRIP water had already a distinct chemical signature before entering the RZ: 80% of the DOC originated from upland riparian wells. This suggests that the chemical role that is associated with RZ's, extends further away from the stream than the traditional fixedwidth buffer management considers. 10 For identification of control points, improving hydrological models, and sustainable forest management practices, a binary approach with little need of local properties can be a very useful tool. However, to understand the underlying mechanisms and the link to the landscape, hydrology of RZ's should be considered non-binary (Klaus and Jackson, 2018). Our LMM's showed that a large part of the variance is explained by the random effects, which contain information regarding the unique properties 15 of individual transects and to a lesser extent the subcatchments. The large variation in non-DRIPs lead to statistically weak contrasts, but this does not mean non-DRIP RZ's are less important. It demonstrated that an important chemical change also occurs in riparian non-DRIP areas, but their complexity overpassed the binary simplifications we have made in this study design. To explain the variance that was accounted for using random factors, it could be of interest to further analyze local landscape characteristics, subsurface soil properties and groundwater level dynamics to decipher whether soil, biology or 20 hydrology define biochemical characteristics throughout the RZ.

Conclusions
Are DOC concentrations in riparian groundwater linked to hydrological pathways in the boreal forest? Yes, based on our findings there is a strong link between the hydrological pathways in the riparian zone, and the DOC concentrations of riparian groundwater. At the confluence of hydrological pathways in the riparian zone, Discrete Riparian Inflow Points (DRIPs), we 25 found groundwater with an organic-rich, stable chemistry, compared to the remaining, drier riparian areas. Combining the organic-rich chemical characteristic and dominant hydrological contributions to headwaters, we propose that DRIPs are control points in the boreal riparian forest for the transport of carbon to small streams. To our knowledge, this study is the first to characterize spatial groundwater chemistry that a priori incorporated hydrological pathways in the study design. 30 further investigate the hydrological activation, and a broader chemical characterization. To understand the mechanisms and processes that link hydrological pathways and groundwater chemistry in boreal forest, we suggest to move towards non-binary approaches incorporating groundwater fluctuations, soil properties and landscape characteristics.

Data availability 5
Data is available upon reasonable request through the first author (stefan.ploum@slu.se). Krycklan data is openly available through the Svartberget database: https://franklin.vfp.slu.se/

Author contributions
LK and HL conceptualized the study design and methodology, supervised data analysis, interpretation and writing process.
AP conducted the well installation, provided field support, reviewed the written text and was involved in discussions. SP was 10 responsible for collection of data and data analysis, figures, interpretation and writing and revision of the manuscript.

Competing interests
HL is a member of the editorial board of this special issue of Hydrology and Earth System Sciences. Table 1 Summary of statistics from LMM models for DOC, pH and EC. The three columns show the response variables DOC, pH, and EC. The upper two rows show the marginal and conditional coefficient of determination (R 2 mar and R 2 con), which explain variance by the fixed effects, and the variance by the fixed and random effects together. For each explanatory variable and the interaction with HP, the p-value and F-statistic is presented. HP differentiates between DRIP and non-DRIPs. POS represents the 5 three positions in along transects being: riparian, transition and upland. TIME represents the three different seasons when sampling has taken place: spring, summer and autumn. Significant codes: p< 0.001 '***', 0.001<p<0.01 '**', 0.01<p<0.05 '*', 0.05<p<0.1 '.', p>0.1 '-'. Explanatory variables with a 'variable1:variable2' represent the interaction between both variables. For the preliminary analysis of DRIP coverage across the Krycklan catchment the following approach was followed: 5

DOC
The stream network was defined by a 10 ha flow initation threshold using a 2 meter DEM. Then a DRIP network was defined using a 2 ha initation threshold. Each point where the 2 ha stream network was incorporated in the 10 ha network, was considered as a DRIP site. The area of the catchment was 62 km. The contributing area of the DRIPs was 35.34 km, which is 57 % of the catchment area. The total length of the stream network was 162.5 km. We considered the total length of both sides of the stream as the riparian zone, which was 325 km. The total length of stream banks where DRIPs flow into the stream 10 network was 20.75 km, when assuming a width of 25 meters for each DRIP (n=830). The total area of DRIPs was 12.8% of the total length of stream banks of the 10 ha stream network.