The paper aims to explain evolution of concentration-discharge patterns in a stream network using a modelling approach. This is an interesting topic, following observations from many systems, where the c-q patterns become homogenised downstream, i.e. from highly variable and positive c-q slopes in first order streams to more linear responses, near chemostatic responses in downstreams. 

The modelling approach adopted here explains one aspect of these previous observations, i.e. changes in curvature. The authors show that in 1^st order streams curvature is larger than in higher order streams which can be explained by hydrological accumulation and homogenisation when moving downstream. Simply speaking, 1^st order streams can have a larger variation in concentration sources compared to bigger streams. And/or activation/deactivation of these sources requires changes in flow discharge that can result in the ‘bent’ c-q relationship or simply speaking different slopes of the relationship for different flows.

The paper is heavy on modelling that can conceal the main findings of the paper, which is that 1) curvature is predominant in 1^st order streams and 2) curvature is better explained by flow characteristics than nitrate uptake velocity. These findings can be explained by higher rate of biological processes in headwater streams (hypothesis investigated in this paper) but there are also other factors that can explain bending of the c-q curves, like stream morphology and flow-stage relationships, activation/deactivation of sources in relation to flow including presence of sewage pollution, drains etc. Thus, I am not convinced that the paper provides the one and only explanation for the observed patterns, rather than provides a plausible explanation for one of the possible explanations. This should be clearly communicated in the paper.

The modelling focus of the paper is, however, very dense and takes precedence over the problem – variations to c-q patterns and their controls. I would suggest ‘moving’ the modelling to the background of the paper and focusing more on the problem. This refocusing would make the paper easier to understand to a non-modelling reader and set it better in the previous research on the topic.

Finally, the original concentration and flow data disappear in the paper convoluted in different models. E.g. linking curvature to other models with inherent uncertainty like Damköhler number or uptake velocity. Showing more raw data in the paper, e.g. providing traditional quantifications of the c-q slopes would be very useful for the reader to link their knowledge of the subject with the new findings of this paper. Also when showing variation in curvature, I would like to know how frequent are concave vs. convex shapes.

Specific comments:

The title could be improved. I am not a big fan of bent c-q, maybe come up with a better term? For example curved c-q relationship as opposed to linear?

Line 13 not clear what you mean by more positive slopes. Be exact.

Lines 13-15 – what about point source pollution impact on low flow concentration?

General: Large number of studies show that high-freq and low-freq c-q relationships are governed by different factors. Please be clear in your paper based on which type of data you derive/base your assumptions on.

This study addresses concentration vs. discharge relationships in streams and rivers. The authors hypothesize that instream uptake will result in “bent” logC vs. logQ relationships, because net instream removal will cause lower concentrations than expected given loadings at low flows.  They apply a new metric (curvature) to quantify this effect using both field data and modeling results from 13 different river networks that range in size and characteristics.  They also do an extensive analysis across parameter space to understand which network characteristics have most influence on curvature and network scale removal.  They find that in stream uptake can indeed lead to more bent logC-logQ relationships (because the curvature parameter is more negative), and that channel hydraulics (the width and depth vs. Q relationships) have the strongest influence.  Uptake velocity (the biological parameter) seems to have less influence, which was surprising.  They suggest that the curvature parameter could be used to quantify network scale removal using only the C vs. Q information, adding a potentially useful tool to understand network scale dynamics.

I think this is overall an interesting analysis and potentially a very useful approach for quantifying network scale uptake.  There are a few things to consider further, emphasize or discuss, and a couple of things that would increase the understandability.

The result depend strongly on the assumption that the relationship between C and Q for loading from the landscape remains linear across seasons (i.e. the parameter b is constant). One of the difficulties getting at broad scale aquatic function is isolating landscape inputs and aquatic processes (inherent in any river network scale analysis).  The constant “b” assumption is what allows inference that the bent C vs. Q relationship results from network-scale nutrient retention. Given that this analysis uses C and Q measured across seasons (as opposed to individual storm events), with seasonality correlated with flow conditions, how likely is that?  That is, the loading C vs. Q relationship will differ between summer and winter, with the former tending to have lower C (e.g. due to higher riparian uptake). Would that also result in bent curves?  I think this is an important consideration, worthy of some discussion.

I agree with one of the findings that the hydraulic dimensions are among the dominant factors when considering network scale removal.  However, in this analysis (if I understand right), a single hydraulic equation is used for width (and depth), i.e. w = Kw * Q ^ a_w, and a single aw is applied over both space and time. However, the hydraulics of rivers are such that the change in width with changing flow at any given site (due to storms) differs from the change in mean flow in the downstream direction (at-a-site vs. downstream hydraulic relationships).  Typically the at-a-site change in w is much lower (~0.1) than in the downstream direction (~0.5) (See Knighton 1998. Fluvial forms and processes: a new perspective).  It appears the best calibrated fit for one of the watersheds (Table C1) was 0.09, closer to the typical at-a-site relationship.  This will greatly affect the pattern of removal within the network (small vs. large rivers) as well as with changing flow.  Note that the constant (Kw) is the width (m) when Q = 1m3/s.  So if you have a low aw, that means large rivers stay relatively narrower and small rivers stay relatively wide (since width doesn’t change much).  The calibrated aw is closer to the at-a-site change (where increasing flow is accommodated mostly by changes in velocity) than the downstream change (where increasing flow is accommodated mostly by changes in width). This may explain why uptake velocity is relatively unimportant (which I was surprised by), and also why water velocity comes out as so important.  It would be worth confirming whether the modeled widths match observations, and reporting the mean width of small headwater rivers (<5km2) and larger rivers (> ~400km2) to evaluate if they are reasonable.

It is interesting and a bit surprising that vf had a relatively small impact.  The authors state that if vf = 0, there is no bending (conservative) – and of course I agree.  But it seems that a low vf would then result in only slight bending, which will only increase as vf increases.   Does this pattern not occur?  The choice of vf in the paper is appropriate for denitrification, but it is on the low side total N uptake (assimilation) which could be 5-10x higher than for denitrification (e.g. Mulholland et al. 2008 found denitrification was ~15% of gross nitrate uptake).  Net assimilation may also be important in watersheds at certain times, particularly during lower flow summers (storing N over medium time scales,  or transforming to PON or DON). Might this ever be a factor in the watershed considered.  Could the Monte Carlo analysis address this possibility by using a higher Vf to determine at what point vf dominates the bending? 

Given that these C vs. Q patterns are based on samples collected over the year there is also the confounding effect of temperature on biological activity.  Denitrification is often represented with Q10 = 2, so winter (cold temperature) reactivity could be much lower.  I know this was not part of the analysis, but given the use of C collected over seasons, it seems important to factor in somehow, at least in the discussion. The temperature effect, correlated with Q, would cause a more rapid shift to saturation with increasing flow (since most of flow change is likely seasonally driven, given the sampling regime).   Should discuss whether this factor is potentially important, why or why not? 

I also had a question about how “bentness” (=curvature) is represented in Figure B1, discussed, and demonstrated.  It would help me a lot (and I assume other readers) if some of the empirical patterns of log C vs. log Q were shown.  Examples for different values of the curvature parameter (end members, the median, and 0) would be helpful.  Especially since one of the conclusions is about the utility of these low frequency empirical data sets (L641) and given that much of the recent literature has used high frequency data to get at C vs. Q relationships. 

Also, I would consider some of the wording regarding “less curvature”.  I initially assumed that meant straighter.  But in fact, “less curvature” meant a more negative curvature parameter, which is actually more bent.  It took me a while to get straight. I

In conceptual figure B1, I think that the bentness as I understand it should show a straight line at high flow parallel to the curvature equal 0 line, but bending down as flows decline.  If the dynamic is saturation, it should approach the slope set by the loading function.  Would it make sense to modify Figure B1 to reflect that (if indeed correct)?  I also think some empirical patterns, showing what the curvature parameters is, would also help increase the intuitiveness of the results.   A demonstration of how curvature is fit would be good in the appendix (to make section 2.1 easier to understand). 

I appreciated the test of the model predictions against observations in the Selke watershed.  The correspondence looks excellent!  But I did not quite understand how the seasonality of concentration emerges give the low removal proportions (I assume this is network scale removal by the entire network), and the fact the loading C vs. Q relationship is flat (b = 0.014).  It seems that loading is fairly constant and removal in Figure 2a is very small (<5% at all times).  So what causes the large drop during summer?  I would add another line that represents the export assuming conservative mixing (Vf = 0).  Also, in Figure 3, add the observed C vs. Q relationship.

What is driving the runoff (water transfer from land to water) variability over time in each watershed?

While the conclusions provide clear and useful summaries, I found the final conclusion seemed underwhelming.  I think more of the implications of these findings could be emphasized, and why they would be useful. Tie back to the big picture of C vs. Q, role of network removal, and management.

Specific Suggestions

Line 116:  should read “log” C-Q

L 135.  Where does the value “402” come from?

L137.  Meaning that at least 10% of the observations come from every season?  Still, less sampled seasons could be underrepresented.  What seasons were most samples collected?

L182.  What does this parameter definition mean?

L196.  Why does the equation have “b+1” rather than just b?

L238.  Explain what PAWN stands for when first introduced.

Table 1.  Kw is not unitless, it has units of the dimension. (it is equivalent to the width at 1 m3/s or whatever units of Q you use).  Same with Kd.

Table 2.  Please add the watershed scale runoff (mm/d) to this table.  It will allow comparison of how the different watersheds function.  Q at the outlet is then just that times the watershed area.  Is median Q the median of all river reaches, or the median at the mouth over time?

Table 3.  Why such small ranges for some of these parameter but not others?

Figure B4.  Define the variables

Table C1.  The parameters for the Selke catchment suggests that inputs of NO3 are relatively chemostatic (fairly low “b”).  This would lead to C vs. Q flattening out at high flows.   It may be helpful to include a “conservative tracer” scenario to each of the catchments, which will be based on the C vs. Q of loading from the landscape.  The divergence (always lower), will indicate bentness.  Consider representing Figure B1 in this way.

L295.  Explain what KSmax means in words and whether high values are better or worse.

L312.  I am not sure that the catchment wide Da adds much to the overall analysis, and could be dropped.

L356-358 and Figure 3.  The comparison of % removed and absolute amount removed within each grid cell is interesting and useful, but not the complete story.  There are many more medium and large river grid cells than headwater grid cells along any nutrient loads flow path.  So cumulative removal by larger rivers likely approaches or maybe even surpassed that of cumulative removal by the headwaters, particularly at high flows (see Wollheim et al. 2006 and 2018).  Consider adding that metric as well. 

L384.  Wouldn’t median over represent low flow periods, rather than total fluxes (since most flows are low, storm flows relatively infrequent).

L416.  It is not clear in the table of watershed characteristics why C1 and C10 have so much higher Lr.perc than the others.    What causes the large variability among watersheds?  Cumulative percent removal should always increase with watershed size.  Are you reporting the median within a watershed?  I think cumulative removal would be a better metric. 

L458.  I have a hard time understanding why catchments results are distinct, when all the parameters are the same. L461 says local loading and uptake differed, but what basis, since all the parameters are the same!  Some of the other explanations in this paragraph are similarly unclear.  It seems the model predictions can be summarized to see ifthe statements are true.

L470.  Is Q higher is some catchments because they are stormier (runoff vs. Q focus). Q integrates watershed size and storminess.

L473.  Is the runoff the same in the small catchments as the large?

L476.  Important point!  What about flow regime (frequency of different runoff events over time).  Are they similar among catchments?

L511.  Replace “Curvature” with “Curvature Parameter” because less curvature is more bent. 

L641.  I think to make this conclusion, you need to include more empirical relationships. 

Figure 3.  Add the observed C vs. Q (fitted relationship, with their R2) as a model test to this figure.  Important to know how close predictions come to observation

Figure 5.  Nice summary of all the correlations, with color coding.

Figure 7.  I fou9nd this figure to be impossible to interpret.  I think more explanation in caption needed.  What are the histograms?  What are the decision values?  Why do variables show up multiple times? Not sure how useful the Cart analysis is based on the discussion here.