Net Radiation Drives Evapotranspiration Dynamics in a Bottomland Hardwood Forest in the Southeastern United States: Insights from Multi-Modeling Approaches

: Evapotranspiration (ET) is a major component of the water budget in Bottomland Hard-wood Forests (BHFs) and is driven by a complex intertwined suite of meteorological variables. The understanding of these interdependencies leading to seasonal variations in ET is crucial in better informing water resource management in the region. We used structural equation modeling and AIC modeling to analyze drivers of ET using Eddy covariance water flux data collected from a BHF located in the Russel Sage Wildlife Management Area (RSWMA). It consists of mature closed-canopy deciduous hardwood trees with an average canopy height of 27 m. A factor analysis was used to characterize the shared variance among drivers, and a path analysis was used to quantify the independent contributions of individual drivers. In our results, ET and net radiation (R n ) showed similar variability patterns with Vapor Pressure Deficit (VPD) and temperature in the spring, summer, and autumn seasons, while they differed in the winter season. The path analysis showed that R n has the strongest influence on ET variations via direct and indirect pathways. In deciduous forests like BHFs, our results suggest that ET is more energy dependent during the growing season (spring and summer) and early non-growing season (autumn) and more temperature dependent during the winter season.


Introduction
In the central and southeastern United States, deciduous forested wetlands located in broad floodplain areas bordering large river systems such as within the Lower Mississippi Alluvial Valley (LMAV) are referred to as Bottomland Hardwood Forests (BHFs) [1].Many ecosystem services provided by BHFs, including water quality regulation, flood control, wildlife habitat, timber production, waste treatment and disturbance regulation, and climate moderation through carbon balance are of global importance [2,3].Like other wetland types, a hydrological regime of alternating wet and dry cycles driven by fluctuating water levels of the associated rivers and groundwater level changes is characteristic of BHFs [4,5].The sustainability of BHFs depends primarily on the longitudinal (upstream to downstream), lateral (river to floodplain to uplands and vice versa), vertical (surface water to groundwater and vice versa), and temporal (seasonal and annual flooding) variability of water availability [6].The primary force controlling the biota in BHFs is the flooding pulse from adjacent water sources, which deposits dissolved nutrients, organic matter, and sediment, also contributing to the formation of young floodplains with each flooding pulse [7].Even small changes in the duration and frequency of water levels can result in a distinct shift in the plant community, as many species are adapted to a certain range of flood tolerance [8].Therefore, a deeper understanding of the water use pattern and surface energy balance in these BHFs is crucial to preserving this dwindling ecosystem along the LMAV.
The study of BHF ecology remains incomplete without the assessment of factors that are integral to the proper functioning of these systems-evapotranspiration (ET), flood regimes, and precipitation.ET has been reported as a major component of BHFs' water balance in several studies, along with LE (heat energy equivalent to ET) dominating the surface heat balance [9,10].Several previous studies have identified biological and climatic drivers of bottomland hardwood ET using chronosequence analyses [9], sapflux measurements [11], and statistical modeling approaches [12] in the southeastern US.These studies report a strong correlation between net radiation (R n ), temperature, and other climatic conditions consistent with the site specificities and variations in ET during different seasonal cycles.The Vapor Pressure Deficit (VPD) has been established as one of the primary drivers of ET and is used increasingly in global simulation studies [13].For example, [14] used maximum likelihood estimation methods to show a complex chain of correlations among ET, VPD, radiation, and temperature in northern high-latitude woodland ecosystems.The complex interrelationships among these atmospheric variables, which are essential for an assessment of the ultimate drivers of variations in ET, can be explored using structural equation modeling (SEM).The use of SEM approaches to diagnose the independent contributions of atmospheric drivers in determining the ET variability from the BHFs largely remains unexplored.
The complexity of interrelationships among multiple variables and their dependencies makes it challenging to quantify the actual contributions of the drivers of ET.SEM, a multivariate statistical modeling technique with factor analysis and path analysis, provides a framework to quantitatively evaluate these interrelationships that need to be untangled to characterize the role and independent contribution of each driver of the variations observed in ET.An SEM analysis carried out with the data measured at multi-site high-latitude regions showed that radiation, temperature, windspeed, and relative humidity (RH) loaded heavily on the first factor during the warm season (May-September) [15,16].Previous studies have applied path analyses to diagnose the drivers of ET in temperate [17], Arctic, and subarctic [15] regions; however, to our knowledge there is not a study reported on the water fluxes in BHFs-where altered hydrologic cycles under global change scenarios are critical and complex.
This paper is a comprehensive assessment of the drivers of variations in bottomland hardwood ET over hourly, daily, and weekly timescales across seasons.The main objectives were (i) to characterize the interrelationships among variables driving ET over these timescales and (ii) to characterize ET dependencies on various factors across different seasons.This study is unique to those mentioned above in various important ways.First, it utilizes an SEM factor analysis and a path analysis to provide a framework for quantitively evaluating the relative importance of drivers of ET variability in this system.Second, this is the first study of the drivers of ET in the BHF of the Russel Sage Wildlife Management Area (RSWMA), a representative of floodplain forests in the entire LMAV.This study will fill the knowledge gap in trying to understand the water flux dynamics of this region.

Study Site
This study was conducted in a BHF in the RSWMA (Figure 1) in northeast Louisiana (32.46 • N, −91.97 • E; elevation 18 m ASL), managed by the Louisiana Department of Wildlife and Fisheries (LDWF).The RSWMA is located within the Bayou Lafourche floodplain and is subjected to annual late winter (December-February) to early spring (March-May) flooding.It currently covers an area of 38,213 acres with mature hardwood stands that vary in age from 90 to 120 years (pers comm Larry Savage, LDWF).The broadleaf deciduous forest canopy consists of the co-dominant canopy species overcup oak (Quercus lyrate Walter) and water hickory (Carya aquatica (F.Michx.)Elliott), along with other canopy species such as green ash (Fraxinus pennsylvanica Marsh.) and sugarberry (Celtis laevigata Willd.) in poorly drained soil [5].In the first bottoms on low ridges, flats, and sloughs, American elm (Ulmus americana L.), sweetgum (Liquidamber styraciflua), winged elm (Ulmus alata Michx.), and red maple (Acer rubrum L.) are abundant.In the newly formed sandbars on river margins, black willow (Salix nigra Marshall), cottonwood (Populus deltoides W. Bartram ex Marshall), river birch (Betula nigra L.), and American sycamore (Platanus occidentalis L.) are prominent.The well-drained bottom ridges are dominated by sweetgum and water oak (Quercus nigra L.), characteristic species of the BHFs of the LMAV.The canopy is relatively flat, with a mean tree height of 27 m.The Leaf Area Index (LAI) of the forested area, derived from Moderate Resolution Imaging Spectroradiometer (MODIS) instrument-based observations, was below 1 throughout the winter season, and on average reached a maximum of 6.53 in July.The soil type is Perry Clay, a fine-textured sediment that has low permeability and a moderate capacity to hold water [18].

Measurements of Above-Canopy Fluxes
The Eddy covariance (EC) technique was used to measure the amount of water vapor that was exchanged between the atmosphere and the BHF ecosystem.The wind components, sonic temperature, and gas concentrations were collected at a 10 Hz frequency using an open-path IRGASON (Campbell Scientific Inc., Logan, UT, USA).The EC system was mounted 12 m above the forest canopy and was directed towards the southwest, the prevailing wind direction at the site, and the fetch of the tower was 2 km.

Measurements of Meteorological, Phenological, and Hydrological Variables
Other measurements relevant to land-plant-atmosphere interactions were carried out from the tower to understand how they affect the hydrodynamics at the study site.The meteorological measurements include air temperature, precipitation, windspeed and direction, barometric pressure, and RH.The VPD was calculated as the difference between saturated and actual vapor pressures at the given temperature, based on the relative humidity and air temperature data.A net radiometer (NR-LITE2, Campbell Scientific Inc., The broadleaf deciduous forest canopy consists of the co-dominant canopy species overcup oak (Quercus lyrate Walter) and water hickory (Carya aquatica (F.Michx.)Elliott), along with other canopy species such as green ash (Fraxinus pennsylvanica Marsh.) and sugarberry (Celtis laevigata Willd.) in poorly drained soil [5].In the first bottoms on low ridges, flats, and sloughs, American elm (Ulmus americana L.), sweetgum (Liquidamber styraciflua), winged elm (Ulmus alata Michx.), and red maple (Acer rubrum L.) are abundant.In the newly formed sandbars on river margins, black willow (Salix nigra Marshall), cottonwood (Populus deltoides W. Bartram ex Marshall), river birch (Betula nigra L.), and American sycamore (Platanus occidentalis L.) are prominent.The well-drained bottom ridges are dominated by sweetgum and water oak (Quercus nigra L.), characteristic species of the BHFs of the LMAV.The canopy is relatively flat, with a mean tree height of 27 m.The Leaf Area Index (LAI) of the forested area, derived from Moderate Resolution Imaging Spectroradiometer (MODIS) instrument-based observations, was below 1 throughout the winter season, and on average reached a maximum of 6.53 in July.The soil type is Perry Clay, a fine-textured sediment that has low permeability and a moderate capacity to hold water [18].

Measurements of Above-Canopy Fluxes
The Eddy covariance (EC) technique was used to measure the amount of water vapor that was exchanged between the atmosphere and the BHF ecosystem.The wind components, sonic temperature, and gas concentrations were collected at a 10 Hz frequency using an open-path IRGASON (Campbell Scientific Inc., Logan, UT, USA).The EC system was mounted 12 m above the forest canopy and was directed towards the southwest, the prevailing wind direction at the site, and the fetch of the tower was 2 km.

Measurements of Meteorological, Phenological, and Hydrological Variables
Other measurements relevant to land-plant-atmosphere interactions were carried out from the tower to understand how they affect the hydrodynamics at the study site.The meteorological measurements include air temperature, precipitation, windspeed and direction, barometric pressure, and RH.The VPD was calculated as the difference between saturated and actual vapor pressures at the given temperature, based on the relative humidity and air temperature data.A net radiometer (NR-LITE2, Campbell Scientific Inc., USA) was used to measure the difference between the incoming and outgoing radiation at the site.A photosynthetically active radiation (PAR) sensor (LI190SB QUANTUM SENSOR, Campbell Scientific Inc., USA) was used to quantify the photosynthetic photon flux density (PPFD).Precipitation was measured by a tipping bucket rain gauge (TE525, Campbell Scientific Inc., USA).

Data Collection, Processing, and Gap-Filling Fluxes
All data were acquired using a solid-state data logger (CR3000, Campbell Scientific Inc., USA).The data were stored on a 2 GB CompactFlash card, which was retrieved monthly.The raw flux data were grouped into individual 30 min files and converted into TOA5 format before processing.These unprocessed data were screened for quality control and gap-filled using LoggerNet, EddyPro, TOVI, and REddyProc R packages [19].To determine the periods of low mixing, which can lead to the underestimation of water and heat fluxes, the frictional velocity (u*) threshold was identified using the Moving Point Test approach in TOVI [20].Further quality control screening was carried out to filter potential data outside the u* threshold value (0.4 m/s).
Data recorded from 1 January 2014 to 31 December 2021 were used in this analysis.There was a major missing data gap due to instrument failure from January 2016 to July 2016 and from November 2018 to December 2020, and these gaps were not included in the analysis.Small gaps (<2 h) due to lower-quality data were gap-filled using the Marginal Distribution Sampling Technique using the R package, REddyProc.Meteorological data from a nearby meteorological station (Monroe Airport-MLU) and NASA's Prediction of Worldwide Energy Resources (POWER) [21] database were used for gap-filling temperature and precipitation data when necessary.The gap-filled data accounted for less than 20% of the whole EC data.The water losses associated with ET were quantified by the conversion of LE values from (W/m 2 ) to mm/day [15].To compare the means and assess the interrelationships among variables, we performed an ANOVA test and the associated post hoc tests and an SEM analysis as required, using R (v4.1.1,R Core Team 2023) [22].
The data were aggregated into hourly, daily, and weekly datasets spanning the growing (March-August) and non-growing (September-February) seasons.Variables with more than 50% of their half-hourly data missing were not included in the SEM analysis.The missing data were then adjusted as described in [23] by dividing each daily and weekly value by the fraction of the data present for that particular unit of time.RH and soil moisture, two of the primary drivers of ET, were not included in our analysis since RH is redundant to the VPD, and soil moisture was not measured in our study.

Structural Equation Modelling (SEM)
As ET is influenced by a suite of atmospheric variables, a powerful multivariate analysis technique was required to identify those strongly driving the variation in the rate of ET.Different empirical approaches [9,11,12,24] and machine learning algorithms [25] have been in practice to characterize the key drivers of ET.However, due to the complexity of a large dataset with interdependent variables, SEM emerges as a suitable multivariate analysis technique to reduce the dimensionality and isolate common variance shared among variables from the residual variance unique to each variable [26].Depending on the types of variables being modeled and their relationships with the predictor variables, several categories of models fall under the suite of SEM.A factor analysis derives the latent constructs from variables that share the most variance with related variables while allowing them to be influenced by the seasonal cycles in the dataset, and a path analysis explains the structural pathways of the interrelationships among variables while testing for underlying causal mechanisms.Together, these analyses allow inferences to be made about the independent contributions of interrelated variables in a dataset.

Factor Analysis
A factor analysis estimates latent variables based on the correlated variations of a dataset (e.g., associations and causal relationships) and can reduce the dimensionality of the dataset, standardize the scale of multiple indicators, and account for the correlations inherent in the dataset [26].The covariance matrix is central to a factor analysis, as described by the various equations in textbooks [26,27].There are two types of factor analyses: exploratory factor analysis (EFA) and confirmatory factor analysis (CFA).In practice, an EFA is often performed to select the useful underlying latent constructs for a CFA when there is little prior knowledge about the latent construct [28].A CFA is applied when the indicator for each latent variable is specified according to the related theories or prior knowledge.However, both types operate under the same basic assumption that "for a set of observed variables there are a suite of underlying factors which explain the interrelationships among the variables" [15].The strongest association with ET is shown by variables that load highly on factors in which ET loads strongly as well.

Path Analysis
A path analysis helps to quantify the direct and indirect relationships among multiple variables while allowing them to covary with other variables in the dataset.It is very powerful in testing and developing the structural hypothesis where variables can influence an outcome directly and indirectly through another variable (e.g., mediation) [29].This specific type of SEM uses a set of exogenous variables (variance is independent of other variables) to predict endogenous variables (variance is dependent on other variables) while allowing the variables to predict each other in the process [15].We used the R-Lavaan package (v4.1.1,R Core Team 2023) to define and run the model to predict ET using R n , VPD, temperature, windspeed, and precipitation.The output yields regression estimates, standard errors, z-values, and p-values for each variable for that specific analysis.The path coefficients (regression coefficients) represent the slope of the linear relationship between the response variable and predictor variables independent of all other variables.The model was run on the hourly, daily, and weekly timescales for a sufficiently large sample size and model performance confidence.Additionally, the results from the factor analysis were used to develop latent constructs, and a structural model was run to characterize the direct and indirect contributions of observed and latent variables in ET variability.
The statistical significance of the differences among path coefficients was tested using ANOVA for hourly, daily, and weekly timescales.For the group of variables with significant results, a Tukey HSD post hoc test was carried out to further examine the pairwise differences in path coefficients.Overall, using these analyses, we were able to categorize the variables/factors having similar variability patterns and quantify the independent contributions of each of those in influencing ET, as variables with large significant path coefficients can be interpreted to have a large influence on ET.

Akaike's Information Criteria (AIC) and Model Selection
An information criterion approach (AIC) was also used to evaluate how well our models fit the data.AIC modeling was used to develop candidate models with different combinations of potentially interacting variables that were also used in the SEM.By calculating and comparing the AIC scores of several candidate models, the model with the lowest AIC value and the highest AIC weight was selected as the best model using the AICmodavg package in R.

Factor Analysis
As a preliminary investigation of the relationships, Figure 2 shows the correlation matrices among all variables for all seasons at the hourly, daily, and weekly timescales.ET was moderately (| r | > 0.2) or highly correlated (| r | > 0.5) with R n , VPD, sensible heat, and temperature at all timescales.Precipitation was weakly correlated with ET at all timescales, with positively increasing coefficients at higher timescales.Interestingly, windspeed showed a weak negative correlation at the hourly timescale, with an increasing negative correlation at the daily and weekly timescales.The strength of the correlation between R n, H, and VPD was dependent on the temporal scale of the comparisons.Fur-ther, temperature and VPD showed a strong correlation with increasing coefficients at larger timescales.
Atmosphere 2024, 15, x FOR PEER REVIEW 6 of 18 ET was moderately (| r | > 0.2) or highly correlated (| r | > 0.5) with Rn, VPD, sensible heat, and temperature at all timescales.Precipitation was weakly correlated with ET at all timescales, with positively increasing coefficients at higher timescales.Interestingly, windspeed showed a weak negative correlation at the hourly timescale, with an increasing negative correlation at the daily and weekly timescales.The strength of the correlation between Rn, H, and VPD was dependent on the temporal scale of the comparisons.Further, temperature and VPD showed a strong correlation with increasing coefficients at larger timescales.
Atmosphere 2024, 15, x FOR PEER REVIEW 7 of 18 As evidenced by the correlation plots above, many of these variables share similar variability patterns and are strongly influenced by each other.This necessitated a method to quantify the independent contributions of each variable to the variability in ET.As described in Section 2.6, a factor analysis was carried out across all seasons and variables shown above (excluding RH) to segregate them into groups showing similar variability patterns.The first factor loadings for these variables across different seasons are shown in Figure 3.As evidenced by the correlation plots above, many of these variables share similar variability patterns and are strongly influenced by each other.This necessitated a method to quantify the independent contributions of each variable to the variability in ET.As described in Section 2.6, a factor analysis was carried out across all seasons and variables shown above (excluding RH) to segregate them into groups showing similar variability patterns.The first factor loadings for these variables across different seasons are shown in Figure 3.As evidenced by the correlation plots above, many of these variables share similar variability patterns and are strongly influenced by each other.This necessitated a method to quantify the independent contributions of each variable to the variability in ET.As described in Section 2.6, a factor analysis was carried out across all seasons and variables shown above (excluding RH) to segregate them into groups showing similar variability patterns.The first factor loadings for these variables across different seasons are shown in Figure 3.  ET, R n , sensible heat flux, temperature, and VPD all loaded highly positive on the first factor for all seasons except for the winter season, where ET and R n loaded comparatively lower than VPD and temperature.Also, sensible heat flux and pressure loaded highly negatively in winter.Windspeed and rain (precipitation) showed low factor loadings, as they become highly variable as the timescale increases and seasons change.
The factor analysis yielded multiple factor loadings, with the first factor capturing approximately 60% of the variability in the dataset.The first and second factor loadings for all variables are shown in score plots in Figure 4 with groupings similar to that in Figure 3. Factor loading pairs were categorized by color for seasons and markers for timescales (hourly, daily, and weekly).Score plots revealed the spread of first and second factor loadings across seasons and timescales.
Atmosphere 2024, 15, x FOR PEER REVIEW 9 of 18 The factor analysis yielded multiple factor loadings, with the first factor capturing approximately 60% of the variability in the dataset.The first and second factor loadings for all variables are shown in score plots in Figure 4 with groupings similar to that in Figure 3. Factor loading pairs were categorized by color for seasons and markers for timescales (hourly, daily, and weekly).Score plots revealed the spread of first and second factor loadings across seasons and timescales.Consistent with Figure 3, ET, Rn, temperature, and VPD loaded moderately or higher on the first factor for spring, summer, and autumn seasons.However, in winter, all three factors loaded moderately or higher in the second factor, while temperature loadings were low.For sensible heat flux and pressure, factor loadings were scattered, while windspeed and rain had loadings more clustered around the plot origin.

Path Analysis
While the factor analysis helped to identify variables with similar variability patterns, the correlations deduced from this analysis alone do not represent the relative importance of one variable in determining ET variability, independent of the effect of other interrelated variables.Hence, a combination of a path analysis with a factor analysis is helpful to demonstrate the direct dependencies of ET on one variable independent of other variables, whereby regression coefficients represent the strength of independent contributions.The direct and indirect dependencies and associations between different variables in terms of physical pathways derived from the path analysis are shown in Figure 5. Sensible heat flux and pressure were not included in the analysis due to poor model fit.Consistent with Figure 3, ET, R n , temperature, and VPD loaded moderately or higher on the first factor for spring, summer, and autumn seasons.However, in winter, all three factors loaded moderately or higher in the second factor, while temperature loadings were low.For sensible heat flux and pressure, factor loadings were scattered, while windspeed and rain had loadings more clustered around the plot origin.

Path Analysis
While the factor analysis helped to identify variables with similar variability patterns, the correlations deduced from this analysis alone do not represent the relative importance of one variable in determining ET variability, independent of the effect of other interrelated variables.Hence, a combination of a path analysis with a factor analysis is helpful to demonstrate the direct dependencies of ET on one variable independent of other variables, whereby regression coefficients represent the strength of independent contributions.The direct and indirect dependencies and associations between different variables in terms of physical pathways derived from the path analysis are shown in Figure 5. Sensible heat flux and pressure were not included in the analysis due to poor model fit.Consistent with Figure 3, ET, Rn, temperature, and VPD loaded moderately or higher on the first factor for spring, summer, and autumn seasons.However, in winter, all three factors loaded moderately or higher in the second factor, while temperature loadings were low.For sensible heat flux and pressure, factor loadings were scattered, while windspeed and rain had loadings more clustered around the plot origin.

Path Analysis
While the factor analysis helped to identify variables with similar variability patterns, the correlations deduced from this analysis alone do not represent the relative importance of one variable in determining ET variability, independent of the effect of other interrelated variables.Hence, a combination of a path analysis with a factor analysis is helpful to demonstrate the direct dependencies of ET on one variable independent of other variables, whereby regression coefficients represent the strength of independent contributions.The direct and indirect dependencies and associations between different variables in terms of physical pathways derived from the path analysis are shown in Figure 5. Sensible heat flux and pressure were not included in the analysis due to poor model fit.5a,b.While correlation and factor loading plots showed Rn having the strongest correlation and shared variability pattern with ET, its direct independent influence on ET was lower.Similarly, Rn not only had a comparable direct effect on ET but also a strong indirect effect via VPD and temperature.VPD showed a strong positive influence on ET, while in reverse, ET had a strong negative impact on VPD.VPD showed the largest impact on ET at the hourly timescale.Air temperature was strongly driven by Rn, which then influenced the VPD, ultimately affecting ET in the process.Windspeed had a negligible negative effect on ET, while precipitation (rain) had a feeble positive impact on ET, although windspeed showed a negative impact on temperature, thereby affecting ET indirectly as well.
In the daily and weekly timescales, these regression coefficients generally increased for all the pathways.For example, Rn emerged as a major direct influencer of ET, while the strengths of indirect linkages through VPD and temperature also increased at larger timescales.The direct effect of rain on ET increased positively, while that of windspeed decreased overall.
The combined results of the path analysis with significant regression coefficients for different variables at hourly, daily, and weekly timescales are shown in Figure 6.Rn stood out, with the highest median regression coefficient for all the timescales.VPD had the second-highest independent contribution to ET variability on the hourly scale, while rain and temperature showed the same at the daily and weekly timescales, respectively.The windspeed showed lesser influence at the hourly timescale, while this increased for the daily and weekly timescales.To test the significance of differences in the regression coefficients of different variables across these timescales, an ANOVA test was performed for all the categories.The results showed significant differences (p < 0.05) in regression coefficients on all timescales for all variables.This was followed by a Tukey HSD post hoc analysis to test the significance for pairwise differences in regression coefficients.Only the regression coefficients for Rn showed significant differences from that of other variables.Thus, it was concluded that Rn is the major direct driver of ET variability in this ecosystem across all timescales.5a,b.While correlation and factor loading plots showed R n having the strongest correlation and shared variability pattern with ET, its direct independent influence on ET was lower.Similarly, R n not only had a comparable direct effect on ET but also a strong indirect effect via VPD and temperature.VPD showed a strong positive influence on ET, while in reverse, ET had a strong negative impact on VPD.VPD showed the largest impact on ET at the hourly timescale.Air temperature was strongly driven by R n , which then influenced the VPD, ultimately affecting ET in the process.Windspeed had a negligible negative effect on ET, while precipitation (rain) had a feeble positive impact on ET, although windspeed showed a negative impact on temperature, thereby affecting ET indirectly as well.
In the daily and weekly timescales, these regression coefficients generally increased for all the pathways.For example, R n emerged as a major direct influencer of ET, while the strengths of indirect linkages through VPD and temperature also increased at larger timescales.The direct effect of rain on ET increased positively, while that of windspeed decreased overall.
The combined results of the path analysis with significant regression coefficients for different variables at hourly, daily, and weekly timescales are shown in Figure 6.R n stood out, with the highest median regression coefficient for all the timescales.VPD had the second-highest independent contribution to ET variability on the hourly scale, while rain and temperature showed the same at the daily and weekly timescales, respectively.The windspeed showed lesser influence at the hourly timescale, while this increased for the daily and weekly timescales.To test the significance of differences in the regression coefficients of different variables across these timescales, an ANOVA test was performed for all the categories.The results showed significant differences (p < 0.05) in regression coefficients on all timescales for all variables.This was followed by a Tukey HSD post hoc analysis to test the significance for pairwise differences in regression coefficients.Only the regression coefficients for R n showed significant differences from that of other variables.Thus, it was concluded that R n is the major direct driver of ET variability in this ecosystem across all timescales.Similar to Figure 6, the path analysis results were compared for all the variabl significant path coefficients across seasons.The distributions of path coefficients f variable during the spring, summer, autumn, and winter seasons are shown in Fi Rn had the highest regression coefficient for the autumn season, followed by summ spring.In winter, the independent contribution of Rn was the lowest in determin variability in ET.However, VPD showed the least contribution to ET variability autumn season, which became maximum in the summer season.This revealed th dependent driving force in summer seasons when there is an abundant amount o vapor in the atmosphere.The effect of temperature was observed to be greater dur Similar to Figure 6, the path analysis results were compared for all the variables with significant path coefficients across seasons.The distributions of path coefficients for each variable during the spring, summer, autumn, and winter seasons are shown in Figure 7. R n had the highest regression coefficient for the autumn season, followed by summer and spring.In winter, the independent contribution of R n was the lowest in determining the variability in ET.However, VPD showed the least contribution to ET variability in the autumn season, which became maximum in the summer season.This revealed the VPD-dependent driving force in summer seasons when there is an abundant amount of water vapor in the atmosphere.The effect of temperature was observed to be greater during the spring and winter seasons compared to the summer and autumn seasons.Windspeed and rain had the highest impact during the summer season.The lesser effect of rain on ET variability during spring and autumn suggested there was sufficient water available for ET.

AIC Results
The results from the AIC model selection suggested that the model with ET as the dependent variable and Rn, VPD, temperature, windspeed, and precipitation independent variables was the best-fit model, given the data.All other models with possible combinations and meaningful interactions were tested and excluded based on the number of parameters (K > 10).The best-fit model, carrying 98% of the cumulative model weight, included the above-mentioned variables with meaningful interaction effects (Table 1).

AIC Results
The results from the AIC model selection suggested that the model with ET as the dependent variable and R n , VPD, temperature, windspeed, and precipitation independent variables was the best-fit model, given the data.All other models with possible combinations and meaningful interactions were tested and excluded based on the number of parameters (K > 10).The best-fit model, carrying 98% of the cumulative model weight, included the above-mentioned variables with meaningful interaction effects (Table 1).The results show that the best model was the interaction model (ET ~Rn *Ta + VPD + WS*P), which corroborates with the findings from the SEM analysis.The 'best' model carried 98% of the cumulative model weight and had the lowest AIC score.The next-best model was more than 2 AIC units higher than the best model (8.25 units) and had only 2% of the cumulative model weight; it was deemed unimportant since it did not meaningfully 'add' to the amount of information already explained by the best-fit model.

Discussion
The SEM approach has been used to diagnose drivers of ET using high-frequency data collected using the EC method for forested systems.The results of this study in a BHF in the RSWMA revealed that the variability in ET is directly influenced by R n during spring, summer, and autumn, primarily vegetatively active seasons.However, during times of vegetative dormancy (i.e., in winter), the variability in ET is largely influenced by VPD and temperature, indirect controls of R n .These results are consistent with the strong seasonal cycle for the variables that gradually increase from winter to summer and gradually diminish from summer to winter.This typical seasonal cycle also suggests that the influence of temperature and VPD on ET is indirectly driven by R n .The greater control of temperature and VPD as drivers of ET during winter is suggestive of temperaturedependent ET, especially when the R n -dependent direct control on ET is minimal.Hence, the direct as well as indirect (through greater control of R n on temperature and VPD) control of R n on ET across all seasons reinforces the role of R n as a primary driver of ET in this forest.
The results from the best-fit model in the AIC modeling support the inference from SEM that there is an interaction between R n and temperature.Furthermore, VPD showed independent influence as a predictor of ET compared to its interaction with R n .This corroborates the increased influence of temperature on ET during vegetatively dormant seasons as opposed to the higher impact of R n on ET during vegetatively active seasons.Although independent, the positive impact of WS on ET has been previously established in many other studies [15,24]; the interaction of WS with other variables like temperature and precipitation is likely to attenuate dryness (latent construct) and augment wetness, subsequently regulating the rate of ET in this system.The best model from AIC suggested that these interactions are highly plausible and significant in determining the dynamics of ET in this BHF.
Consistent with the findings from the data collected in a BHF in the RSWMA from 2014 to 2021, which showed R n as the major direct and indirect driver of ET variability across different times of the year, research reports from most other forest types have characterized R n as a primary driver of ET.For example, [15] carried out an SEM analysis in various boreal, tundra, and permafrost ecosystems of high-latitude regions to demonstrate that R n is the major driver of ET variability, albeit having a smaller independent contribution due to its control on other variables.On the contrary, a similar path analysis conducted in a mid-latitude agricultural site in northern China reported that R n had the largest direct independent contribution to ET [17].Brown [10] suggested that the increase in the amount of R n received by a BHF in Missouri resulted in higher ET.Its vegetation composition, however, was silver maple (Acer saccharinum), eastern cottonwood (Populus deltoides), boxelder (Acer negundo), and sycamore (Platanus occidentalis), which is mostly different from that found around the US-ULM tower in the RSWMA.Mackay [30] reported R n as a major driver in upland hardwood growth forests and VPD as a major driver in wetland ecosystems during the vegetatively active time of the year in northern Wisconsin.The dominant hardwood vegetation composition at this site was sugar maple (Acer saccharum), basswood (Tilia americana L.), and green ash (Fraxinus pennsylvanica Marsh), mostly different from the hardwood community found in the RSWMA.However, a much higher effect of precipitation has been reported in water-limited ecosystems such as BHFs and seasonal cycles of canopy greenness in energy-limited ecosystems in higher latitudes using path analyses [31].Similarly, in humid boreal regions, VPD and R n were characterized as major drivers of sap flow and thus transpiration during the growing season, as well as during drought [32].
In closed-canopy deciduous BHFs like the one in this study, R n controls the variability of ET through two different pathways: first is the direct pathway, in which throughout the growing period (spring to summer) and early autumn, R n directly promotes transpiration, which contributes about 80-90% of the total ET, as shown by [24]; second is the indirect pathway, in which when there are no leaves in the vegetatively dormant period, the direct impact of R n is somehow attenuated, and the R n influences ET variability indirectly via temperature and VPD.This is also supported by the consistent variability pattern shown by R n and ET in seasonal plots (except winter) in Figure 3, and as reported by others [33,34].
Seasonally, WS and precipitation play a critical role during summer when the atmospheric humidity is higher compared to other seasons, albeit a small one.The increase in WS and precipitation positively influences wetness (one of the latent constructs), subsequently affecting ET negatively in the process (path coefficients = 0.70 and −0.26 for dryness-ET and wetness-ET, respectively).However, the negative effect of WS on ET, as seen in Figure 5a, could be due to its significant negative impact on temperature, which in turn has a strong influence on ET through the VPD.As observed in the structural model in Figure 5b, WS contributes significantly as a major driver of atmospheric wetness, a latent construct with a significant negative effect on ET.In a study by Lobos-Roco [35] in the Atacama Desert ecosystem of Chile, it was shown that strong winds in the afternoon enhance mechanical turbulence and increase evaporation.A similar path analysis conducted in a mid-latitude agricultural site in China found similar results, with WS having the least direct and indirect effects on ET [17].On the other hand, precipitation has a minor positive impact as a driver of ET, as seen in the correlograms (Figure 2).However, there is an increasing negative impact of precipitation, as seen in Figure 3, and it has a significant negative contribution to latent construct wetness; hence, over larger timescales the impact is suggestive of a seasonal reduction in ET, as precipitation is observed for longer periods, especially during growing season.From similar research carried out at high-latitude regions, Thunberg [15] reported similar control of precipitation as a driver of ET and suggested it as a potentially relevant driver of ET in mid-latitude regions such as these BHFs.These findings of seasonal relations of ET with meteorological variables are consistent with those from similar studies in Canadian forest ecosystems [36].This also strengthens the conclusion that precipitation contributes positively to short-term enhancement in ET and negatively in the long run.These conclusions can be further strengthened by the simultaneous measurement of soil water content and heat flux on the site, one of the limitations of this study.
Since temperature and VPD are largely controlled by R n , the importance of these variables as drivers of ET is more complex to understand, multifaceted, and largely dependent on the direct and indirect influence of the timing, duration, and intensity of solar energy in association with the seasonal phenological characteristics.The results from the factor analysis suggested that their importance as drivers of ET becomes more prominent only during vegetatively dormant seasons, when the direct control of R n on ET remains lower.However, these thermal variables share communality in their variability patterns and have a greater impact on ET by increasing the dryness (one of the latent constructs in the structural model) of the atmosphere, as observed in Figure 5b.The ambient temperature not only positively influences the VPD but also negatively impacts R n .This could be due to a lower retention of incoming radiation as the canopy becomes saturated with heat and a higher loss of longwave radiation from the canopy as the temperature of the canopy increases.For example, with the increase in R n before noon, all variables, ET, temperature, and VPD, increase consistently until the canopy becomes saturated with heat, thereafter leading to a decrease in R n and ET in the afternoon while the temperature and VPD increase further.This has implications for the stomatal regulation of water loss, GPP, and canopy temperature regulation in this forested system, as also suggested by [18].The independent contributions of several other phenological and hydro-meteorological variables including temperature and VPD need to be further investigated to better understand the key role played by these variables as drivers of ET in this forested system.

Conclusions
During this study, the EC method was used to measure the water fluxes from 2014 to 2021 in the BHF of the RSWMA in northeast Louisiana.The results show that R n is a major direct driver of ET during the vegetatively active season and an indirect driver of ET via temperature and VPD during the vegetatively dormant season, implying a crucial role played by phenological changes in the process.The timing and duration of WS and precipitation play a significant role as drivers of ET in different seasons, albeit with smaller independent contributions compared to other variables.This research needs to be expanded further to include other types of BHFs in the LMAV region to create a holistic understanding of water use patterns that can have implications for developing better-informed strategies for natural resource management in the whole region.This is even important in light of the findings of this research, which suggest that unusual changes in weather patterns leading to unprecedented anomalies in temperature, VPD, and precipitation patterns could alter the water use dynamics, exacerbating the deviations from the typical water cycle in this system.The advancement of machine learning and robust data analytical methods and their ability to process multiple variables simultaneously to disentangle the complex relationships among variables has the potential to accomplish this goal, the results of which have implications for more accurate modeling and forecasting of water use dynamics for informed natural resource management in the region.In light of climate change, the ecosystem models studying the potential impact of unusual climatic conditions on water fluxes can help better plan and prepare various mitigating efforts for the future of these temperate forests.

Figure 1 .
Figure 1.(a) Location map showing the study area, the Russel Sage Wildlife Management Area in northeast Louisiana, with the position of the US-ULM tower location indicated by the arrow tip, and (b) the study site shown flooded, as is typical during the late winter and early spring (Photo: J.B.).

Figure 1 .
Figure 1.(a) Location map showing the study area, the Russel Sage Wildlife Management Area in northeast Louisiana, with the position of the US-ULM tower location indicated by the arrow tip, and (b) the study site shown flooded, as is typical during the late winter and early spring (Photo: J.B.).

18 Figure 3 .
Figure 3. Factor loadings for evapotranspiration (ET), net radiation (Rn), sensible heat flux (H), temperature of air (Ta), Vapor Pressure Deficit (VPD), windspeed (WS), precipitation (P), and pressure (Pa) on the first pattern of factor analysis for spring, summer, autumn, and winter seasons.ET, Rn, sensible heat flux, temperature, and VPD all loaded highly positive on the first factor for all seasons except for the winter season, where ET and Rn loaded comparatively lower than VPD and temperature.Also, sensible heat flux and pressure loaded highly negatively in winter.Windspeed and rain (precipitation) showed low factor loadings, as they become highly variable as the timescale increases and seasons change.

Figure 4 .
Figure 4. Score plots of the first and second factor loadings by various timescales across different seasons.Variables include evapotranspiration (ET), net radiation (R n ), sensible heat flux (H), Vapor Pressure Deficit (VPD), temperature of air (Ta), pressure (Pa), windspeed (WS), and precipitation (P).

Figure 5 .
Figure 5.The (a) path diagram and (b) SEM of all the observed variables (ET, Rn, Ta, VPD, WS, and P) and latent constructs (dryness and wetness) within the path analysis with positive (black) and negative (red) path coefficients for hourly timescale across all seasons.The dashed line indicates direct contribution.The * sign indicates significance level and n.s.indicates non-significance.The standardized solutions for the best-fit path model and structural model [χ 2 = 752 (df = 4, n = 42471), p < 0.0001], Comparative Fit Index ((CFI) = 0.98), and Tucker-Lewis Index ((TLI) = 0.96, RMSEA = 0.06, SRMR = 0.02) are shown in Figure5a,b.While correlation and factor loading plots showed Rn having the strongest correlation and shared variability pattern with ET, its direct independent influence on ET was lower.Similarly, Rn not only had a comparable direct effect on ET but also a strong indirect effect via VPD and temperature.VPD showed a strong positive influence on ET, while in reverse, ET had a strong negative impact on VPD.VPD showed the largest impact on ET at the hourly timescale.Air temperature was strongly driven by Rn, which then influenced the VPD, ultimately affecting ET in the process.Windspeed had a negligible negative effect on ET, while precipitation (rain) had a feeble positive impact on ET, although windspeed showed a negative impact on temperature, thereby affecting ET indirectly as well.In the daily and weekly timescales, these regression coefficients generally increased for all the pathways.For example, Rn emerged as a major direct influencer of ET, while the strengths of indirect linkages through VPD and temperature also increased at larger timescales.The direct effect of rain on ET increased positively, while that of windspeed decreased overall.The combined results of the path analysis with significant regression coefficients for different variables at hourly, daily, and weekly timescales are shown in Figure6.Rn stood out, with the highest median regression coefficient for all the timescales.VPD had the second-highest independent contribution to ET variability on the hourly scale, while rain and temperature showed the same at the daily and weekly timescales, respectively.The windspeed showed lesser influence at the hourly timescale, while this increased for the daily and weekly timescales.To test the significance of differences in the regression coefficients of different variables across these timescales, an ANOVA test was performed for all the categories.The results showed significant differences (p < 0.05) in regression coefficients on all timescales for all variables.This was followed by a Tukey HSD post hoc analysis to test the significance for pairwise differences in regression coefficients.Only the regression coefficients for Rn showed significant differences from that of other variables.Thus, it was concluded that Rn is the major direct driver of ET variability in this ecosystem across all timescales.

Figure 5 .
Figure 5.The (a) path diagram and (b) SEM of all the observed variables (ET, R n , T a , VPD, WS, and P) and latent constructs (dryness and wetness) within the path analysis with positive (black) and negative (red) path coefficients for hourly timescale across all seasons.The dashed line indicates direct contribution.The * sign indicates significance level and n.s.indicates non-significance.The standardized solutions for the best-fit path model and structural model [χ 2 = 752 (df = 4, n = 42471), p < 0.0001], Comparative Fit Index ((CFI) = 0.98), and Tucker-Lewis Index ((TLI) = 0.96, RMSEA = 0.06, SRMR = 0.02) are shown in Figure5a,b.While correlation and factor loading plots showed R n having the strongest correlation and shared variability pattern with ET, its direct independent influence on ET was lower.Similarly, R n not only had a comparable direct effect on ET but also a strong indirect effect via VPD and temperature.VPD showed a strong positive influence on ET, while in reverse, ET had a strong negative impact on VPD.VPD showed the largest impact on ET at the hourly timescale.Air temperature was strongly driven by R n , which then influenced the VPD, ultimately affecting ET in the process.Windspeed had a negligible negative effect on ET, while precipitation (rain) had a feeble positive impact on ET, although windspeed showed a negative impact on temperature, thereby affecting ET indirectly as well.In the daily and weekly timescales, these regression coefficients generally increased for all the pathways.For example, R n emerged as a major direct influencer of ET, while the strengths of indirect linkages through VPD and temperature also increased at larger timescales.The direct effect of rain on ET increased positively, while that of windspeed decreased overall.The combined results of the path analysis with significant regression coefficients for different variables at hourly, daily, and weekly timescales are shown in Figure6.R n stood out, with the highest median regression coefficient for all the timescales.VPD had the second-highest independent contribution to ET variability on the hourly scale, while rain and temperature showed the same at the daily and weekly timescales, respectively.The windspeed showed lesser influence at the hourly timescale, while this increased for the daily and weekly timescales.To test the significance of differences in the regression coefficients of different variables across these timescales, an ANOVA test was performed for all the categories.The results showed significant differences (p < 0.05) in regression coefficients on all timescales for all variables.This was followed by a Tukey HSD post hoc analysis to test the significance for pairwise differences in regression coefficients.Only the regression coefficients for R n showed significant differences from that of other variables.Thus, it was concluded that R n is the major direct driver of ET variability in this ecosystem across all timescales.

Figure 6 .
Figure 6.Distributions of regression coefficients from all significant results of path analysis hourly, daily, and weekly timescales with sample size (n = 12) for each variable.Black so represent median values and boxes represent interquartile range.

Figure 6 .
Figure 6.Distributions of regression coefficients from all significant results of path analysis SEM at hourly, daily, and weekly timescales with sample size (n = 12) for each variable.Black solid lines represent median values and boxes represent interquartile range.

Figure 7 .
Figure 7. Distributions of regression coefficients from all significant results of path analysis SEM across different seasons with sample size (n = 12) for each variable.Black solid lines represent median values and boxes represent interquartile range.

Figure 7 .
Figure 7. Distributions of regression coefficients from all significant results of path analysis SEM across different seasons with sample size (n = 12) for each variable.Black solid lines represent median values and boxes represent interquartile range.

Table 1 .
Results of AIC model selection with shortlisted candidate models with the number of parameters (K), Akaike's Information Criteria (AIC), delta AIC, AIC weights, and log-likelihood (LL) values.

Table 1 .
Results of AIC model selection with shortlisted candidate models with the number of parameters (K), Akaike's Information Criteria (AIC), delta AIC, AIC weights, and log-likelihood (LL) values.