Observational partitioning of water and CO₂ fluxes at National Ecological Observatory Network (NEON) sites: a 5-year dataset of soil and plant components for spatial and temporal analysis

Long-term time series of transpiration, evaporation, plant net photosynthesis, and soil respiration are essential for addressing numerous research questions related to ecosystem functioning. However, quantifying these fluxes is challenging due to the lack of reliable and direct measurement techniques, which has left gaps in the understanding of their temporal cycles and spatial variability. To help address this open challenge, we generated a dataset of these four components by implementing five (conventional and novel) approaches to partition total evapotranspiration (ET) and CO₂ fluxes into plant and soil fluxes across 47 National Ecological Observatory Network (NEON) sites. The final dataset (https://doi.org/10.5281/zenodo.12191876; Zahn and Bou-Zeid, 2024) spans a 5-year period and covers various ecosystems, including forests, grasslands, and agricultural terrain. This is the first comprehensive dataset covering such a wide spatial and temporal distribution. Overall, we observed good agreement across most methods for ET components, increasing confidence in these estimates. Partitioning of CO₂ components, on the other hand, was found to be less robust and more dependent on prior knowledge of water use efficiency. This highlights some limitations of these present methods that we discuss, emphasizing the broader challenge posed by the lack of an accurate reference method to validate against. Despite these limitations, this dataset has several potential applications, especially in addressing critical questions regarding the response of ecosystems to extreme weather events, which are expected to become more severe and frequent with climate change.

1 Introduction

Plant transpiration and photosynthesis are important components of the global water and CO₂ cycles. At the ecosystem level, these quantities reflect the complex interactions within the soil–vegetation–atmosphere continuum. Given the complexity of these interactions and the challenges associated with measuring or modeling leaf–atmosphere exchanges, many open questions remain. For instance, the competing effects of rising CO₂ and temperature on plant water use efficiency (Mengis et al., 2015; Kirschbaum and McMillan, 2018; Dusenge et al., 2019; Baslam et al., 2020; Wang et al., 2022) are still being debated, together with their response and resilience to more frequent and intense droughts and changes in soil moisture conditions (Maxwell and Condon, 2016; Lesk et al., 2021). In addition to the challenges of investigating such nonlinear processes, one of the main barriers to advancing research in this area has been the lack of long-term observations and/or reliable estimates of plant transpiration and CO₂ assimilation over various types of vegetation and soils.

Eddy covariance (EC) is a reliable approach to continuously monitor evapotranspiration (ET) and net CO₂ (F_c) fluxes across ecosystems in a standardized manner. However, separate measurement of the individual ecosystem components, i.e., plant transpiration (T) and soil evaporation (E), together with plant gross primary production (GPP) and ecosystem respiration (R_eco), is still difficult. Traditional measurement techniques, such as sap flow and lysimeters, often lack the temporal and spatial resolutions necessary for capturing ET flux components comprehensively (Kool et al., 2014; Stoy et al., 2019). Other techniques, such as soil chambers (Kool et al., 2014), stable isotopes (Good et al., 2014), and carbonyl sulfide as a tracer for plant fluxes (Wohlfahrt et al., 2011), are not yet deployable at large spatial and temporal scales. Additionally, all of these methods carry their own uncertainties, can be logistically challenging to implement, and are often costly to scale up. Partitioning of CO₂ fluxes into GPP and R_eco can also be done based on nighttime CO₂ flux extrapolation (Reichstein et al., 2005) and light-response curves (Lasslop et al., 2010); this method has been widely implemented and is part of FLUXNET products (Pastorello et al., 2020). Nonetheless, it cannot partition ET components. As a result, there is a scarcity of datasets representing the temporal and spatial variabilities of both ET and CO₂ flux components across a broad range of ecosystems and weather conditions.

An alternative approach to obtain long-term time series of flux components is partitioning of evapotranspiration and net CO₂ exchanges into their respective plant and soil components based strictly on the analysis of their turbulent perturbation signals. This approach is the focus of this paper, which applies five recently developed and tested partitioning methods to disentangle the contributions of transpiration, evaporation, net photosynthesis, and soil respiration to the total ecosystem fluxes measured across all the EC towers of the National Science Foundation's National Ecological Observatory Network (NEON). The methods employed in this study have been rigorously evaluated in previous studies that combined experimental data from EC sites and virtual EC data from large-eddy simulations. Although, as with all partitioning approaches, the validation is hindered by a lack of a universally accepted accurate standard method, previous research has highlighted conditions under which these turbulence methods agree and are more reliable, together with conditions under which their partitioning estimates may have larger errors. In the dataset we produce in this work, we consolidate the results from five partitioning methods and calculate their ensemble average, providing a measure of uncertainty for each flux period at each site.

With a total of 47 sites spanning diverse ecosystems in the United States, this dataset can be an invaluable resource for future investigation of the spatial variability and ecological and meteoclimatic drivers shaping flux dynamics across different landscapes. It has the potential to contribute to the advancement of our understanding of ecosystem processes and their responses to current and future climates, in turn better informing ecosystem management and climate change mitigation and adaptation practices.

2 Study area and input data

NEON is a comprehensive initiative designed to collect long-term open-access environmental observations in the continental United States (Metzger et al., 2019). It operates a network of field sites strategically located to capture the diversity of ecosystems and environmental conditions present in the country. One of the main features of NEON is the standardized collection of data across sites, including biological, hydrological, and meteorological observations. As part of its effort to understand land–atmosphere exchanges of energy, water, and CO₂, NEON has a total of 47 eddy-covariance towers spread across 24 states (see Table 1), where most sites are located in deciduous and evergreen forests (Fig. 1).

Figure 1 Locations of the 47 eddy-covariance towers maintained by NEON. Note that Alaska, Hawaii, and Puerto Rico are not scaled proportionally to the continental US area.

Table 1 Description of the 47 sites included in the present paper. In addition to the NEON site identification (ID) and state, the table lists the geographical coordinates, the ratio between the eddy-covariance measurement height and the mean canopy height ($z/h$), and the predominant vegetation type around the tower following Fiorella et al. (2021). The last two columns represent the availability of ET fluxes during the daytime hours (when solar radiation > 10 W m²) in the 5 years of record and the availability of transpiration estimates, based on the partitioning approaches used here, in the same period. The range represents the availability according to the different methods (least available to most available). Availability is computed after gap-filling as described in Sect. 3.4.2.

For the purpose of implementing the partitioning methods, we downloaded instantaneous measurements from the eddy-covariance towers from 2019 to 2023 for all 47 locations. All sites were equipped with the same gas analyzer (LICOR LI7200) and sonic anemometers (Campbell Scientific, CSAT-3 3-D) collecting data at a frequency of 20 Hz. These time series were processed and prepared to be used by the five partitioning methods, as described in the next section. In addition to raw data, we also downloaded the processed eddy-covariance data (product level 4) for the same period (NEON, 2024a, b), which include evapotranspiration and net ecosystem CO₂ fluxes as processed by eddy4R (Metzger et al., 2017). Fluxes computed by eddy4R include all standard corrections, such as high-frequency spectral correction, some of which are only applied to the final turbulent fluxes and not the instantaneous time series. Thus, we use the total turbulent fluxes as quantified by NEON to rescale all the flux components that we computed (T, E, R, and P), ensuring that they are consistent with the corrected fluxes reported by the network.

As the last step of postprocessing, a gap-filling algorithm was implemented to maximize the availability of flux-partitioning estimates. Additional meteorological variables were downloaded using the R package neonUtilities (Lunch et al., 2024) and used as features to implement the Extreme Gradient Boosting (XGBoost) algorithm. The selected variables (inputs) included air temperature (NEON, 2024g, h), incoming solar radiation (NEON, 2024i, j), relative humidity (NEON, 2024e, f), photosynthetically active radiation (NEON, 2024c, d), and wind speed (NEON, 2024k, l). All quantities were measured at the top of the tower. The sensible heat flux and friction velocity, included in the bundled data products, were also used. Due to the presence of long gaps and/or noisy measurements, precipitation and soil moisture were not included as features.

3 Methods

Below, we describe the methods used to partition the data and the processing steps followed to generate the final dataset. A diagram summarizing all the preprocessing and postprocessing steps is shown in Fig. 2. We implemented five methods to partition evapotranspiration into soil evaporation and plant transpiration as well as net CO₂ flux into soil respiration and net plant photosynthesis. It is important to note that these quantities differ from GPP and ecosystem respiration (see Wohlfahrt and Gu, 2015, for a detailed discussion of the different definitions of CO₂ flux components). Since these methods are based on principles of turbulent transport, they do not provide the information needed to separate soil respiration from plant respiration. Instead, the framework effectively partitions below-canopy fluxes, which include respiration and evaporation from soil, roots, and possibly the lower parts of tree trunks, and above-canopy fluxes, which encompass plant transpiration and net photosynthesis. Thus, none of these approaches allows for the distinction between canopy and ground respiration. Another limitation of these models is that they cannot distinguish between transpiration and evaporation of intercepted water on the leaves. Nonetheless, since EC data are often lost through poor data quality during and immediately following rainfall periods, many such periods are likely removed from the time series already. With these definitions in mind, we refer to net photosynthesis simply as photosynthesis (P) and plant transpiration simply as transpiration (T) in the following sections, while evaporation (E) and respiration (R) are considered to originate solely from the ground or soil levels. We advise future users to carefully note the specific definitions of each component, as this is essential for making meaningful comparisons across different partitioning products.

Figure 2 Diagram representing all steps of data processing used to generate the final dataset.

3.1 Partitioning approaches

The partitioning methods implemented in this study are based on the idea of transport similarity between CO₂ and water vapor. It assumes that turbulent structures, or “eddies”, simultaneously transport evaporation and respiration from the soil towards the sensor, where the respective concentrations of these “event samples” are measured. Similarly, turbulence also transports air parcels from the canopy, where net photosynthesis and plant transpiration take place. By separating the signature of soil and canopy events in the time series, these methods try to infer the respective magnitude of each flux component. Note that, while both evapotranspiration components are sources of water vapor (i.e., they are positive fluxes), CO₂ fluxes are a combination of a source (soil respiration) and a sink (photosynthesis). This feature of the CO₂ flux is what creates this framework, since it works as a tracer that identifies the origin of water vapor. However, another consequence of this feature is that, while the components E and T are bounded by E+T, the magnitudes of R>0 and P<0 are not constrained by $F_{\mathrm{c}}=R+P$. As a result, partitioning results for CO₂ components tend to be more uncertain. Nonetheless, as will be explained in the next section, two of the partitioning methods implemented here require the water use efficiency, $W=P/T$, as an input. While the need for a parameterization for W is a disadvantage for the method, the inclusion of W adds a bound to P and R.

In total, five partitioning methods were implemented: flux-variance similarity (FVS) (Scanlon and Sahu, 2008; Scanlon and Kustas, 2010; Scanlon et al., 2019), modified relaxed eddy accumulation (MREA) (Thomas et al., 2008), conditional eddy covariance (CEC) (Zahn et al., 2022), modified CEC that also utilizes water use efficiency (CECw) (Zahn et al., 2024), and conditional eddy accumulation (CEA) (Zahn et al., 2024). The advantage of these models is that they mostly rely on high-frequency eddy-covariance observations, requiring no or few input parameters. While MREA, CEC, and CEA are based solely on statistics computed from the high-frequency data, FVS and CECw additionally require water use efficiency as an input.

Three of the methods, FVS, MREA, and CEC, have been investigated extensively in previous studies (Thomas et al., 2008; Klosterhalfen et al., 2019; Scanlon et al., 2019; Zahn et al., 2022). More recently, these three methods were tested by Zahn et al. (2024) using large-eddy simulations; they also formulated and tested the CECw and CEA approaches to probe the assumptions and performance of all the methods. Their simulation-based analyses indicated that all the methods result in reliable results for a range of flux combinations that are expected in real ecosystems. While that study was not able to pinpoint conditions under which the accuracy of the results in field experiments can be guaranteed, it delineated flux combinations that may result in higher uncertainties, it quantified the impact of errors in W on the results of FVS and CECw, and it identified the most consequential assumptions in the various methods. The broad recommendations of that study are as follows:

1. The measurement height z, i.e., the height where the EC system is placed, should be as close as possible to the mean canopy height, h, to better distinguish between the soil and plant signals. However, while $z/h<\mathrm{3}$ is a good recommendation, we note that the performance also depends on other factors such as canopy density and the specific method, as some approaches were found to perform well even outside this region.

2. All flux components from the soil (evaporation and respiration) and from the vegetation (transpiration and photosynthesis) should be non-negligible. This criterion arises from the water vapor–CO₂ coupling, which needs the CO₂ signal to determine the origin of the water vapor.

3. The (anti)correlation between CO₂ and water vapor, ρ_c,q, should not be perfect, i.e., $\left|{\mathit{\rho }}_{\mathrm{c},q}\right|<\mathrm{1}$. Perfect correlation or anticorrelation indicates that turbulence fully mixed the scalars from the soil and from the canopy. In this case, the framework of all five partitioning models is no longer valid. Most often, high correlation is seen when measurements are too far from the canopy top or when one of the flux components is negligible. Thus, this condition should be satisfied if points 1 and 2 are addressed first.

Overall, by comparing the outputs of the five models and identifying periods or sites where they consistently agree, we can increase our confidence in the results obtained, even without a “true” value for comparison. In addition, the range of predictions from the various methods can be used to quantify the uncertainty in the overall partitioning.

3.2 Data preprocessing

We downloaded and processed high-frequency eddy-covariance turbulent exchange (ECTE) data using our in-house Python routines. The following variables were extracted from the original datasets for each site: the mixing ratios of CO₂ and water vapor (r_c and r_q); the three components of the velocity field in the x, y, and z directions (u, v, and w); the sonic temperature (T_s); and the atmospheric pressure (P). In addition, sensor flags describing the quality of the measurements and signal were used for quality control (see Table 2 for a description of the variables used). Our algorithm was implemented following the preprocessing of daily files as follows:

1. Instantaneous measurements that were assigned 1 (poor quality) by any of the flags described in Table 2 were discarded.

2. CO₂ and H₂O mixing ratio observations were discarded if the signal strength indicator was smaller than 0.7.

3. Air density (kg m⁻³) and mass concentrations of CO₂ and H₂O (mg m⁻³ and g m⁻³, respectively) were computed.

4. Daily files with less than 50 % raw high-frequency data were discarded, with very few cases falling into this category. Files that passed this test were split into half-hour blocks of 36 000 points each (at 20 Hz measurements).

5. Instantaneous observations were checked for physical plausibility, such as positive scalar concentrations.

6. Outliers were removed using a despiking algorithm (Zahn et al., 2016).

7. The velocity field was rotated using double-coordinate rotation.

8. CO₂ and H₂O mass concentrations, c and q, were lag-corrected, since the LI7200 is a closed-path sensor with tubing, using the cross-correlation technique (Rebmann et al., 2012).

9. Turbulent fluctuations were computed for all the quantities by removing the linear trends for all the respective variables.

The outcomes of these steps are 30 min time series comprising cleaned turbulence data. In the next section we describe how these data were used to partition ET and F_c fluxes following the five methods described in Sect. 3.1.

Table 2 Variables used for data preprocessing.

3.3 Implementation of the partitioning methods

Since all partitioning methods are based on turbulence transport and similarity, their implementation requires the computation of turbulent statistics such as covariances and correlation coefficients. Details of their implementation can be found in Zahn et al. (2022, 2024). These statistics are computed for each 30 min block. Using these turbulence statistics, we first implement the three methods that do not require extra inputs, i.e., CEC, CEA, and MREA. For each partitioning method, we thus obtain T_part, E_part, P_part, and R_part, where the subscript “part” represents any of the three methods. For instance, for CEC we have T_CEC, E_CEC, P_CEC, and R_CEC.

The implementation of FVS and CECw additionally requires an estimate of the water use efficiency W. Note that, in this context, W is defined as the ratio between net photosynthesis and transpiration, which again differs from other definitions. This quantity was also estimated for each 30 min block using the following expression:

$$\begin{array}{}\text{(1)}& W=\mathrm{0.65}{\displaystyle \frac{\stackrel{\mathrm{‾}}{c_{\mathrm{c}}}-\stackrel{\mathrm{‾}}{c_{\mathrm{s}}}}{\stackrel{\mathrm{‾}}{q_{\mathrm{c}}}-\stackrel{\mathrm{‾}}{q_{\mathrm{s}}}}},\end{array}$$

where $\stackrel{\mathrm{‾}}{q_{\mathrm{c}}}$ and $\stackrel{\mathrm{‾}}{c_{\mathrm{c}}}$ are the H₂O and CO₂ atmospheric mean concentrations near the canopy, while $\stackrel{\mathrm{‾}}{q_{\mathrm{s}}}$ and $\stackrel{\mathrm{‾}}{c_{\mathrm{s}}}$ are the mean intercellular concentrations. The near-canopy concentrations ($\stackrel{\mathrm{‾}}{q_{\mathrm{c}}}$ and $\stackrel{\mathrm{‾}}{c_{\mathrm{c}}}$) are obtained using the logarithmic profile and mean concentrations from EC data, while q_s is calculated under the assumptions of stomatal saturation (at air temperature) and a well-coupled leaf. $\stackrel{\mathrm{‾}}{c_{\mathrm{s}}}$, on the other hand, is the most challenging variable to obtain. Following previous studies (Wagle et al., 2021; Zahn et al., 2022), we implemented five models for $\stackrel{\mathrm{‾}}{c_{\mathrm{s}}}$, giving a total of five estimates of W. These parameterizations are summarized in Table 3. For each 30 min period, FVS and CECw were thus implemented five times, one for each parameterization of W. As a result, each flux component is estimated five times; for instance, following the nomenclature presented in Table 3, after implementing FVS we obtain $T_{{\mathrm{FVS}}_{\mathrm{CC}}}$, $T_{{\mathrm{FVS}}_{\mathrm{CR}}}$, $T_{{\mathrm{FVS}}_{\mathrm{lD}}}$, $T_{{\mathrm{FVS}}_{\mathrm{sD}}}$, and $T_{{\mathrm{FVS}}_{\mathrm{OPT}}}$. These five estimates are then averaged to $T_{{\mathrm{FVS}}_{\mathrm{MEAN}}}$ when at least one of the five outputs is available. The same averaging is done to obtain $E_{{\mathrm{FVS}}_{\mathrm{MEAN}}}$, $P_{{\mathrm{FVS}}_{\mathrm{MEAN}}}$, and $R_{{\mathrm{FVS}}_{\mathrm{MEAN}}}$. Likewise, the same procedure is applied to CECw, which also takes W as an input. Our final dataset includes the outputs from all five models as well as the averaged values, as described in the next section.

Scanlon et al. (2019)

Table 3 Descriptions of all five parameterizations adopted for W. More details can be found in Wagle et al. (2021) and Zahn et al. (2022). All five parameterizations were used to partition fluxes by FVS and CECw, thus resulting in five partitioning results for each partitioning approach.

3.4 Data postprocessing

3.4.1 Data cleaning and flux rescaling

The last section described how all five partitioning methods were implemented to obtain estimates of soil and canopy fluxes at a 30 min time interval. In postprocessing, we first eliminated unrealistic high flux magnitudes and counter-gradient fluxes by discarding periods when any of the following conditions were identified:

• T_part<0 W m⁻² or E_part<0 W m⁻²;

• P_part>0 mg m⁻² s⁻¹ or R_part<0 mg m⁻² s⁻¹;

• $T_{\mathrm{part}}+E_{\mathrm{part}}>\mathrm{1000}$ W m⁻²; and

• $|P_{\mathrm{part}}+R_{\mathrm{part}}|>\mathrm{10}$ mg m⁻² s⁻¹.

As previously described, all the partitioning methods relied on high-frequency time series for turbulent statistics and flux computation. Consequently, the resulting flux components lacked certain eddy-covariance corrections commonly applied to the final fluxes, including high-frequency spectral adjustments. To address this, we rescaled our fluxes based on those obtained from the NEON bundled data processed by eddy4R. We first computed correction factors for the ET and F_c fluxes as

$$\begin{array}{}\text{(2)}& {\mathrm{cor}}_{\mathrm{ET}}={\displaystyle \frac{{\mathrm{ET}}_{\mathrm{NEON}}}{{\mathrm{ET}}_{\mathrm{part}}}},\end{array}$$

$$\begin{array}{}\text{(3)}& {\mathrm{cor}}_{{\mathrm{F}}_{\mathrm{c}}}={\displaystyle \frac{F_{{\mathrm{c}}_{\mathrm{NEON}}}}{F_{{\mathrm{c}}_{\mathrm{part}}}}},\end{array}$$

where ET_NEON and $F_{{\mathrm{c}}_{\mathrm{NEON}}}$ are the evapotranspiration and net CO₂ fluxes downloaded directly from NEON, while ET_part and $F_{{\mathrm{c}}_{\mathrm{part}}}$ are the fluxes computed by our routines using the high-frequency data following the preprocessing described in Sect. 3.2. Note that the total fluxes, ${\mathrm{ET}}_{\mathrm{part}}=T_{\mathrm{part}}+E_{\mathrm{part}}$ and $F_{{\mathrm{c}}_{\mathrm{part}}}=P_{\mathrm{part}}+R_{\mathrm{part}}$, are the same for all five partitioning methods. The corrections were then applied as

$$\begin{array}{}\text{(4)}& {\displaystyle }& {\displaystyle }{T}_{\mathrm{part},\mathrm{cor}}={\mathrm{cor}}_{\mathrm{ET}}{T}_{\mathrm{part}}\mathrm{and}{E}_{\mathrm{part},\mathrm{cor}}={\mathrm{cor}}_{\mathrm{ET}}{E}_{\mathrm{part}},\text{(5)}& {\displaystyle }& {\displaystyle }{P}_{\mathrm{part},\mathrm{cor}}={\mathrm{cor}}_{{\mathrm{F}}_{\mathrm{c}}}{P}_{\mathrm{part}}\mathrm{and}{R}_{\mathrm{part},\mathrm{cor}}={\mathrm{cor}}_{{\mathrm{F}}_{\mathrm{c}}}{R}_{\mathrm{part}}.\end{array}$$

By rescaling our fluxes based on this comparison, we ensure that our fluxes are consistent with NEON products, mitigating potential discrepancies arising from variations in data processing methodologies. For simplicity, we drop the subscript “cor” in the remainder of the text and in the data files.

Additionally, we utilized NEON's bundled eddy-covariance products to assess the quality of both evapotranspiration and CO₂ flux data. Although our initial data preprocessing removed numerous periods with poor quality, we refrained from conducting further tests, such as assessing data stationarity, integral turbulence characteristics, and spectral analysis. Consequently, we integrated the flags as outlined in Metzger et al. (2022) into our final dataset. These flags are classified into two categories: 0, indicating compliance with all the quality criteria and rendering the data valid and reliable for analysis, and 1, indicating failure to meet one or more quality tests and requiring cautious utilization or possible exclusion. Retaining the original flags in our dataset facilitates the identification of periods characterized by high-quality data. Nonetheless, we opted not to automatically discard data based solely on these flags, leaving such decisions to the discretion of subsequent users.

3.4.2 Gap-filling

As a consequence of the inherent challenges associated with eddy-covariance measurements – such as sensor malfunction, data transmission, or even external factors such as rain and dust collection – flux time series are often incomplete. For our dataset, Table 1 shows the availability of ET fluxes during the 5-year period for each site. Considering only daytime periods (selected based on incoming solar radiation greater than 10 W m²), this reveals that some sites have available fluxes only 34 % of the time (BARR), while others have daytime ET estimates for up to 77 % of the 5-year record (KONA).

Due to mathematical constraints and/or invalidity of physical assumptions, partitioning methods are not guaranteed to converge to a solution for every 30 min time series, leading to further gaps in the record of flux components. As a consequence, the dataset generated as described in the previous sections contains missing periods, with the length of these gaps depending on the site and method. Table 4 illustrates that the methods with the lowest convergence rates – defined as the fraction of solutions found per available period – range from least to most solutions found as follows: FVS, CECw, MREA, CEA, and CEC. The lowest rate of convergence for the FVS method (on average 50 %) is not surprising and has been discussed in previous studies (Wagle et al., 2021; Zahn et al., 2022, 2024). Additionally, the low rate of solutions for FVS and CECw when combined with W_OPT is due to the model for this water use efficiency itself not always converging (i.e., W_OPT was not available).

Table 4 Success rate of each partitioning approach in finding a solution, expressed as a percentage of the available half-hour periods (i.e., only considering cases when the time series was available and the partitioning methods were implemented). The table shows the average and standard deviation across all the sites before and after gap-filling. For FVS and CECw, the convergence using each water use efficiency parameterization is also shown.

To increase the availability of periods with valid partitioning estimates for each method, we implemented XGBoost (Chen and Guestrin, 2016) for regression to predict transpiration. By training the algorithm to partition ET into E and T based on selected environmental features, the goal is to apply the trained algorithm to estimate T only when ET is available but a given partitioning method did not yield a valid estimate. That is, an XGBoost model is developed for each partitioning approach independently and then applied to fill the gaps in the time series of that approach. Additional machine learning algorithms, including neural networks and k-means clustering, were investigated but yielded similar or inferior results compared to XGBoost. Previous studies focusing on gap-filling eddy-covariance data also observed superior performance by XGBoost (Huang and Hsieh, 2020; Irvin et al., 2021). These findings are corroborated by research indicating that tree-based models outperform deep learning for tabular data used in regression tasks (Grinsztajn et al., 2022). Thus, the selected algorithm is robust and a very suitable choice for the problem at hand.

Only daytime conditions were considered, a constraint imposed by selecting periods when the incoming solar radiation was >10 W m². In addition, we ensured that only cases with positive ET were included. For regression purposes, we further filtered out periods where one or more input features were unavailable. As a result of these criteria, not every segment of the time series was suitable for gap-filling. Nonetheless, this approach maximizes the data availability of T and E, in particular during daytime, expanding the use of this dataset for future studies. Note that gap-filling was not applied to the F_c components. This decision was made due to the noisier nature of CO₂ fluxes. Furthermore, as discussed later, estimates for respiration and photosynthesis are less reliable across methods (larger variance between methods), and thus gap-filling for these components was not deemed as robust as with ET components.

To maximize the periods with available data features, we only selected meteorological variables with short periods of missing data and high quality. Following tests on different meteorological variables, the following features were retained: ET, air temperature (T_AIR), photosynthetically active radiation (PAR), relative humidity (RH), vapor pressure deficit (VPD), mean wind speed (U), sensible heat flux (H), and friction velocity (u_*). In addition, daily and seasonal cycles were represented as

$$\begin{array}{}\text{(6)}& {\displaystyle }& {\displaystyle }{X}_{\mathrm{diurnal}}={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\sin \left({\displaystyle \frac{\mathrm{2}\mathit{\pi }}{\mathrm{24}}}\mathrm{hour}-{\displaystyle \frac{\mathit{\pi }}{\mathrm{2}}}\right)+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}},\text{(7)}& {\displaystyle }& {\displaystyle }{X}_{\mathrm{yearly}}={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\sin \left({\displaystyle \frac{\mathrm{2}\mathit{\pi }}{\mathrm{365}}}\mathrm{julian}-{\displaystyle \frac{\mathit{\pi }}{\mathrm{2}}}\right)+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}},\end{array}$$

where X_diurnal and X_yearly represent the diurnal and yearly cycles corresponding to the hour of the day (in local time) and the Julian day, respectively. While there are other variables that could be relevant for this problem, they were excluded due to poor data quality or availability. Nonetheless, using only the listed variables resulted in excellent performance metrics.

For each site and partitioning method – including the different combinations of the FVS, CECw, and water use efficiency models taken individually – the dataset was divided into training (60 %), validation (20 %), and testing (20 %) sets. Across the different cases, the size of the training set varied between 1700 and 12 000 points. To monitor the performance of the algorithm, we used the coefficient of determination, R², defined as

$$\begin{array}{}\text{(8)}& R^{\mathrm{2}}=\mathrm{1}-{\displaystyle \frac{\sum (y-\widehat{y}{)}^{\mathrm{2}}}{\sum (y-\stackrel{\mathrm{‾}}{y}{)}^{\mathrm{2}}}},\end{array}$$

where y represents the observed values with average $\stackrel{\mathrm{‾}}{y}$, while $\widehat{y}$ represents the predictions. This metric indicates how well the independent variables explain the variance in the dependent variable. Thus, the upper bound $R^{\mathrm{2}}=+\mathrm{1}$ is expected only if the selected features explain 100 % of the variance.

After experimenting with different hyperparameters, the same values for the learning rate (0.01), maximum depth of a tree (5), and number of trees (1000) were used to train all the sites and methods individually. Finally, gap-filling was only performed for a particular case, when its respective R² on the test set was greater than 0.7. A few cases had an R² below this threshold, which was found to correspond most often to datasets combining FVS and W_CC or W_CR over shrubland or tundra sites, combinations which resulted in the lowest number of valid partitioning estimates (smallest training sets). Excluding these few exceptions, the majority of the cases featured R²>0.8 (and bias < 1 W m⁻²) and were deemed reliable for gap-filling. The final metrics for all the models and sites can be found in the Supplement. All auxiliary data used to implement the gap-filling algorithm are also available in the final dataset.

As expected, the algorithm identified ET as the most relevant variable across all the cases (approximately 80 %). Note that, in this context, we do not use ET to explain the variability in E or T in the sense of an environmental forcing. Instead, ET is used to teach the model how to partition the fluxes, effectively functioning as a machine-learning-based partitioning method. In Sect. 6, we briefly explore feature importance for T and $T/\mathrm{ET}$ in the context of environmental forcings, where only environmental variables are included as features (excluding ET).

While the variability in feature importance was small across the vegetation types, a larger variability was observed between the methods. Notably, FVS combined with the water use efficiency models W_CC and W_CR often attributed a lower importance to ET (as low as 27 %) and a higher importance to VPD (as high as 25 %). Previous work (Wagle et al., 2023) has shown that these two water use efficiency models are very sensitive to the choice of parameters and greatly overestimate W compared to the other models, especially early in the morning. While investigating the possible causes of this discrepancy across water use efficiency models is beyond the scope of this paper, deeper analyses will be enabled by this dataset.

Following the algorithm evaluation, transpiration was estimated for the missing data periods and defined as $T_{\mathrm{part}}^{\mathrm{xgb}}$, where the superscript “xgb” refers to the gap-filled variable using XGBoost. Only half-hour periods with $\mathrm{0}<T_{\mathrm{part}}^{\mathrm{xgb}}<\mathrm{ET}$ were gap-filled. Evaporation was then estimated as $E_{\mathrm{part}}^{\mathrm{xgb}}=\mathrm{ET}-T_{\mathrm{part}}^{\mathrm{xgb}}$. For FVS and CECw, we also computed the average across their different water use efficiency outputs, corresponding to $T_{{\mathrm{FVS}}_{\mathrm{MEAN}}}^{\mathrm{xgb}}$ and $T_{{\mathrm{CECw}}_{\mathrm{MEAN}}}^{\mathrm{xgb}}$.

The convergence rate to valid solutions after gap-filling is shown in the last column of Table 4. Most notably, gap-filling increased the presence of estimates from the FVS method for 80 % to 90 % of the periods when partitioning was applied. To further illustrate the completeness of the flux-partitioning record in the 5-year dataset, Table 1 shows the availability of flux components as a fraction of the total number of half-hour periods in the record. Overall, all the methods cover a similar temporal distribution of flux partitioning and are potential candidates for ensemble averaging.

4 Description of the final dataset

The final dataset is available for download at https://doi.org/10.5281/zenodo.12191876 (Zahn and Bou-Zeid, 2024). It is organized into different folders for each site, with each site containing a .csv file for each method. This format is selected to be user-friendly and accessible in various programming languages and software packages. For FVS and CECw, in addition to their ensemble averages for different water use efficiency parameterizations, we include each individual output. This allows future users to compare the impact of water use efficiency with the partitioning output.

All the results are in UTC, where this represents the start time of measurements for a half-hour period (i.e., starting at minute “00” or “30”). Each file contains all fluxes in continuous time stamps of 30 min. In cases where a specific variable was unavailable, the corresponding cell in the .csv file was left blank. A Python script showing how to read the final data files is included with the dataset.

In each file, the following variables are included:

• Tpart, transpiration [W m⁻²];

• Tpart_xgb, transpiration after gap-filling [W m⁻²];

• Ppart, plant net photosynthesis [mg${}_{{\mathrm{CO}}_{\mathrm{2}}}$ m⁻² s⁻¹];

• Epart, below-canopy evaporation [W m⁻²];

• Epart_xgb, below-canopy evaporation after gap-filling [W m⁻²];

• Rpart, below-canopy respiration [mg${}_{{\mathrm{CO}}_{\mathrm{2}}}$ m⁻² s⁻¹];

• Fc, net CO₂ flux [mg${}_{{\mathrm{CO}}_{\mathrm{2}}}$ m⁻² s⁻¹];

• ET, evapotranspiration [W m⁻²];

• qfFinalCO2, final NEON quality flag for the CO₂ flux (0: Fc passed all quality assurance tests; 1: Fc was flagged by one or more tests); and

• qfFinalH2O, final NEON quality flag for the ET flux (0: ET passed all quality assurance tests; 1: ET was flagged by one or more quality assurance tests).

For the flux components, “part” represents any of the methods. In addition, a separate file, WUEinput.csv, contains the water use efficiency values from all five parameterizations, defined as follows:

• const_ppm [kg${}_{{\mathrm{CO}}_{\mathrm{2}}}$(kg${}_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}{)}^{-\mathrm{1}}$];

• const_ratio [kg${}_{{\mathrm{CO}}_{\mathrm{2}}}$(kg${}_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}{)}^{-\mathrm{1}}$];

• linear [kg${}_{{\mathrm{CO}}_{\mathrm{2}}}$(kg${}_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}{)}^{-\mathrm{1}}$];

• sqrt [kg${}_{{\mathrm{CO}}_{\mathrm{2}}}$(kg${}_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}{)}^{-\mathrm{1}}$]; and

• opt [kg${}_{{\mathrm{CO}}_{\mathrm{2}}}$(kg${}_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}{)}^{-\mathrm{1}}$].

Finally, a file named MeteorologicalVariables.csv is also included, which contains all the auxiliary meteorological variables used for gap-filling and feature importance analysis. This file is synchronized with the partitioning data file to contain the same records. The following variables are available:

• SoilMoisture, soil moisture measured at the second level below the surface (depth varies by site) [m³ m⁻³];

• AirTemperature, air temperature [°C];

• SolarRadIn, incoming solar radiation measured at the top of the tower [W m⁻²];

• RH, relative humidity (%);

• Rain, accumulated precipitation at the top of the tower during the 30 min period [mm];

• PAR, photosynthetically active radiation measured at the top of the tower [µmolec. m⁻² s⁻¹];

• fluxTemp, sensible heat flux at the top of the tower [W m⁻²];

• ustar, friction velocity at the top of the tower [m s⁻¹];

• WindSpeed, mean wind speed at the top of the tower [m s⁻¹]; and

• vpd, vapor pressure deficit [kPa].

5 Comparison of the partitioning approaches

A comparison of $T/\mathrm{ET}$ across all the seasons for different types of vegetation is presented in Fig. 3. The averages were computed over half-hour periods obtained between 06:00 and 18:00 LT (local time). For all the vegetation types, the average $T/\mathrm{ET}$ ratios range from 0.5 in winter to 0.7 in summer. Except for the results over the tundra locations, there is a notable similarity across the five methods, particularly during spring and fall. In contrast, the dissimilarity across the methods for tundra is possibly linked to the heterogeneity of these sites, which breaks down the assumptions of the partitioning methods more often. A surprising result is that good agreement across the methods was observed regardless of the measurement height. Even above the agricultural sites, where the ratio $z/h$ can be as large as 30, good agreement was observed for many locations.

Given the underlying uncertainties with the partitioning algorithms, one route for verifying the reliability of partitioning estimates is to compare different methods. In Table 5, we compare the different transpiration estimates obtained by CEC, FVS, CEA, and CECw, where the slope α of the linear regression is shown. The results for MREA were similar to those for CEC and are not shown. To help interpret these results, they were separated into four groups: agreement in the 10 % range, i.e., $\mathrm{0.9}<\mathit{\alpha }<\mathrm{1.1}$; slopes in the range $\mathrm{0.8}<\mathit{\alpha }<\mathrm{1.2}$; slopes in the range $\mathrm{0.7}<\mathit{\alpha }<\mathrm{1.3}$; and disagreement greater than 30 %, i.e., α<0.7 or α>1.3.

Table 5 Slopes of the linear regression for transpiration for all the methods across all the sites. Only half-hour periods between 06:00 and 18:00 LT, without gap-filling, were included. For FVS and CECw, their ensemble average across the water use efficiency options was used. The sites were grouped according to the dominant vegetation type representing the tower footprint. To facilitate visual inspection of the results, slopes α in the 10 % range ($\mathrm{0.9}<\mathit{\alpha }<\mathrm{1.1}$) are shown in bold, while slopes in the 20 % range ($\mathrm{0.8}<\mathit{\alpha }<\mathrm{1.2}$) appear in italic.

Figure 3 $T/\mathrm{ET}$ ratios averaged over all the half-hour periods across seasons and sites with similar vegetation types. Only half-hour periods during the day (06:00 and 18:00 LT) were included, and only when solutions for all the methods were available. For FVS and CECw, the ensemble averages of outputs with different water use efficiencies were used. For each season or vegetation type, the black dashed line indicates the average across all the methods.

Overall, T_FVS and T_CEA are in very good agreement, with differences often within the 10 % range across various sites. Only at two sites (PUUM and BARR) were the slopes found to be greater than 1.3 (1.55 and 1.39, respectively). Similarly good agreement is seen for the pairs T_FVS and T_CECw as well as T_CEA and T_CECw. The close agreement for these three methods was also observed in their average diurnal cycles, in particular between 09:00 and 15:00 LT (Fig. A1), and even under conditions when $z/h\gg \mathrm{3}$. Nonetheless, note that the diurnal cycles for FVS exhibit higher values early in the morning and late in the afternoon. Previous studies have highlighted significant variability among the five W options, with differences that can span 2 orders of magnitude during the early hours (Wagle et al., 2023), which inevitably influences the output of the FVS method. Therefore, without additional information on the actual water use efficiency, we cannot explain this behavior with certainty.

The good agreement between FVS and CEA is particularly encouraging given that the CEA method does not require water use efficiency as an input, using only turbulence information. Moreover, we showed that the outputs of FVS using W_CR and W_CC were the ones that differed most from CEA, the remaining three FVS outputs being more similar; W_OPT produced the estimates closest to CEA (Fig. A2). Finally, transpiration estimates obtained by CEC tend to disagree more often with other methods (similar trends were observed for MREA; data not shown). Nonetheless, increased agreement is observed when the data are separated by seasons, as shown in Figs. 3 and A1.

The main challenge of flux partitioning is verifying the accuracy of these results, so we cannot initially assume that CEC's results are inherently inaccurate or incorrect. However, previous numerical experiments (Zahn et al., 2024) found that CEC and MREA are more restrictive than the other approaches, being stricter in terms of measurement height and canopy density. Thus, it is possible that the denser canopy during summer will directly affect these assumptions, decreasing these two methods' performance. They were also more prone to larger errors when the ratios $T/E$ and $P/R$ were more dissimilar. Unfortunately, a metric to determine the “goodness” of partitioning is not feasible since in experiments the true answer is never known. Only estimates of these partitioned fluxes are available since even the sap flow, isotope, and flux chamber methods have their limitations.

Given this uncertainty, a prudent recommendation is to implement the ensemble averages of T and E for all the methods in future studies. Ideally, the results should converge regardless of the partitioning methods used. However, depending on the specific NEON site or analyses of interest, it may be advantageous to prioritize the FVS, CEA, and CECw methods under certain circumstances given their close agreement, potentially adding either CEC or MREA since these two methods give similar results and should thus be treated as one member of the ensemble. To this end, the slopes shown in Table 5 offer a reference.

Table 6 Slopes of the linear regression for net primary productivity (P) for all the methods across all the sites. Only half-hour periods between 06:00 and 18:00 LT were included. For FVS and CECw, their ensemble averages across the water use efficiency options were used. The sites were grouped according to the dominant vegetation type. To facilitate visual inspection of the results, slopes α in the 10 % range ($\mathrm{0.9}<\mathit{\alpha }<\mathrm{1.1}$) are shown in bold, while slopes in the 20 % range ($\mathrm{0.8}<\mathit{\alpha }<\mathrm{1.2}$) appear in italic.

In contrast to ET flux partitioning, the estimated CO₂ flux components showed more disagreement across the methods (Table 6). Generally, methods that do not require water use efficiency (CEC, CEA, and MREA) underestimated P and R compared to FVS and CECw, where most slopes indicated less than 10 % agreement. Our previous findings using large-eddy simulations to test all the methods (Zahn et al., 2024) found larger errors for CO₂ components with CEC, CEA, and MREA. We found that, while CO₂ is essential as a tracer for H₂O in these methods' formulations, it cannot be reliably partitioned between soil and canopy components given the lack of an upper bound on F_c components, which combine a positive (respiration) flux and a negative (photosynthesis) flux. In contrast, FVS and CECw limit the magnitudes of these components by imposing $W=P/T$, thus adding a constraint to the magnitude of $F_{\mathrm{c}}=R+P$. Note, however, that P (and R) will be subjected to uncertainties present in the water use efficiency estimate. Overall, these findings indicate that CO₂ flux partitioning is more challenging, and care must be taken when using and interpreting these results. At present, for CO₂ flux components, we thus recommend taking the ensemble average of CECw and FVS (across all five W estimates), potentially adding CEA as a third member of the ensemble for the selected sites, in particular above deciduous forests following Table 6.

6 Exploring research opportunities with a flux-partitioning dataset

The final dataset will enable hydrologists and ecologists to investigate a variety of research questions related to flux partitioning. Considering the dataset's breadth of coverage, its spatial and temporal distributions enable comparison across diverse conditions, encompassing various forest structures, water availability, climate and weather patterns, and more. Due to remaining gaps in the final flux time series, evaluation of long-term budgets over the period of 5 years is not recommended without additional work to complete the temporal coverage. Nonetheless, budgets over smaller periods and for specific NEON sites are possible and can be explored. In this case, for nighttime periods, when most partitioning results are unavailable, researchers might fill gaps by assuming zero plant fluxes (e.g., for T and P), E=ET, and R=F_c.

One relevant research question regarding ecosystem functioning revolves around the main drivers of plant and soil fluxes. To briefly explore this question, we implemented XGBoost to estimate the importance of a few environmental variables for T and $T/\mathrm{ET}$, where we first averaged the output from all five partitioning models (excluding gap-filled periods). Note that this analysis is different from the procedure used for gap-filling, where ET was used as a predictor. Here, we simply investigate the importance of each environmental variable as a driver of plant fluxes, and thus ET is not used as an input. The selected features were air temperature, photosynthetically active radiation, relative humidity, vapor pressure deficit, mean wind speed, sensible heat flux, and friction velocity.

Figure 4 Importance of seven environmental variables as drivers of transpiration (top-row panels) and $T/\mathrm{ET}$ (bottom-row panels) across the vegetation types. From left to right, the legend includes air temperature (T_AIR), photosynthetically active radiation (PAR), relative humidity (RH), vapor pressure deficit (VPD), mean wind velocity (U), sensible heat flux (H), and friction velocity (u_*). The pie charts display the average over all sites with the same vegetation type, where bars represent the standard deviation of feature importance across the sites ( %). As a reference, the outer circle represents a standard deviation of 20 %.

The importance of each environmental variable is shown in Fig. 4 for T (top row) and $T/\mathrm{ET}$ (bottom row) averaged across sites with the same vegetation type. The standard deviation was also computed across the sites and is represented by the bars in the figure (the outer circle represents a standard deviation of 20 %). Notably, these seven environmental variables can explain up to 77 % (R²=0.77 over the agricultural sites) of the variability in transpiration. With the exception of the tundra sites, where sensible heat was the most important variable, photosynthetically active radiation, followed by air temperature, was the most relevant feature. Intersite variability was also observed: for instance, the standard deviation of PAR and T_AIR was as large as 20 % for evergreen and mixed forests. The variability across the methods is not included in this analysis since the ensemble average was used; nonetheless, similar to our discussion regarding gap-filling, we observed that FVS combined with W_CC and W_CR often differed from the other methods (results not shown). This discrepancy is also left to future studies.

The ratio $T/\mathrm{ET}$, on the other hand, is more complex and dependent on additional variables not included in this analysis. Across the vegetation types, these seven features were able to explain as little as 23 % (shrublands) or no more than 48 % (deciduous forests). While relevant environmental variables such as soil moisture and accumulated rainfall were not included in our analyses due to availability and/or data quality, these variables can be included for selected sites and/or smaller temporal records. By exploring different environmental forcings, this dataset presents an opportunity to address specific inquiries, such as the impact of heat waves, droughts, and other extreme events on plant transpiration. Research endeavors of this nature are vital, as they contribute to a deeper comprehension of how present and future climate conditions affect local water and energy budgets, together with plant water use efficiency (Hatfield and Dold, 2019).

The link between fluxes and canopy structure features another set of critical scientific questions that can be addressed using this dataset. Previous work (Wang et al., 2014; Wei et al., 2017) compiled results from several case studies and summarized the correlation between $T/\mathrm{ET}$ and leaf area index (LAI), which has been turned into a simple regression model often used in other studies (Nelson et al., 2020). Future investigations could expand this study and utilize NEON products or satellite imagery to derive seasonal variations in LAI, thereby enhancing our understanding of how vegetation dynamics influence transpiration. By combining the most important meteorological features and vegetation characteristics, simpler and more accurate partitioning models can be derived and applied more broadly.

In addition to observational studies, modeling or hybrid studies (Rafi et al., 2019; Kozii et al., 2020) can also benefit from the dataset compiled here. This approach can be explored to develop models for transpiration partitioning or modeling (Bright et al., 2022), where observations can validate the model or help determine relevant parameters, such as plant conductance. For instance, Schreiner‐McGraw et al. (2022) used transpiration data from FVS to optimize parameters for ParFlow.CLM, a hydrological model coupling groundwater flow to a land surface model (Maxwell and Miller, 2005). This methodology is particularly useful for evaluating the performance of land surface models, such as NOAH-MP (Niu et al., 2011; Li et al., 2021) and the Community Land Model (Lawrence et al., 2007), thereby improving assessments of their impact on weather prediction and climate projections (Berg and Sheffield, 2019; Dong et al., 2022).

Future studies can also expand this dataset by including additional partitioning models, including the products derived from the CO₂ partitioning algorithm proposed by Reichstein et al. (2005) and Lasslop et al. (2010). In terms of ET partitioning, for instance, different methods have been proposed that take GPP as an input (Perez-Priego et al., 2018; Zhou et al., 2016; Nelson et al., 2018); such models were compared by Nelson et al. (2020) across FLUXNET sites. By comparing different algorithms, we can further explore their uncertainties and focus on model improvement. Finally, as more data become available, other options can be used to train machine learning algorithms, focusing on gap-filling.

7 Code and data availability

The dataset is available at https://doi.org/10.5281/zenodo.12191876 (Zahn and Bou-Zeid, 2024). In addition to all the flux components, it contains the auxiliary meteorological inputs used to implement the Extreme Gradient Boosting algorithm for gap-filling and feature importance analysis. The scripts used to implement all five partitioning methods can be found at https://doi.org/10.5281/zenodo.13244979 (Zahn, 2024).

8 Conclusions

This paper described the first dataset of flux partitioning across all NEON eddy covariance sites from 2019 to 2023. Estimates of transpiration, evaporation, soil respiration, and plant net CO₂ assimilation were obtained using five partitioning approaches. These methods had undergone extensive testing in both experimental and numerical settings in prior studies (Zahn et al., 2022, 2024), aiding in understanding their limitations and reliability.

Overall, ET components exhibited notable consistency across the methods, which was particularly evident in FVS, CEA, and CECw, thereby enhancing confidence in these estimates. In contrast, the CO₂ components are more challenging, and care must be taken when interpreting these estimates. The final dataset, available at https://doi.org/10.5281/zenodo.12191876 (Zahn and Bou-Zeid, 2024), can be used to address numerous research questions. We recommend utilizing an ensemble average of the methods available for a particular site. However, based on the comparison presented in this study, a combination of the methods might be adopted for specific analyses or NEON sites. In any case, the user should keep in mind the following considerations:

1. None of the methods can estimate plant respiration and gross photosynthesis separately: only an estimate of their difference and plant net photosynthesis (i.e., plant net fluxes) is computed. Similarly, below-canopy respiration – which we assume to be mostly from the soil – also includes respiration from roots and litter.

2. Evaporation from canopy-intercepted water is neglected.

3. While nighttime fluxes are also included, we do not recommend the inclusion of nighttime periods, since CO₂ assimilation is not expected at night for C₃ and C₄ plants, thus breaking the assumption of CO₂–H₂O coupling.

4. Uncertainties in W will be transferred to the partitioning estimates for FVS and CECw.

5. Given the differences in the method frameworks and their assumptions, the data gaps vary. Thus, selecting the periods of interest and the methods to average must be done with care to avoid biases.

This dataset opens up the many research pathways that we outlined above and that can be the basis for other derived datasets. The partitioning approaches can also be applied to other periods of NEON observations or implemented in the NEON workflow so that the partitioned fluxes are part of the standard NEON outputs.

Appendix A: Diurnal cycle of $T/\mathrm{ET}$

Figure A1 Diurnal cycle of the partitioning ratio $T/\mathrm{ET}$ for all the methods across season and vegetation type. Half-hourly estimates during the daytime hours (06:00 to 18:00 LT) are only included when all the methods have valid solutions. For FVS and CECw, the ensemble average across the water use efficiency options was used.

Figure A2 Diurnal cycle of the partitioning ratio $T/\mathrm{ET}$ comparing the FVS output from different water use efficiency models (and their ensemble average) and the CEA method (dashed line). Half-hourly estimates during the daytime hours (06:00 to 18:00 LT) are only included when all the methods have valid solutions.