# Partitioning net carbon dioxide fluxes into photosynthesis and respiration using neural networks

Gianluca Tramontana, Mirco Migliavacca, Martin Jung, Markus Reichstein, Trevor Keenan, Gustau Camps-valls, Jérôme Ogée, Jochem Verrelst, Dario Papale

► To cite this version:

Gianluca Tramontana, Mirco Migliavacca, Martin Jung, Markus Reichstein, Trevor Keenan, et al. Partitioning net carbon dioxide fluxes into photosynthesis and respiration using neural networks. Global Change Biology, 2020, 26 (9), pp.5235 - 5253. 10.1111/gcb.15203. hal-03173753

HAL Id: hal-03173753

https://hal.inrae.fr/hal-03173753

Submitted on 18 Mar 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Partitioning net carbon dioxide fluxes into photosynthesis and respiration using neural networks

Gianluca Tramontana<sup>1,2 ID</sup> | Mirco Migliavacca<sup>3</sup> | Martin Jung<sup>3</sup> | Markus Reichstein<sup>3</sup> |
Trevor F. Keenan<sup>4,5</sup> | Gustav Camps-Valls<sup>2</sup> | Jerome Ogee<sup>6 ID</sup> | Jochem Verrelst<sup>2</sup> |
Dario Papale<sup>1,7</sup>

## Correspondence

Gianluca Tramontana, Image Processing Laboratory (IPL), Universitat de València, Parc Científic Universitat de València, Calle Catedrático José Beltrán, 2, 46980 Paterna (València), Spain.

Email: gianluca.tramontana@uv.es; g.tramontana@unitus.it

## Funding information

European Commission, Grant/Award Number: 640176, 755617, 647423 and 703944; NASA, Grant/Award Number: NNH17AE861; the EU's Horizon 2020 Programme, Grant/Award Number: GA 820852; the EU's Horizon 2020 Programme, Grant/Award Number: GA 703944

## Abstract

The eddy covariance (EC) technique is used to measure the net ecosystem exchange (NEE) of CO₂ between ecosystems and the atmosphere, offering a unique opportunity to study ecosystem responses to climate change. NEE is the difference between the total CO₂ release due to all respiration processes (RECO), and the gross carbon uptake by photosynthesis (GPP). These two gross CO₂ fluxes are derived from EC measurements by applying partitioning methods that rely on physiologically based functional relationships with a limited number of environmental drivers. However, the partitioning methods applied in the global FLUXNET network of EC observations do not account for the multiple co-acting factors that modulate GPP and RECO flux dynamics. To overcome this limitation, we developed a hybrid data-driven approach based on combined neural networks (NN<sub>C-part</sub>). NN<sub>C-part</sub> incorporates process knowledge by introducing a photosynthetic response based on the light-use efficiency (LUE) concept, and uses a comprehensive dataset of soil and micrometeorological variables as fluxes drivers. We applied the method to 36 sites from the FLUXNET2015 dataset and found a high consistency in the results with those derived from other standard partitioning methods for both GPP (R<fnref n="2" /> > .94) and RECO (R<fnref n="2" /> > .8). High consistency was also found for (a) the diurnal and seasonal patterns of fluxes and (b) the ecosystem functional responses. NN<sub>C-part</sub> performed more realistic than the traditional methods for predicting additional patterns of gross CO₂ fluxes, such as: (a) the GPP response to VPD, (b) direct effects of air temperature on GPP dynamics, (c) hysteresis in the diel cycle of gross CO₂ fluxes, (d) the sensitivity of LUE to the diffuse to direct radiation ratio, and (e) the post rain respiration pulse after a long dry period. In conclusion, NN<sub>C-part</sub> is a valid data-driven approach to provide GPP and RECO estimates and complementary to the existing partitioning methods.

## KEYWORDS

carbon dioxide fluxes partitioning, ecosystem respiration (RECO), eddy covariance, gross primary production (GPP), machine learning, net ecosystem exchange, neural network

<footnote n="1"><sup>1</sup>DIBAF, Department for Innovation in Biological, Agro-food and Forestry Systems, University of Tuscia, Viterbo, Italy</footnote>

<footnote n="2"><sup>2</sup>Image Processing Laboratory (IPL), Parc Científic Universitat de València, Universitat de València, Paterna, Spain</footnote>

<footnote n="3"><sup>3</sup>Max Planck Institute for Biogeochemistry, Jena, Germany</footnote>

<footnote n="4"><sup>4</sup>Department of Environmental Science, Policy and Management, UC Berkeley, Berkeley, CA, USA</footnote>

<footnote n="5"><sup>5</sup>Earth and Environmental Sciences Area, Lawrence Berkeley National Lab, Berkeley, CA, USA</footnote>

<footnote n="6"><sup>6</sup>INRAE, UMR 1391 ISPA, Villenaved'Ornon, France</footnote>

<footnote n="7"><sup>7</sup>Euro-Mediterranean Center on Climate Change (CMCC), Viterbo, Italy</footnote>
# 1 | INTRODUCTION

The eddy covariance (EC) technique offers a unique opportunity for monitoring carbon and energy exchanges between land ecosystems and the atmosphere (Baldocchi, 2003) allowing near-continuous measurements integrated at the ecosystem scale. The number of study sites equipped with EC systems has increased over the years (Baldocchi, 2014; Chu, Baldocchi, John, Wolf, & Reichstein, 2017; Pastorello et al., 2017), and now we estimate that more than 1,400 globally distributed sites (Chu et al., 2017) are monitoring the most representative land ecosystems in different climate conditions. Most of the EC study sites are organized in regional networks such ICOS, AmeriFlux, and AsiaFlux, and contribute to the global FLUXNET network (Baldocchi, 2014).

The EC method allows for the measurement of the net ecosystem exchange (NEE) which is the difference between two larger flux components: gross primary production (GPP) and ecosystem respiration (RECO). GPP is the gross amount of carbon uptake by photosynthesis from vegetation while RECO is the total carbon efflux by the respiration processes of all organisms. Estimating GPP and RECO is a key step to better understand the underlying mechanisms constraining ecosystem function. Moreover, GPP and RECO estimates from EC are useful for modeling, supporting process-based model parameterization and validation, data assimilation, plant trait retrieval by model inversion (Dutta, Schimel, Sun, van der Tol, & Frankenberg, 2019; Pacheco-Labrador et al., 2019), upscaling (Jung et al., 2020; Tramontana et al., 2016), as well as photosynthesis estimates based on remote sensing (e.g., Arnone et al., 2008; Verrelst et al., 2016; Zhang et al., 2016).

The EC technique does not directly measure GPP and RECO, and for this reason, numerical methods (termed partitioning methods) have been proposed for estimating GPP and RECO from NEE measurements (e.g., Desai et al., 2008; Keenan et al., 2019; Lasslop et al., 2010; Reichstein et al., 2005; Sulman, Tyler Roman, Scanlon, Wang, & Novick, 2016). The most widely used approaches are based on the use of NEE measurements for fitting simple physiologically based nonlinear relationships to estimate GPP and RECO using few meteorological drivers. These functional relationships are in general either light response functions linking global incoming radiation and GPP (Gilmanov et al., 2003), or respiration functions based on temperature (Reichstein et al., 2005), or also combinations of the two approaches (Keenan et al., 2019; Lasslop et al., 2010). The simple implementation and the robustness of the results (Lasslop et al., 2010) have led to their adoption as standard processing tools in the FLUXNET community (Pastorello et al., 2020, accepted).

However, these partitioning methods rely on important assumptions about the flux dynamics and their relationship with environmental drivers. Importantly, the assumed functional relations used can be similar to those applied in ecosystem models that make use of these data for their validation, creating a risk of circularity. In addition, although the functional relationships used in standard partitioning approaches are known to be valid at the organ level (where

they can be measured), their direct application at the ecosystem spatial scale should be evaluated carefully (Medlyn, 1998; Medlyn et al., 2017; Musavi et al., 2016). Furthermore, to guarantee their wide applicability, partitioning methods use a limited number of drivers (i.e., air temperature, vapor pressure deficit, and global radiation). However, fluxes dynamics are also potentially affected by many important environmental factors, such as soil moisture, soil temperature, or the ratio of diffuse to direct radiation, among others, that are generally not considered in the partitioning methods (Lasslop et al., 2012; Wohlfahrt & Galvagno, 2017). These limitations are partially compensated by the use of short temporal moving windows for parameters estimation, which takes more slowly changing factors indirectly into consideration (such as phenology, water, or substrate availability), but does not handle fast ecosystem responses well, like for example respiration pulses following rain events (Williams, Hanan, Scholes, & Kutsch, 2009).

Some studies have attempted to solve the limitation of FLUXNET standard partitioning methods by developing more comprehensive approaches. For instance, Scanlon and Kustas (2010) coupled CO₂ and H₂O fluxes dynamics in the flux variance similarity approach for simultaneously partitioning carbon and water fluxes. The method, although interesting, requires canopy scale estimates of water-use efficiency (WUE), which introduces assumptions and uncertainty. More recently, the estimation of gross CO₂ fluxes from NEE has been approached using the EC method in combination with parallel measurements of a trace gas such as carbonyl sulfide (COS; Commane et al., 2015) or ¹³C isotopes (Ogée, Peylin, et al., 2003; Oikawa et al., 2017; Wehr et al., 2016; Wehr & Saleska, 2015) with the aim to disentangle the photosynthesis signals from respiration in daytime NEE measurements. Both methods are promising and are starting to be applied in the field; however, they currently require extensive and expensive instrumentation and the uncertainty in the results is still large (Oikawa et al., 2017; Whelan et al., 2018), limiting, for now, their application to a restricted number of study sites. These methods, however, are all subject to assumptions that could affect the resulting estimates of GPP and RECO. A possible alternative approach, based on existing measurements and not subject to the limitation of the standard methods, could be provided by machine learning methods.

Machine learning methods are generally used to model underlying complex relationships linking a set of predictors with one or more outputs. The core characteristics of machine learning methods are that predictors are not prescribed and the relationships between input (drivers) and output (fluxes) are inferred from the data. Several studies have reported the capability of machine learning to reproduce complex patterns in ecological studies, and in particular in relation to EC measurements (e.g., see Moffat, Beckstein, Churkina, Mund, & Heimann, 2010; Moffat et al., 2007; Papale & Valentini, 2003; Reichstein et al., 2019; Tramontana et al., 2016).

In this study, we develop and test an empirical machine learning-based approach for retrieving GPP and RECO from NEE. The main motivation was to evaluate if a purely empirical approach could provide estimates of the two components without any predefined
relationships and with the flexibility to use multiple sets of meteorological and biological forcings.

In this experiment, we used artificial neural networks (ANNs), which have already been tested as a partitioning tool in a cross-site intercomparison among several partitioning methods (Desai et al., 2008; Oikawa et al., 2017), without appreciable differences with respect to the other methods. In those studies, the ANNs were trained using only nighttime measurements and then applied to extrapolate RECO during daytime. Here, we revisited the ANN implementation by proposing a new scheme in which the ANN is constrained using the measured net CO₂ fluxes to directly estimate the two components GPP and RECO.

Our aim was to design a machine learning-based method flexible enough to: (a) use a large set of drivers, accounting for the complexity of the GPP and RECO responses, (b) ensure general applicability of the method to allow ANN-based partitioning in the context of FLUXNET thanks to the use of generally measured variables. In order to evaluate the results, we compared the GPP and RECO produced by the ANN with the estimates obtained by two partitioning methods used in FLUXNET, by analyzing: (a) the consistency of the estimates at yearly, seasonal, daily, and (half) hourly time steps; (b) the dynamics of the seasonal and diurnal cycles of the estimated gross CO₂ fluxes; (c) the functional relationships between micrometeorological inputs and fluxes in output as reproduced by the ANN; and (d) the effectiveness of the proposed method to predict additional patterns of gross CO₂ fluxes not (or less) accounted for by the standard partitioning methods.

# 2 | MATERIALS AND METHODS

## 2.1 | ANN algorithm

ANNs are nonlinear and nonparametric methods for regression and classification that artificially emulate the functioning of a biological brain (Haykin, 1999). The base unit of an ANN is the “neuron” where the numeric information in the input is weighted, condensed, and transformed (by a linear/nonlinear transfer/activation function) to be transferred to other neurons. Neurons are organized in layers that are interconnected: the outputs of m neurons in one layer are the inputs for n neurons of the next layer and the signals are transferred through connections associated with multiplicative weights. The learning (or training) process of an ANN consists of adjusting the weights of the network on the basis of specific examples provided as input (supervised training).

In this experiment, we developed a new ANN architecture designed to emulate<fnref n="1" /> the ecosystems processes driving NEE. This customized neural network (hereafter NN<sub>C-part</sub>) is based on the concept

that NEE, measured by the EC system, is the difference between RECO and GPP (Equation 1):

$$NEE = RECO - GPP, \qquad (1)$$

and that these two fluxes have individual drivers and dependencies. The uniqueness of our approach consists of imposing physical constraints in the proposed neural network, on the basis of the known properties of RECO and GPP. The overall structure is based on two subnetworks (Figure 1): one for retrieving RECO (sub neural network, SNN<sub>RECO</sub>) and the other for retrieving GPP (SNN<sub>GPP</sub>). The output of these two subnetworks (S<sub>RECO</sub> and S<sub>GPP</sub> for SNN<sub>RECO</sub> and SNN<sub>GPP</sub>, respectively) is combined in the last neuron of the overall structure, where NEE is calculated and compared with the measurements in order to optimize the network's weights (network training). The output signals of both subnetworks are constrained to be always positive, except for SNN<sub>GPP</sub>, which is constrained to be 0 during nighttime, as photosynthesis requires light. However, RECO and GPP have opposite signs; therefore, the connection weight of RECO is fixed to be positive (w<sub>RECO</sub> = 1) and that of GPP to be negative (w<sub>GPP</sub> = -1), mirroring the sign convention of the measured NEE. Since the last transfer function is unbounded and linear, and the bias is fixed to 0, the last node effectively reproduces the equation of NEE as:

$$NEE_{ANN} = RECO_{ANN} + GPP_{ANN}, \qquad (2)$$

where

$$RECO_{ANN} = S_{RECO} * W_{RECO}, \qquad (3)$$

and

$$GPP_{ANN} = S_{GPP} * W_{GPP}, \qquad (4)$$

The two subnetworks have different input drivers but a similar structure including two hidden layers: the first hidden layer has n hidden neurons (different for the two fluxes) and a hyperbolic tangent sigmoid activation function (“Tansig”), which improved the performance of network optimization; the second has only one neuron and a log-sigmoidal activation function (“Logsig”) that constrains the output to positive values, allowing the implementation of Equation (2) in the last node. The SNN<sub>GPP</sub> used to simulate GPP, has also a third layer where the measured incoming shortwave radiation enters as input in a node where a product is applied. This last node has a positive linear transfer function (“Poslin”); thus, the output from the SNN<sub>GPP</sub> is 0 at nighttime and positive during the day. With this structure, the last node of the SNN<sub>GPP</sub> mimics the light-use efficiency (LUE) approach: the output coming from the two previous hidden layers is an LUE that is then multiplied by the available shortwave radiation to estimate GPP. In this way, the drivers are used to estimate an instantaneous proxy of the LUE which takes into consideration seasonal and diurnal variability of carbon uptake, including photosynthesis saturation at high light, while shortwave radiation mainly defines the magnitude of

<footnote n="1">¹The use of this term is not incidental: surrogate modeling, also known as emulation, is about developing statistical models that learn to mimic costly physical models, such as radiative transfer or climate models, using a representative dataset of simulation/forward runs. ML algorithms learn such parameter-observation relation with high accuracy, and after training can be used for forward simulation very efficiently (Bastos and O'Hagan, 2009; Camps-Valls et al., 2019; Rivera et al., 2015; Vicent, et al 2018).</footnote>
<figure><img src="image_2.png" /><figcaption>**FIGURE 1** The scheme of the customized neural network applied for retrieving gross primary production (GPP) and ecosystem respiration (RECO) from net ecosystem exchange (NEE) measurements. The two subnetworks ($S_{\text{GPP}}$ left, in green, $S_{\text{RECO}}$ right, in brown) are connected in the last step to estimate NEE. Inputs for RECO are air and soil temperature; soil water content; wind speed and wind direction; the day of the year (sine and cosine values of angular day of year); the average value of nighttime NEE. Input GPP are, in addition to the shortwave incoming radiation used in the product, air temperature; vapor pressure deficit; soil water content; potential and actual shortwave incoming radiation; wind speed and wind direction; a proxy of the mean seasonal GPP dynamic derived from the nighttime and daytime averaged NEE. See Section 2.3 of the main text for details</figcaption></figure>

the flux. This decreases the weight given to the incoming shortwave radiation in relation to other driving variables. Up to this step, the output of the two subnetworks ($S_{\text{RECO}}$ and $S_{\text{GPP}}$) have positive signs; the sign convention of the NEE equation is finally restored in the last neuron of the overall structure by associating a positively signed weight to $S_{\text{RECO}}$ and a negatively signed weight to $S_{\text{GPP}}$.

## 2.2 | Dataset used

For the purpose of this experiment, we used data from the FLUXNET2015 dataset (http://fluxnet.fluxdata.org) and, as an additional comparison, a synthetic dataset generated by the Multi-layer Simulator of the Interactions between a vegetation Canopy and the Atmosphere (MuSICA) process-based model (Ogée, Brunet, Loustau, Berbigier, & Delzon, 2003, see section S4).

FLUXNET2015 includes meteorological and EC measurements that were quality checked and processed with standard tools (Papale et al., 2006; Pastorello et al., 2020) and provided with per-variable quality flags. For more details, see the webpage http://fluxnet.fluxd ata.org/data/fluxnet2015-dataset/.

Data used for training and validation of the neural network were taken from the "FULLSET" "TIER 1" collection. Among the sites available in this collection, we selected a subset of 36 study sites (listed in Table 1) on the basis of the data quality and data availability

in order to ensure the best conditions for the partitioning methods comparison. Site-years were selected if the percentage of meteorological gap-filled data was less than 20% and the measured NEE covered at least 10% of both daytime and nighttime periods.

Among the variables stored in the FLUXNET2015 dataset, we used the GPP and RECO obtained with the two standard methods: the nighttime method (NT) from Reichstein et al., 2005, and the daytime method (DT) from Lasslop et al., 2010. In the NT method, the Arrhenius-type temperature-response curve of respiration as modeled by the Lloyd and Taylor equation (Lloyd & Taylor, 1994), driven by air temperature, is used for estimating RECO. This method makes use of the nighttime (when global incoming radiation <20 W/m²) NEE observations (assumed as representative of RECO given the absence of photosynthesis at night) to fit the Lloyd and Taylor equation. There are two parameters to fit in the Lloyd and Taylor equation: the activation energy (E0) and the respiration rate at the reference temperature ($R_{\text{ref}}$). In the NT method, E0 is estimated at an annual scale by calculating E0 values every 15 days and then averaging the three with smaller uncertainty. Once E0 is fixed, the $R_{\text{ref}}$ parameter is estimated using short-term moving windows (4 days). On the basis of Equation (1), GPP is finally calculated by subtracting NEE from the daytime RECO extrapolated by the fitted model.

The DT method combines a rectangular hyperbola light response curve (Falge et al., 2001) for estimating GPP and the Lloyd and Taylor equation for estimating RECO (as in the NT method). In the
<table><caption>TABLE 1 List of the subset of FLUXNET2015 study sites used in this experiment for net ecosystem exchange partitioning. The study sites used also for NN<sub>C-part</sub> validation (see section S1) are marked with the symbol (*) in the "Validation" column. Selected study sites represent the following International Geosphere-Biosphere Programme (IGBP) vegetation classes: Evergreen Needleleaf Forest (ENF), Deciduous Broadleaf Forest (DBF), Mixed Forest, Evergreen Broadleaf Forest (EBF), Grassland (GRA), Cropland (CRO), Woody SAvanna (WSA), Closed SHrubland (CSH), and Open SHrubland (OSH). For accurate coordinates, please refer to the FLUXNET2015 website</caption><thead><tr><th>ID</th><th>Site code</th><th>IGBP</th><th>Lat</th><th>Lon</th><th>Validation</th><th>Reference</th></tr></thead><tbody><tr><td>1</td><td>AU-Cpr</td><td>SAV</td><td>-34,00</td><td>140,59</td><td>*</td><td>Meyer, Kondrlovà, and Koerber (2015)</td></tr><tr><td>2</td><td>AU-DaP</td><td>GRA</td><td>-14,06</td><td>131,32</td><td></td><td>Beringer et al. (2011)</td></tr><tr><td>3</td><td>AU-Dry</td><td>SAV</td><td>-15,26</td><td>132,37</td><td></td><td>Cernusak, Hutley, Beringer, Holtum, and Turner (2011)</td></tr><tr><td>4</td><td>AU-How</td><td>WSA</td><td>-12,49</td><td>131,15</td><td></td><td>Beringer, Hutley, Tapper, and Cernusak (2007)</td></tr><tr><td>5</td><td>AU-Stp</td><td>GRA</td><td>-17,15</td><td>133,35</td><td>*</td><td>Beringer et al. (2011)</td></tr><tr><td>6</td><td>BE-Lon</td><td>CRO</td><td>50,55</td><td>4,75</td><td></td><td>Moureaux, Debacq, Bodson, Heinesch, and Aubinet (2006)</td></tr><tr><td>7</td><td>BE-Vie</td><td>MF</td><td>50,31</td><td>6,00</td><td></td><td>Aubinet et al. (2001)</td></tr><tr><td>8</td><td>CA-Qfo</td><td>ENF</td><td>49,69</td><td>-74,34</td><td>*</td><td>Giasson, Coursolle, and Margolis (2006)</td></tr><tr><td>9</td><td>DE-Geb</td><td>CRO</td><td>51,10</td><td>10,91</td><td>*</td><td>Anthoni et al. (2004)</td></tr><tr><td>10</td><td>DE-Gri</td><td>GRA</td><td>50,95</td><td>13,51</td><td>*</td><td>Prescher, Grünwald, and Bernhofer (2010)</td></tr><tr><td>11</td><td>DE-Kli</td><td>CRO</td><td>50,89</td><td>13,52</td><td>*</td><td>Prescher et al. (2010)</td></tr><tr><td>12</td><td>DE-Obe</td><td>ENF</td><td>50,79</td><td>13,72</td><td>*</td><td></td></tr><tr><td>13</td><td>DE-Tha</td><td>ENF</td><td>50,96</td><td>13,57</td><td>*</td><td>Grünwald and Bernhofer (2007)</td></tr><tr><td>14</td><td>DK-Sor</td><td>DBF</td><td>55,49</td><td>11,64</td><td>*</td><td>Pilegaard, Ibrom, Courtney, Hummelshøj, and Jensen (2011)</td></tr><tr><td>15</td><td>FI-Hyy</td><td>ENF</td><td>61,85</td><td>24,29</td><td>*</td><td>Suni et al. (2003)</td></tr><tr><td>16</td><td>FR-LBr</td><td>ENF</td><td>44,72</td><td>-0,77</td><td></td><td>Berbigier, Bonnefond, and Mellmann (2001)</td></tr><tr><td>17</td><td>GF-Guy</td><td>EBF</td><td>5,28</td><td>-52,92</td><td></td><td>Bonal et al. (2008)</td></tr><tr><td>18</td><td>IT-BCi</td><td>CRO</td><td>40,52</td><td>14,96</td><td></td><td>Vitale, Di Tommasi, D'Urso, and Magliulo (2016)</td></tr><tr><td>19</td><td>IT-Cp2</td><td>EBF</td><td>41,70</td><td>12,36</td><td>*</td><td>Fares, Savi, Muller, Matteucci, and Paoletti (2014)</td></tr><tr><td>20</td><td>IT-Cpz</td><td>EBF</td><td>41,71</td><td>12,38</td><td>*</td><td>Garbulsky, Penuelas, Papale, and Filella (2008)</td></tr><tr><td>21</td><td>IT-MBo</td><td>GRA</td><td>46,01</td><td>11,05</td><td>*</td><td>Marcolla et al. (2011)</td></tr><tr><td>22</td><td>IT-Noe</td><td>CSH</td><td>40,61</td><td>8,15</td><td></td><td>Papale, Black, et al. (2015)</td></tr><tr><td>23</td><td>IT-Ro1</td><td>DBF</td><td>42,41</td><td>11,93</td><td></td><td>Rey et al. (2002)</td></tr><tr><td>24</td><td>IT-SRo</td><td>ENF</td><td>43,73</td><td>10,28</td><td>*</td><td>Chiesi et al. (2005)</td></tr><tr><td>25</td><td>NL-Loo</td><td>ENF</td><td>52,17</td><td>5,74</td><td>*</td><td>Moors (2012)</td></tr><tr><td>26</td><td>RU-Fyo</td><td>ENF</td><td>56,46</td><td>32,92</td><td>*</td><td>Kurbatova, Li, Varlagin, Xiao, and Vygodskaya (2008)</td></tr><tr><td>27</td><td>US-ARM</td><td>CRO</td><td>36,61</td><td>-97,49</td><td></td><td>Fischer, Billesbach, Riley, Berry, and Torn (2007)</td></tr><tr><td>28</td><td>US-GLE</td><td>ENF</td><td>41,37</td><td>-106,24</td><td>*</td><td>Frank, Massman, Ewers, Huckaby, and Negrón (2014)</td></tr><tr><td>29</td><td>US-MMS</td><td>DBF</td><td>39,32</td><td>-86,41</td><td>*</td><td>Dragoni et al. (2011)</td></tr><tr><td>30</td><td>US-NR1</td><td>ENF</td><td>40,03</td><td>-105,55</td><td>*</td><td>Monson et al. (2002)</td></tr><tr><td>31</td><td>US-SRG</td><td>GRA</td><td>31,79</td><td>-110,83</td><td>*</td><td>Scott, Biederman, Hamerlynck, and Barron-Gafford (2015)</td></tr><tr><td>32</td><td>US-SRM</td><td>WSA</td><td>31,82</td><td>-110,87</td><td>*</td><td>Scott, Jenerette, Potts, and Huxman (2009)</td></tr><tr><td>33</td><td>US-UMB</td><td>DBF</td><td>45,56</td><td>-84,71</td><td>*</td><td>Gough et al. (2013)</td></tr><tr><td>34</td><td>US-Whs</td><td>OSH</td><td>31,74</td><td>-110,05</td><td>*</td><td>Scott et al. (2015)</td></tr><tr><td>35</td><td>US-Wkg</td><td>GRA</td><td>31,74</td><td>-109,94</td><td></td><td>Scott, Hamerlynck, Jenerette, Moran, and Barron-Gafford (2010)</td></tr><tr><td>36</td><td>ZA-Kru</td><td>SAV</td><td>-25,02</td><td>31,50</td><td></td><td>Archibald et al. (2009)</td></tr></tbody></table>

DT method, nighttime NEE is used only to estimate the E0 parameter of the Lloyd and Taylor equation, while the remaining unknown parameters, including R<sub>ref</sub>, are fitted using daytime measurements. In DT, the light response curve is driven by the incoming shortwave

radiation (SW_IN), but it is adjusted for the effect of stomata closure due to the atmospheric evaporative demand using VPD as an additional driver of GPP (Lasslop et al., 2010). Similar to NT: (a) the Lloyd and Taylor model is driven by air temperature, (b) the model
parameters are fitted on short-term moving time windows to account for slow changing factors.

## 2.3 | Input variables and data preparation

The variable used as target values in the NN<sub>C-part</sub> training was the half hourly NEE (µmol CO<sub>2</sub> m<sup>-2</sup> s<sup>-1</sup>) measured with the EC technique. In particular, we used the NEE_CUT_USTAR50 variable (Pastorello et al., 2020). We used a comprehensive subset of micrometeorological variables measured at the flux towers as candidate drivers for NN<sub>C-part</sub> and additional variables derived from NEE, with the aim to explain variability of carbon fluxes at seasonal, daily, and half hourly resolutions. In particular, in terms of meteorological drivers, we used the SW_IN (W/m<sup>2</sup>), VPD (kPa), air temperature (TA [°C]), soil temperature (TS [°C]), and soil water content (SWC [%]) in the upper soil layer (depth function of the site, see Pastorello et al., 2020), wind speed (WS [m/s]), and wind direction (WD [degrees]). Each one of these variables is involved to different degrees in carbon exchange processes or in NEE measurements. SW_IN carries the photosynthetically active radiation that is a key variable for the light-dependent reactions of photosynthesis, while leaves regulate stomatal conductance in response to VPD (Lasslop et al., 2010); temperature has a key role in chemical reactions of biological processes; thus, it is involved both in photosynthesis and in the RECO (Falge et al., 2001). SW_IN and VPD are used in the DT method to estimate photosynthesis while both DT and NT use TA to estimate RECO only. Wind-related variables (WS and WD) affect the footprint of flux measurements (Arriga et al., 2017; Kljun, Calanca, Rotach, & W. and Schmid, H. P., 2015) which is important in particular if WD changes systematically (e.g., nighttime vs. daytime, morning vs. afternoon) and if the surrounding land cover is heterogeneous. The Julian day of the year (DOY) and the potential incoming radiation (SW_IN_POT [W/m<sup>2</sup>]) were also used in order to provide information about seasonality and length of the day. The SW_IN_POT was also aggregated at daily resolution to better track the seasonality of light conditions; the first derivatives of half hourly and daily SW_IN_POT were also calculated to add specific information about the seasonal and diurnal dynamics of light (Bodesheim, Jung, Gans, Mahecha, & Reichstein, 2018; Papale, Black, et al., 2015; Papale, Migliavacca, et al., 2015).

Finally, since plant photosynthesis and RECO change seasonally also due to substrate availability and management (in case of managed sites), daily average of nighttime NEE (NEE<sub>NIGHT</sub>) and a proxy of GPP (GPP<sub>prox</sub>) were calculated and used in input. In this case, gap-filled (Reichstein et al., 2005) half hourly NEE was used in order to have a continuous time series. Proxy of GPP was calculated by using NEE<sub>NIGHT</sub>, and the average of daytime NEE (NEE<sub>DAY</sub>) as follow:

$$GPP_{prox} = (NEE_{DAY} - NEE_{NIGHT}) * k, \quad (5)$$

where k is the fraction of daytime hours for each day.

Some of the drivers are periodic, which means that extreme values have similar meanings (e.g., the DOY 0 and 365 or the WD 0 and

360). To take this into consideration, we used their sine and cosine transformations instead of the original variable's value.

The two subnetworks that estimate RECO and GPP use a different set of drivers, selected on the basis of the expected role and tests on performances and that are listed in Table 2. In particular, SNN<sub>RECO</sub> used TA, TS, SWC, WS, NEE<sub>NIGHT</sub>, and sine and cosine of WD and DOY while the drivers for the SNN<sub>GPP</sub> were SW_IN, SW_IN_POT (half hourly and daily and their first derivatives), VPD, TA, SWC, WS, GPP<sub>prox</sub>, and sine and cosine values of WD. Note that SW_IN is then used both as driver and as single input in the GPP subnetwork for the product with the LUE (see Figure 1).

## 2.4 | ANN training

Because our goal is to partition NEE into its two fluxes components, the NN<sub>C-part</sub> was trained at the site level and year-by-year; we used only high quality measured half hourly NEE (see Section 2.2) as training target. The training of the NN<sub>C-part</sub> was carried out using both nighttime and daytime measurements of NEE.

As commonly used in ANN training and application, inputs and outputs were normalized in the range +1/-1 by Equation (6):

$$X_{norm} = 2 * \left[ \left( \frac{X - X_{min}}{X_{max} - X_{min}} \right) - 0.5 \right], \quad (6)$$

where X<sub>min</sub> and X<sub>max</sub> are estimated as the maximum of the absolute value of X then negative signed in the case of X<sub>min</sub> and positive signed in the case of X<sub>max</sub>; this preserved the zero value in both normalized

<table><caption>TABLE 2 List of the variables used as drivers to estimate gross primary production (GPP) and ecosystem respiration (RECO)</caption><thead><tr><th>Variable</th><th>Time resolution</th><th>GPP</th><th>RECO</th></tr></thead><tbody><tr><td rowspan="6">Measured</td><td rowspan="6">Half hourly</td><td>SW_IN<sup>a</sup></td><td>TA</td></tr><tr><td>SW_IN_POT</td><td>TS</td></tr><tr><td>VPD</td><td>SWC</td></tr><tr><td>TA</td><td>WS</td></tr><tr><td>SWC</td><td></td></tr><tr><td>WS</td><td></td></tr><tr><td rowspan="8">Calculated</td><td rowspan="4">Half hourly</td><td>SW_IN_POT</td><td>Cos WD</td></tr><tr><td>First derivative</td><td></td></tr><tr><td>Cos WD</td><td>Sin WD</td></tr><tr><td>Sin WD</td><td></td></tr><tr><td rowspan="4">Daily</td><td>Daily average</td><td>Sin(DOY)</td></tr><tr><td>SW_IN_POT</td><td></td></tr><tr><td>Daily average</td><td>Cos(DOY)</td></tr><tr><td>SW_IN_POT first derivative</td><td></td></tr><tr><td></td><td>GPP prox</td><td>Net ecosystem exchange night</td></tr></tbody></table>

<footnote><sup>a</sup>SW_IN is then used both as driver and as single input in the GPP subnetwork for the product with the LUE.</footnote>
<figure><img src="image_3.png" /></figure>

and original values without affecting the linearity of relationships between the original and the transformed time series.

Available records with all the inputs and the measured output were randomly split 25 times in sets for training (60% of the observations and, in general, at least 2,000 examples), test (20%) to avoid overfitting, and validation (20%) to assess the performances. This results in 25 datasets used for the ANN training. The ANN training used the mean square error as the cost function and the Levenberg-Marquardt backpropagation algorithm.

Five ANNs with different numbers of neurons in the first hidden layer of both subnetworks were trained: 10–13, 11–14, 12–15, 13–16, and 14–17 for SNN<sub>RECO</sub> and SNN<sub>GPP</sub>, respectively. For each one of the 25 training datasets, the five networks with different structures were trained five times, changing the initial weights randomly, and then selecting the one with the best performance on the validation set. In case the sign of RECO and GPP did not respect the expected convention, the weights were reinitialized and the training was repeated. At the end of the process, there were 125 trained ANNs for each site-year (25 training datasets × 5 structures). To select the best among them, which were then used for the application, the ANN with the highest performances on the validation set was selected for each training dataset (obtaining 25 selected ANNs), and then, the five ANNs with the best model efficiency in predicting the overall NEE time series were finally selected and applied. RECO and GPP were calculated by averaging the output of the two subnetworks SNN<sub>GPP</sub> and SNN<sub>RECO</sub> from the five selected ANNs. The capacity of the ANN to simulate NEE was also analyzed (see Data S1, section S1).

## 2.5 | Statistical analysis and evaluation

Results of this experiment were evaluated by comparing GPP and RECO from NN<sub>C-part</sub> with the NT- and DT-based partitioning methods in the FLUXNET2015 collection. As complementary information, we repeated the same analysis for the gross CO₂ fluxes derived from the synthetic NEE dataset provided by the MUSICA model (hereafter, we explicitly refer the synthetic CO₂ fluxes modeled with MUSICA as GPP<sub>MUSICA</sub>, RECO<sub>MUSICA</sub>, and NEE<sub>MUSICA</sub>). The GPP and RECO estimates from the different methods were compared at half hourly, daily, and yearly time resolutions and were further evaluated about the consistency among methods for the seasonal cycle and its daily anomalies. The seasonal cycle was estimated by averaging the fluxes over 5 days while the daily anomalies were estimated by subtracting the seasonal cycle from each daily value. Finally, in order to compare the seasonal cycles across sites-years having different flux magnitudes, the annual average daily value was subtracted from each seasonal cycle. In addition to the seasonal cycles, the mean diurnal cycles were calculated and compared. To verify that the different drivers used in the NN<sub>C-part</sub> had the role that is ecologically meaningful, functional relationships between micrometeorological drivers and retrieved gross CO₂ fluxes were analyzed. In particular, the functional relationship between SW_IN and GPP, VPD and GPP, and

TA and both GPP and RECO were considered and analyzed. This was carried out in two different ways: (a) we simulated the flux responses in the NN<sub>C-part</sub> varying artificially one driver at time, while the other drivers were kept fixed to the average condition of the season at midday; and (b) we compared the seasonal variations of gross CO₂ fluxes with respect to the meteorological forces for both fluxes calculated using the NN<sub>C-part</sub> results and those obtained using the NT and DT methods. The responses were obtained by averaging the values of 10 equal intervals step of the fluxes and drivers for 2 months.

The results of the machine learning approach were finally evaluated by looking at three specific responses of gross CO₂ fluxes that are not explicitly accounted for by the DT and NT methods but can be indirectly considered in the NN<sub>C-part</sub> method: (a) the instantaneous LUE responses to the diffuse to direct radiation ratio, (b) the hysteresis of the diel cycles of predicted RECO and GPP with respect to TA and SW_IN, and (c) the respiration pulse after a rapid increase in SWC occurring frequently in dry environments. Finally, the NN<sub>C-part</sub> algorithm implemented here allows for an analysis of the role and importance of the different drivers (see Data S1).

# 3 | RESULTS

## 3.1 | Cross consistency between NN<sub>C-part</sub> and standard partitioning methods (DT and NT)

The half hourly GPP and RECO retrieved by NN<sub>C-part</sub> were consistent with the output obtained by the NT and DT methods, with a slightly higher correlation (R² and lower RMSE) with the DT method for GPP estimation and with the NT method for RECO, in both cases higher than the correlation between DT and NT methods (Table 3). The mean squared error of partitioned gross CO₂ fluxes between NN<sub>C-part</sub> and the NT/DT methods was <1.53 µmol CO₂ m⁻² s⁻¹ on average and comparable to the average value of the estimated random uncertainty of NEE ("NEE_CUT_USTAR50_RANDOM," variable provided by FLUXNET 2015). Looking at the bias, the magnitude of both GPP and RECO retrieved by NN<sub>C-part</sub> was slightly lower than the DT and NT methods, with the predicted values by NN<sub>C-part</sub> closer to the DT method (though the spread among methods was very small). The average values of bias among methods were low on average (<27 g C m⁻² year⁻¹) and comparatively lower with respect to the NEE random uncertainty reported in the FLUXNET2015 database (<10% of the reported "NEE_CUT_USTAR50_RANDOM"). The agreement between NN<sub>C-part</sub> and the two FLUXNET standard methods was significantly lower in the case of RECO suggesting that the differences occurred in the diurnal and seasonal cycles of RECO between NN<sub>C-part</sub> and both DT and NT methods were comparatively larger than GPP.

High correlation and consistency among NN<sub>C-part</sub> and the DT and NT methods in both GPP and RECO were found for daily, seasonal, and yearly values and also for the daily average anomalies with respect to the seasonal values (Figure 2). Daily, growing season, and yearly average values from NN<sub>C-part</sub> were highly correlated with the two standard methods (R² > .83), in particular for GPP (R² > .97).
<table><caption>TABLE 3 Cross consistency among the retrieved gross primary production (GPP) and ecosystem respiration (RECO) by $NN_{C-part}$, daytime method (DT), and nighttime method (NT) methods at half hourly time step, also after removing the mean daily value of the fluxes. Statistics reported in the table are the median of the study sites values, in brackets the 25th and 75th percentile</caption><thead><tr><th rowspan="2">Gross CO<sub>2</sub> flux</th><th rowspan="2">Variability</th><th rowspan="2">Comparison</th><th colspan="3">Statistics</th></tr><tr><th>R<sup>2</sup></th><th>RMSE (µmol CO<sub>2</sub> m<sup>-2</sup> s<sup>-1</sup>)</th><th>Bias (µmol CO<sub>2</sub> m<sup>-2</sup> s<sup>-1</sup>)</th></tr></thead><tbody><tr><td rowspan="6">GPP</td><td rowspan="3">Overall</td><td>NN<sub>C-part</sub> versus DT</td><td>.96 (.93/.97)</td><td>1.22 (0.76/1.52)</td><td>-0.042 (-0.17/0.039)</td></tr><tr><td>NN<sub>C-part</sub> versus NT</td><td>.94 (.89/.95)</td><td>1.53 (0.89/2.02)</td><td>-0.068 (-0.18/0.025)</td></tr><tr><td>DT versus NT</td><td>.90 (.84/.93)</td><td>1.84 (1.20/2.43)</td><td>0.020 (-0.16/0.14)</td></tr><tr><td rowspan="3">Removing mean daily value</td><td>NN<sub>C-part</sub> versus DT</td><td>.95 (.93/.97)</td><td>0.93 (0.59/1.11)</td><td></td></tr><tr><td>NN<sub>C-part</sub> versus NT</td><td>.90 (.86/.92)</td><td>1.34 (0.82/1.86)</td><td></td></tr><tr><td>DT versus NT</td><td>.86 (.81/.89)</td><td>1.55 (1.02/2.13)</td><td></td></tr><tr><td rowspan="6">RECO</td><td rowspan="3">Overall</td><td>NN<sub>C-part</sub> versus DT</td><td>.73 (.54/.83)</td><td>0.98 (0.70/1.29)</td><td>-0.023 (-0.19/0.12)</td></tr><tr><td>NN<sub>C-part</sub> versus NT</td><td>.87 (.72/.92)</td><td>0.68 (0.47/0.83)</td><td>-0.072 (-0.16/0.022)</td></tr><tr><td>DT versus NT</td><td>.78 (.65/.86)</td><td>0.95 (0.64/1.21)</td><td>-0.0045 (-0.21/0.16)</td></tr><tr><td rowspan="3">Removing mean daily value</td><td>NN<sub>C-part</sub> versus DT</td><td>.21 (.09/.32)</td><td>0.39 (0.28/0.51)</td><td></td></tr><tr><td>NN<sub>C-part</sub> versus NT</td><td>.24 (.11/.38)</td><td>0.35 (0.26/0.47)</td><td></td></tr><tr><td>DT versus NT</td><td>.68 (.62/.76)</td><td>0.20 (0.14/0.26)</td><td></td></tr></tbody></table>

<figure><img src="image_4.png" /><figcaption>FIGURE 2 Scatter density plot showing the cross-consistency between partitioned gross primary production and ecosystem respiration by $NN_{C-part}$ (x-axis) and the ones by the daytime method and nighttime method (y-axis) aggregated at daily time step (a, e, i, m), yearly (b, f, j, n), and looking for the seasonal cycle (c, g, k, o) and daily anomalies (d, h, l, p)</figcaption></figure>

In terms of bias, results by $NN_{C-part}$ were more consistent with the DT results confirming the findings from the subdaily comparison. High consistency was also found in the case of GPP daily anomalies ($R^2 > .88$ and RMSE < 0.52 µmol CO$_2$ m$^{-2}$ s$^{-1}$), but it decreased in the case of RECO (Figure 2).

This general consistency between $NN_{C-part}$ and standard methods (DT and NT) was confirmed by comparing the gross CO$_2$ fluxes derived by training the methods on synthetic NEE ($NEE_{MUSICA}$, see section S4.1 in Data S1). The main findings from this comparison

were that: (a) $NN_{C-part}$ and standard methods accurately retrieved both $GPP_{MUSICA}$ and $RECO_{MUSICA}$ at daily, yearly, seasonal timescales, and also the daily anomalies (see Figures S5 and S6); (b) $NN_{C-part}$ slightly outperformed the standard methods to estimate $GPP_{MUSICA}$ and $RECO_{MUSICA}$, particularly in the case of daily anomalies; (c) all methods slightly underestimated both $GPP_{MUSICA}$ and $RECO_{MUSICA}$ (bias ranging between -0.2 and -0.45 µmol CO$_2$ m$^{-2}$ s$^{-1}$) with the largest discrepancy found in the case of DT method (for more details, see section S4.1, Figures S5 and S6).
<figure><img src="image_5.png" /></figure>

## 3.2 | Consistency among the seasonal and diurnal dynamics of partitioned gross CO₂ fluxes

### 3.2.1 | Seasonal cycle

The seasonal pattern of the retrieved gross CO₂ fluxes by NN<sub>C-part</sub> closely matched that of the DT and NT methods, with the only exception being that predicted RECO by NN<sub>C-part</sub> was slightly lower in late spring/early summer (Figure 3, left). The predicted seasonal patterns were also similar among the represented International Geosphere-Biosphere Programme vegetation classes (Figure 3, right). Results showed a large degree of agreement among the three methods, in particular between NN<sub>C-part</sub> and NT.

The seasonal dynamic of gross CO₂ fluxes derived from NEE<sub>MUSICA</sub> by NN<sub>C-part</sub>, NT, and DT agreed consistently. However, all methods underestimated both GPP<sub>MUSICA</sub> and RECO<sub>MUSICA</sub> (in particular the DT method) during the growing season (between April and July, see section S4.2 and Figure S7 in Data S1).

### 3.2.2 | Mean diurnal cycle of GPP

We found a good agreement between the pattern of the mean diurnal cycle of GPP predicted by NN<sub>C-part</sub> and FLUXNET standard methods (Figure 4). The GPP diurnal dynamic followed an asymmetrically bell-shaped curve with the maximum value reached

before midday (Figure 4, top panels). Analyzing the differences more in depth, the mean diurnal cycle of GPP predicted by NN<sub>C-part</sub> systematically diverged from the one predicted by the DT standard method, following a pattern that was consistent across the seasons (Figure 4, bottom panels). In particular, the GPP estimated by NN<sub>C-part</sub> was lower during the first hours of the morning (just after the dawn) and in the late afternoon; instead, the GPP estimated by NN<sub>C-part</sub> was slightly higher than DT from the morning until midday, where NN<sub>C-part</sub> estimated higher GPP values. A similar pattern appeared when comparing the two standard methods (difference calculated as NT-DT) although shifted toward positive values, indicating higher values of GPP predicted by the NT method compared to the NN<sub>C-part</sub>. The GPP derived from NEE<sub>MUSICA</sub> confirmed the differences we highlighted between DT and the other methods (see Section 4.3; Figure S8).

### 3.2.3 | Mean diurnal cycle of RECO

In comparison to GPP, the diurnal cycle of RECO modeled by NN<sub>C-part</sub> showed larger differences with respect to the DT and NT methods, in particular regarding the magnitude of the flux. All the methods showed an increased respiration when moving toward midday, in particular during the growing season (Figure 5, top plots); nevertheless, the NN<sub>C-part</sub> estimated values were lower than predictions from the other methods. Conversely, the differences between the two

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

**FIGURE 3** The mean seasonal cycles of gross primary production (GPP; a) and ecosystem respiration (RECO; c) obtained by NN<sub>C-part</sub>, daytime method, and nighttime method. The consistency among methods was evaluated by the determination coefficients (R²) among the seasonal patterns of GPP (b) and RECO (d) per International Geosphere-Biosphere Programme (IGBP) vegetation class (the number in brackets after each IGBP vegetation class stands for the number of sites in each category). Only FLUXNET2015 study sites at northern hemisphere (latitude > +15°N) were used for that analysis. The vegetation classes are the same as Table 1
<figure><img src="image_10.png" /></figure>

FIGURE 4 The trends of the mean diurnal cycle of gross primary production (GPP; a–f) retrieved by $NN_{C-part}$ (black line), daytime method (DT; red line), and nighttime method (NT; green line). The differences between diurnal cycle of retrieved GPP are also shown (g–l): $NN_{C-part}$-DT (red line), $NN_{C-part}$-NT (green line), and NT-DT (dashed blue line). Only FLUXNET2015 study sites at northern hemisphere (latitude > +15°N), having at least 2 years of data, were used for the analysis

<figure><img src="image_11.png" /></figure>

FIGURE 5 The trends of the mean diurnal cycle of ecosystem respiration (RECO; a–f) retrieved by $NN_{C-part}$ (black line), daytime method (DT; gray line), and nighttime method (NT; dashed gray line). The differences between diurnal cycle of retrieved RECO (g–l) are also shown: $NN_{C-part}$-DT (continuous black line), $NN_{C-part}$-NT (dashed black line), and NT-DT (gray line). Only FLUXNET2015 study sites at northern hemisphere (latitude > +15°N), having at least 2 years of data, were used for that analysis

standard methods were negligible with a slightly higher RECO by NT in comparison to DT (Figure 5, lower plots, calculated as NT-DT).

This finding was not fully mirrored in the patterns of RECO<sub>MUSICA</sub> (see section S4.3 and Figure S9). In that case, $NN_{C-part}$ predicted higher

RECO in comparison to both DT and NT, but more consistent with RECO<sub>MUSICA</sub>. Conversely, $NN_{C-part}$ predicted lower RECO than the NT method (and more consistent with DT) during the morning; the differences between DT and NT fitted on NEE<sub>MUSICA</sub> were also enhanced.
<figure><img src="image_12.png" /></figure>

These differences increased during the nighttime hours, while the predicted RECO by the two standard methods became closer during the daytime hours (particularly when air temperature reached its maximum value). It is interesting to note that the peak of the maximum RECO by NN<sub>C-part</sub> was shifted compared to the DT and NT methods and closer to the one of RECO<sub>MUSICA</sub>, suggesting a better capacity to reproduce diurnal patterns when multiple drivers were involved (Figure S9).

## 3.3 | Functional relationships between partitioned fluxes and meteorological drivers

The functional response obtained by artificially varying the drivers in NN<sub>C-part</sub> was consistent with the current knowledge. The response of GPP to SW_IN was as expected, with an increased photosynthesis due to the increased light that resulted as saturated at high radiation values (Figure 6a-d). The response to VPD clearly showed the effect of VPD limitation on photosynthesis (Körner, 1995) that was evident also at relatively low values of VPD (Figure 6e-h).

The response of photosynthesis to air temperature showed a strong and positive response at the start of the growing season that became flat in the summer period (Figure 6i-l). The functional relationships between RECO and TA followed the expected increase of RECO with TA (Figure 6m-p) except for the dry ecosystems (red points in Figure 6) where TA was not the main driving factor. These general patterns of relationships were consistent across the seasons and in the different vegetation types, despite the expected variability due to the ecosystem-specific properties.

When applied using the original measurements, the functional ecosystem response to the micrometeorological forces was largely consistent with the ones obtained using the DT and NT standard methods, although affected by multiple co-acting and confounding factors (Figure 7). The shape of the relations was slightly different in the four periods analyzed, but it was important to recall that these were average responses from sites in different climatic conditions. It could be noted that an almost perfect match was present in the case of GPP at this broad scale; in the case of RECO, the NN<sub>C-part</sub> showed in general lower respiration fluxes at high temperatures (so mainly daytime) as opposed to the other methods, although the pattern of response was basically the same. This could be an indication of a water resources limitation during the summer period and at higher temperatures, that could be detected by the NN<sub>C-part</sub> since it also used SWC as a driver. In fact, the functional relationships between RECO and SWC retrieved by NN<sub>C-part</sub> highlighted a clear direct relationship between SWC and RECO in drought conditions, in particular for non-forested ecosystems (Figure S16).

## 3.4 | Additional ecological patterns simulated by NN<sub>C-part</sub>

### 3.4.1 | NN<sub>C-part</sub> instantaneous LUE and the ratio between diffuse to direct radiation

The instantaneous LUE estimated by NN<sub>C-part</sub> showed an increasing trend with respect to the increased diffuse to direct radiation ratio (here estimated by the proxy 1-SW_IN/SW_IN_POT; Figure 8a-d).

<figure><img src="image_13.png" /><figcaption>**FIGURE 6** NN<sub>C-part</sub> predicted responses of gross primary production as function of SW_IN ([W/m²], a-d), VPD, ([kPa], e-h), Ta ([°C], i-l), and of ecosystem respiration as function of Ta ([°C], m-p). The functional responses are simulated in the study sites at the northern hemisphere (latitude > +15°N) and in fixed seasonal conditions. The average site-specific responses (mean values of at least 3 years of data) are reported as colored lines, in particular: forested sites are in green, grassland and cropland sites are in black, and dry sites are reported in red. The singular site-year response are dashed light gray lines</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>FIGURE 7 Univariate functional ecosystem relationships between gross CO₂ fluxes retrieved by NN<sub>C-part</sub> (black lines) and FLUXNET standard methods (green and red lines for nighttime method and daytime method, respectively): gross primary production (GPP; µmol CO₂ m⁻² s⁻¹) and SW_IN (W/m²; a–d), GPP (µmol CO₂ m⁻² s⁻¹) and VPD (kPa; e–h), ecosystem respiration (µmol CO₂ m⁻² s⁻¹) and TA (°C; i–l). These patterns are derived from study sites located at the northern hemisphere (latitude > +15°N)</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>FIGURE 8 Sensitivity of the light use efficiency (LUE) of NN<sub>C-part</sub> to the diffuse to direct radiation ratio. The trend of LUE by NN<sub>C-part</sub> with respect to the diffuse radiation ratio is reported in the panels (a–d; here LUE is normalized for comparison purpose). The difference between LUE (µmol CO₂/W) as estimated by NN<sub>C-part</sub> and that from the daytime method (DT) standard partitioning method is also reported (e–h); the statistics were aggregated per diffuse to direct radiation class (here estimated by the proxy 1-SW_IN/SW_POT and reported as percentage). We reported the following statistics of x (x = LUE(NN<sub>C-part</sub>) – LUE(DT)): the median values (bar), the mean value (*), the interquartile range (by brackets). Only LUE at midday is used in this figure</figcaption></figure>
<figure><img src="image_16.png" /></figure>

When compared against the LUE estimated by the DT method, $NN_{C-part}$ showed a moderate but consistent higher sensitivity to the fraction of diffuse to direct light with a higher LUE when the diffuse to direct radiation ratio was higher than 40% (Figure 8e-h).

### 3.4.2 | Hysteresis of the diel cycle of gross CO₂ fluxes

The diurnal cycle of gross CO₂ fluxes exhibited hysteresis, which we analyzed through the diurnal cycle of GPP with respect to

$SW_{IN}$ (Figure 9a-c) and the diurnal cycle of RECO with respect to TA (Figure 9e-g). We generally observed a clockwise cycle for GPP, with higher values of GPP in the morning in comparison to the afternoon. All the methods showed this pattern although it was more evident in the $NN_{C-part}$ and NT methods. We provided numerical estimates by calculating the integral area included in the hysteresis pattern (larger area indicates larger hysteresis) that was reported in Figure 9d.

For RECO, only $NN_{C-part}$ detected an appreciable hysteresis in the diurnal cycle (Figure 9h); this is expected because both the DT and NT methods use the same TA-dependent and invariant

<figure><img src="image_17.png" /></figure>

**FIGURE 9** Hysteresis of the gross primary production (GPP) diel cycle predicted by $NN_{C-part}$ with respect to $SW_{IN}$ in three sampled study sites (a-c). Three different International Geosphere-Biosphere Programme classes are reported: DBF (a), ENF (b), and GRA (c), while the period of interest was July-August. (d) Distribution of the GPP hysteresis integral area (then normalized) estimated by $NN_{C-part}$, nighttime method (NT), and daytime method (DT) partitioning methods. The hysteresis of the diel cycle of ecosystem respiration (RECO; here log-transformed) with respect to the TA dynamic, in three sampled study sites, is reported (e-g). The distribution of the RECO hysteresis integral area (then normalized) estimated by $NN_{C-part}$, NT, and DT partitioning methods is reported (h)

**FIGURE 10** The sensitivity of half-hourly ecosystem respiration (RECO) predicted by $NN_{C-part}$ to the SWC "pulse". The investigated study site was US-SRG. Fluxes are reported in panel (a) and more specifically we reported: RECO predicted by $NN_{C-part}$ (blue line), nighttime method (green) and daytime method (red) and the half-hourly net ecosystem exchange (gray). Meteorological variables are reported in panel (b); more specifically we reported the SWC dynamics (blue line) and the rain events (P [mm]; red line)

<figure><img src="image_18.png" /></figure>
relationships for the whole diel cycle. In general, the hysteresis of RECO is realized as a counterclockwise cycle, with the value of RECO, for the same TA, higher during the afternoon and nighttime hours and lower during the morning. However, the counterclockwise cycle observed in the RECO patterns was not always observed, and in some cases, the opposite hysteresis cycle was found.

### 3.4.3 | RECO pulses due to SWC variations

The evaluation of the capacity to correctly reproduce rapid responses of RECO (respiration pulse) to changes in SWC (e.g., a post rain event after a prolonged drought period in arid environments) is in general difficult due to the sporadic nature of the events and the noise associated with the measurements. We evaluated the performances of the three models in a specific event that was clearly visible in the NEE time series (US-SRG site, Figure 10). The analysis of the predicted RECO shows that $NN_{C-part}$ was much more sensitive to the fast response of respiration pulses. The NT and DT methods, that did not have SWC as a direct input, showed a slow changing pattern with the peak of respiration significantly delayed compared to the measurements.

## 4 | DISCUSSION

### 4.1 | Consistency of $NN_{C-part}$ functional response with theory on plant physiology and with the DT and NT standard methods

The retrieved GPP and RECO as obtained from $NN_{C-part}$ were consistent with the estimates by the DT and NT methods both when trained directly on EC measurements and using synthetic data with NEE$_{MUSICA}$ (Figure 2; Figures S3 and S4). This was valid also for the relationship between micrometeorological drivers and gross CO$_2$ fluxes calculated by $NN_{C-part}$ and those implemented in the DT and NT partitioning methods (Figures 6 and 7).

The GPP response to light used by $NN_{C-part}$ is curvilinear and consistent with the pattern of that used in DT method (despite some small but systematic divergences). The pattern of the GPP-VPD relationship captured by $NN_{C-part}$ is consistent with current knowledge about the protective (physiological) mechanisms of stomata closure carried out by plants in response to the increase in atmospheric evaporative demands (Körner, 1995). There is also a general consistency of the GPP-VPD dynamics retrieved by $NN_{C-part}$ with those from the other methods (Figure 7) although the relationships include also the effect of other covarying factors, in particular incoming radiation. This stomatal regulation of photosynthesis is explicitly inserted in the formulation only in the DT approach. Interestingly, the $NN_{C-part}$ shows an effect of VPD also at very low values (Figures 6e-h and 7e-h), lower than the 10 hPa used as threshold in the DT method. The patterns of

the TA-RECO relationships are also consistent across methods and in agreement with the expected trend due to the biochemical reactions involved in the organic matter respiration processes that find their optimal conditions at high (but not extreme) temperatures, if not limited by water availability.

The largest mismatch we found between $NN_{C-part}$ and the two standard methods of partitioning (NT and DT) is in late spring/early summer, with lower fluxes of daytime RECO predicted by $NN_{C-part}$ compared to both standard methods. When evaluated on the synthetic dataset, the mismatch between $NN_{C-part}$ and the standard methods was significantly lower except for the anomalies, which were better reproduced by $NN_{C-part}$. In addition, the estimates from $NN_{C-part}$ better matched the diurnal dynamics of RECO$_{MUSICA}$ than the estimates from the standard methods (Figure S9); conversely the method that showed the highest mismatch with respect to RECO$_{MUSICA}$ was the DT method. All the methods underestimated both the reference GPP$_{MUSICA}$ and RECO$_{MUSICA}$ (despite the differences were very tiny) without any relevant effect on NEE. This could be related to a certain limitation of these methods to infer the complex relationships of the MUSICA model with the reduced set of drivers used in this experiment, but we have also to consider that both $NN_{C-part}$ and standard methods of partitioning were trained on noisy NEE signals (see section S4 in Data S1 for details) while both GPP$_{MUSICA}$ and RECO$_{MUSICA}$ used for reference are noise-free.

Recent papers highlighted possible limitations of the standard methods for partitioning stemming from the fact that they do not consider the inhibition of daytime leaf respiration (the Kok effect; Keenan et al., 2019; Wehr et al., 2016) that produces an overestimation of daytime RECO and GPP in the NT method (Wehr et al., 2016) and an underestimation of nighttime respiration (see Keenan et al., 2019) in DT. Although some of the differences between the RECO estimated by $NN_{C-part}$ and the other methods (see Figures 3 and 5 and Figure S9) are consistent with experimental data focused on this topic (e.g., Keenan et al., 2019; Wehr et al., 2016), it is not possible to demonstrate that the $NN_{C-part}$ method, as implemented in this experiment, is able to reproduce the light inhibition of leaf respiration. Moreover, the consistency between the DT and NT partitioned fluxes used as reference in this experiment, and the pattern emerged by other modeling experiences (see e.g., Jung et al., 2020, where DT and NT partitioned fluxes were globally upscaled) suggests that the Kok effect could have a minor relevance on biases of RECO and GPP estimates compared to other sources of NEE measurement uncertainties.

In terms of ecological responses, the $NN_{C-part}$ outputs showed a clear and interesting response to the diffuse/direct light with an increased LUE in diffuse conditions and with hysteresis of the diel cycle of GPP and RECO that are only partially shown by the standard methods. The latter showed also strong limitation in reproducing the respiration pulse while $NN_{C-part}$ showed this capacity, although tested only in one case. In summary, these findings highlight the strength of this approach that is able to reproduce these patterns even without being trained specifically for this.
<figure><img src="image_19.png" /></figure>

## 4.2 | The advantages of the proposed NN<sub>C-part</sub> approach

One relevant difference between NN<sub>C-part</sub> and the standard partitioning methods is the absence of prescribed relationships between drivers and fluxes. In fact, only a few constraints have been set in the NN<sub>C-part</sub> structure (RECO and GPP signals must be positive, weights in the last node have to be consistent with the sign convention of NEE) and the relationship between inputs and outputs is set on the basis of the data without any assumption (as common in all the machine learning methods). This property could be a great advantage when used to study a phenomenon that does not respond with the same invariant pattern. The synthetic datasets used in this experiment are generated using functional relationships which are more complex than the ones implemented in the FLUXNET standard methods. The better matching by NN<sub>C-part</sub> highlights the potential of this algorithm to capture additional patterns of carbon fluxes dynamics. As a direct consequence, this led to some systematic differences observed in the GPP and RECO patterns compared to standard methods. For instance, the shape of the GPP diurnal cycle predicted by the DT method is almost symmetrically bell shaped: GPP reached a maximum around midday and only VPD, used as a downregulation function, when higher than 10 hPa, can modify this pattern that is otherwise dictated by incoming radiation. This is also visible in the modeled NEE, by comparing NN<sub>C-part</sub> and DT outputs with measured data (Figure S10). Conversely, the pattern of GPP obtained by the NT method is less affected by prescribed relationships because GPP is calculated as the difference between RECO and the measured NEE. The diurnal dynamic of GPP estimated by NN<sub>C-part</sub> is more similar to the NT estimates, which suggests that the formulation of the DT approach could be improved, for example, by using the non-rectangular light response function (instead of the rectangular used in DT, see Gilmanov et al., 2003) or a different threshold (or a different function) for the VPD downregulation effect. Although these differences did not lead to large mismatches in terms of magnitude and/or seasonal dynamics of GPP, they can have an effect on the estimates of ecosystem functional properties such as LUE, WUE, or other important physiological parameters (Keenan et al., 2019; Reichstein, Bahn, Mahecha, Kattge, & Baldocchi, 2014).

Machine learning algorithms, such as the one implemented in this experiment, are also more effective at extracting relevant features from the comprehensive set of drivers used as input. For example, NN<sub>C-part</sub> indirectly derived information on diffuse radiation from the measured and potential incoming radiation used as drivers to scale the LUE, as consistently found in other studies (Alton, North, & Los, 2007). In addition, NN<sub>C-part</sub> directly uses additional input to improve the estimates of GPP. For example: (a) TA shows a role in spring and then a less strong effect possibly due to acclimation (Figure 6i-l); (b) SWC could play a role in downregulating the GPP response to VPD, with the effect of the atmospheric evaporative demand being stronger in the case of reduced soil water availability from soil and effective water stress condition.

In the case of RECO, both DT and NT used the same modeling approach, based on the Lloyd and Taylor model driven by the air temperature, resulting in high consistency between their results. However, it is known that soil temperature is also an important driver of RECO in particular to define its temporal (diurnal) dynamic (Lasslop et al., 2012; Wohlfahrt & Galvagno, 2017) because it is directly linked to soil respiration, which is a large contribution to the total RECO (Misson et al., 2007). The NN<sub>C-part</sub> method uses soil temperature (in addition to air temperature) as a direct driver of RECO and the analysis of the functional responses highlighted an important and direct effect (Figure S3). By looking at the analysis of RECO driver importance (see section S8 and Figure S20 in Data S1), it seems that TS and TA had, in general, a similar importance, but their relative weights changed by site. For this reason, it could be difficult to use both TA and TS in the same Lloyd and Taylor relationship without a priori information about the relative importance of the two variables. The multiple drivers approach implemented in NN<sub>C-part</sub> where both TA and TS are used as drivers of RECO without any a priori assumptions also allows for the detection of hysteresis in the diurnal cycle of RECO (see Wohlfahrt & Galvagno, 2017 for details on the topic).

Another important driver of RECO that is not explicitly accounted for in the standard methods was SWC that is particularly important in dry sites (see section S8). The DT and NT methods use a dynamic parameter ($R_{ref}$) estimated using moving windows to indirectly consider the slow changes in water availability (Reichstein et al., 2005). Soil water content, however, also has instantaneous effects in some ecosystems (e.g., the respiration pulses after short rain events, see Jarvis et al., 2007) that need the direct use of SWC as input, like in the NN<sub>C-part</sub> method, to be effectively modeled.

Other slow changing factors affecting the fluxes, such as physiological state, substrate availability, and management or other disturbance events, that are considered indirectly through $R_{ref}$ in the standard methods, are represented in the NN<sub>C-part</sub> approach through the use of averaged NEE-derived quantities as drivers.

Finally, the standard methods do not account for EC footprint variability, which can significantly affect the magnitude of the measured fluxes if the fetch is small and not fully homogeneous (see section S8). This could be considered, for example, by fitting the models per wind sectors and wind speed classes (or using footprint models), but the amount of available data could become critically low in less represented conditions. In the NN<sub>C-part</sub>, we used wind variables (direction, speed) as input in order to indirectly take into consideration footprint variability. Wind-related variables appeared to be important in a few study sites (see Figures S15, S16, and S17). By comparing the ANN trained by including/excluding wind variables, we detected an effect on performance (increase of RMSE), which highlights an important effect on the instantaneous estimated flux values. However, the role of wind variables is reduced in the yearly budget calculation of gross CO<sub>2</sub> fluxes (~6 g C m<sup>-2</sup> year<sup>-1</sup>, average value across study site) even though differences were detected site by site (site-specific differences ranged between -169 and +122 g C m<sup>-2</sup> year<sup>-1</sup>). In summary, the
implemented machine learning approach gives the possibility to use variables as drivers which would be difficult to consider in process-based approaches, which require that the role and effects of each single variable being known and correctly encoded (at the ecosystem scale).

## 4.3 | Uncertainty and limitation of the proposed method

Despite the encouraging results, there are also limitations and uncertainty sources to consider when a pure machine learning approach is used. The absence of prescribed relationships and the few constraints fixed made the method robust against incorrect or incomplete formulations used in the standard methods but also led to higher sensitivity to the uncertainty in the data used in the training phase of machine learning, particularly when both NEE and the meteorological variables used as drivers are affected by long gaps. In the NN<sub>C-part</sub>, the weights of the ANN are optimized against NEE; thus, the uncertainty of NEE could have a significant impact on their estimation. This could be particularly relevant for nighttime NEE measurements that can be affected by higher uncertainty (e.g., advection, large footprint) that, even if filtered, lead to long gaps. In these situations, a model based on ecological responses is in general more robust, and needs less data for a proper calibration than a data-driven approach.

Another limitation is related to the availability of the large set of input data that this approach requires in order to fully exploit its advantages, which are often not fully available at all FLUXNET sites. In this situation, the model can be used with a reduced set of drivers, but this could also lead to higher uncertainty and lower performance in specific sites, while models based on ecological processes like the NT and DT may be more robust because the response pattern is already prescribed.

## 5 | CONCLUSIONS

In this study, we propose a machine learning method for NEE partitioning (NN<sub>C-part</sub>) using an ANN approach with a tailored structure to simultaneously retrieve GPP and RECO. From a methodological point of view, this approach is an example of simple hybrid modeling (Reichstein et al., 2019), combining a neural net with a "theoretical" ecosystem equation based on LUE.

The NN<sub>C-part</sub> is designed to use as input a comprehensive set of drivers selected among the most commonly measured variables at EC sites. This property of NN<sub>C-part</sub> (and of machine learning in general) can be exploited in future studies, for example, using as input other variables such as diffuse radiation, COS, ¹³C isotopes, or SIF. Additionally, NN<sub>C-part</sub> does not make use of prescribed relationships between inputs and outputs and only a few constraints are introduced. Results in terms of the consistency of the gross CO₂ estimated fluxes by NN<sub>C-part</sub> with the estimates from FLUXNET standard

methods were good and the proposed model is able to reproduce magnitudes, patterns, and functional responses. However, the systematic differences emerged by the cross-consistency analysis and the ecological patterns captured by NN<sub>C-part</sub> suggest that the implemented ANN reproduces fluxes dynamics more realistic than standard methods of partitioning. Being that said, further investigation are required in order to clarify the impact of the additional drivers used in NN<sub>C-part</sub>, or possible missing functions and relationships in the standard methods.

Altogether, this new method provides GPP and RECO estimates that are based on purely data-driven empirical relationships without any assumption of the driver-output relations. This product would be of high interest for process model parameterization and validation as it avoids unwanted circularity, and for this reason, we propose it as a complement to the standard processing in large networks such ICOS, AmeriFlux, or FLUXNET.

## ACKNOWLEDGEMENTS

G.T. and D.P. conceived the experiment, G.T. implemented the method and did the analysis, J.O. provided the simulated data from MUSICA model, all the co-authors contributed in the preparation of the study. The authors thank Carlo Trotta for its help in the processing of the MUSICA data using the FLUXNET standard partitioning methods. This work used EC data acquired and shared by the FLUXNET community, including these networks: AmeriFlux, AfriFlux, AsiaFlux, CarboAfrica, CarboEuropeIP, Carboltaly, CarboMont, ChinaFlux, Fluxnet-Canada, GreenGrass, ICOS, KoFlux, LBA, NECC, OzFlux-TERN, TCOS-Siberia, and USCCC. The ERA-Interim reanalysis data are provided by ECMWF and processed by LSCE. The FLUXNET EC data processing and harmonization was carried out by the European Fluxes Database Cluster, AmeriFlux Management Project, and Fluxdata project of FLUXNET, with the support of CDIAC and ICOS Ecosystem Thematic Center, and the OzFlux, ChinaFlux, and AsiaFlux offices. G.T. was supported by the projects H2020 BACI (grant agreement 640176) and by the European Research Council (ERC) under the ERC-2017-STG SENTIFLEX project (grant agreement 755617). J.V. would like to acknowledge the support from the European Research Council (ERC) under the ERC-2017-STG SENTIFLEX project (grant agreement 755617). G.C.V. would like to acknowledge the support from the European Research Council (ERC) under the ERC Consolidator Grant 2014 project SEDAL (grant agreement 647423). T.F.K. was supported by the NASA Terrestrial Ecology Program IDS Award NNH17AE86I. M. J. thanks the support of the e-shape project funded by the EU's Horizon 2020 Programme (GA 820852). D.P. thanks for the support the RINGO project funded by the EU's Horizon 2020 Programme (GA 703944) and the e-shape project funded by the EU's Horizon 2020 Programme (GA 820852).

## DATA AVAILABILITY STATEMENT

All the data used in the paper are included in the FLUXNET2015 collection available at https://fluxnet.fluxdata.org/ and shared under the CC4BY license. The codes for preprocessing and preparing data for partitioning, train the algorithm, select the best ANN, and produce the outputs are also shared under the BSD-Clause3 open license in GitHub at the address https://github.com/icos-etc/Partitioning_ANN.
ORCID

Gianluca Tramontana ID https://orcid.org/0000-0003-4014-7958

Jerome Ogee ID https://orcid.org/0000-0002-3365-8584

REFERENCES

Alton, P. B., North, P. R., & Los, S. O. (2007). The impact of diffuse sunlight on canopy light-use efficiency, gross photosynthetic product and net ecosystem exchange in three forest biomes. *Global Change Biology*, 13, 776-787. https://doi.org/10.1111/j.1365-2486.2007.01316.x

Anthoni, P. M., Knohl, A., Rebmann, C., Freibauer, A., Mund, M., Ziegler, W., ... Schulze, E.-D. (2004). Forest and agricultural land-use-dependent CO₂ exchange in Thuringia, Germany. *Global Change Biology*, 10(12), 2005-2019. https://doi.org/10.1111/j.1365-2486.2004.00863.x

Archibald, S. A., Kirton, A., van der Merwe, M. R., Scholes, R. J., Williams, C. A., & Hanan, N. (2009). Drivers of inter-annual variability in Net Ecosystem Exchange in a semi-arid savanna ecosystem, South Africa. *Biogeosciences*, 6(2), 251-266. https://doi.org/10.5194/bg-6-251-2009

Arnone III, J. A., Verburg, P. S. J., Johnson, D. W., Larsen, J. D., Jasoni, R. L., Lucchesi, A. J., ... Schimel, D. S. (2008). Prolonged suppression of ecosystem carbon dioxide uptake after an anomalously warm year. *Nature*, 455(7211), 383-386. https://doi.org/10.1038/nature07296

Arriga, N., Rannik, Ü., Aubinet, M., Carrara, A., Vesala, T., & Papale, D. (2017). Experimental validation of footprint models for eddy covariance CO₂ flux measurements above grassland by means of natural and artificial tracers. *Agricultural and Forest Meteorology*, 242, 75-84. https://doi.org/10.1016/j.agrformet.2017.04.006

Aubinet, M., Chermanne, B., Vandenhaute, M., Longdoz, B., Yernaux, M., & Laitat, E. (2001). Long term carbon dioxide exchange above a mixed forest in the Belgian Ardennes. *Agricultural and Forest Meteorology*, 108(4), 293-315. https://doi.org/10.1016/S0168-1923(01)00244-1

Baldocchi, D. (2003). Assessing the eddy covariance technique for evaluating carbon dioxide exchange rates of ecosystems: Past, present and future. *Global Change Biology*, 9, 479-492. https://doi.org/10.1046/j.1365-2486.2003.00629.x

Baldocchi, D. (2014). Measuring fluxes of trace gases and energy between ecosystems and the atmosphere - The state and future of the eddy covariance method. *Global Change Biology*, 20(12), 3600-3609. https://doi.org/10.1111/gcb.12649

Bastos, L. S., & O'Hagan, A. (2009). Diagnostics for Gaussian process emulators. *Technometrics*, 51(4), 425-438. https://doi.org/10.1198/TECH.2009.08019

Berbigier, P., Bonnefond, J.-M., & Mellmann, P. (2001). CO₂ and water vapour fluxes for 2 years above Euroflux forest site. *Agricultural and Forest Meteorology*, 108(3), 183-197. https://doi.org/10.1016/S0168-1923(01)00240-4

Beringer, J., Hutley, L. B., Hacker, J. M., Neininger, B., & Paw U, K. T. (2011). Patterns and processes of carbon, water and energy cycles across northern Australian landscapes: From point to region. *Agricultural and Forest Meteorology*, 151(11), 1409-1416. https://doi.org/10.1016/j.agrformet.2011.05.003

Beringer, J., Hutley, L. B., Tapper, N. J., & Cernusak, L. A. (2007). Savanna fires and their impact on net ecosystem productivity in North Australia. *Global Change Biology*, 13(5), 990-1004. https://doi.org/10.1111/j.1365-2486.2007.01334.x

Bodesheim, P., Jung, M., Gans, F., Mahecha, M. D., & Reichstein, M. (2018). Upscaled diurnal cycles of land-atmosphere fluxes: A new global half-hourly data product. *Earth System Science Data*, 10, 1327-1365. https://doi.org/10.5194/essd-10-1327-2018

Bonal, D., Bosc, A., Ponton, S., Goret, J.-Y., Burban, B., Gross, P., ... Granier, A. (2008). Impact of severe dry season on net ecosystem exchange in the Neotropical rainforest of French Guiana. *Global Change Biology*, 14(8), 1917-1933. https://doi.org/10.1111/j.1365-2486.2008.01610.x

Camps-Valls, G., Sejdinovic, D., Runge, J., & Reichstein, M. (2019). A perspective on Gaussian processes for Earth observation. *National Science Review*, 6(4), 616-618. https://doi.org/10.1093/nsr/nwz028

Cernusak, L. A., Hutley, L. B., Beringer, J., Holtum, J. A. M., & Turner, B. L. (2011). Photosynthetic physiology of eucalypts along a sub-continental rainfall gradient in northern Australia. *Agricultural and Forest Meteorology*, 151(11), 1462-1470. https://doi.org/10.1016/j.agrformet.2011.01.006

Chiesi, M., Maselli, F., Bindi, M., Fibbi, L., Cherubini, P., Arlotta, E., ... Seufert, G. (2005). Modelling carbon budget of Mediterranean forests using ground and remote sensing measurements. *Agricultural and Forest Meteorology*, 135(1-4), 22-34. https://doi.org/10.1016/j.agrformet.2005.09.011

Chu, H., Baldocchi, D. D., John, R., Wolf, S., & Reichstein, M. (2017). Fluxes all of the time? A primer on the temporal representativeness of FLUXNET. *Journal of Geophysical Research Biogeosciences*, 122, 289-307. https://doi.org/10.1002/2016JG003576

Commane, R., Meredith, L. K., Baker, I. T., Berry, J. A., Munger, J. W., Montzka, S. A., ... Wofsy, S. C. (2015). Seasonal fluxes of carbonyl sulfide in a mid-latitude forest. *Proceedings of the National Academy of Science of the United States of America*, 112(46), 14162-14167. https://doi.org/10.1073/pnas.1504131112

Desai, A. R., Richardson, A. D., Moffat, A. M., Kattge, J., Hollinger, D. Y., Barr, A., ... Stauch, V. J. (2008). Cross-site evaluation of eddy covariance GPP and RE decomposition techniques. *Agricultural and Forest Meteorology*, 148, 821-838. https://doi.org/10.1016/j.agrformet.2007.11.012

Dragoni, D., Schmid, H. P., Wayson, C. A., Potter, H., Grimmond, C. S. B., & Randolph, J. C. (2011). Evidence of increased net ecosystem productivity associated with a longer vegetated season in a deciduous forest in south-central Indiana, USA. *Global Change Biology*, 17, 886-897. https://doi.org/10.1111/j.1365-2486.2010.02281.x

Dutta, D., Schimel, D. S., Sun, Y., van der Tol, C., & Frankenberg, C. (2019). Optimal inverse estimation of ecosystem parameters from observations of carbon and energy fluxes. *Biogeosciences*, 16, 77-103. https://doi.org/10.5194/bg-16-77-2019

Falge, E., Baldocchi, D., Olson, R., Anthoni, P., Aubinet, M., Bernhofer, C., ... Wofsy, S. (2001). Gap filling strategies for defensible annual sums of net ecosystem exchange. *Agricultural and Forest Meteorology*, 107(1), 43-69. https://doi.org/10.1016/S0168-1923(00)00235-5

Fares, S., Savi, F., Muller, J., Matteucci, G., & Paoletti, E. (2014). Simultaneous measurements of above and below canopy ozone fluxes help partitioning ozone deposition between its various sinks in a Mediterranean Oak Forest. *Agricultural and Forest Meteorology*, 198-199, 181-191. https://doi.org/10.1016/j.agrformet.2014.08.014

Fischer, M. L., Billesbach, D. P., Riley, W. J., Berry, J. A., & Torn, M. S. (2007). Spatiotemporal variations in growing season exchanges of CO₂, H₂O, and sensible heat in agricultural fields of the southern great plains. *Earth Interactions*, 11, 1-21. https://doi.org/10.1175/EI231.1

Frank, J. M., Massman, W. J., Ewers, B. E., Huckaby, L. S., & Negrón, J. F. (2014). Ecosystem CO₂/H₂O fluxes are explained by hydraulically limited gas exchange during tree mortality from spruce bark beetles. *Journal of Geophysical Research. Biogeoscience*, 119(6), 2013JG002597. https://doi.org/10.1002/2013JG002597

Garbulsky, M. F., Penuelas, J., Papale, D., & Filella, I. (2008). Remote estimation of carbon dioxide uptake by a Mediterranean forest. *Global Change Biology*, 14, 2860-2867. https://doi.org/10.1111/j.1365-2486.2008.01684.x

Giasson, M.-A., Coursolle, C., & Margolis, H. A. (2006). Ecosystem-level CO₂ fluxes from a boreal cutover in eastern Canada before and after scarification. *Agricultural and Forest Meteorology*, 1, 23-40. https://doi.org/10.1016/j.agrformet.2006.08.001

Gilmanov, T. G., Verma, S. B., Sims, P. L., Meyers, T. P., Bradford, J. A., Burba, G. G., & Suyker, A. E. (2003). Gross primary production and light
response parameters of four Southern Plains ecosystems estimated using long-term CO₂-flux tower measurements. *Global Biogeochemical Cycles*, 17(2), 1071. https://doi.org/10.1029/2002GB002023

Gough, C. M., Hardiman, B. S., Nave, L. E., Bohrer, G., Maurer, K. D., Vogel, C. S., ... Curtis, P. S. (2013). Sustained carbon uptake and storage following moderate disturbance in a Great Lakes forest. *Ecological Applications*, 23(5), 1202–1215. https://doi.org/10.1890/12-1554.1

Grünwald, T., & Bernhofer, C. (2007). A decade of carbon, water and energy flux measurements of an old spruce forest at the Anchor Station Tharandt. *Tellus B Chemical and Physical Meteorology*, 59(3), 387–396. https://doi.org/10.1111/j.1600-0889.2007.00259.x

Haykin, S. (1999). *Neural networks – A comprehensive foundation* (2nd ed.). Prentice Hall.

Jarvis, P., Rey, A., Petsikos, C., Wingate, L., Rayment, M., Pereira, J., ... Valentini, R. (2007). Drying and wetting of Mediterranean soils stimulates decomposition and carbon dioxide emission: The “Birch Effect”. *Tree Physiology*, 27, 929–940. https://doi.org/10.1093/treephys/27.7.929

Jung, M., Schwalm, C., Migliavacca, M., Walther, S., Camps-Valls, G., Koirala, S., ... Reichstein, M. (2020). Scaling carbon fluxes from eddy covariance sites to globe: Synthesis and evaluation of the FLUXCOM approach. *Biogeosciences*, 17, 1343–1365. https://doi.org/10.5194/bg-17-1343-2020

Keenan, T. F., Migliavacca, M., Papale, D., Baldocchi, D., Reichstein, M., Torn, M., & Wutzler, T. (2019). Widespread inhibition of daytime ecosystem respiration. *Nature Ecology & Evolution*, 3, 407–415. https://doi.org/10.1038/s41559-019-0809-2

Kljun, N., Calanca, P., Rotach, M. W., & Schmid, H. P. (2015). A simple two-dimensional parameterisation for Flux Footprint Prediction (FFP). *Geoscientific Model Development*, 8, 3695–3713. https://doi.org/10.5194/gmd-8-3695-2015

Körner, C. (1995). Leaf diffusive conductances in the major vegetation types of the globe. In E. D. Schulze & M. M. Caldwell (Eds.), *Ecophysiology of photosynthesis, springer study edition* (Vol. 100, pp. 463–490), Berlin, Heidelberg: Springer. https://doi.org/10.1007/978-3-642-79354-7_22

Kurbatova, J., Li, C., Varlagin, A., Xiao, X., & Vygodskaya, N. (2008). Modeling carbon dynamics in two adjacent spruce forests with different soil conditions in Russia. *Biogeosciences*, 5(4), 969–980. https://doi.org/10.5194/bg-5-969-2008

Lasslop, G., Migliavacca, M., Bohrer, G., Reichstein, M., Bahn, M., Ibrom, A., ... Cescatti, A. (2012). On the choice of the driving temperature for eddy-covariance carbon dioxide flux partitioning. *Biogeosciences*, 9, 5243–5259. https://doi.org/10.5194/bg-9-5243-2012

Lasslop, G., Reichstein, M., Papale, D., Richardson, A. D., Arneth, A., Barr, A., ... Wohlfahrt, G. (2010). Separation of net ecosystem exchange into assimilation and respiration using a light response curve approach: Critical issues and global evaluation. *Global Change Biology*, 16, 187–208. https://doi.org/10.1111/j.1365-2486.2009.02041.x

Lloyd, J., & Taylor, J. A. (1994). On the temperature dependence of soil respiration. *Functional Ecology*, 8, 315–323. https://doi.org/10.2307/2389824

Marcolla, B., Cescatti, A., Manca, G., Zorer, R., Cavagna, M., Fiora, A., ... Zampedri, R. (2011). Climatic controls and ecosystem responses drive the inter-annual variability of the net ecosystem exchange of an alpine meadow. *Agricultural and Forest Meteorology*, 151(9), 1233–1243. https://doi.org/10.1016/j.agrformet.2011.04.015

Medlyn, B. E. (1998). Physiological basis of the light use efficiency model. *Tree Physiology*, 18(3), 167–176. https://doi.org/10.1093/treephys/18.3.167

Medlyn, B. E., De Kauwe, M. G., Lin, Y.-S., Knauer, J., Duursma, R. A., Williams, C. A., ... Wingate, L. (2017). How do leaf and ecosystem measures of water-use efficiency compare? *New Phytologist*, 216, 758–770. https://doi.org/10.1111/nph.14626

Meyer, W. S., Kondrlová, E., & Koerber, G. R. (2015). Evaporation of perennial semi-arid woodland in southeastern Australia is adapted for

irregular but common dry periods. *Hydrological Processes*, 29(17), 3714–3726. https://doi.org/10.1002/hyp.10467

Misson, L., Baldocchi, D. D., Black, T. A., Blanken, P. D., Brunet, Y., Curiel Yuste, J., ... Goldstein, A. H. (2007). Partitioning forest carbon fluxes with overstory and understory eddy-covariance measurements: A synthesis based on FLUXNET data. *Agricultural and Forest Meteorology*, 144, 14–31. https://doi.org/10.1016/j.agrformet.2007.01.006

Moffat, A. M., Beckstein, C., Churkina, G., Mund, M., & Heimann, M. (2010). Characterization of ecosystem responses to climatic controls using artificial neural networks. *Global Change Biology*, 16(10), 2737–2749. https://doi.org/10.1111/j.1365-2486.2010.02171.x

Moffat, A. M., Papale, D., Reichstein, M., Hollinger, D. Y., Richardson, A. D., Barr, A. G., ... Stauch, V. J. (2007). Comprehensive comparison of gap-filling techniques for eddy covariance net carbon fluxes. *Agricultural and Forest Meteorology*, 147(3–4), 209–232. https://doi.org/10.1016/j.agrformet.2007.08.011

Monson, R. K., Turnipseed, A. A., Sparks, J. P., Harley, P. C., Scott-Denton, L. E., Sparks, K. L., & Huxman, T. E. (2002). Carbon sequestration in a high-elevation, subalpine forest. *Global Change Biology*, 8(5), 459–478. https://doi.org/10.1046/j.1365-2486.2002.00480.x

Moors, E. J. (2012). *Water use of forests in the Netherlands*. PhD thesis, Vrije Universiteit Amsterdam, The Netherlands.

Moureaux, C., Debacq, A., Bodson, B., Heinesch, B., & Aubinet, M. (2006). Annual net ecosystem carbon exchange by a sugar beet crop. *Agricultural and Forest Meteorology*, 139(1–2), 25–39. https://doi.org/10.1016/j.agrformet.2006.05.009

Musavi, T., Migliavacca, M., van de Weg, M. J., Kattge, J., Wohlfahrt, G., van Bodegom, P. M., ... Mahecha, M. D. (2016). Potential and limitations of inferring ecosystem photosynthetic capacity from leaf functional traits. *Ecology and Evolution*, 6, 7352–7366. https://doi.org/10.1002/ece3.2479

Ogée, J., Brunet, Y., Loustau, D., Berbigier, P., & Delzon, S. (2003). MuSICA, a CO₂ water and energy multilayer, multileaf pine forest model: Evaluation from hourly to yearly time scales and sensitivity analysis. *Global Change Biology*, 9, 697–717. https://doi.org/10.1046/j.1365-2486.2003.00628.x

Ogée, J., Peylin, P., Ciais, P., Bariac, T., Brunet, Y., Berbigier, P., ... Bonnefond, J. M. (2003). Partitioning net ecosystem carbon exchange into net assimilation and respiration using ¹³CO₂ measurements: A cost effective sampling strategy. *Global Biogeochemical Cycles*, 17(2), 1070. https://doi.org/10.1029/2002GB001995

Oikawa, P. Y., Sturtevant, C., Knox, S. H., Verfaillie, J., Huang, Y. W., & Baldocchi, D. D. (2017). Revisiting the partitioning of net ecosystem exchange of CO₂ into photosynthesis and respiration with simultaneous flux measurements of ¹³CO₂ and CO₂, soil respiration and a biophysical model, CANVEG. *Agricultural and Forest Meteorology*, 234–235, 149–163. https://doi.org/10.1016/j.agrformet.2016.12.016

Pacheco-Labrador, J., Perez-Priego, O., El-Madany, T. S., Julitta, T., Rossini, M., Guan, J., ... Migliavacca, M. (2019). Multiple-constraint inversion of SCOPE. Evaluating the potential of GPP and SIF for the retrieval of plant functional traits. *Remote Sensing of Environment*, 234, 111362. https://doi.org/10.1016/j.rse.2019.111362

Papale, D., Black, T. A., Carvalhais, N., Cescatti, A., Chen, J., Jung, M., ... Raduly, B. (2015). Effect of spatial sampling from European flux towers for estimating carbon and water fluxes with artificial neural networks. *Journal of Geophysical Research - Biogeosciences*, 120, 1941–1957. https://doi.org/10.1002/2015JG002997

Papale, D., Migliavacca, M., Cremonese, E., Cescatti, A., Alberti, G., Balzarolo, M., ... Valentini, R. (2015). Carbon, water and energy fluxes of terrestrial ecosystems in Italy. In R. Valentini & F. Miglietta (Eds.), *The greenhouse gas balance of Italy: An insight on managed and natural terrestrial ecosystems* (pp. 11–45). Berlin, Heidelberg: Springer. ISBN: 9783642324239. https://doi.org/10.1007/978-3-642-32424-6_2

Papale, D., Reichstein, M., Aubinet, M., Canfora, E., Bernhofer, C., Kutsch, W., ... Yakir, D. (2006). Towards a standardized processing of Net Ecosystem Exchange measured with eddy covariance technique:
Algorithms and uncertainty estimation. *Biogeosciences*, 3, 571–583.
https://doi.org/10.5194/bg-3-571-2006

Papale, D., & Valentini, R. (2003). A new assessment of European forests carbon exchanges by eddy fluxes and artificial neural network spatialization. *Global Change Biology*, 9(4), 525–535. https://doi.org/10.1046/j.1365-2486.2003.00609.x

Pastorello, G., Trotta, C., Canfora, E., Chu, H., Christianson, D., Cheah, Y.-W. (2020). The FLUXNET2015 dataset and the ONE Fluxprocessing pipeline for eddy covariance data. Scientific Data (in press).

Pastorello, G. Z., Papale, D., Chu, H., Trotta, C., Agarwal, D. A., Canfora, E., ... Torn, M. S. (2017). A new data set to keep a sharper eye on land-air exchanges. *Eos*, 98. https://doi.org/10.1029/2017EO071597

Pilegaard, K., Ibmro, A., Courtney, M. S., Hummelshøj, P., & Jensen, N. O. (2011). Increasing net CO₂ uptake by a Danish beech forest during the period from 1996 to 2009. *Agricultural and Forest Meteorology*, 151(7), 934–946. https://doi.org/10.1016/j.agrformet.2011.02.013

Prescher, A.-K., Grünwald, T., & Bernhofer, C. (2010). Land use regulates carbon budgets in eastern Germany: From NEE to NBP. *Agricultural and Forest Meteorology*, 150, 1016–1025. https://doi.org/10.1016/j.agrformet.2010.03.008

Reichstein, M., Bahn, M., Mahecha, M. D., Kattge, J., & Baldocchi, D. D. (2014). Linking plant and ecosystem functional biogeography. *Proceedings of the National Academy of Sciences of the United States of America*, 111, 13697–13702. https://doi.org/10.1073/pnas.1216065111

Reichstein, M., Camps-Valls, G., Stevens, B., Jung, M., Denzler, J., Carvalhais, N., & Prabhat, M. (2019). Deep learning and process understanding for data-driven Earth system science. *Nature*, 566(7743), 195–204. https://doi.org/10.1038/s41586-019-0912-1

Reichstein, M., Falge, E., Baldocchi, D., Papale, D., Aubinet, M., Berbigier, P., ... Valentini, R. (2005). On the separation of net ecosystem exchange into assimilation and ecosystem respiration: Review and improved algorithm. *Global Change Biology*, 11, 1424–1439. https://doi.org/10.1111/j.1365-2486.2005.001002.x

Rey, A., Pegoraro, E., Tedeschi, V., De Parri, I., Jarvis, P. G., & Valentini, R. (2002). Annual variation in soil respiration and its components in a coppice oak forest in Central Italy. *Global Change Biology*, 8(9), 851–866. https://doi.org/10.1046/j.1365-2486.2002.00521.x

Rivera, J. P., Verrelst, J., Gómez-Dans, J., Muñoz-Mari, J., Moreno, J., & Camps-Valls, G. (2015). An emulator toolbox to approximate radiative transfer models with statistical learning. *Remote Sensing*, 7, 9347–9370. https://doi.org/10.3390/rs70709347

Scanlon, T. M., & Kustas, W. (2010). Partitioning carbon dioxide and water vapor fluxes using correlation analysis. *Agricultural and Forest Meteorology*, 150, 89–99. https://doi.org/10.1016/j.agrformet.2009.09.005

Scott, R. L., Biederman, J. A., Hamerlynck, E. P., & Barron-Gafford, G. A. (2015). The carbon balance pivot point of southwestern U.S. semiarid ecosystems: Insights from the 21st century drought. *Journal of Geophysical Research. Biogeoscience*, 120(12), 2612–2624. https://doi.org/10.1002/2015JG003181

Scott, R. L., Hamerlynck, E. P., Jenerette, G. D., Moran, M. S., & Barron-Gafford, G. A. (2010). Carbon dioxide exchange in a semidesert grassland through drought-induced vegetation change. *Journal of Geophysical Research. Biogeoscience*, 115(G3), G03026. https://doi.org/10.1029/2010JG001348

Scott, R. L., Jenerette, G. D., Potts, D. L., & Huxman, T. E. (2009). Effects of seasonal drought on net carbon dioxide exchange from a woody-plant-encroached semiarid grassland. *Journal of Geophysical Research. Biogeoscience*, 114(G4), G04004. https://doi.org/10.1029/2008JG000900

Sulman, B. N., Tyler Roman, D. T., Scanlon, T. M., Wang, L., & Novick, K. A. (2016). Comparing methods for partitioning a decade of carbon dioxide and water vapor fluxes in a temperate forest. *Agricultural and Forest Meteorology*, 226–227, 229–245. https://doi.org/10.1016/j.agrformet.2016.06.002

Suni, T., Rinne, J., Reissell, A., Altimir, N., Keronen, P., Rannik, Ü., ... Vesala, T. (2003). Long-term measurements of surface fluxes above a Scots pine forest in Hyytiälä, southern Finland, 1996–2001. *Boreal Environment Research*, 8, 287–301.

Tramontana, G., Jung, M., Schwalm, C. R., Ichii, K., Camps-Valls, G., Ráduly, B., ... Papale, D. (2016). Predicting carbon dioxide and energy fluxes across global FLUXNET sites with regression algorithms. *Biogeosciences*, 13, 4291–4313. https://doi.org/10.5194/bg-13-4291-2016

Verrelst, J., Van der Tol, C., Magnani, F., Sabater, N., Rivera, J. P., Mohammed, G., & Moreno, J. (2016). Evaluating the predictive power of sun-induced chlorophyll fluorescence to estimate net photosynthesis of vegetation canopies: A SCOPE modeling study. *Remote Sensing of Environment*, 176, 139–151. https://doi.org/10.1016/j.rse.2016.01.018

Vicent, J., Verrelst, J., Rivera-Caicedo, J. P., Sabater, N., Muñoz-Mari, J., Camps-Valls, G., & Moreno, J. (2018). Emulation as an accurate alternative to interpolation in sampling radiative transfer codes. *IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing*, 11(12), 4918–4931. https://doi.org/10.1109/JSTARS.2018.2875330

Vitale, L., Di Tommasi, P., D'Urso, G., & Magliulo, V. (2016). The response of ecosystem carbon fluxes to LAI and environmental drivers in a maize crop grown in two contrasting seasons. *International Journal of Biometeorology*, 60(3), 411–420. https://doi.org/10.1007/s00484-015-1038-2

Wehr, R., Munger, J. W., McManus, J. B., Nelson, D. D., Zahniser, M. S., Davidson, E. A., ... Saleska, S. R. (2016). Seasonality of temperate forest photosynthesis and daytime respiration. *Nature*, 534, 680–683. https://doi.org/10.1038/nature17966

Wehr, R., & Saleska, S. R. (2015). An improved isotopic method for partitioning net ecosystem-atmosphere CO₂ exchange. *Agricultural and Forest Meteorology*, 214–215, 515–531. https://doi.org/10.1016/j.agrformet.2015.09.009

Whelan, M. E., Lennartz, S. T., Gimeno, T. E., Wehr, R., Wohlfahrt, G., Wang, Y., ... Campbell, J. E. (2018). Reviews and syntheses: Carbonyl sulfide as a multi-scale tracer for carbon and water cycles. *Biogeosciences*, 15, 3625–3657. https://doi.org/10.5194/bg-15-3625-2018

Williams, C. A., Hanan, N., Scholes, R. J., & Kutsch, W. (2009). Complexity in water and carbon dioxide fluxes following rain pulses in an African savanna. *Oecologia*, 161, 469–480. https://doi.org/10.1007/s00442-009-1405-y

Wohlfahrt, G., & Galvagno, M. (2017). Revisiting the choice of the driving temperature for eddy covariance CO₂ flux partitioning. *Agricultural and Forest Meteorology*, 237–238, 135–142. https://doi.org/10.1016/j.agrformet.2017.02.012

Zhang, Y., Guanter, L., Berry, J., Van der Tol, C., Yang, X., Tang, J., & Zhang, F. (2016). Model-based analysis of the relationship between sun-induced chlorophyll fluorescence and gross primary production for remote sensing applications. *Remote Sensing of Environment*, 187, 145–155. https://doi.org/10.1016/j.rse.2016.10.016

## SUPPORTING INFORMATION

Additional supporting information may be found online in the Supporting Information section.

**How to cite this article:** Tramontana G, Migliavacca M, Jung M, et al. Partitioning net carbon dioxide fluxes into photosynthesis and respiration using neural networks. *Glob Change Biol.* 2020;26:5235–5253. https://doi.org/10.1111/gcb.15203