Technical Note: Linking soil -- and stream-water chemistry

Technical Note: Linking soil – and stream-water chemistry based on a riparian flow-concentration integration model J. Seibert, T. Grabs, S. Köhler, H. Laudon, M. Winterdahl, and K. Bishop Dept. of Geography, University of Zurich, Winterthurerstr. 190, 8057 Zurich, Switzerland Dept. of Physical Geography and Quaternary Geology, Stockholm University, 10691 Stockholm, Sweden Dept. of Aquatic Sciences and Assessment, Swedish University of Agricultural Sciences (SLU), Box 7050, 75007 Uppsala, Sweden Dept. of Forest Ecology and Management, Swedish University of Agricultural Sciences (SLU), Umeå, Sweden


Introduction
In the effort to understand how stream water is influenced by catchment inputs, the riparian zone (RZ) has been identified as a key part of the catchment, especially when considering the short term dynamics of water chemistry which are of great ecological importance (Fiebig et al., 1990;Gregory et al., 1991;Hill, 2000;Buffam et al., 2001;Smart et al., 2001).As a simple starting point, the RZ can be thought of as a vertical array of soil solute sources that behave like chemostats.The temporal variation of flow pathways through the riparian zone associated with changing groundwater levels connects different combinations of these solute sources in the RZ soil profile to the Figures

Back Close
Full stream.The dynamics of stream water chemistry thus reflect the chemical "fingerprint" of the combination of chemostat-like soil solute sources in the RZ that connect to the stream at a given flow rate.The RZ is of special importance for total organic carbon and related constituents, since TOC concentrations increase markedly when passing through the RZ, relative to upslope concentrations (Bishop et al., 1990;Hinton et al., 1998;Cory et al., 2007;Laudon et al., 2007;K öhler et al., 2009;Sanderman et al., 2009).The RZ is a dominant control in first order catchments, which in turn are crucial for our understanding of water chemistry dynamics in larger aquatic systems (Bishop et al., 2008).Hillslopes, i.e., the upslope area draining through the riparian zone, are the most extensive landscape units in nearly all catchments and provide both the majority of the water and most solutes in terms of fluxes.In the last meters before reaching the stream this matter passes through the RZ.It is here where important features of temporal variation in stream water chemistry are determined (Hooper et al., 1998).The riparian zone, with its often distinctive soil properties, therefore plays a crucial role in catchment hydrology and biogeochemistry.While hillslope flows to the RZ over the course of millennia have shaped the chemical and hydrological structure of the RZ, the RZ controls the short-term stream-water-chemistry on an episodic and inter-annual time scale.
One example of where the RZ is called on to explain the complexity of soil-surface water interactions is the so called double paradox of runoff generation (Kirchner, 2003).
During hydrological events there is often a rapid mobilization of old, or pre-event water as indicated by the isotopic composition of stream water.This "old water" often exhibits a highly temporally variable chemistry for some constituents, while variations are damped for others.In other words, there is a strong variability in the hydrochemical response to hydrological episodes associated with flowrates changing by an order of magnitude while the dominant source of streamflow remains "old water".
We previously proposed a perceptual model that builds on the interaction of the vertical profiles of lateral water fluxes and soil solution chemistry in riparian soils to explain the double paradox (Bishop et al., 2004).In this paper we describe and test this Riparian Profile Flow-Concentration Integration Model (RIM) in more detail.While the juxtaposition of lateral water fluxes and soil water chemistry has been used to explain different aspects of the way in which the RZ imparts a characteristic fingerprint on the temporal patterns in stream water chemistry, we focus here mainly on TOC.
For stream water TOC this RZ fingerprint has been quantified by K öhler et al. (2009).In their study, long-term riparian patterns of soil solution hydrochemistry within 4 m of the stream captured a large part of the observed variation of stream TOC especially at low to medium flow conditions.As a consequence these stream TOC patterns reflect the small fraction of riparian soils in the catchment much more than the extensive upland soils that cover a major part of the catchment.Similar results have been obtained by others.For example, Dosskey and Bertsch (1994) quantified the amount of carbon originating from upland soils covering 94% of the catchment area to contribute less than 10% while the riparian zone covering only 6% of the area contributed the remaining flux.
There are several ways to describe how streamflow and stream water chemistry are generated.These formulations often involve a number of discrete reservoirs, or "boxes".Our concept of the RZ function builds on continuity in the vertical patterns of riparian soil solution chemistry and the distribution of lateral fluxes across that soil profile.While upslope conditions are of course important for the characteristics of the RZ, RIM does not aim at describing these long-term controls.Rather we focus on the linkage between the characteristics of the RZ and what is seen in the streamflow at the time-scale of hours to years.Describing this parsimoniously is essential for making robust models.
To test and develop this concept further, a mathematical formulation is needed.We therefore present a quantitative formulation of this riparian flow-concentration integration for the RIM model as well as an explicit discussion of the assumptions in the model.We also present a way to test the model by backward calculation of soil solute concentration profiles needed to produce observed stream water concentration.These calculations are then compared to observed soil solution concentration profiles from the riparian zone.This allows us to demonstrate quantitatively that the juxtaposition of water fluxes and soil water chemistry is a plausible mechanism for observed stream water concentration dynamics of some important solutes.

Site description
The model development has been carried out in the 13 ha V ästrab äcken catchment.This catchment is a tributary of the 50 ha Svartberget-Ny änget catchment (monitored since 1980) which is located within the Krycklan basin in Northern Sweden (Fig. 1).The catchment has a mean annual air temperature of ∼1 • C, and mean annual precipitation of 600 mm with an average runoff of 325 mm (Ottosson L öfvenius et al., 2003;K öhler et al., 2008).Vegetation is mainly Norway spruce (Picea abies) in wet areas and Scots pine (Pinus sylvestis) in drier soils away from the stream channel.Podzols have developed on glacial till that extends to gneissic bedrock 5-10 m below the ground surface.Peat deposits are found in riparian areas close to the stream, and are up to 1 m deep.Discharge was computed on an hourly basis from water level measurements (using a pressure transducer connected to a Campbell Scientific data logger) at a thin plate, 90 • V-notch weir at the outlet of the catchment.The rating curve was derived using manual, instantaneous discharge measurements (50 measurements in the range from 0.1 to 10 mm day −1 ).The stream sampling program is based on weekly samples and more frequent sampling during some episodes, in particular, spring flood.For the soil measurements three soil profiles are located along a 25 m transect (Fig. 1b).This so called S-transect was established in 1992 and has been sampled approximately monthly since 1996 when suction lysimeters were installed (Laudon et al., 2004;Cory et al., 2007;Petrone et al., 2007) topography, to follow the assumed lateral flow paths of the groundwater towards the stream.The riparian soil profile closest to the stream, S04, is dominated by organic material with a transition from organic to organic-rich mineral soil at 30 cm depth.The organic enrichment continues to a depth of 60 cm.The upslope site, S22, is located in a typical podzolic soil with a 10 to 15 cm organic layer overlying the mineral soil.
Site S12 is between the riparian and the upslope site in an organic rich mineral soil with a transition from organic to organic-rich mineral soil at 20 cm depth.The organic enrichment continues to a depth of 50 cm.Details about flow pathways, soil properties, water content and temperature, as well as the role of soil frost in the hydrology of the transect can be found in previous publications (Nyberg et al., 2001;St ähli et al., 2001).

Riparian profile flow-concentration integration model (RIM)
In many soils physical and chemical properties, including lateral transmissivity and soil water concentration, vary systematically with depth.With the riparian profile flowconcentration integration model (RIM) these variations are used to estimate the stream water concentration (or load).The idea behind the model is that the variations of water and solute fluxes with depth for a representative soil profile can be used to explain the variations of stream chemistry with varying runoff (Fig. 2).The flow-concentration integral is basically the juxtaposition of lateral water flux, q[L 2 T −1 ], and soil-water chemistry.The lateral mass flux of a constituent at a certain depth z [L] is the water flux at this depth, q(z), times the concentration of the soil water solution at depth z, c soilwater (z) [M L −3 ].The total mass flow rate of this constituent, i.e., the constituent load, L[M T −1 ], in the runoff is then the integration of these lateral mass fluxes over depth from a certain base level (bedrock), z 0 , to the groundwater table, z 1 (Eq.1).Lateral fluxes above the groundwater table and below z  and the upslope recharge along the stream network, these values for streamflow and load can be assumed to represent the values at the catchment outlet.
The variation of water flux and concentration with depth can be described by different functions.In this study we used exponential functions based on the observed variations of flow as a function of groundwater depth (parameters a and b) and soil solution (parameters c 0 and f ) as a function of soil depth z.Stream constituent loads can then be calculated using Eq. ( 2).
Several parameters in this equation can be estimated independently from each other.The higher the flow rate, the more the soil profile contributes to the runoff and, thus is identifiable from the stream water chemistry based on RIM.Introduction

Conclusions References
Tables Figures

Back Close
Full Screen / Esc Printer-friendly Version Interactive Discussion

Assumptions of the riparian flow-concentration integration
For the lateral groundwater-flow computations it is assumed that the direction of soil water fluxes is parallel to the ground surface and lateral flow only occurs in saturated soils below the water table, as well as that a unique relationship exists between groundwater level and streamflow.This means that all the additional lateral flow associated with ground water table rise occurs in the newly saturated soil layer.Matrix flow or Darcian flow with a constant hydraulic gradient can be one way to fulfil this assumption.The large variations in lateral water fluxes with changes in the groundwater level would then be explained by the decrease of hydraulic conductivities with depth.This transmissivity feedback mechanism has been suggested as a runoff generation mechanism in till soils, where there is extensive evidence of saturated hydraulic conductivity increasing upwards in the soil profile (Lundin, 1982;Rodhe, 1989;Bishop et al., 1990;Nyberg, 1995;Laudon et al., 2004).However, flow does not have to be Darcian to fulfil the above assumption.As long as there is no bypass flow and any additional water flow always occurs in the soil layer just above the previous groundwater table, the computation of the lateral fluxes is valid.A vertical gradient in saturated flow is a useful visualization appropriate to till soils with a strong vertical profile in saturated hydraulic conductivity, but not a requirement in future applications of the model.
For the concentration computations it is assumed that the soil water chemical signal is imprinted in the riparian zone instantaneously, or at least at a rate that is faster than the rate at which the water traverses the riparian zone, before entering the stream.The solute concentration of the soil water entering the riparian zone from upslope does not have any influence on the concentration of solutes of the water flowing out of the riparian zone.In other word, the RZ functions as chemostat.As a consequence, the predicted load will be a function of flow only, i.e., the resulting concentration will be the same as long as the flow is the same, and the flow will be the same as long as the groundwater level is the same.This assumption can easily be relaxed by letting the concentration soil profile vary temporally.It is also assumed that the solutes are Introduction

Conclusions References
Tables Figures

Back Close
Full transported by advection only to the stream and that diffusion is negligible.

Analytical solution of the riparian integration
The integral expression in Eq. ( 1) can be solved analytically for certain functions used to describe the variations with depth.Below we show this for the exponential functions (Eq.2).Note that the existence of an analytical solution is convenient but not necessary since Eq. ( 1) can also be solved numerically.
The key for the analytical solution is to integrate over streamflow, which is possible since streamflow is directly related to the riparian groundwater table position.Substituting the profile depth z by streamflow (Eqs.3 and 4) and rewriting Eq. ( 2) results in a simple power law where the constituent load is directly related to streamflow Q (Eqs.5 and 9).
The constituent load is then related directly to runoff and soil parameters through a simple power function (Eq.9): It is interesting to note that this type of power-law function is well-established as an empirical way to describe the relationship between discharge and load (or concentration); sometimes also referred to as the rating-curve method.Originally developed for estimating stream sediment loads (Campbell and Bauder, 1940;Cohn, 1995) it is also frequently applied to calculate nutrient loads including phosphorus and DOC (Crawford, 1996;Smith et al., 1996;Cooper and Watts, 2002).Here we demonstrate how this relationship can be derived based on an integration of the lateral fluxes and constituent concentrations in the soil water solution at different depths through the riparian zone.One difference between the numerical integration and the analytical solution are the integration boundaries.While we integrate from the soil surface to a certain depth for the numerical solution, the integration continues to infinite depth in the analytical solution.A lower limit could of course be introduced in the analytical solution but in order to simplify the equations we choose to integrate to minus infinity.In testing of the RIM described below we used a lower limit of 1 m, and with reasonable parameter settings there was virtually no difference between the numerical and the analytical solution.In other words, the difference in the lower boundary had no effect because of the assumed exponential decrease of the lateral water fluxes and concentrations with depth.Introduction

Conclusions References
Tables Figures

Back Close
Full Screen / Esc Printer-friendly Version Interactive Discussion

Testing the riparian profile integration model
We tested the riparian flow-concentration integration approach using observed data from the V ästrab äcken catchment, both streamflow and chemistry at the outlet of the catchment, and the instrumented transect within the catchment.Here we assumed that this transect is a fair representation of the whole catchment.This assumption can be motivated by vegetation and soil morphology of the transect which was similar to other sites within the catchment.Depths of the boundaries between mineral soil, organicrich mineral soil and peat soil were within a few centimetres of the median depth for these interfaces measured at 30 transects on the east side of the stream (Bishop et al., 2004).In addition, a tracer investigation of water sources and flow pathways using stable isotope measurements (Laudon et al., 2004) demonstrated that the proportion of event and pre-event water leaving the transect during spring flood was similar to that of the entire catchment.
Using groundwater level observations from the riparian site at the S-transect we established an exponential relationship between catchment runoff and groundwater levels (Fig. 3).This approach is supported by previous studies showing that riparian groundwater levels are in phase with runoff in the conditions similar to those found in the study catchment (Seibert et al., 2003).The differentiation of this groundwater level -streamflow relationship with depth yielded the parameters (a and b) for estimating the water flux at a certain depth under the assumption of an exponential transmissivity profile.
For 132 instances with stream concentration observations collected between 1993 and 2001 (K öhler et al., 2008), we calculated the corresponding soil water concentration profiles for TOC according to Eq. ( 2).All measurement periods had a streamflow corresponding to a specific runoff of more than 0.5 mm/d.The soil solution concentration at 1 m depth was set to 5 mg/l based on the stream concentration observed at low flow periods of around 0.15 mm/d, which generates less than 2% of the total annual runoff in that catchment.The only parameter varying between the different measurement periods was the shape factor, f .After fitting f , we compared the backwards calculated Introduction

Conclusions References
Tables Figures

Back Close
Full Screen / Esc

Printer-friendly Version
Interactive Discussion soil concentration profiles with profiles that had been observed on 14 occasions.For comparison we used the profiles computed for a window of up to 10 days from the day at which the soil solute concentration measurements were taken.This meant that we compared each of the observed 14 soil concentration profiles with those estimated based on one to six different stream water concentration measurements (35 in total).
We also tested the riparian integration model for three other constituents, namely Ca, Mg and Cl.Additionally, here an exponential function was adjusted to the concentration profile and the following values were used for the concentration at 1 m depth: 3 mg/l (Ca), 1 mg/l (Mg) and 1 mg/l (Cl), again based on low flow concentrations.Four instances with measured stream and soil water concentrations were available for this test.
For TOC we investigated how the backwards calculated soil profiles changed seasonally.For this we compiled distributions of the shape factors, f , for the different months based on the 132 computations of the TOC soil profiles based on different observed streamflow concentrations.

Results
The backward calculated soil profiles of TOC for the stream-and soil-water concentration observations agreed well with the observed soil concentration profiles for TOC on all 14 instances with soil concentration observations (Fig. 4).To evaluate the estimated soil concentration profiles we computed the RMSE based on the differences between observations and estimations for the different depths for each instance.The median RMSE for the different instances was 3.7 mg/l, with 80% of the RMSE values ranging from 1.4 to 11.6 mg/l.For comparison one could estimate the soil profile by one constant value corresponding to the observed stream water concentration.The median RMSE of this benchmark was 7.8 mg/l and the RIM soil profile was better than the benchmark on 80% of the measurement periods.
For the other constituents the backwards calculated concentration profiles also Introduction

Conclusions References
Tables Figures

Back Close
Full generally agreed with the observed profiles for the four instances where observations were available (Fig. 5).For Ca, the median RMSE was 0.47 mg/l, for Mg 0.09 mg/l and for Cl 0.04 mg/l.These values corresponded to about 20% of the mean value for Ca, 15% for Mg, 5% for Cl and were only slightly higher or close to the analytical precision that may be expected for those solutes.
For locations further away from the stream the backward calculated TOC concentration profiles did not agree with the observations (not shown here).For S22 (Fig. 1), for instance, the median RMSE for TOC was 17.2 mg/l.
There was a considerable variation in the estimated shape factor values resulting in a large variation of the computed soil profiles (Fig. 6a).Grouping the f -values according to the month of the observation revealed some seasonal pattern (Fig. 6b).In general, higher values for f were obtained for flow situations during summer and fall compared to spring.

Discussion
The RIM approach allowed us to infer soil concentration profiles from observed stream concentrations.This is similar to the approach of inferring rainfall and evaporation based on streamflow fluctuations by "doing hydrology backward" (Kirchner, 2009).The backward calculated profiles of several solutes agreed, in general, well with the observed soil solute concentration profile in the RZ.This failure to falsify the RIM model is not sufficient for the model and its assumptions to be valid, since there might be other ways to explain the observed relationships.Our results do, however, demonstrate that the suggested qualitative explanation (Bishop et al., 2004) is quantitatively plausible and can be predicted in a model that captures major features of the riparian zone architecture regarding both hydrology and soil solution chemistry using just four parameters.The backward calculation of soil concentration profiles is especially valuable because stream water chemistry is easier to measure than soil solute profiles.In this study, we observed a seasonal pattern of the shape factor values that indicate increasing TOC concentrations in the soil through spring and summer.End-member mixing analysis (EMMA) (Hooper et al., 1990) is a commonly used technique for interpreting stream water chemistry variations and the importance of different distinct sources of streamflow (Hagedorn et al., 2000;Katsuyama and Ohte, 2002;Inamdar andMitchell, 2006, 2007).The RIM approach differs from EMMA by allowing for a continuous transition between waters of different chemical composition such as TOC-rich soil water in the upper soil layers to TOC-poor water in deeper layers.
We suggest that the RIM approach can provide a framework for investigating the temporal variations of soil solute concentrations, and allows this to be done with a much higher temporal resolution than would be possible with soil solution measurements alone.As mentioned before, there is a long tradition of using power-laws in hydrology to describe the relationship between streamflow and load (or concentration) (Campbell and Bauder, 1940;Cohn, 1995;Crawford, 1996;Smith et al., 1996;Cooper and Watts, 2002;Godsey et al., 2009).The analytical solution of flow-concentration integration presented here provides a physically plausible explanation for the load function as a power function for certain constituents such as TOC, Ca, Mg, and Cl.The analytical derivation is similar to the one presented by Godsey et al. (2009) but differs in that it relies solely on soil water concentrations in the riparian zone as explanatory endmembers rather than on the ensemble of reactive soil surface areas along an entire hillslope.
The power function, which we derived based on the assumed exponential profiles for flow and concentration, represents a special case of the more general RIM concept.The derivation also does not constitute a true proof of concept as other processes might lead to the similar power-law-shaped relations.This is illustrated by the fact that power-law-shaped relations observed between streamflow and stream sediment loads have often been reported (Campbell and Bauder, 1940;Cohn, 1995) but are obviously not related to riparian soil concentration profiles.In addition, flow and concentration profiles could be represented by other combinations of functions than two exponential functions, which would lead to differently shaped streamflow-load relations.However, Introduction

Conclusions References
Tables Figures

Back Close
Full applied to several solutes in this study, the power-law function offers an appealingly simple and yet physically plausible explanation based on observed stream and soil water concentrations.The RIM approach did not hold for linking stream water chemistry to soil solutions profiles of TOC further away from the stream.This demonstrates how, at least for some solutes, it is the "fingerprint" of the RZ that determines the stream water chemistry dynamics, rather than the much more extensive upslope soils.
The implicit assumptions of spatial and temporal invariance in the very basic type of model presented here may not be valid in different types of catchments or landscapes.Subagyono et al. (2005), for instance, provided for a Japanese catchment evidence that the not the riparian zone not always resets the hillslope-water chemical signature.In alpine catchment the spatial variation of snow accumulation and melt might be an important control on stream chemistry (Boyer et al., 1997).
With regards to the assumption of spatial uniformity, we assume that stream recharge conditions as well as the hydraulic and chemical properties of riparian soils in the V ästrab äcken catchments act as one homogenous hydromorphic unit, and that this unit is adequately represented by the study transect.If a catchment has a less homogeneous topography or even other types of landscape elements, such as wetlands or lakes, a disaggregation of the different landscape elements or hydromorphic units would be appropriate.In this way also the expected variation of the chemical fingerprint of the riparian soils with varying vegetation (forest) type and age structure in the RZ could be considered.The continuous changes in concentrations of solutes in the soil with depth, which we assumed, might not always be valid, especially not in strongly structured soils where more abrupt changes are observed.
With regards to changes in time, the only temporally varying feature of the model presented is the groundwater level and associated flows.With time-invariant parameters in the RIM any hysteresis in the flow-concentration relationship are, thus, neglected, although these can be important for chemistry variations during events (Evans and Davies, 1998;Chanat et al., 2002;Rice et al., 2004;Hood et al., 2006).One may Introduction

Conclusions References
Tables Figures

Back Close
Full anticipate that the hydraulic properties of the RZ may also change with time as a response to natural long-term or short term climate variability.In the current formulation of the RIM, hysteresis effects from antecedent soil wetness on the hydraulic properties of the riparian zone are considered negligible.Changes in subsurface flow patterns may occur during prolonged periods of wetness or very dry spells.Finally we neglect any overland flow component that would lead to dilution effects over time.
There may also be temporal changes in the soil solution profile in the soil.The grouping of the back-calculated soil solution profile shape parameter, f , by month demonstrated that there was a seasonal variation in the properties of the soil solution profile.This will provide guidance on how to extend the chemically static RIM approach presented here in the direction of a dynamic model with time-variant soil solution concentrations.Mechanisms that might influence the temporal variation of the concentration profile parameters f and c D in Eq. ( 2) could potentially include periods of prolonged soil wetness, elevated soil temperature, combinations of both, or seasonality of plant growth that influences cation uptake and/or organic matter quality (Christ and David, 1996), and temperature dependant soil weathering processes.Such processes and their influence on soil concentration profiles could be included into the model by using time-variant formulations for f and c D and describing these profile parameters as functions of the relevant variables such as temperature, season, antecedent flow or degree-days.Also plant growth dynamics could be accounted for in this way.For TOC simple formulations such as those suggested by Boyer et al. (1996) might be used.

Concluding remarks
We have presented a simple model to explain stream water chemistry variations based on the fingerprint of the riparian zone.Based on backward calculated soil concentration profiles the concept was demonstrated to be a possible mechanism for several dissolved constituents.We have previously used this riparian flow-concentration integration approach in a forward way to quantify soil-stream linkages for a variety of Figures
In this paper we presented this useful approach as the RIM model in a mathematically formalized way which also shows how this approach ties into power-law representations of catchment outputs.This mathematical model facilitates testing the concept that the architecture of how flow paths and soil properties consistently interact in the riparian zone is a useful framework for explaining stream water chemistry dynamics.
The simple formulation of RIM presented here should be seen as a starting point for identifying spatial and temporal variability in catchment properties controlling hydrological and chemical processes in both streams and riparian soils.
0 are neglected.Under the assumption that all water passes the RZ before reaching the stream the water-flux integral equals the streamflow Q [L 3 T −1 ].The concentration in the runoff, c runoff [M L −3 ], can then be computed as the mass flux integral divided by the water-flux integral (Eq.1).If we can assume uniform conditions for the riparian zone The parameters related to the variation of flow with depth (a and b) can be estimated by the differentiation of the functional relationship between groundwater levels and streamflow.We further assume that we can fix the concentration of the soil solution at a certain depth, d , to a value c d that can be related to stream water concentrations at low flow when only the deepest flow paths in the soil profile are transmitting lateral fluxes and thus contributing to the total streamflow and constituent loads.The value of c d allows us to calculate c 0 (c 0 = c d e f d ).This leaves us with one free parameter, f , for any instance where we have simultaneous observations of streamflow and concentration at flow rates above low flow.In other words we can back-calculate the soil concentration which would be needed to produce the observed stream water concentration.
Screen / EscPrinter-friendly Version Interactive Discussion

Fig. 2 .Figure 3 .
Fig. 2. Schematic figure of the riparian profile flow-concentration integration model (RIM) which explains stream water concentration as the juxtaposition of lateral water flux and soil-water chemistry.