A perturbed biogeochemistry model ensemble evaluated against in situ and satellite observations

Anugerahanti, Prima; Roy, Shovonlal; Haines, Keith

doi:https://doi.org/10.5194/bg-15-6685-2018

Articles | Volume 15, issue 21

https://doi.org/10.5194/bg-15-6685-2018

© Author(s) 2018. This work is distributed under
the Creative Commons Attribution 4.0 License.

https://doi.org/10.5194/bg-15-6685-2018

© Author(s) 2018. This work is distributed under
the Creative Commons Attribution 4.0 License.

Articles | Volume 15, issue 21

Research article

|

12 Nov 2018

Research article |

| 12 Nov 2018

A perturbed biogeochemistry model ensemble evaluated against in situ and satellite observations

Prima Anugerahanti, Shovonlal Roy, and Keith Haines

Download

Final revised paper (published on 12 Nov 2018)
Supplement to the final revised paper
Preprint (discussion started on 05 Apr 2018)
Supplement to the preprint

Interactive discussion

Status: closed

AC: Author comment | RC: Referee comment | SC: Short comment | EC: Editor comment

- Printer-friendly version

- Supplement

RC1: 'Review of Anugerahanti et al', Anonymous Referee #1, 19 Apr 2018
- AC1: 'Authors' Response', Prima Anugerahanti, 07 Jun 2018
RC2: 'Review', Anonymous Referee #2, 24 Apr 2018
- AC2: 'Authors' Response', Prima Anugerahanti, 07 Jun 2018
RC3: 'Review', Anonymous Referee #3, 14 May 2018
- AC3: 'Authors' Response', Prima Anugerahanti, 07 Jun 2018

Peer-review completion

AR: Author's response | RR: Referee report | ED: Editor decision

ED: Reconsider after major revisions (11 Jun 2018) by Katja Fennel

AR by Prima Anugerahanti on behalf of the Authors (28 Jun 2018) Author's response Manuscript

ED: Referee Nomination & Report Request started (29 Jun 2018) by Katja Fennel

RR by Anonymous Referee #1 (12 Jul 2018)

RR by Anonymous Referee #2 (23 Jul 2018)

Suggestions for revision or reasons for rejection

The updated manuscript shows some improvement in places and the authors have added some useful analyses. However, some of my main criticisms of the manuscript still apply.

general comments

In my first review I wrote that "The results of the default run can have knock-on effects on the ensemble: in multiple parts of the manuscript the authors note that when there is a large bias between the model (ensemble) and the observation, that the ensemble spread is too low when really other model aspects may be to blame for the bias. In other words, problems with the parametrization, the physical model, or the 1D nature of the model cannot be explained by structural uncertainty in the biogeochemical model." To make it more explicit: if the ensemble does not envelop the observations, this does not necessarily imply that the ensemble spread is too low. A bad model (physical model, unresolved BGC processes) with a large bias may be to blame. This is now touched on in the discussion section but ignored before, and phrases like "the observation is outside the ensemble range ... making the ensemble spread too narrow" appear too often, especially in the phytoplankton phenology section 3.4 (blooms are often controlled by physical processes).

With regard to my earlier comment (reviewer 2, general comment (2)): The plot in the response shows little difference between the regular and the log-transformed fit. This is surprising to me. The linear fit produces a somewhat good match to the sigmoidal function only for high phytoplankton concentrations. This can be easily explained as the linear scale emphasizes high phytoplankton values. For the log-scale fit, what I suggested was to fit the linear function across the range of phytoplankton concentrations shown in the figure using the same log-scaling. I would expect a much lower mortality value in this case, much closer to 0.05 than the 0.085 that were obtained. This may yield significantly different results. The best idea would probably be to let the data determine the fit, i.e. minimize the misfit for the values occurring in the model output (or the observations if these are available). In any case, since many of the results hinge on which functional form is used, the authors need to include more detail on how exactly the fit was obtained. This should include the range of phytoplankton concentrations considered in the fit.

At the risk of sounding pedantic: there are still many sentences in the manuscript that are not clearly formulated and that can be interpreted in different ways. Readers familiar with the topic will likely know what is meant but others will not. I have pointed out some of these instances below but I would recommend that the authors go through the manuscript again carefully.

specific comments:

p1

l3: "We assessed the impact of structural sensitivity ... by modelling the chlorophyll and nitrogen concentrations at five different oceanographic stations spanning three different regimes": This appears to imply (incorrectly) that structural sensitivity is linked to the model location.

l3: "nitrogen": Do you mean "nitrate"? This comment applies to later mentions of nitrogen as well.

l6: Agreed, though there are other applications for BGC models than climate change assessments.

p2

l24: Mention some causes of mortality here.

p3

l14 "demonstrated ... that linear density-dependent mortality produces the biggest difference in diatoms". As opposed to other phytoplankton types? This is not clear.

p4

l29: "The difference in shape of the optimised functional forms are more obvious at low nutrient concentrations." Likely due to the way the fitting is performed, see general comment (1).

p5

l3: "graze on small phytoplankton, non-diatoms, and detritus": This is confusing, as small phytoplankton and non-diatoms are the same variable.

l13: "hyperbolic": at this point it is not clear what the hyperbolic refers to, please use "hyperbolic" and "sigmoidal" when G_1 and G_2 are introduced initially, in line 6.

l28: Here it is very important to note if the area below the function is computed in log-space or not and over what range of concentrations. Figure 1 makes it look like a log transformation was used and concentrations between 0.001 and 10 mmol N m^-3 were considered in the fitting.

p6

l7: "as" -> "is"

l30: "mortalities" -> "mortality formulations"

p7

l18: There are several 1D models now, do they all use the same depth levels (at every station)?

l21: The station abbreviations should be mentioned earlier when the stations are first mentioned (section 2.5).

p8

l3: I don't quite understand the reasoning behind using the same concentration at each depth level. Are the initial conditions quickly forgotten? Mention this. A short model spinup could have been useful.

l17: "for to" something is missing here.

p10

l9: How does the ensemble range not cover the whole ensemble and what exactly is meant by "if we only allow one process function at a time to change"?

p13

l23: "phenology metrics" -> "observed phenology metrics" or "phenology metrics obtained from observations"

l27: How can a range be late; rephrase?

l28 and again l30: "For initiation...": add "bloom"

p17

l23: This long sentence mentions "active prey selection" twice.

Fig 1: Should it be "nitrate" in the axis label in (a)? Same comment applies to other figures.

Fig 3: missing units in (b) and (e)

Fig 7: "(d)" is used twice

Fig 9: The order in which panels are labeled is confusing and inconsistent. It lead to a problem in the description as well. White line is in panel (a) not (b).

Fig 10: Use same y-axis across all panels.

Fig 11: "(e)-(f)" should be "(e)-(h)"

Hide

RR by Anonymous Referee #3 (04 Aug 2018)

Suggestions for revision or reasons for rejection

As already noted in my first review, I really appreciate the very thorough and exhaustive attempt to investigate the sensitivity of a biogeochemical model that will be, at a later stage, used for projections of climate change and other applied tasks. This paper could provide important information on global biogeochemical model sensitivity with respect to structural uncertainty, even for users not applying this particular model. However, so far, and for this purpose, I find the discussion of results somehow incomplete. I think discussing the outcome of these experiments before the background of earlier studies that also deal with global model uncertainty and performance would enhance the paper's impact, and put this work into a broader context.

For example, I am not convinced that the results indicate a larger importance of structural sensitivity compared to uncertainties related to physical or parametric uncertainty. To my impression, so far there is little (with respect to the metrics applied) evidence that mean or median perform much better than the default model (see below, point B). I think this should be discussed more, e.g. before the background of the study by Kwiatkowski et al. (2014), who showed that "No model is shown to consistently outperform all other models across all metrics." when comparing six models (among them MEDUSA) against surface tracers such as DIN and Chl. Likewise, Galbraith et al. (2016) or Kriest (2017) found that models of varying complexity performed quite similar with respect to metrics such as RMSE etc. The importance of parametric uncertainty was addressed in several studies, even at a global scale. There are also many studies that deal with the uncertainty due to physics (starting with the "classic" study by Najjar et al., 2007), and show a large impact of circulation or phenomena at smaller (meso) scales.

Secondly, different equations, as applied in this study, necessarily result in different nutrient or food affinities of plankton, as briefly discussed by the authors; therefore, to some extent the structural sensitivity also includes some parametric sensitivity. Again, this might also be worthwhile discussing (see below, point A, and specific comments for p15).

I would thus recommend an overhaul of the discussion section, as already suggested by Reviewer #1, to give this paper a wider impact and "to address the big picture".

(A) Authors response to my comment 2 on fitting procedure of different structural forms

"However, we are trying to capture the whole range of nutrient and phytoplankton at all the different region, [...]"

As shown in Figures 3, 9 and 14, the full range of observed and simulated nutrient concentration is between ~0-10 mmol/m3, and the range of Chl ~0-1 (PAP, Cariaco), 0-0.3 (ALOHA and BATS) and ~0-6 (L4), i.e. more than an order of magnitude lower than used for curve fitting. I would suggest to give the range over which the functions were fitted in the paper (method section).

Further, in their response the authors state that

"Suppose we are optimising the nutrient uptake on the similar range of station BATS and ALOHA (with maximum nitrogen and phytoplankton concentration of 5 mmol N m- 3, shown on Figure 3, although at stations like Cariaco, PAP, and L4, we may see nitrogen larger than 5 mmol N m-3), the functions still deviate at low nitrogen and phytoplankton concentration."

This is very important information, and I would suggest to mention this in the paper, together with a full description on how the fitting was carried out. (Was the range between 0-100 discretized, and if so, into how many bins? Were the discrete values equally distributed over the range?)

I still think this issue (fitting different functional forms, the side effects on affinities and rates in different oceanic regimes) might merit more in depth discussion and detailed presentation.

B) I find the analysis of model performance for the different structural types (ensemble mean or median vs default) somewhat biased and incomplete: For example, Table 3 shows, that out of the 51 criteria listed in that table for r, RMSE and Bias for DIN and Chl (profile, surface, integrated), the default model is the best model w.r.t. to the bias and also r, compared to ensemble mean or median. It also outperforms mean and median at BATS, and, all three criteria combined, is as good as the ensemble mean at ALOHA. Also, for RMSE some of the values do not seem to provide a clear evidence that the ensemble mean or median is better: for example, at ALOHA the RMSE for DIN and Chl (profile or surface) is very similar among the different models.

Specific comments:

p1, line 6-7: "that describe the key biogeochemical processes" - I would suggest to rephrase this to "that describe some key biogeochemical processes", as, at a larger scale, phytoplankton light affinity, recycling and sinking organic matter may also be regarded as key processes.

p1, line 13: "Changing mortality" - which mortality: both zooplankton and phytoplankton?

p1, line 14-15: "The RMSEs between in situ observations and the ensemble mean and median are mostly reduced compared to the default model output." - See above, point B. Why not address the bias or correlation coefficient here, and in the discussion?

p2, Introduction: This is probably just a minor point, but usually global models are NPZD models, with detritus (and/or DOP) being an important component in the recycling and vertical redistristribution of nutrients.

p2, line 17ff: "Moreover, in order to investigate the effect of global climate change and anthropogenic activities in the ocean, marine biogeochemical models are now being embedded into earth system models." - Simple marine biogeochemical models were first embedded into GCMs in the 90ies, and into ESMs possibly around 2005 (e.g., UVic, Schmittner et al., 2005). Perhaps it would be better to rephrase this as "[...] activities in the ocean, MORE COMPLEX marine biogeochemical models are now being embedded into earth system models."

p3, line 14: "A few studies have investigated the effects of biogeochemical process formulations, e.g. Yool et al. (2011) has demonstrated [...]" - Perhaps better "A few studies have investigated the effects of biogeochemical process formulations. For example, Yool et al. (2011) have demonstrated [...]" ?

p4, section 2.1: I suggest to give range of nutrient concentrations over which the functions were fitted, and some details of the fitting procedure here (see A).

p5, line 5 "which are high quality food sources" - In the model? In reality? Is this of importance here, and if so, how is it embedded in the model?

p5, line 12: "These functions both become constant at a maximum grazing rate." - Perhaps better: "These functions approach a maximum grazing rate at high concentrations of prey."

p5, line 25: "as shown on Fig. 1(b)." - Fig 1(c)?

p5, line 27ff "integrated over the range of prey density" - Give range here (see p4, section 2.1)

p6, line 5: "From a previous 3-D MEDUSA run, in the oligotrophic regions show" - skip "in".

p6, line 9: "MEDUSA also contains both slow and fast detritus sinking factors." - What does this mean? Sinking speeds? Remineralisation length scales? Perhaps better "MEDUSA also parameterises slow and fast sinking detritus" ?

p6, line 14f: " We chose a lower sinking rate of 0.1 m day to prevent depletion of state variables particularly at the shallower stations." - Does this mean sinking organic matter is buried at the sea floor? If so, I would mention it here.

p7, line 18f: "Our 1D model uses these same 63 depth levels vertical resolution." - The same as the global (MEDUSA) model? Then, on p8, line 9: "The model is simulated at 37 depth levels, [...]" - I am a bit confused - How many levels did the 1D setup have - 63 or 37?

p7, line 18f: "we use the integrated nitrogen over 200 m (integrated nitrogen / depth)" - is this meant to be "dissolved inorganic nitrogen averaged over the upper 200 m"? I would suggest to refer to "averaged" when appropriate (also in some of the figure captions); Additionally, if nutrients are mean it should be "nitrate" or DIN or "dissolved inorganic nitrogen", because "nitrogen" alone can also include the organic forms.

p 9, line 25: "Chlorophyll and nitrogen profile 10 year means"? What does this mean? "Observed mean profiles of chlorophyll and DIN"?

p15, lines 15ff: " In order to maintain phenomenological similarity, these functions are calibrated using non-linear least squares, while keeping the maximum process rates fixed." - To me it seems as if the phenomenological similarity is maintained with respect to the rates integrated over the nutrient and Chlorophyll range considered; i.e. between 0-100 mmol/m3 or mg/m3. This could downweigh the affinity at very low nutrient or chlorophyll concentrations. (See above, A.) Thus, phenomenological similarity is maintained with regard to a certain criterion. As stated by the authors in their reply, fitting the curves over a narrower range results [0-5] in almost the same curves, which would strengthen the above statement of phenomenological similarity. I would suggest to add a few sentences about this here.

p15, lines 19ff: "Applying structural sensitivity in the 1-D framework has also allowed a large parameter space of concurrent variations to be explored for several different oceanographic regions, and with minimal computational cost." - How is the parameter space explored, particularly when considering the above statement about phenomenological similarity?

p15, line 24: "have" -> has

p15, line 25: "This is expected as at low concentrations, using the G2 function would graze more phytoplankton, as shown on Fig. 1(b)." - Again, here the effect of the fitting procedure for different functions shows up: perhaps adjusting the k's of the different function such that they become similar at very low concentrations might have resulted in similar results at the oligotrophic stations BATS and ALOHA? Thus, the effect of different functional forms seems to be mingled with the effects of parametric uncertainty. Of course, this cannot be avoided, but I think I think it could be discussed more.

p16, lines 19ff: "At most of the stations, the ensemble mean produced lower RMSE compared to the default run, suggesting that the structural ensemble with a wide range of predictions covering the in situ observations, is likely to produce a mean field closer to the observation, than a single-structure model." - I am still not fully convinced about this - see my above comment B.

p18, line 8: "widely used in the community" - In which community?

References:

Galbraith et al. (2016), Complex functionality with minimal computation: Promise and pitfalls of reduced-tracer ocean biogeochemistry models. Journal of Advances in Modeling Earth Systems, doi:10.1002/2015MS000463

Kriest (2017), Calibration of a simple and a complex model of global marine biogeochemistry, Biogeosciences, doi:10.5194/bg-14-4965-2017

Najjar et al. (2007), Impact of circulation on export production, dissolved organic matter, and dissolved oxygen in the ocean: Results from Phase II of the Ocean Carbon-cycle Model Intercomparison Project (OCMIP-2), Global Biogeochemical Cycles, doi:10.1029/2006GB002857, 2007

Schmittner et al. (2005), A global model of the marine ecosystem for long-term simulations: Sensitivity to ocean mixing, buoyancy forcing, particle sinking, and dissolved organic matter cycling, Global Biogeochemical Cycles, doi:10.1029/2004GB002283

Hide

ED: Reconsider after major revisions (08 Aug 2018) by Katja Fennel

AR by Prima Anugerahanti on behalf of the Authors (30 Aug 2018) Author's response Manuscript

ED: Referee Nomination & Report Request started (12 Sep 2018) by Katja Fennel

RR by Anonymous Referee #2 (28 Sep 2018)

Suggestions for revision or reasons for rejection

The updated manuscript has improved in many aspects. I have no more issues with methods and results that are presented but the way they are presented needs improvement before the manuscript can be published.

general comments:

I found the results sections were difficult to read and not well structured. For example, section 3.1 starts out by describing the data. The second paragraph covers model DIN seasonality but in the last sentence also starts the model-data comparison. In the third paragraph, it is initially not clear if the chlorophyll data is described or the ensemble. The fourth paragraph then launches the model data comparison and mentions the same result from the second paragraph again (NRR=1.25 for DIN). This back and forth is confusing to the reader. I suggest to reorganize each results section: start out by describing the ensemble characteristics (which are the focus here) and then compare them to the data.

Section 4 is very dense and difficult to read. Again there is a lack of structure and the first few paragraphs each cover a range of different aspects without much consistency. Some run-on sentences are almost impossible to follow (for example, p18, l15). The manuscript would greatly benefit from tidying up this summary section and a clear focus on the key results.

In order to improve readability, I'd suggest to rename the functions that are used throughout the manuscript, so it is easier for the reader to keep track. For example, new names for the uptake functions could be:
U_1: U_hyp
U_2: U_sig
U_3: U_exp
U_4: U_tri
The same naming scheme can easily be applied to the other functions.

specific comments:

p2

l6: "NPZ-D": "NPZD" is much more commonly used.

l31: Change "parameters" to "biological parameters".

p3

l23: Anderson et al. (2015) is cited twice in the same sentence.

p4

l10: Here the manuscript shifts briefly from present into past tense. Use on consistently.

p6

l11: "has" -> "have"

l18: "study" -> "studies"

l18: 0.1 is not lower than zero, rephrase to something like "towards the lower end of the range of literature values"

l26: "the 3-D models' predictive skill" -> "the predictive skill of 3-D models"

l33: "ensemble model" -> "model ensemble"

p7

l2: "This provides a total number of 128 ..." This sentence should probably precede the previous one. Altogether, the ensemble generation is needs a bit more explanation.

l9: What is "in situ satellite SST"?

l15: "there is no" rephrase to "is zero"

p8

l12: This paragraph is very difficult to understand. How does a simulation from 6-1200m minimize the computational cost? Is the maximum depth of the model at Cariaco also at 500m? What does "The level thickness is similar to that in 63 depth levels." mean? Do the 2 stations with reduced depth have fewer layers? You do not need to answer these questions in your reply, please rephrase the paragraph so that the reader knows.

l19: "We use ..." The sentence does not say what the metrics are used for.

p9

l9: I don't this paragraph is very helpful. Instead of explaining in what order the metrics will be used (the order is not followed anyway: in 3.1, NRR is mentioned before RMSE), it could explain the order the following subsections will be in.

p10

l11: "range": Is this the ensemble range? Make this explicit.

l11: "If we only allow one process function at a time to change": I still don't think the new explanation helps much to understand what is going on. The "at a time" seems to imply that something is changing over time. I would suggest: "by generating an ensemble in which the functional form of only one of the processes is modified"

l13: "covering 84% of the all ensemble members": what does this mean, please explain better.

p11

l5: "A DCM occurs when lower chlorophyll is detected at the surface": Not sure if a description of DCM is needed but this one is more confusing than helpful.

p13

l3: "the in situ" add "data"

l30: "At most stations, the observed phenology metrics are covered by the ensemble range.": "covered by" does not make it clear what you are trying to say, I would suggest to rephrase using "fall within".

p15

l26: "Through this approach, we provide a perturbed biology ensemble conditioned upon process structural uncertainties." This is a very dense sentence

p16

l2: "Even though fitting..." This sentence is difficult to understand, please rephrase.

Hide

ED: Publish subject to minor revisions (review by editor) (03 Oct 2018) by Katja Fennel

AR by Prima Anugerahanti on behalf of the Authors (16 Oct 2018) Author's response Manuscript

ED: Publish as is (23 Oct 2018) by Katja Fennel

AR by Prima Anugerahanti on behalf of the Authors (26 Oct 2018)

Download

Article (6053 KB)
Full-text XML

Short summary

Minor changes in the biogeochemical model equations lead to major dynamical changes. We assessed this structural sensitivity for the MEDUSA biogeochemical model on chlorophyll and nitrogen concentrations at five oceanographic stations over 10 years, using 1-D ensembles generated by combining different process equations. The ensemble performed better than the default model in most of the stations, suggesting that our approach is useful for generating a probabilistic biogeochemical ensemble model.