Uncertainty analysis of gross primary production partitioned from net ecosystem exchange measurements

Gross primary production (GPP) can be separated from flux tower measurements of net ecosystem exchange (NEE) of CO2. This is used increasingly to validate processbased simulators and remote-sensing-derived estimates of simulated GPP at various time steps. Proper validation includes the uncertainty associated with this separation. In this study, uncertainty assessment was done in a Bayesian framework. It was applied to data from the Speulderbos forest site, The Netherlands. We estimated the uncertainty in GPP at half-hourly time steps, using a non-rectangular hyperbola (NRH) model for its separation from the flux tower measurements. The NRH model provides a robust empirical relationship between radiation and GPP. It includes the degree of curvature of the light response curve, radiation and temperature. Parameters of the NRH model were fitted to the measured NEE data for every 10-day period during the growing season (April to October) in 2009. We defined the prior distribution of each NRH parameter and used Markov chain Monte Carlo (MCMC) simulation to estimate the uncertainty in the separated GPP from the posterior distribution at half-hourly time steps. This time series also allowed us to estimate the uncertainty at daily time steps. We compared the informative with the non-informative prior distributions of the NRH parameters and found that both choices produced similar posterior distributions of GPP. This will provide relevant and important information for the validation of process-based simulators in the future. Furthermore, the obtained posterior distributions of NEE and the NRH parameters are of interest for a range of applications.

of net ecosystem exchange (NEE).They use a Bayesian approach, moving from the regression analysis of rectangular hyperbola fitting daytime data.The argument treated is within the scope of Biogeosciences and the computational instrument they developed is promising.Nevertheless, at this current stage the study suffer from several limitations, not only in the presentational form but also in the substance.In fact, several sources of uncertainty exist in the partitioning the GPP from eddy covariance data.GPP is not directly measured and must be extrapolated from available NEE.Both extrapolation approaches, from night-time or daytime data, can suffer from systematic errors.At least, the authors should acknowledge that a) day respiration can be significantly different from night respiration in reason of the different processes occurring at leaf level (photorespiration or dark respiration).Lower respiration values are expected during the day, see Sun et al. (2015), although compensatory effect could occur (Reichstein et al, 2005); b) the shape of the light response curve measured by eddy covariance can be significantly biased by an inadequate quantification of the storage contribution of the NEE flux, particularly if measurements are taken above high canopies like in the present study.I'm unsure that the authors can quantify these potential sources of bias using the data they have, but at least they should clearly state that they analysed only a component of the possible uncertainty sources.Overall, the approach used by Beer et al., 2010, still seems more solid.
AR 1: Indeed these are important issues.We used a temperature dependent respiration term that is equal for day and night.The fact that we have not parameterized respiration separately for day and night (or separately for vegetation and soil) can be considered as limitations of our model.Indeed the respiration of the vegetation may depend on other factors such as irradiance (Sun et al., 2015).Some other terms are also not included our respiration estimates: photorespiration (because it is nearly proportional to GPP), and respiration terms of which the produced CO2 remains in the trees (Teskey et al., 2008).Our GPP estimates are uncorrected for these particular respiration terms.It is technically possible to extend the model by including more parameters and statistically test whether the model fit with additional parameters performs C8220 better, but we doubt whether this will improve the estimates of GPP, given the limitations of our data: We do not have independent measurements of respiration terms (by means of gas chambers).The night time storage can indeed affect the diurnal shape of the NEE measurements, in particular when CO2 builds up below the sensor height during stable nights which is released when the surface layer becomes unstable after sunrise.Even though some of the stable night time data have been excluded by the quality filtering due to insufficient turbulence, we cannot exclude that the diurnal cycle is affected by storage effects.In the revision, we will acknowledge these limitations of our study.Nevertheless, we believe that the method to establish credible intervals (metrics of uncertainty) in the parameter estimates for the non-rectangular hyperbola (NRH) model with temperature dependent respiration term is useful.For most sites this is still the best we can do, and by providing an algorithm at least the uncertainties of the model parameters can be estimated.The approach of Beer et al. (2010) is an extensive study with a wide scope.They partitioned GPP (at FLUXNET sites) from NEE both using a rectangular hyperbola (RH) light-response curve (Lasslop et al., 2010) and conventional night-time data based approach (Reichstein et al., 2005).They used these partitioned GPP to calibrate five highly diverse diagnostic models for GPP to produce the distribution of global GPP.In our study, we tried to understand better the uncertainty in partitioning GPP using NRH light-response curve.Future research can build on our findings and extend our approach in order to estimate uncertainty across various flux tower sites with different GPP characteristics.We have already highlighted the importance of our findings in the manuscript (P13988 L21 to L 23).We will clarify all the above points in our revised manuscript.
Comments of Referee #1 on the structure of the paper: RC 2: Much of the text is used to present and define the computational approach used.If this could perfectly fit for a Journal like Geoscientific Model Development, it could become probably excessively complicated for a wide audience like that of Biogeosciences.Even more importantly, the sections defined as '4.Results' and '5.Dis-C8221 cussion', are similarly presenting some of the results and partially discuss them.In the revision, I recommend to decide if the paper will have a discussion section separated from the presentation of the results, since what is done in the current text version is confusing.Overall, I recommend this study for publication, but only after a thorough work of revision.
AR 2: We understand that computational approach used in this study may be complicated to some of the audience of Biogeosciences.We have tried to explain the flux partitioning in a Bayesian framework in a clear way and we will give this further attention in the revised manuscript.We hope the audience of Biogeoscience will not have any problem to understand our computational approach.We will follow the advice of the referee to combine the Results and Discussion section in the revised manuscript.
Specific comments of Referee #1: RC 3: P13972 L23 'ecosystem scientists' use positive values of GPP.To my knowledge, GPP is always positive since it represents a production, and production was represented with positive values much earlier than Brahmagupta, in the 7th century, described negative numbers.If we use the micrometeorological approach, the most correct term will be probably gross ecosystem exchange (GEE), which is defined as negative since it represents the quantity of CO2 which enters in the ecosystem.
AR 3: We have also used positive sign for GPP.To avoid the confusion that only ecosystem scientists use positive value for GPP, we will remove the term "ecosystem scientist" in line P13971 L23 in the revised manuscript.We will also rephrase P13971 L23 to L25 to clarify that GPP is always positive as it represents production.
RC 4: P13972L16-18 Units are missing in the passage from 'Pa' to 'GPP'.In fact, NEE is generally expressed in terms of micromoles m-2 s-1, and GPP with the same units or in terms of g m-2 d-1.Here, in table 1, both quantities are expressed in terms of g CO2 m-2 s-1, so the proposed conversion factor of 12/44 shouldn't be present.In any case, it converts for instance micromoles of CO2 into micrograms of C.

C8222
AR 4: The conversion of unit from µmoles CO2 m-2 s-1 to g C m-2 s-1 requires the multiplication of factor by 12/44.Originally, we had the unit of NEE in µmoles CO2 m-2 s-1 and we converted it to mg CO2 m-2 s-1 (This unit appears in Table 1 for NEE and also for Pa) that does not require the multiplication by factor 12/44.We have expressed GPP in g C m-2 s-1, so we multiplied Pa by the factor 12/44 in order to account only carbon (C) in the unit of GPP.We agree that this is not clearly expressed in the manuscript.We will give clear explanation of unit and its conversion in the revised manuscript.
RC 5: P13974 L13-14 'NEE data were corrected for storage of CO2 in the air between the sensor and the ground'.The storage is a relevant component or in the mixing ratio conservation equation (e.g, Kowalski and Argueso, 2011) and hence in the correct computation of NEE.Although the storage terms tend to cancel out when producing annual sums, it is well established that they can asymmetrically influence the apparent light response curve, presenting opposite values in the evening and in the morning and a significant hysteresis.It has been clearly established that a convenient number of measurement points have to be established when measurements are dove above forests (Yang et al., 2007), but there is no mention about the instrumentation used.
More information is needed!AR 5: We regret that we have made a mistake in the data description.We have not carried out a storage correction of the CO2 flux since profile measurements were unavailable.We had installed a second Eddy Covariance system at 35 m height, and there was no significant difference in magnitude of the flux and shape of the diurnal cycle compared to the height of 46 m.However, we did not have flux or concentration measurements below this height, thus within the canopy, the place where storage effects will be most significant.The storage effects on daily or annual NEE may be limited, but we agree that it most likely affects the light response curves.The morningafternoon hysteresis effect will contribute to higher scatter in the model fit, and thus to more uncertainty in the retrieved parameters and the GPP estimates.We will discuss C8223 this in the revised manuscript.
RC 6: P13982 'credible interval spanned zero'.After this sentence there is one paragraph of discussion.What the authors mean exactly?Was the mean zero and the distribution of values partially above it, or there was a bias?AR 6: A Bayesian prediction gives a probability distribution, this can be summarized by the median and the 95% credible interval as a quantification of uncertainty.Hence the actual value is likely to be somewhere in this interval and not necessarily at the median.The predicted median is negative during the night; however, the credible interval spans zero, implying that the actual value could be positive.Clearly this is not physically possible, but is an artefact of the statistical approach and highlights that we are indeed uncertain about these predictions.We have explained the possible reasons in P13983 from L1 to L12 and will clarify this in the revised manuscript.
AR 7: OK, we will provide the reference here in the revised manuscript.
RC 8: P13978: 'the photosynthetic capacity. ..is reached when the photosynthesis is Rubisco limited varies among different tree species'.At canopy level, the photosynthetic capacity depends also on the structure of the canopy, when multiple leaf layers are present Amax increases.Please read carefully the cited paper of Ruimy et al., 1995.AR 8: Non-rectangular hyperbola (NRH) model used in this study requires Amax at canopy level.We agree with the referee that the photosynthesis capacity at the canopy level also depends on the structure of the canopy (i.e.arrangement of the canopy leaves) and the area of leaves available to absorb photons.Both are determined by leaf area index (Ruimy et al., 1995).It should be noted that we have also used leaf area index to determine canopy Amax using a model that incorporates a radiative transfer scheme and a vertically declining needle-level Amax.We will rephrase the sentence C8224 in P13977 L25 to L27 to give clear explanation of canopy level Amax in the revised manuscript.
RC 9: P13978 L15.'in the literature were 0.0097. ..', again, please be specific on the sources and also add units (is mg CO2 m-2s-1 valid for all the numbers reported?).
AR 9: We have already provided the sources in L16 at the end of the sentence.The unit of mg CO2 m-2s-1 is valid for all numbers reported.
RC 10: P13980 L 13-14: 'but was short enough that we could observe temporal change between the 10-day blocks'.This is strange, possibly the authors were meaning the contrary (long enough)?
AR 10: When selecting the block length there was a trade-off between selecting a block that was LONG enough to contain sufficient measurements for accurate modelling yet SHORT enough to allow us to observe changes between the blocks.Thus the temporal change is observed between consecutive blocks, not within a block.We will clarify this in the revised manuscripts.
RC 11: P13982 L6-7: 'The chains were thoroughly interdigitating, indicating that the the Markov chains had mixed and converged. ..'Besides the repeated article, I cannot understand.In any case, I recommend to avoid lab jargon.
AR 11: We have explained in P13973 L20 to L27 that the stationary distribution of the Markov chains is the posterior distribution of parameters.The stationary distribution can be assessed graphically (Fig. S3 in the supplementary file) when the Markov chains interlock with each other (in other words, chains thoroughly interdigitate) showing the proper mixing of the chains.We will clarify this in the revised manuscript.AR 14: At the end of figure 4 caption, we have mentioned that information about NRH parameters (showed as symbol in the Y axis of the plots) is given in Table 1 that shows the meaning and the unit of the symbol of NRH parameters.
RC 15: Figure 5: What are the units in the Y axes?I add that in this set of images, a magnifier is needed to distinguish what is reported along the axes, at least for many readers including myself.
AR 15: We have written "as Fig. 4" at the start of the figure 5 caption.This means that our response AR 14 is also valid here for the unit.We will redraw the figure 4 and 5 to make the values along the axis more visible in the revised manuscript.
Minor/language remarks of referee 2: RC 16: Page (P) 1397 Line (L) 21 'a non linear empirical models': please check the consistency between article and noun.P1397 L26-27: 'a single optimized values', same as above.
AR 16: We will correct this in the revised manuscript.
RC 17: P13974 L8 Cambell->Campbell AR 17: We will correct this in the revised manuscript.
RC 18: P13982: Lay->Laid AR 18: We will correct this in the revised manuscript.
RC 19: P13984 L24-25.'In order to undertake a Bayesian analysis it is necessary to specify the prior distributions on the NRH parameters.' Written in this way, it seems C8226 that the use of NRH parameters is a general rule in Bayesian analysis, but it is not.
AR 19: We agree with the referee.We will rephrase this sentence in the revised manuscript.
RC 20: 13985 L8-9: 'This wide variation in Amax was chosen as the non-informative priors led to spikes in the value of Amax in the posterior (Fig. 5e).' Remove 'was' and possibly reformulate the sentence.
AR 20: We will rephrase this sentence in the revised manuscript.

RC 12 :
Figure 1: Please define what Y axis represents.AR 12: Y axis represents the density of the distribution of NRH parameter.We will add this in the caption of figure 1 in the revised manuscript.C8225 RC 13: Figure 3: What are the frequency units in the Y axes?AR 13: In the histogram, the Y axis has no unit.Y axis represents the frequency or the number of GPP points.RC 14: Figure 4: What are the units in the Y axes?