The effect of rainfall amount and timing on annual transpiration in a grazed savanna grassland

The role of precipitation (P) variability with respect to evapotranspiration (ET) and its two components, transpiration (T) and evaporation (E), from savannas continues to draw significant research interest given its relevance to a number of ecohydrological applications. Our study reports on 6 years of measured ET and estimated T and E from a grazed savanna grassland at Welgegund, South Africa. Annual P varied significantly with respect to amount (508 to 672 mm yr⁻¹), with dry years characterized by infrequent early-season rainfall. T was determined using annual water-use efficiency and gross primary production estimates derived from eddy-covariance measurements of latent heat flux and net ecosystem CO₂ exchange rates. The computed annual T for the 4 wet years with frequent early wet-season rainfall was nearly constant, 326±19 mm yr⁻¹ ($T/\mathrm{ET}=\mathrm{0.51}$), but was lower and more variable between the 2 dry years (255 and 154 mm yr⁻¹, respectively). Annual T and $T/\mathrm{ET}$ were linearly related to the early wet-season storm frequency. The constancy of annual T during wet years is explained by the moderate water stress of C₄ grasses as well as trees' ability to use water from deeper layers. During extreme drought, grasses respond to water availability with a dieback–regrowth pattern, reducing leaf area and transpiration and, thus, increasing the proportion of transpiration contributed by trees. The works suggest that the early-season P distribution explains the interannual variability in T, which should be considered when managing grazing and fodder production in these grasslands.

1 Introduction

Similar to other semiarid areas, wooded grasslands in central South Africa deliver essential ecosystem services such as grazing land and fodder (Bengtsson et al., 2019). In such semiarid zones, evapotranspiration (ET) approximately matches annual precipitation (P≈500 mm yr⁻¹; Zhang et al., 2001). The transpiration (T) component accounts for water loss from the leaf stomata of the sparse tree component, seasonal grasses, and the minor forb component. The evaporation (E) component is large following rain events, as intercepted water and near-surface soil water evaporate; the latter may continue over periods longer than a week (Perez-Priego et al., 2018). The partitioning of ET between E and T may affect the net radiation (R_n) and surface temperature on short timescales (sub-daily). However, the processes that increase the proportion (and amount) of water used in T, facilitating greater carbon uptake and subsequent fodder production for cattle, occur over timescales of weeks or longer. Given the link between T and carbon uptake from the atmosphere, there is growing interest in how ET is partitioned into E and T in semiarid ecosystems (Merbold et al., 2009; Sankaran et al., 2004; Scanlon et al., 2002, 2005; Scholes and Archer, 1997; Scott and Biederman, 2017; Volder et al., 2013; Williams and Albertson, 2004; Xu et al., 2015; Yu and D'Odorico, 2015). The aim here is to explore this partitioning of ET using a long-term data set of measured fluxes of energy, water, carbon dioxide (CO₂), and vegetation activity from a grazed wooded grassland. The focus is restricted to processes operating over daily to seasonal timescales, commensurate with controls over the annual partitioning of P into T, and the resulting carbon uptake in gross primary production (GPP). These longer timescales are of interest in the valuation of ecosystem productivity and their services when assessing climatic shifts (Godde et al., 2020). The results presented here on the partitioning of ET must be viewed as necessary but insufficient for developing best practices for the management of grazing or fodder production.

The contrasting vegetation layers of wooded grasslands have distinct seasonal dynamics of leaf area and physiological activity. In southern African drylands, Vachellia erioloba (camel thorn tree) is an important woody component of grasslands (Barnes et al., 1997); this species is a deep-rooted semi-deciduous tree with a low leaf turnover rate, resulting in minor leaf area changes. Furthermore, V. erioloba has been shown to absorb 37 % of its water below a depth of 1 m, partially decoupling its physiological activity from recent precipitation and shallow soil water content (Beyer et al., 2018). Little interannual variation in tree water use has been shown in many semiarid ecosystems (Do et al., 2008; Hutley et al., 2001; Montaldo et al., 2020). However, pronounced seasonality in tree transpiration of semi-deciduous and deciduous savanna tree species has been observed in South Africa, even in riparian forest (Scott-Shaw and Everson, 2019). Perennial C₄ grass species in the area have shallower rooting systems and are physiologically responsive to intermittent rainfall events (Sankaran, 2019). The grass dependence on the temporal distribution of P has been demonstrated by a positive correlation between rainfall frequency and productivity at the field scale (Swemmer et al., 2007) and by a positive relationship between wet-season rainfall frequency and grass cover in the grassland ecosystems of sub-Saharan Africa (P≤630 mm yr⁻¹; D'Onofrio et al., 2019). Compared with C₃ trees, the C₄ grass layer has CO₂ concentrated in the bundle sheath that enables greater light- and water-use efficiencies of CO₂ uptake in the warmer inter-canopy spaces (Ripley et al., 2010). Due to their ability to regulate intercellular CO₂ concentrations, C₄ grasses have a higher photosynthetic ability per unit leaf area that can be sustained even in moderate-water-stress situations (Taylor et al., 2014). However, for the same volume of soil, a C₄ grass is an intensive and fast user of soil water when compared with C₃ trees. Indeed, due to their shallow rooting depth, severe droughts may alter both their water-use efficiency per unit leaf area and their leaf area dynamics. Thus, our study objective is not only to partition measured ET into T and E but also to quantify the effect of environmental variables on the seasonality of the grass activity.

Three methods that link T to GPP are used to estimate monthly $T/\mathrm{ET}$. These methods were chosen because previous applications have shown some success in partitioning ET into E and T when applied to multi-site data sets. These three methods provide an estimate of ecosystem-scale T, albeit with differing assumptions and uncertainty (Stoy et al., 2019). Previous method comparisons have shown that the linear regression method with optimality assumption produces lower $T/\mathrm{ET}$ estimates than the machine learning approach to ET partitioning for grassland and savanna ecosystems (Nelson et al., 2020; Scott et al., 2021). Comparing these methods allows for the selection of the most suitable partitioning scheme for water-limited ecosystems in general and savannas in particular. Furthermore, an agreement between the methods lends confidence to the estimates of $T/\mathrm{ET}$ and the drivers of T (e.g., precipitation). Disagreements between the methods may also identify potential uncertainties in the hydroclimatic or land-cover conditions explored here. Hence, a corollary goal is to understand the opportunities and limitations of these methods when combined with eddy-covariance-measured ET, such as those supplied by FLUXNET (Baldocchi et al., 2001).

The main question addressed here is how T and $T/\mathrm{ET}$ vary with P at monthly and annual timescales in a grazed savanna grassland ecosystem. Available MODIS vegetation indexes and Landsat 8 grass and tree normalized difference vegetation index (NDVI) information allow for the quantification of vegetation dynamics at the site. Soil-moisture-based grass transpiration and ecosystem-scale T offer a new perspective on the relation between water fluxes at seasonal and annual timescales as well as new ways to examine the role of grass and trees in water budgets. Our study objectives are (i) to quantify the variation in annual P, ET, and T; (ii) identify the main drivers of the annual and monthly T and $T/\mathrm{ET}$; and (iii) relate the growth dynamics of tree and grass components to the hydrological balance.

2 Materials and methods

2.1 Site description

The Welgegund measurement site is located in a grazed savanna grassland in South Africa (26^∘34^′10^′′ S, 26^∘56^′2^′′ E; 1480 m a.s.l.), as shown in Fig. 1. The research site is part of a large-scale commercial ranch with an annual cattle head count of 1300±300. The cattle grazing area is approximately 6000 ha.

Figure 1 (a) Photograph of the measurement station and surrounding vegetation taken on 22 August 2016. (b) Daytime flux footprint estimated with measured flow statistics and a 2D footprint model (Kljun et al., 2015). Footprint contour lines (black) are shown in 10 % increments from 10 % to 90 %. The station is located at the center, and the green cross indicates the farm well. The Landsat 8 NDVI was calculated from the area marked by the blue square for grass and by the green square for trees.

The area experiences two seasonal periods: a warm rainy season from October to April and a cool dry season from May to September. The 16-year mean annual rainfall determined at a nearby weather station (town of Potchefstroom) was 540 mm yr⁻¹ ± 112 mm yr⁻¹ (Räsänen et al., 2017). The soil around the site is loamy sand in the top 1 m. Although the water table depth is not known, the farm well has a continuous water supply at 30 m below the surface (Fig. 1).

The vegetation in the area is an open thornveld. Eragrostis trichophora, Panicum maximum, and Setaria sphacelata are the dominant perennial C₄ grass species. The mean maximum grass height across sampling plots was 0.1 m in 2011 (Räsänen et al., 2017). Tree cover is 15 %, and the dominant tree species is V. erioloba, with other less-prominent species such as Celtis africana and Searsia pyroides. Dicoma tomentosa, Hermannia depressa, Pentzia globosa, and Selago densiflora are the dominant forb species. Details on the site and vegetation cover may be found elsewhere (Jaars et al., 2016; Räsänen et al., 2017).

2.2 Measurements

Atmospheric aerosols, trace gases, and meteorological variables were measured continuously at the site (Beukes et al., 2015; Petäjä et al., 2013). The measurements directly related to energy fluxes and water balance are briefly described. The eddy-covariance (EC) system consisted of a triaxial sonic anemometer (METEK USA-1) and a LI-COR (LI-7000) closed-path infrared gas analyzer, which were positioned 9 m a.g.l. (meters above ground level). The sampling frequency of the EC system was 10 Hz. The gas analyzer was calibrated every month with a high-precision CO₂ span gas using synthetic air with CO₂ < 0.5 ppm as a reference gas. The meteorological measurements included atmospheric temperature and pressure, mean wind speed and direction, and mean air relative humidity. The radiation measurements were made using Kipp & Zonen PAR-lite sensors, CMP-3 pyranometers, and an NR Lite2 net radiometer positioned 3 m a.g.l. with a field of view at the grass level. These sensors measure photosynthetically active radiation (PAR), direct and reflected global radiation, and net radiation. The soil surface heat flux was measured with a Hukseflux HFP01 heat flux plate at 5 cm below the soil surface. The meteorological variables were sampled every 1 min (radiation every 10 s), and 15 min averages were then recorded.

Precipitation was measured 1.5 m a.g.l. with two tipping bucket gauges (Vaisala QMR102 and Casella). Most of the precipitation values reported here are from the Vaisala gauge, with the Casella gauge only used to gap-fill missing values in the Vaisala gauge time series. No significant differences were observed between the two sensors. Between December 2011 and February 2012, the measured rainfall was underestimated due to the high intensity of the rainfall, and it was corrected using nearby weather station measurements (Sect. S1, Figs. S1 and S2 in the Supplement). Wind-induced underestimation is a well-known problem with pointwise rainfall measurements. Thus, the measured precipitation was corrected by multiplying it by 1.094. This corresponds to the 9.4 % bias that was determined for the Casella rain gauge at the height of 0.5 m at a measurement site with a similar mean wind speed (5 m s⁻¹) and annual rainfall (P=700–1000 mm yr⁻¹) (Pollock et al., 2018).

Stored soil water changes were determined using two separate soil moisture profiles. The measurements of individual soil moisture sensors at depths of 0.05, 0.2, and 0.5 m (Delta-T ML2) were converted to a single average soil moisture value using the weights of 125, 225, and 200 mm, respectively. These soil moisture measurements covered the complete experiment period from September 2010 to August 2016. Starting from March 2012 onwards, a Delta-T PR$\mathrm{2}/\mathrm{6}$ probe was installed to record soil moisture at depths of 0.1, 0.2, 0.3, 0.4, 0.6, and 1.0 m. This profile measurement was converted to stored soil water using the weights 150, 100, 100, 150, 300, and 200 mm, respectively, and it was used to calculate the annual change in soil water storage (ΔΘ_1 m) over the entire 1 m soil column.

The site was visited once or twice a week during the 6-year period to check the status of the sensors and correct errors if necessary. Measurement records were used to identify anomalies, outliers, or erroneous measurement periods. Further details on the site and EC measurements are presented elsewhere (Aurela et al., 2009; Räsänen et al., 2017). The annual energy balance closure was also computed, which varied from 0.75 to 0.85. This lack of closure is comparable to those reported for FLUXNET sites (Stoy et al., 2013; Wilson et al., 2002). Given the heterogeneity in vegetation cover and that EC measurements sense a different footprint from the footprint representing the difference between net radiation and soil heat flux, no Bowen ratio adjustments were performed to force an energy balance closure.

2.3 Flux calculation and gap filling

The details on the turbulent flux calculations are presented in Räsänen et al. (2017). Briefly, the turbulent fluxes were calculated as 30 min block averages after double rotation and by applying the Webb–Pearman–Leuning (WPL) density correction (Webb et al., 1980). The low-frequency flux correction was performed according to Moore (1986), and high-frequency losses were corrected using empirical transfer functions determined using sensible heat flux as a reference scalar. The sensible and latent heat flux values were discarded when the measured friction velocity u_* was below 0.28 m s⁻¹, which was deemed as a state of low turbulence mixing. The steady-state test of Foken and Wichura (1996) was used to screen the latent heat flux data for nonstationary conditions within each 30 min averaging period. The data were discarded if the relative nonstationarity defined by this test exceeded a threshold, which was set to 30 % and 100 % for the data used for gap filling and final analysis, respectively. Latent heat fluxes were checked for an acceptable H₂O concentration range and variance to detect anomalous spikes due to condensation or rainfall. Heat flux values were filtered for outliers by considering values for each month of all of the measurement years and removing outliers using an adjusted boxplot (Hubert and Vandervieren, 2008). The steady-state check resulted in less than 30 % filtered fluxes, which were gap filled using marginal distribution sampling (MDS) from the REddyProc package (Reichstein et al., 2005). Both daytime and nighttime fluxes were gap filled using this approach, given the significant role that nighttime evaporation and respiration play in the water and carbon balances. Gap filling of nighttime evaporation is of significance at Welgegund, as most of the rainfall occurs in the late afternoon and early evening. The meteorological parameters were also gap filled using the MDS approach (Reichstein et al., 2005).

The flux footprint was estimated using the daytime measured flow statistics for the 6-year period and a standard 2D footprint model (Kljun et al., 2015), which are presented in Fig. 1. These calculations suggest that 80 % of the ET fluxes originate from the homogeneous thornveld.

The EC-inferred GPP was used to derive the water-use efficiencies to partition measured ET into T and E. The measured net ecosystem CO₂ exchange (NEE) was partitioned into GPP and ecosystem respiration using nighttime mean respiration values. These values were assumed to be the same for daytime respiration, and GPP was determined as the difference between NEE and daytime ecosystem respiration. Nighttime mean respiration was used instead of the exponential temperature function, as only 2 % of the fitting windows had a linear or exponential relation between EC-based ecosystem respiration and soil temperature. The difference between the mean monthly transpiration from these two methods was small, with transpiration from the exponential temperature function being 4 % higher than transpiration from the nighttime mean method (Fig. S3). The GPP fit parameters and the nighttime mean respiration were calculated in a moving data window that was defined for each day with an initial length of 6 d. The moving window was expanded by up to 20 d if necessary, to include at least 50 measurement points. The measured NEE had one large 25 d gap in September 2013, and the fit parameters were linearly interpolated for this gap. The preprocessing of NEE was performed with the same filters as the heat fluxes, as discussed in Räsänen et al. (2017).

The potential ET (PET) was calculated using the Priestley–Taylor formulation given by Priestley and Taylor (1972):

$$\begin{array}{}\text{(1)}& \mathrm{PET}={\mathit{\alpha }}_{\mathrm{PT}}{\displaystyle \frac{\mathrm{\Delta }}{\mathrm{\Delta }+{\mathit{\gamma }}_{\mathrm{p}}}}\left(R_{\mathrm{n}}-G\right),\end{array}$$

where α_PT=1.26 is the Priestley–Taylor coefficient, $\mathrm{\Delta }=\mathrm{d}{e}^{*}/\mathrm{d}{T}_{\mathrm{a}}$ (Pa K⁻¹), e^* (Pa) is the saturation vapor pressure given by the Clausius–Clapeyron equation and evaluated at the measured air temperature T_a (K), R_n (W m⁻²) is net radiation, G (W m⁻²) is soil heat flux, and γ_p (Pa K⁻¹) is the psychrometric constant. The energy balance closure (EBC) slope was estimated for each year by regressing all measured 30 min values of R_n–G against the sum of the measured latent and sensible heat fluxes for the same period.

2.4 Uncertainty of annual ET estimates

The friction velocity (u^*) threshold was estimated using a bootstrapping technique from 200 artificial replicates of the data set (Wutzler et al., 2018). The mean u^* estimate value for the whole data set was 0.28 m s⁻¹, with heat flux and NEE values being discarded when u^* was lower than this limit. The 5th, 50th, and 95th percentiles of the estimates were 0.27, 0.29, and 0.32 m s⁻¹, respectively. The data set was u^* filtered and gap filled with these three u^* limits. The annual u^* uncertainty range was calculated for each k year as follows:

$$\begin{array}{}\text{(2)}& E_{u^{*},k}={\displaystyle \frac{{\mathrm{ET}}_{\max ,k}-{\mathrm{ET}}_{\min ,k}}{{\mathrm{ET}}_{\mathrm{median},k}}}{\mathrm{ET}}_k,\end{array}$$

where $E_{u^{*},k}$ is the u^* uncertainty for year k, and ET_k is the evapotranspiration for year k.

The MDS gap-filling algorithm estimates random error for each 30 min value based on the standard deviation of the observed latent heat flux with similar meteorological conditions in a moving window. The annual random error was estimated as the root-mean-squared error:

$$\begin{array}{}\text{(3)}& E_{\mathrm{rand},k}=\sqrt{\sum _{i=\mathrm{1}}^{n_k}{\mathit{\sigma }}_i^{\mathrm{2}}},\end{array}$$

where n_k is the number of 30 min periods in year k, and σ_i is the standard deviation of latent heat flux from the MDS gap-filling algorithm. The total uncertainty of the annual ET was calculated by adding the random error and u^* uncertainty in quadrature to yield

$$\begin{array}{}\text{(4)}& E_{\mathrm{tot},k}=\sqrt{E_{u^{*},k}^{\mathrm{2}}+E_{\mathrm{rand},k}^{\mathrm{2}}}.\end{array}$$

2.5 Rainfall interception

The total rainfall interception (I_t) was not measured; instead, it was estimated by modeling grass, litter, and tree interception. The interception was estimated for each storm using discrete rainfall events separated by at least 1 h. The grass interception for one storm event was calculated using an expression derived for crops (Moene and Van Dam, 2014):

$$\begin{array}{}\text{(5)}& I_{\mathrm{grass}}=a\mathrm{LAI}\left(\mathrm{1}-{\displaystyle \frac{\mathrm{1}}{\mathrm{1}+\frac{c_{\mathrm{g}}{P}_{\mathrm{g}}}{a\mathrm{LAI}}}}\right),\end{array}$$

where a is the scale parameter; c_g is the grass cover fraction; LAI is the leaf area index, estimated here from satellite (Sect. 2.7); and P_g is the rainfall amount per storm. The scale parameter was set to 0.5 mm per event, which corresponds to a maximal 1 mm interception loss for LAI = 2. The grass cover fraction was estimated using LAI, as follows:

$$\begin{array}{}\text{(6)}& c_{\mathrm{g}}=\mathrm{1}-e^{-k\mathrm{LAI}},\end{array}$$

where the extinction coefficient (k) is set to 0.4. Tree interception was estimated using the revised model for a sparse canopy (Gash et al., 1995). The model assumes that rainfall events consist of wetting, saturation, and drying phases. The interception for small events that do not saturate the canopy was estimated separately from large storms that saturate the canopy. The rainfall to fill the canopy storage is

$$\begin{array}{}\text{(7)}& P_{\mathrm{g}}^{\prime }=-S_{\mathrm{c}}\left({\displaystyle \frac{\stackrel{\mathrm{‾}}{R}}{\stackrel{\mathrm{‾}}{E}}}\right)\ln \left(\mathrm{1}-\left({\displaystyle \frac{\stackrel{\mathrm{‾}}{E}}{\stackrel{\mathrm{‾}}{R}}}\right)\right),\end{array}$$

where $S_{\mathrm{c}}=S/c_{\mathrm{t}}$ is the canopy storage capacity per unit cover, $\stackrel{\mathrm{‾}}{R}$ is the mean rainfall, and $\stackrel{\mathrm{‾}}{E}$ is the mean evaporation rate during a storm. The measured ET was used to calculate the mean evaporation rate for each event. The tree cover fraction c_t was set to a constant 0.15, and storage capacity S was set to 1.07 mm, corresponding to a measured value for Acacia mearnsii (Bulcock and Jewitt, 2012). The total tree interception is then determined as follows:

$$\begin{array}{}\text{(8)}& I_{\mathrm{tree}}=c_{\mathrm{t}}\sum _{i=\mathrm{1}}^m{P}_{\mathrm{g},i}+\sum _{j=\mathrm{1}}^n\left[c_{\mathrm{t}}{P}_{\mathrm{g}}^{\prime }+c_{\mathrm{t}}{\displaystyle \frac{\stackrel{\mathrm{‾}}{E}}{\stackrel{\mathrm{‾}}{R}}}\left(P_{\mathrm{g},j}-P_{\mathrm{g}}^{\prime }\right)\right].\end{array}$$

The first sum accounts for the m small events that do not saturate the canopy, and the second sum accounts for the n large events. The litter interception was assumed to be 1 mm per rainfall event (Scholes and Walker, 1993), and it was multiplied by c_t.

2.6 Partitioning ET

Prior to presenting the three ET partitioning approaches, the link between GPP and T is reviewed. From definitions, the flux-based water-use efficiency (WUE) is as follows:

$$\begin{array}{}\text{(9)}& \mathrm{WUE}={\displaystyle \frac{\mathrm{GPP}}{T}}\propto {\displaystyle \frac{c_{\mathrm{a}}}{\mathrm{VPD}}}\left(\mathrm{1}-c_{\mathrm{i}}/c_{\mathrm{a}}\right),\end{array}$$

where c_i and c_a are the respective intercellular and ambient atmospheric CO₂ concentrations, and VPD is the vapor pressure deficit. Based on stomatal optimization theories that maximize carbon gain for a given amount of water loss in the rooting system per unit leaf area, the ratio of CO₂ concentrations ($\mathrm{1}-c_{\mathrm{i}}/c_{\mathrm{a}}$) is proportional to $\sqrt{\mathrm{VPD}}$, as demonstrated in several studies reviewed elsewhere (Hari et al., 2000; Katul et al., 2009, 2010). Combining these theories with the definition of WUE (Eq. 9) makes T proportional to GPP×VPD^0.5 provided that c_a does not vary appreciably. The proportionality constant in this expression ($T\propto \mathrm{GPP}\times {\mathrm{VPD}}^{\mathrm{0.5}}$) is linked to the so-called marginal water-use efficiency (or the Lagrange multiplier in optimal stomatal control theories), which differs from the intrinsic water-use efficiency $\mathrm{iWUE}=(\mathrm{1}-c_{\mathrm{i}}/c_{\mathrm{a}})c_{\mathrm{a}}$. It must be externally supplied or determined from EC measurements during conditions when T approximately equals ET. When this proportionality constant is known, an EC-based GPP estimate (together with VPD) can be used to infer T and, via subtraction from ET, produce an estimate of E.

Table 1 Summary of the ET partitioning methods applied to the 6 years of data.

Three approaches were used to divide ET into E and T using the abovementioned link between GPP and T (Table 1). The first method was presented by Berkelhammer et al. (2016) (hereafter, B16). Here, it was applied to each year individually to allow for the large interannual variation in vegetation phenology. The method assumes that ET is linearly related to GPP×VPD^0.5 only when T is the dominant term in ET. Moreover, the $T/\mathrm{ET}$ ratio is assumed to approach unity intermittently. To estimate the $T/\mathrm{ET}$ value for each 30 min period, the product GPP×VPD^0.5 was plotted against ET for each year, and the minimum value of ET was then selected as the 5th percentile for each equal-sized GPP×VPD^0.5 bin. The bin was defined by discretizing the 30 min GPP×VPD^0.5 values into 50 bins, with each containing the same number of measurements but encompassing different value ranges, for reasons provided elsewhere (Berkelhammer et al., 2016). A linear regression line of these bins defines the ET value for which T dominates ET. Any value falling below the line is considered to have $T/\mathrm{ET}=\mathrm{1}$. For points above the regression line, $T/\mathrm{ET}$ is defined as the ratio between the minimum ET that represents T and the observed ET:

$$\begin{array}{}\text{(10)}& {\displaystyle \frac{T}{\mathrm{ET}}}={\displaystyle \frac{{\min }_{\mathrm{GPP}}\Vert \mathrm{ET}\Vert }{{\mathrm{ET}}_{\mathrm{flux}}}},\end{array}$$

where ${\min }_{\mathrm{GPP}}\Vert \mathrm{ET}\Vert$ is the minimum ET value, and ET_flux is the observed ET value. The calculation of 30 min $T/\mathrm{ET}$ values for 1 year is illustrated in Fig. S4. The daily $T/\mathrm{ET}$ values were calculated by taking the mean of these 30 min $T/\mathrm{ET}$ values, and it was used to calculate daily T and E using the measured ET. The regression slope and intercept of the T=ET line are related to the inverse of water-use efficiency for each year based on the 30 min data. The 30 min data points used for the $T/\mathrm{ET}$ estimation were also filtered with additional quality criteria, i.e., only data points with measured ET, positive GPP, and R_n were used (see Zhou et al., 2016). However, rainy days were included in the estimation to capture the rainfall interception events measured by the EC system while maintaining the data-stationarity filter. Shortly after rain, water droplets remaining on the sonic anemometer transducers can block the detection of sound waves emitted and received, leading to anomalous vertical velocity and friction velocity measurements for these 30 min runs. However, as the sonic anemometer transducers are inclined, smooth, and have small surface area, they dry out faster than the leaves, thereby allowing the EC system to operate shortly after each rainfall event. At an annual scale, the estimated I_t was used to calculate soil evaporation (E_s) by subtracting I_t from E.

As previously mentioned, two other methods were also used to estimate T in order to identify the method most appropriate to water-limited ecosystems. The second approach entailed fitting the T=ET line using quantile regression with zero intercept for each year; the slopes of these fitted lines are termed uWUEp (Zhou et al., 2016) (hereafter, Z16). The apparent uWUE slopes (uWUEa) were defined for each day separately by fitting the 30 min ET values to GPP×VPD^0.5 values using linear regression with zero intercept. The daily $T/\mathrm{ET}$ value is the ratio of the uWUEa slope and uWUEp of each year. The difference between the B16 method and the uWUE method is primarily in the process of fitting the T=ET line.

The third approach is a random forest regressor that first isolates the most likely periods when T is equal to ET, after which it trains on GPP and T relations during these periods to infer T from measured GPP (Nelson et al., 2018) (hereafter, N18). A summary of these methods and their requirements is featured in Table 1 for convenience. In addition, an ET partitioning method that does not assume equality between T and ET during any periods was tested, but it could not be applied at this site because the monthly multiyear correlations between ET and GPP were not significant (Fig. S5) (Scott and Biederman, 2017).

2.7 Stage-2 soil evaporation

The estimated daily E was assessed using the stage-2 soil evaporation theory each year during the early dry season. During stage-2 conditions of soil evaporation, the evaporation is controlled by the soil moisture and soil physical properties (desorptivity) (Brutsaert and Chen, 1995; Hu and Lei, 2021). After the rainfall event, the cumulative daily E can be expressed as $D_{\mathrm{e}}{t}_{\mathrm{d}}^{\mathrm{1}/\mathrm{2}}$, where D_e is the soil desorptivity to be determined, and t_d is the dry-down duration in days. By regressing cumulative daily E inferred from the aforementioned partitioning methods upon $\sqrt{t_{\mathrm{d}}}$ for a single dry-down period, the D_e can be computed and compared to literature values. The rainfall events were chosen from the end of April onwards, and all rainfall events were higher than 10 mm. The dry-down periods varied from 12 to 30 d. The expected range of D_e based on several experiments is about 3 to 6 mm d${}^{-\mathrm{1}/\mathrm{2}}$ for sandy soils (Brutsaert and Chen, 1995).

2.8 Grass transpiration

The grass transpiration (T_g) was estimated from step-shaped diurnal changes in soil moisture for the last 4 years (Jackisch et al., 2020). The method estimates root water uptake from the daily change in root zone soil moisture, and it was applied to each soil moisture measurement from 10 to 60 cm. Total daily grass transpiration is the sum of daily changes in each soil moisture layer. The algorithm considers only periods when soil moisture measurement has characteristic root water uptake decline and excludes periods of percolation (Jackisch et al., 2020). Twelve percent of the estimated T_g values exceeded daily ET and were replaced with the daily ET values. The daily grass transpiration was set to zero during the dry season, defined by a MODIS LAI value of less than 0.3. The portion of soil evaporation was removed from the grass transpiration estimate by removing a constant 0.28 mm d⁻¹, which is the mean daily dry-season evaporation from the B16 method, from the estimated daily grass transpiration. This constant value was used instead of the daily E from ET partitioning methods because the EC-estimated daily E values were unrealistically high compared with the estimated grass transpiration. This is because the soil moisture measurements may not fully capture the soil evaporation of the EC footprint.

2.9 Satellite data

Changes in vegetation cover were quantified using the monthly average of the MODIS 16 d enhanced vegetation index (EVI) with a 250 m spatial resolution (MOD13Q1, Collection 6) (Didan, 2015). The monthly average of MODIS 8 d LAI (MOD15A2H, Collection 6) with a 500 m spatial resolution was used to relate monthly $T/\mathrm{ET}$ to LAI, comparing estimated $T/\mathrm{ET}$ to variations in vegetation phenology. The EVI signal is a ratio of spectral bands, whereas the LAI has corrected units of foliage area per ground area. For the last 3 years, Landsat 8 L2 16 d NDVI with a 30 m spatial resolution was used to determine separate grass and tree NDVI for the EC footprint (Fig. 1). Cloud- and cloud-shadow-affected Landsat pixels were removed using automated cloud cover identification (Braaten et al., 2015).

2.10 Rainy-season timing and green-up dates

Rainy-season length (T_wet) was estimated based on a climatological threshold of 5 % of the mean annual rainfall (Guan et al., 2014). The start of the rainy season was defined as the day when cumulative rainfall of the hydrological year (September to August) reached the threshold value of 27 mm, which was based on the long-term mean annual rainfall (540 mm yr⁻¹). Similarly, the end of the rainy season was estimated as the first day when cumulative rainfall, starting backward from the end of the hydrological year (August), reached the same threshold value. Early wet-season (September to November) precipitation was characterized by estimating the mean daily rainfall statistics using the daily mean precipitation amount (α) and daily mean storm frequency (λ). The daily mean precipitation amount was calculated as the mean precipitation of rainy days, and the mean storm frequency was calculated as the inverse of the mean time between rainy days.

The tree green-up date, estimated from the raw 16 d EVI time series (Archibald and Scholes, 2007), was defined as the day when the EVI signal was higher than the moving average of the previous four time steps at the beginning of the hydrological year, which is the time when the EVI time series experiences a sudden increase.

3 Results

Before addressing the study objectives, we first present the variability in precipitation and LAI (or EVI). Next, the outcomes of the three ET partitioning methods summarized in Table 1 are featured. The tree and grass dynamics are outlined analyzing the monthly variability in grass and tree NDVI, T, and T_g. Finally, the annual water balance components are presented. Hereafter, hydrological years are defined as the time period from September to August and are referred to by the year in which they began.

Figure 2 Time series of daily precipitation, soil moisture at 10 and 60 cm depth, EC-based evapotranspiration, GPP, VPD, EVI, and 7 d average $T/\mathrm{ET}$ for the B16, Z16, and N18 methods. The GPP and VPD lines show the 7 d running mean. Rainy-season length T_wet (days), daily mean precipitation amount (λ), and storm frequency (α) for the early wet season (September to November) are indicated for each year in panel (a).

3.1 Site meteorology and ET partitioning

The early-season rainfall was frequent in all years except 2011 and 2015 (Fig. 2a). The tree green-up days and start of the rainy period were not linearly related (R²=0.03, p=0.753; Table S1 in the Supplement): the earliest tree green-up date occurred in 2011, 72 d earlier than the start of the rainy season. The year 2015 was an extreme drought year in South Africa, and it was characterized by the lowest early-season rainfall frequency and by a rainy-season length that was nearly twice as long as in other years (Fig. 2a). Water entered the deeper soil layers in wet years 2012 and 2013 but not in the wet year 2014 nor in the drought year 2015 (Fig. 2b). There was a 2-week dry spell at the end of January 2011 and another dry spell at the end of November 2015, which are visible through the low measured ET and soil moisture values (Fig. 2b, c). Grass experienced dieback and regrowth in 2015, as shown by decreased GPP during the early wet season (Fig. 2d). This is also seen in the EVI trend, which decreases after initial increase (Fig. 2e). This period was also characterized by high VPD.

The daily $T/\mathrm{ET}$ estimated using the N18 method was consistently higher than the $T/\mathrm{ET}$ from the B16 and Z16 methods (Fig. 2f). The difference between N18 and the other methods was most significant during the end of wet seasons and during the drought year. The $T/\mathrm{ET}$ estimated from Z16 was higher than $T/\mathrm{ET}$ from B16 during the dry seasons. In the B16 method, the annually fitted line between the variable GPP×VPD^0.5 and ET established the empirical link between GPP and T. The bin values of the variable GPP×VPD^0.5 were linearly related to the 5th percentile of measured ET (Fig. S6), with the largest scatter occurring during the drought year (R²=0.85, in 2015). For all years, the mean surface soil moisture during T=ET instances was 0.1 m³ m⁻³ or less (Table S1). The annual slope of the T=ET was related to the rainy-season length, with the year 2011 falling below the 95 % confidence interval of the mean (Fig. 3a). The slope represents the T=ET values, and only 67 % of those values were from the rainy season in 2011, as opposed to 75 %–84 % in the other years (Table S1). The greater slope value in 2011 means lower water-use efficiency. The slope and intercept of the T=ET line were also linearly related (Fig. 3b). Thus, when most T=ET values are observed during the rainy season, it is possible to estimate the annual T=ET line based only on the rainy-season length for the B16 method. The annually fitted uWUEp for the Z16 method was related to the sum of R_n during the wet season (Fig. 3c).

Figure 3 (a) The relationship between the slope of the T=ET line (B16 method) and rainy-season length. (b) The relationship between the y intercept and the slope of the T=ET line for the B16 method. (c) The relationship between annual uWUEp (Z16 method) and wet-season R_n. Dashed lines demarcate the 95 % confidence interval.

Figure 4 The relationship between cumulative daily soil evaporation and (day)${}^{\mathrm{1}/\mathrm{2}}$ for the stage-2 evaporation events for the B16, Z16, and N18 methods. The slope of the regression line is soil desorptivity. Titles indicate the start date of each dry-down event.

The early dry-season daily evaporation was assessed according to the stage-2 theory of soil evaporation. The regression between cumulative daily E and $\sqrt{t_{\mathrm{d}}}$ is linear for all methods (Fig. 4). The derived D_e values are the highest for the B16 method, ranging from 2.23 to 4.20 mm d${}^{-\mathrm{1}/\mathrm{2}}$. The largest difference between the methods occurred in the late wet season in 2015, with D_e values of 2.91, 1.92, and 1.08 for the B16, Z16, and N18 methods, respectively (Fig. 4f). The N18 method has the lowest D_e values except in the middle of the dry season in 2015 (Fig. 4g). The estimated D_e values were linearly related to the first day air temperature for the B16 and Z16 methods (Fig. S7). Overall, the D_e values from the B16 method most closely match the reported range of D_e from other studies of sandy soils (Brutsaert and Chen, 1995; Hu and Lei, 2021).

Figure 5 (a) Time series of monthly $T/\mathrm{ET}$ estimated with the B16, Z16, and N18 methods. (b) Time series of monthly T and T_g. (c) Time series of 16 d MODIS EVI and Landsat 8 NDVI for grasses and trees. The dotted vertical line indicates the start of the hydrological year (1 September).

3.2 Monthly tree and grass transpiration

A comparison of the three different monthly $T/\mathrm{ET}$ estimates shows that $T/\mathrm{ET}$ according to the N18 method is consistently higher than the $T/\mathrm{ET}$ from the other two methods during the wet seasons (Fig. 5a). The B16 and Z16 methods have similar $T/\mathrm{ET}$ seasonality during the wet seasons, whereas the N18 $T/\mathrm{ET}$ is higher than the other methods during the late wet seasons. The monthly $T/\mathrm{ET}$ from N18 is nearly constant during the wet season in 2011, whereas $T/\mathrm{ET}$ from B16 and Z16 shows a decreasing trend after the early-season peak. The largest difference between $T/\mathrm{ET}$ estimated with the Z16 and B16 methods occurred from March to June 2015 (Fig. 5a). During that period, according to all methods, the monthly GPP decreased while $T/\mathrm{ET}$ increased. The B16 soil desorption was closest to published values during this period (Fig. 4f), suggesting that Z16 and N18 overestimated $T/\mathrm{ET}$ during the late wet season of 2015.

Estimates of monthly T were similar based on the B16 and Z16 methods, but they were consistently higher based on the N18 method (Fig. 5b). The rainy season began 48 d later in 2014 than in 2012 (Table S1). This delay is consequently reflected in the monthly course of T (B16 method) in 2014, which lagged behind that of 2012 until January (Fig. 5b). The dry spell in 2011 is clearly shown by reduced T and ET during this period. The soil-moisture-based estimate of grass transpiration shows a similar trend to estimated T in 2012 and 2013, whereas the T_g trend is less variable than the T trend in 2014 and 2015 (Fig. 5b). In 2014, T_g was nearly constant during the wet season despite increasing grass NDVI (Fig. 5b, c).

Both the tree and grass peak NDVI values are lower during the drought year compared with the wet years 2013 and 2014 (Fig. 5c). The grass dieback is reflected in low grass NDVI after the initial increase during the early drought year. In addition, the T_g decreases during this period. The second peak in grass NDVI is similar in magnitude to the first, suggesting modest grass growth after the dieback. The grass NDVI and T_g decrease during the middle and late wet season, whereas the tree NDVI and T stay nearly constant. This means that the tree contribution to the total transpiration increases from the middle of the wet season onwards in 2015.

The monthly $T/\mathrm{ET}$ to LAI relation was scattered for all methods (Fig. 6a). The $T/\mathrm{ET}$ was higher for the Z16 method than for B16 at low LAI values. The VPD response of monthly $\mathrm{GPP}/T$ was most nonlinear for the B16 method and least nonlinear for the N18 method (Fig. 6b). The monthly LAI–T relation was linear for all methods and was similar between the B16 and Z16 methods (Fig. 6c). The T_g was increasing for higher LAI values in 2012 and 2013, whereas it was nearly constant at 15 mm per month for a wide range of LAI values in 2014 and 2015 (Fig. 6d).

Figure 6 The relationship between monthly $T/\mathrm{ET}$ and MODIS LAI. The black line is the relationship for shrub and grass ecosystems that is based on measured LAI values (Wei et al., 2017). The equal-width bins are fitted using the following equation: a⋅(LAI)^b. (b) The relationship between monthly $\mathrm{GPP}/T$ and monthly mean VPD. The fit line is $a\cdot \ln \left(\mathrm{VPD}\right)+b$. (c) The relationship between monthly T and MODIS LAI. (d) The relationship between monthly T_g and MODIS LAI, where colors indicate data from the 4 different years.

The monthly GPP and T were linearly related (R²=0.97, p<0.001; Fig. S8), allowing for the estimate of an effective (constant) ecosystem water-use efficiency using a zero-intercept regression. The constant water-use efficiency was 2.83, 2.78, and 2.29 g C per kg H₂O for B16, Z16, and N18 methods, respectively (Fig. S8).

3.3 Interannual variation

Annual P was close to or above the mean annual rainfall of 591 mm yr⁻¹ (adjusted here for the undercatch) for every year except in 2015, which was an extreme drought year in South Africa (Table 2). The EC-measured annual ET was close to annual P for all years (Table 2). The annual P−ET ranged from −31 to 45 mm yr⁻¹, and it was inversely related to the annual maximum EVI (R²=0.87, p=0.007). The annual change in soil water storage was small (1 to 14 mm yr⁻¹) compared with the variation in other water balance components, and it was unrelated to the annual P−ET (R²=0.45, p=0.332). The frequent evening and nighttime precipitation resulted in nighttime evapotranspiration (ET_N), which varied from 58 to 85 mm yr⁻¹ (12 % of annual ET). The annual P−ET would be positive for all years if ET_N were assumed to be zero. The annual estimated rainfall interception ranged from 69 to 97 mm yr⁻¹ and was linearly related to ET_N (R²=0.75, p=0.025).

Table 2 Annual sum of water balance components for each hydrological year (September to August). The total uncertainty (Eq. 4) is indicated for ET after the “±” sign. ET_N is the annual nighttime evapotranspiration. Transpiration and evaporation estimated using the B16 method are shown. The tree transpiration (T_tree) was calculated by subtracting T_g from T. The PET was determined from Eq. (1). The EBC slope stands for the slope of the energy balance closure, with the ordinate defined by measured R_n-G and the abscissa defined by the sum of the measured latent and sensible heat fluxes.

The estimated annual $T/\mathrm{ET}$ ratio varied from 0.33 to 0.56 for the B16 method (Table 2). The annual $T/\mathrm{ET}$ from Z16 was similar to that from the B16 method (with a 0.04 higher 6-year mean), whereas the N18 $T/\mathrm{ET}$ mean was 0.15 higher than the B16 mean (Table S2). The $T/\mathrm{ET}$ in 2011 was similar to other wet years for the N18 method but reduced compared with other wet years for the B16 and Z16 methods (Tables 2, S2). The annual $T/\mathrm{ET}$ and T were linearly related to the early wet-season storm frequency for the B16 and Z16 methods (Fig. 7), while the relation of early-season P, middle wet-season P, annual P, or mean annual EVI with annual $T/\mathrm{ET}$ and T were more scattered (Figs. S9, S10, S11). For the N18 method, the annual T was linearly related to mid-wet-season P (Fig. S10), and annual T and $T/\mathrm{ET}$ were linearly related to mean EVI (Fig. S11).

Figure 7 The relationship between early wet-season storm frequency and annual $T/\mathrm{ET}$ and T for the B16, Z16, and N18 methods. Dashed lines demarcate the 95 % confidence intervals.

Annual transpiration was nearly constant at 326±19 mm yr⁻¹ for the 4 years with frequent early-season rainfall (Table 2). The annual T was similar for the wet years of 2012 and 2014, despite 97 mm lower T_g and the late start of the rainy season in 2014. The lower T in 2011 is explained by the infrequent early-season rainfall and the 2-week dry spell in January. The average dry-season transpiration was 9 mm (per quarter) for the B16 method (Table S3), which suggests a minimum tree transpiration of 36 mm yr⁻¹. The annual T−T_g representing tree transpiration ranged from 88 to 160 mm yr⁻¹ for the 3 wet years (Table 2). The annual T_g was similar for 2014 and 2015, despite a significantly lower grass NDVI in 2015 (Fig. 5c). Therefore, during the drought year, T_g might be an overestimate due to the contribution from soil evaporation, which would explain the low value (23 mm yr⁻¹) of tree transpiration.

4 Discussion

At annual timescales, P was approximately equal to annual ET, consistent with other studies from sites with similar annual rainfall (Gwate et al., 2018; Scholes and Walker, 1993). Annual T was nearly constant (326±19 mm yr⁻¹) during the 4 years with frequent early wet-season precipitation (Table 2), as has been found in different types of forest ecosystems (Oishi et al., 2010; Tor-ngern et al., 2017; Ward et al., 2018). However, it was lower in years with infrequent early wet-season rainfall, producing intermittent dry spells. The variation in annual T and $T/\mathrm{ET}$ was explained by the variation in the early wet-season storm frequency (Fig. 7). The monthly T was linearly related to LAI, whereas $T/\mathrm{ET}$ had an expected nonlinear relation to LAI (Fig. 6a, c). The constancy of annual T reflects the moderate water stress of C₄ grasses shown by similar T_g during the 2 wet years (Table 2). The annual T was not reduced during the year with a late rainy-season start and lower T_g, showing that trees can use the water that is not transpired by the grass layer (Fig. 5b). However, during the extreme drought year, the peak grass and tree NDVI values were lower than during the wet years (Fig. 5c). The T_g was decreasing and grass regrowth was modest after the grass dieback from the middle to late wet season, while the tree NDVI was constant, showing a greater contribution of tree transpiration during the late wet season (Fig. 5b, c). The rainfall timing control on annual T as well as the grass dieback and regrowth dynamics show that the grasses' early-season development foreshadows the interannual variation in T and $T/\mathrm{ET}$, which reflects the fact that the temporal distribution of P is a good predictor of the productivity at these Highveld grasslands (Swemmer et al., 2007).

4.1 Transpiration

The small interannual variation in transpiration is likely due to moderate water stress of the C₄ grass layer as well as trees' ability to use water from shallow and deeper soil layers, depending on water availability. The photosynthesis reduction in C₄ grass is more related to non-stomatal limitations compared with C₃ grass, which is predominantly limited by stomatal control (Ripley et al., 2010). Under South African field conditions over one growing season, the C₄ grass layer tended to maintain a constant difference between the predawn and midday leaf water potential, with similar transpiration at rain-fed and irrigated pot trials (Taylor et al., 2014). The annual T_g at Welgegund was similar during the 2 wet years, but it was lower during the year with a late start to the rainfall. During the 2015 drought year, the grass cover underwent a dieback–regrowth cycle, tree peak NDVI was less than during wet years, and the annual T was reduced by 51 % from the mean of wet years (Fig. 5b, c). The tree contribution to the total transpiration increased from the middle of the wet season onwards during the drought year, as T_g decreased while T increased (Fig. 5b). A similar grass growth pattern was also observed in Kruger National Park, where the grass biomass decreased and vast areas were barren during the drought, but the grasses quickly recovered once the rains returned (Wigley-Coetsee and Staver, 2020).

The trees at Welgegund likely rely more on topsoil water when available, switching to a greater reliance on deeper layers during dry seasons and droughts. A similar tree water-use pattern has been observed using stable isotopes in the South African savanna (Priyadarshini et al., 2016). Welgegund is located at the wet end of the distribution range of the dominant tree species V. erioloba. The estimated radiocarbon age of these trees is approximately 20 years (Steenkamp et al., 2008). The roots of V. erioloba are deep and reportedly extend to a depth of up to 60 m (Jennings, 1974); in one study, the roots absorbed 37 % of the transpired water below a depth of 1 m (Beyer et al., 2018). In addition, the horizontal extent of the roots of this species can exceed 20 m (Wang et al., 2007). When soil moisture profile measurements at the site were installed 15 m away from the nearest tree, tree roots were observed 0.4 m from the surface and deeper. The mean annual tree transpiration was 127 mm yr⁻¹ (40 % of annual T) for the 3 wet years (Table 2). This estimate is similar to the tree transpiration of 126 mm yr⁻¹ at a site in South Africa (Nylsvley) (P=586 mm yr⁻¹) with 30 % tree cover and shallower tree roots (measurements and modeling; Scholes and Walker, 1993). At a savanna site in this region (P=241 mm yr⁻¹), V. erioloba (5 m tall) had an annual/dry-season transpiration ratio of 6.5 ($\mathrm{3.9}/\mathrm{0.6}$ mm d⁻¹; Tfwala et al., 2019). Multiplying this ratio by the 9 mm dry-season tree transpiration in our study results in 59 mm yr⁻¹ annual tree transpiration. This shows that the dry-season-based tree transpiration estimate is lower than the T−T_g estimated that combines the eddy-covariance estimate of T and soil-moisture-based T_g. The T−T_g tree transpiration estimate includes any error made in the T_g estimation. The overestimation of T_g likely explains the low tree transpiration during the drought year. Removing the constant soil evaporation from the T_g estimate may not be adequate during the drought year, leading to a soil evaporation contribution to T_g. The T_g estimation algorithm also has lower performance in dryer soils (Jackisch et al., 2020). The estimated long-term water-use efficiency was 2.83 and 2.78 g C per kg H₂O for the respective B16 and Z16 methods. This is somewhat higher than the wet-season value of 2.4 g C per kg H₂O (Z16 method) at a C₄ grassland site in southeastern Arizona, USA (Scott et al., 2021), and the field-scale long-term grass community value of 2.15 g C per kg H₂O for the aforementioned savanna with shallower-rooted trees and 30 % tree cover in Nylsvley (Scholes and Walker, 1993).

During water-stressed years, the partitioning of the tree and grass contribution to LAI and $T/\mathrm{ET}$ may be needed to derive meaningful relations at the monthly scale. The drought year was characterized by different grass and tree NDVI trends (Fig. 5c) as well as nearly constant grass transpiration for a wide range of LAI values (Fig. 6d). New remote sensing products may be able to separate these contributions, as shown by a recent study that successfully separated tree and grass leaf area using a canopy height model, Sentinel vegetation indexes (10 m spatial resolution), and a Sentinel radar band during the 2015 drought in Kruger National Park (Urban et al., 2018). The effect of dieback–regrowth on annual transpiration is also interesting, as a stochastic model based on measured precipitation statistics with explicit bare soil, grass, and tree cover showed that vegetation dynamics had little effect on annual transpiration (Williams and Albertson, 2005).

Water availability for grass is the dominant factor in transpiration at Welgegund. The annual $T/\mathrm{ET}$ range (0.33 to 0.56) at Welgegund is slightly wider than the mean annual $T/\mathrm{ET}$ range of 0.35 to 0.46 at a C₄ grassland site in southeastern Arizona, USA (P=317 mm yr⁻¹), where the annual $T/\mathrm{ET}$ correlates with annual P and mean LAI (Scott et al., 2021). At Welgegund, the early-season rainfall frequency was a better predictor of annual T than annual P or mean EVI (Figs. 7, S10, S11). These relations might be due to the heavy grazing at the site, which limits peak EVI and emphasizes the early-season grass development. The analysis of the different C₄ grass species during the 2014–2016 South African drought suggests that their bundle sheath morphology explains the differences in drought tolerance (Wigley-Coetsee and Staver, 2020). Therefore, it is difficult to generalize whether the invariance in annual transpiration during the wet years would hold for sites with higher grass LAI or a different grass species composition.

4.2 Uncertainty

The 6-year mean $\mathrm{ET}/P$ ratio was 1.0 at Welgegund (Table 2), slightly higher than a long-term $\mathrm{ET}/P$ ratio of 0.96 at a C₄ grassland in Arizona, USA (Scott and Biederman, 2019). The lowest annual P−ET was −31 mm yr⁻¹, which is more negative than the estimated annual ET uncertainty but less than the uncertainty related to ET_N gap filling (Table 2). Due to frequent afternoon and nighttime precipitation, the ET_N was 12 % of the annual ET. The ET_N values here may appear high, but they are commensurate with reported values for forested ecosystems (Novick et al., 2009) in regions with higher precipitation and LAI. The gap-filled ET_N may be an overestimate because only 30 % of the values were measured, and these values were determined during high wind speeds ($u_{*}>\mathrm{0.28}$ m s⁻¹).

The ratio of annual ET uncertainty to annual ET was 1.3 %, which is lower than the 5 % to 9 % range reported from eddy-covariance ET measurements from a cultivated area in Benin (Mamadou et al., 2016). The difference can be ascribed to different error terms used in the uncertainty estimation. The ET in the middle of the dry season ranged from 45 to 68 mm (per quarter; middle dry-season monthly value multiplied by 3) at this cultivated site in Benin, which has isolated trees (height < 10 m) and bare soil during the dry season (Mamadou et al., 2014, 2016). This is higher than the 29 to 52 mm (per quarter) range measured in our study. These differences may be attributed to the relatively shallow water table (a depth of 3 m during the dry season) and the higher annual precipitation (P=1200 mm yr⁻¹) at the Beninese site.

4.3 ET partitioning methods

The B16 and Z16 transpiration estimates were more similar and closer to reported grassland $T/\mathrm{ET}$ values than the N18 estimate (Fig. 5a). The N18 estimate was also higher than the Z16 estimate at the C₄ grassland site in southeastern Arizona, USA (Scott et al., 2021). The annual maximum of monthly $T/\mathrm{ET}$ ranged from 0.57 to 0.67 for the 4 wet years (B16 and Z16 methods), which is similar to the maximum value of 0.60 at a C₄ grassland site in southeastern Arizona, USA, estimated using the ET partitioning method that does not assume equality between T and ET (Scott and Biederman, 2017). The assessment of soil evaporation according to stage-2 theory (Fig. 4) showed that the B16-estimated D_e most closely matched the reported range of D_e (3 to 6 mm d${}^{-\mathrm{1}/\mathrm{2}}$) from other studies of sandy soils (Brutsaert and Chen, 1995; Hu and Lei, 2021). The largest difference in D_e between the B16 and Z16 methods occurred in the late wet season in 2015, when the Z16 $T/\mathrm{ET}$ was deemed an overestimate based on the low D_e value of 1.92 mm d${}^{-\mathrm{1}/\mathrm{2}}$. The estimated D_e values were linearly related to the air temperature on the first day for the B16 and Z16 methods (Fig. S7). Similar dependence has been observed under laboratory conditions with the full wetting of sandy soil columns (Ben Neriah et al., 2014). For the Z16 method, a one-to-one T=ET line is fitted using quantile regression for each year combined with the intercept forced through zero, whereas the T=ET line is fitted using linear regression over the GPP×VPD^0.5 bins for the B16 method. Regarding the N18 method, T was likely an overestimate for several reasons. The estimated D_e values were much lower than the literature values, and the annual $T/\mathrm{ET}$ values had much higher range (0.54 to 0.65) than other methods and studies (Scholes and Walker, 1993; Scott et al., 2021). In addition, the tree transpiration estimated using N18 (by subtracting estimated grass transpiration from N18 estimated annual T) results in an average annual tree transpiration of 193 mm over the 3 wet years, which is much higher than values reported in the literature. The N18 algorithm does not use measured soil moisture in the training period; instead, it uses the P and ET water balance, which may explain the small interannual variance in the maximum $T/\mathrm{ET}$ values (Nelson et al., 2018).

The low surface soil moisture values during T=ET periods and their concentration during the rainy season give assurance that the annual fitted T=ET lines correspond to periods when T equals ET (Table S1). The annual T=ET line could be predicted using the rainy-season length, except in 2011, which experienced the earliest green-up of trees and the highest number of T=ET moments outside of the rainy season (Fig. 4a, Table S1). The water balance analysis, focused on monthly and annual timescales using the ET partition methods, has shown good agreement with independent estimates (Berkelhammer et al., 2016; Zhou et al., 2018). Berkelhammer et al. (2016) showed that a 3 d running mean of the 30 min $T/\mathrm{ET}$ estimates reduced the root-mean-square difference between the Berkelhammer method (B16) and the isotopic estimate of $T/\mathrm{ET}$ to ≤0.2. Therefore, the random error of the monthly means of the 30 min $T/\mathrm{ET}$ estimates in this study can be assumed to be small. For a Mediterranean tree–grass savanna, $T/\mathrm{ET}$ was shown to rarely exceed 0.8 (Perez-Priego et al., 2018). In contrast to the Mediterranean site, the Welgegund site has sandy soil, deep-rooted trees, and no clay horizon close to the soil surface. More importantly, the mean surface soil moisture was 0.1 m³ m⁻³ or below for the 30 min runs when T=ET at Welgegund. This low soil moisture resulted in small diffusion-limited soil evaporation and, thus, periods when T=ET. This is another independent confirmation of the partitioning of ET into E and T (even at such short timescales).

5 Conclusion

The measurements reported here show that annual transpiration is nearly constant during years with frequent early-season rainfall but can be lower because the C₄ grass cover reacts to dry spells. Trees at the site are likely able to utilize water that is not used by the grass layer as well as water from deeper layers during extreme drought. Our work highlights precipitation control over T and the annual variation in the T–LAI relationship. These results can be used to assess the water resources and fodder production of grassland grazing systems. Although further work is required to determine the generality of these conclusions to other savanna systems, the methodologies developed and tested here can be employed when investigating a wide range of arid and semiarid ecosystems experiencing water shortages in times of drought.