Spatial and Temporal Variability of Throughfall among Oak and Co-Occurring Non-Oak Tree Species in an Upland Hardwood Forest

Canopy throughfall comprises the largest portion of net precipitation that is delivered to the forest floor. This water flux is highly variable across space and time and is influenced by species composition, canopy foliage, stand structure, and storm meteorological characteristics. In upland forests throughout the central hardwoods region of the Eastern United States, a compositional shift is occurring from oak-hickory to more mesic, shade-tolerant species such as red maple, sweetgum, and winged elm. To better understand the impacts of this shift on throughfall flux and the hydrologic budget, we monitored throughfall for one year in Northern Mississippi under the crowns of midstory and overstory oak (post oak and southern red oak) and non-oak species (hickory, red maple, and winged elm). In general, oak had more throughfall than co-occurring non-oak species in both canopy levels. In the overstory during the leaf-off canopy phase, white oak had relatively higher throughfall partitioning (standardized z-score = 0.54) compared to all other species (z-score = −0.02) (p = 0.004), while in the leaf-on canopy phase, red maple had relatively lower throughfall (z-score = −0.36) partitioning compared to all other species (z-score = 0.11). In the midstory, red maple was the only species to exhibit a difference in throughfall between canopy phases, with much lower throughfall in the leaf-off compared to the leaf-on canopy phase (z-score = −0.30 vs. 0.202, p = 0.039). Additionally, throughfall under oak crowns was less variable than under non-oak crowns. These results provide evidence that the spatial and temporal distribution of throughfall inputs under oak crowns are different than non-oak species, likely due to differences in crown architecture (i.e., depth and density). As oak dominance diminishes in these forests, it is possible that the portion of rainfall diverted to throughfall may decrease as well. The net impacts to watershed hydrology are still unknown, but these results provide one mechanism by which the distribution of water resources may be affected.


Introduction
In upland hardwood forests throughout the Eastern United States, there has been a shift in dominance from oak (Quercus spp.) to shade-tolerant species such as maples (both Acer rubrum and A. saccharum), ashes (Fraxinus spp.), and elms (Ulmus spp.) [1][2][3]. This shift may be caused by "mesophication" a hypothesized feedback loop by which systematic exclusion of the natural fire regime has allowed fire-sensitive, and often shade-tolerant, species (aka "mesophytes") to establish and create cooler, moister, and darker understory conditions that are unfavorable for oak regeneration and survival [4]. These compositional shifts may impact the hydrologic budget at the watershed scale. For example, mesophytes have higher rates of transpiration compared to co-occurring oaks, such that Table 1. Overstory and midstory species (abbreviations given in parentheses) selected for study with mean ± standard error of diameter at breast height (DBH), basal area (BA), crown area (CA), and density. Throughfall gauges were constructed from 3.8 L polyethylene collectors fitted with a 12.5 cm diameter funnel. Underneath each overstory tree, one throughfall gauge was deployed adjacent to the tree bole and a second throughfall gauge was deployed at distance midway between the tree bole and crown edge, termed "mid-crown". In the midstory, due to the small size of tree crowns, a single gauge was deployed at the mid-crown point. Throughfall volume was measured at 4-week intervals, after which the gauges were emptied and placed back in the same position.

Species
Monthly throughfall volume was converted into a depth equivalent through the equation: where D is the depth equivalent in cm, V is the volume of throughfall in cm 3 , and A is the orifice area of the throughfall gauge in cm 2 . In order to compare across sampling periods and rainfall amounts, throughfall depth was standardized to a z-score through the equation: where TF z is the throughfall z-score for each gauge during the sampling interval, µ D is the mean depth of all throughfall gauges for a sampling interval, and σ D is the standard deviation of all throughfall gauges for a sampling interval.
Linear mixed models were used to explore the impact of species, canopy leaf phase, and collector position on TF z . In the overstory, species, canopy leaf phase, and collector positions and their interactions were fixed effects, while individual trees within a species were a random effect. We chose to use a mixed model in order to incorporate individual trees and their repeated measurement as random effects. In this way, we are able to explain variation across species while the model accounts for variation within a species. Models were run using the lmer function in R and fit using the restricted maximum-likelihood criterion [30]. Least-Squares Means post-hoc pairwise comparisons were performed when statistical differences in sample means were detected (α = 0.05) using the emmeans package in R [31]. In the midstory, only species and canopy leaf phase were set as fixed effects because there was only one gauge per tree. Prior to model development, TF z was tested to meet the assumptions of normality and homogeneity of variance. In the overstory, TF z was not normally distributed, so it was cube root transformed. Transformed TF z values were used for statistical analysis but un-transformed TF z values are presented in figures for ease of interpretation.
To assess throughfall variability, the coefficient of variation (CV) was calculated with the equation: where σ is the standard deviation of throughfall depth and µ is the mean throughfall depth for each species during each sampling interval, and by collector location for overstory species. Linear models were used to explore the impact of throughfall depth, species, leaf canopy phase, and collector position on CV.

General Throughfall Trends
From winter 2017/18 through winter 2018/19, throughfall was monitored approximately monthly for a total of 12 sampling periods. Total rainfall over the period was 203.3 cm and total throughfall, averaged across species and gauge location, was 181.3 cm (89.2%) ( Table 2). During the sampling period November 2017 through February 2019, more throughfall than rainfall was recorded at the study site, which could be due to a highly localized storm cell that impacted the study site but not the meteorological station. Linear mixed models were developed separately for overstory and midstory species to test the effects of canopy leaf phase and gauge position on throughfall distribution. In the overstory, TF z differed between species (p = 0.025) and as an interaction between species and canopy leaf phase (p = 0.006) (Table 3). Overall, white oaks had higher TF z than other species (Figure 1), which means that partitioning of throughfall by white oaks was greater than the average across both the leaf-off and leaf-on canopy phases. Specifically, in the leaf-off canopy phase, white oak had relatively higher throughfall partitioning compared to all other species. In the leaf-on canopy phase, white oak had more similar throughfall partitioning while, red maple had relatively lower throughfall partitioning compared to all other species (Figure 1). All species except southern red oak exhibited a decrease in TF z from the leaf-on canopy phase to the leaf-off canopy phase, with the largest decrease in white oak. The primary interaction of species and leaf phase was between southern red oak and white oak, whereby TF z was different between the two species in the leaf-off canopy phase (p = 0.004) but not during the leaf-on canopy phase (p = 0.923). Among species, mid-canopy throughfall gauges had marginally greater standardized throughfall than gauges located next to tree boles (mid-canopy = 0.152 vs. Bole = 0.035, p = 0.104). Linear mixed models were developed separately for overstory and midstory species to test the effects of canopy leaf phase and gauge position on throughfall distribution. In the overstory, TFz differed between species (p = 0.025) and as an interaction between species and canopy leaf phase (p = 0.006) (Table 3). Overall, white oaks had higher TFz than other species (Figure 1), which means that partitioning of throughfall by white oaks was greater than the average across both the leaf-off and leaf-on canopy phases. Specifically, in the leaf-off canopy phase, white oak had relatively higher throughfall partitioning compared to all other species. In the leaf-on canopy phase, white oak had more similar throughfall partitioning while, red maple had relatively lower throughfall partitioning compared to all other species ( Figure 1). All species except southern red oak exhibited a decrease in TFz from the leaf-on canopy phase to the leaf-off canopy phase, with the largest decrease in white oak. The primary interaction of species and leaf phase was between southern red oak and white oak, whereby TFz was different between the two species in the leaf-off canopy phase (p = 0.004) but not during the leaf-on canopy phase (p = 0.923). Among species, mid-canopy throughfall gauges had marginally greater standardized throughfall than gauges located next to tree boles (mid-canopy = 0.152 vs. Bole = 0.035, p = 0.104).  In the midstory, TFz was not different among species (p = 0.233) and was only marginally different due to canopy leaf phase (p = 0.077), with clear differences exhibited as an interaction between species and leaf phase (p = 0.052). In the leaf-off canopy phase, TFz was lower (−0.225) compared to the leaf-on canopy phase (−0.078). Red maple was the only species to exhibit a difference In the midstory, TF z was not different among species (p = 0.233) and was only marginally different due to canopy leaf phase (p = 0.077), with clear differences exhibited as an interaction between species and leaf phase (p = 0.052). In the leaf-off canopy phase, TF z was lower (−0.225) compared to the leaf-on canopy phase (−0.078). Red maple was the only species to exhibit a difference in throughfall partitioning between canopy phases, with much lower throughfall in the leaf-off canopy phase compared to the leaf-on (−0.303 vs. 0.202, p = 0.039) (Figure 2). in throughfall partitioning between canopy phases, with much lower throughfall in the leaf-off canopy phase compared to the leaf-on (−0.303 vs. 0.202, p = 0.039) (Figure 2).

Species-Level Variability
For overstory species, the coefficient of variation in throughfall was lower in the leaf-off canopy phase than in the leaf-on canopy phase (11.0 vs 16.5) (p = 0.002) (Figure 3). Red maple and white oak had the highest overall CVs (15.5 and 15.1, respectively), followed by hickory (12.5) and southern red oak (11.9), although there was no difference in the regression equations among species (p = 0.120) ( Figure 3). In the leaf-on canopy phase, throughfall depth was a moderate predictor of CV (r 2 = 0.11, p = 0.068, y = −0.35x + 20.88) such that as throughfall increased, the variability among gauges decreased. However, throughfall depth was not a strong predictor of CV in the leaf-off canopy phase (r 2 = 0.00, p = 0.466).

Species-Level Variability
For overstory species, the coefficient of variation in throughfall was lower in the leaf-off canopy phase than in the leaf-on canopy phase (11.0 vs. 16.5) (p = 0.002) (Figure 3). Red maple and white oak had the highest overall CVs (15.5 and 15.1, respectively), followed by hickory (12.5) and southern red oak (11.9), although there was no difference in the regression equations among species (p = 0.120) (Figure 3). In the leaf-on canopy phase, throughfall depth was a moderate predictor of CV (r 2 = 0.11, p = 0.068, y = −0.35x + 20.88) such that as throughfall increased, the variability among gauges decreased. However, throughfall depth was not a strong predictor of CV in the leaf-off canopy phase (r 2 = 0.00, p = 0.466).

Species-Level Variability
For overstory species, the coefficient of variation in throughfall was lower in the leaf-off canopy phase than in the leaf-on canopy phase (11.0 vs 16.5) (p = 0.002) (Figure 3). Red maple and white oak had the highest overall CVs (15.5 and 15.1, respectively), followed by hickory (12.5) and southern red oak (11.9), although there was no difference in the regression equations among species (p = 0.120) (Figure 3). In the leaf-on canopy phase, throughfall depth was a moderate predictor of CV (r 2 = 0.11, p = 0.068, y = −0.35x + 20.88) such that as throughfall increased, the variability among gauges decreased. However, throughfall depth was not a strong predictor of CV in the leaf-off canopy phase (r 2 = 0.00, p = 0.466).   In the midstory, red maple and hickory had the greatest CV (15.8 and 14.6, respectively), followed by white oak (13.0), winged elm (12.4), and southern red oak (7.8). Similar to the overstory, CV was greater in the leaf-on canopy phase (14.3) compared to the leaf-off canopy phase (11.3). In the leaf-off canopy phase, there was essentially no relationship between throughfall quantity and throughfall CV (r 2 = 0.00, p = 0.717) (Figure 4). In the leaf-on canopy phase, all species exhibited a decrease in CV with increasing throughfall quantity (r 2 = 0.28, p = 0.060, y = −0.04x + 4.67). Interestingly, the CV of southern red oak increased with increasing throughfall, leading to a predictive model significantly different than the other species (p = 0.012).
Geosciences 2019, 9, 405 7 of 13 In the midstory, red maple and hickory had the greatest CV (15.8 and 14.6, respectively), followed by white oak (13.0), winged elm (12.4), and southern red oak (7.8). Similar to the overstory, CV was greater in the leaf-on canopy phase (14.3) compared to the leaf-off canopy phase (11.3). In the leaf-off canopy phase, there was essentially no relationship between throughfall quantity and throughfall CV (r 2 = 0.00, p = 0.717) (Figure 4). In the leaf-on canopy phase, all species exhibited a decrease in CV with increasing throughfall quantity (r 2 = 0.28, p = 0.060, y = −0.04x + 4.67). Interestingly, the CV of southern red oak increased with increasing throughfall, leading to a predictive model significantly different than the other species (p = 0.012).

Gauge-Level Variability
During the leaf-off canopy phase, throughfall gauges underneath crowns of white oaks consistently had higher standardized throughfall quantities ( Figure 5A,B). For gauges located next to tree boles, six out of eight of the highest median z-scores were white oaks, while red oaks and hickories often had z-scores near zero ( Figure 5A). For gauges underneath the mid-crown point, five of the top eight gauges were underneath white oak crowns ( Figure 5B). Gauges underneath hickory crowns ranked lower at the mid-crown point compared to next to tree boles. During the leaf-on canopy phase, non-oak species had the lowest z-scores in four out of five gauges located next to tree boles ( Figure 5C) and at mid-crown points ( Figure 5D).

Gauge-Level Variability
During the leaf-off canopy phase, throughfall gauges underneath crowns of white oaks consistently had higher standardized throughfall quantities ( Figure 5A,B). For gauges located next to tree boles, six out of eight of the highest median z-scores were white oaks, while red oaks and hickories often had z-scores near zero ( Figure 5A). For gauges underneath the mid-crown point, five of the top eight gauges were underneath white oak crowns ( Figure 5B). Gauges underneath hickory crowns ranked lower at the mid-crown point compared to next to tree boles. During the leaf-on canopy phase, non-oak species had the lowest z-scores in four out of five gauges located next to tree boles ( Figure 5C) and at mid-crown points ( Figure 5D). In the midstory, southern red oak throughfall gauges had the highest z-scores in the leaf-off canopy phase but not in the leaf-on canopy phase ( Figure 6). In the leaf-off canopy phase, winged elm and hickory dominated the number of gauges with low z-scores ( Figure 6A). In the leaf-on canopy phase, seven out of 10 red maple gauges had higher than average z-scores ( Figure 6B).  In the midstory, southern red oak throughfall gauges had the highest z-scores in the leaf-off canopy phase but not in the leaf-on canopy phase ( Figure 6). In the leaf-off canopy phase, winged elm and hickory dominated the number of gauges with low z-scores ( Figure 6A). In the leaf-on canopy phase, seven out of 10 red maple gauges had higher than average z-scores ( Figure 6B). In the midstory, southern red oak throughfall gauges had the highest z-scores in the leaf-off canopy phase but not in the leaf-on canopy phase ( Figure 6). In the leaf-off canopy phase, winged elm and hickory dominated the number of gauges with low z-scores ( Figure 6A). In the leaf-on canopy phase, seven out of 10 red maple gauges had higher than average z-scores ( Figure 6B).

Discussion
At our study site, southern red oak had the largest diameter at breast height and occupied the most basal area, but had the smallest crown area. In contrast, red maple trees had the smallest basal area but the second highest crown area, illustrating how stem size is not necessarily an indicator of crown size. As a result, the larger crowns of red maple lead to greater rainfall interception and less throughfall compared to a co-occurring oak species of the same size but with a smaller crown. Forest canopies are spatially and temporally heterogeneous for many reasons, but some inherent structure can be attributed to differences in function and ecological niches. For example, pioneer species, which are shade-intolerant, allocate resources to grow tall and fast, resulting in shallower crowns with fewer branches. In contrast, climax species, which are moderate to very shade-tolerant, grow slower and are able to photosynthesize at lower light levels, maintaining multiple layers of branches and foliage, resulting in denser and deeper crowns. In evergreen species, this tradeoff between growth strategies in relation to shade tolerance was observed with vertical branching most common among pioneer species (i.e., shade-intolerant species) and horizontal branching most common among climax species (i.e., shade-tolerant species) [32]. As such, the structural traits of species reflect differences in functional traits associated with light acquisition and photosynthesis [33,34]. Consequently, these structural traits may also impact the redistribution of rainwater by the forest canopy.
With regards to canopy water partitioning, the architecture of tree crowns has direct impacts on the number of layers to intercept rainwater. In tropical forest ecosystems where species diversity is highest, differences in DBH and crown size were observed and were most pronounced among midstory trees [35]. At our study site in a temperate forest, the crown to basal area ratio was 2.5-to 25-times greater for red maple overstory trees compared to those of white oak and southern red oak, respectively. In the midstory, red maple and winged elm had crown to basal area ratios up to 11 times greater than either oak species. Thus, both canopy layers displayed differences in structure among species, but these traits were stronger in the overstory.
In this study, the rank of standardized throughfall was different between gauges located at the bole of the tree trunk compared to the mid-crown point underneath overstory trees ( Figure 5). Although there were no consistent trends between gauge location or the interaction with species (p > 0.050). However, overall, the variability of throughfall at this study was within the range (<50%) of CV observed in other systems [15,36,37]. Compared to these studies, we also observed greater variability (i.e., the range of boxplots in Figure 5 and CVs in Figures 3 and 4) during the leaf-on canopy phase compared to the leaf-off canopy phase, especially for the mid-crown gauges. More variability when foliage is present would be expected, as the number of intercepting surfaces increases with new foliage, whereby the spatial complexity of branching architecture combined with that of foliage distribution leads to even higher throughfall variability.
Keim et al. [38] considered the spatial scale over which throughfall collectors were correlated, and found that during the leaf-on canopy phase, the structure of individual tree crowns resulted in correlation between gauges, but in the leafless canopy, that correlation weakened significantly [38]. When foliage is present, the architecture of the tree crown may be more pronounced, leading to stronger spatial patterns. Differences between species in leaf surface characteristics may also influence retention of water by the tree crown. Not only do oaks have smaller crowns, but they have waxier leaf surfaces that help shed intercepted rainfall [10,39], both characteristics leading to greater throughfall than non-oaks. Species that are more shade-tolerant also tend to have larger leaves to capitalize on light interception and photosynthesis [40], which would also increase rainfall interception and decrease throughfall.
Compounding these interactions is also the myriad of storm conditions (i.e., rainfall intensity, wind speed and intensity, synoptic seasonality [41][42][43]), adding to the complexity of an already spatially complex process. In this study, throughfall was collected over 4-week sampling intervals, so any fine-scale temporal influence that species have on throughfall partitioning was not evaluated. However, it is likely that more frequent sampling intervals would help elucidate differences between species.
For example, the timing of leaf life cycles in a deciduous forest from leaf emergence to full canopy occupancy to senescence, termed phenoseasons [44], determine the quantity of radiation making it to the forest floor. At the same time, the quantity of rainfall will also be impacted by phenoseasons. Growth strategies and physiology determine leaf life span [45,46], so it is reasonable to hypothesize that these attributes differ between species and will impact throughfall. Additionally, coatings on leaf surfaces change throughout the growing season, such that leaf hydrophobicity decreases as leaves age, increasing interception storage capacity, thus reducing throughfall [47]. Had this study been able to sample at shorter-time intervals, differences in phenoseasons beyond leaf-on vs. leaf-off canopy phase may have been expressed more strongly and led to periods during the year when throughfall was even more different among species.
The forest hydrologic budget is a balance of inputs via throughfall and stemflow (i.e., net precipitation) and outputs via interception, evaporation, transpiration (i.e., evapotranspiration) and streamflow. This study considered how structural traits vary among species to influence throughfall specifically. Here, we have presented preliminary evidence that suggests throughfall may be reduced during mesophication. However, net precipitation also includes stemflow and the smoother and thinner bark of mesophytes has been shown to substantially increase stemflow production compared to co-occurring oaks [20,21]. Evidence also suggests that mesophytes transpire more water than co-occurring oaks [5]. The impacts of mesophication on forest hydrologic cycles is still unresolved but these few studies provide evidence that watershed yield is decreasing as a result of changing species composition.
Mesophication and changing species composition and structure may have implications beyond forest hydrology. Throughfall is enriched with nutrients that accumulate in forest canopies during antecedent dry periods which are subsequently washed off, and also leached from plant surfaces during rainfall [48][49][50]. Species-level differences in throughfall have been observed across a range of forest functional types [51][52][53], and recent work has linked canopy-derived nutrient inputs with soil microbiological communities [54][55][56][57]. For example, presence of canopy epiphytes increased nutrient fluxes in throughfall and led to differences in abundance of ammonia-oxidizing soil bacteria, compared to soils underneath trees where epiphytes were not present [54].

Conclusions
A shift away from oak dominance throughout upland hardwood forests in the Eastern United States is widespread, and the magnitude and direction of ecosystem functions are uncertain. This paper provides preliminary evidence that changing tree species composition will impact the redistribution of throughfall water, the largest canopy-derived water flux, to the forest floor. As red maple and other mesophytes increase in dominance, throughfall volume will likely decline and throughfall variability will likely increase, both as a result of denser tree crowns associated with these shade-tolerant species. Whether the net impact of denser tree crowns increases soil moisture (by providing more shade and less evaporation) or decreases soil moisture (by reducing throughfall inputs) is yet to be determined. However, these differences in throughfall redistribution could have effects on the net result of mesophication by impacting both growing conditions and forest floor flammability. Furthermore, these differences in water fluxes will be compounded with differences in nutrient fluxes among species, which may alter other biogeophysical components, including soil biogeochemistry.