¹Instituto Nacional de Investigaciones Forestales, Agrícolas y Pecuarias, Centro Nacional de Investigación Disciplinaria en Conservación y Mejoramiento de Ecosistemas Forestales. México.

²Instituto Nacional de Investigaciones Forestales, Agrícolas y Pecuarias. Sitio Experimental Tlaxcala. México.

³Instituto Nacional de Investigaciones Forestales, Agrícolas y Pecuarias. Campo Experimental Valle de México. México.

In forest ecosystems, throughfall (TF) is the primary source of water entering the soil and is determined both by the characteristics of rainfall events and by the canopy structure. The influence of forest management and structural attributes on TF was evaluated in a temperate forest. Circular sampling plots of 1 000 m² were established in stands managed using the Silvicultural Development Method, the Mexican Method for Irregular Forest Management, and a reference stand (REF). Over a three-year period, total rainfall (TR) and daily precipitation (TF) were recorded using rain gauges during the rainy season. The data were analyzed using generalized linear mixed models, incorporating the TR, canopy cover, and silvicultural treatment, as well as the nested hierarchical structure of the data and the temporal dependence between measurements. The TR had a greater influence on the TF, with an average increase of 1.8 % for every additional millimeter of rainfall; canopy cover had a negative effect, reducing the TF by approximately 1.04 % for every percentage point increase in cover. Compared to REF, thinning (T) showed a significant reduction of approximately 12 %. The model showed good predictive performance, with a Correlation coefficient (r) of 0.98 between observed and predicted TF values and Mean absolute error of 1.16 mm. The results indicate that forest management alters the canopy structure, which in turn impacts rainfall distribution. Long-term studies are required to guide adaptive management.

Keywords: Water balance, rainfall interception, rainfall allocation, total rainfall, hydrological ecosystem services, silvicultural treatments.

La precipitación directa (Pd) es la principal entrada de agua al suelo forestal y está determinada por las características de la lluvia y la estructura del dosel. La influencia del manejo forestal y de los atributos estructurales sobre la Pd se evaluó en un bosque templado. El establecimiento de unidades de muestreo circulares de 1 000 m² se realizó en rodales manejados con el Método de Desarrollo Silvícola, el Método Mexicano de Ordenación de Bosques Irregulares y un rodal de Referencia (REF). Durante tres años se registró, en la época de lluvias, la precipitación total (Pt) y la Pd con pluviómetros. Los datos se analizaron mediante modelos lineales mixtos generalizados, incorporando la Pt, cobertura del dosel y tratamiento silvícola, así como la estructura jerárquica anidada de los datos y la dependencia temporal entre mediciones. La Pt tuvo mayor influencia sobre la Pd, con un aumento promedio de 1.8 % por milímetro adicional de lluvia; la cobertura del dosel presentó un efecto negativo, reduciendo la Pd en aproximadamente 1.04 % por unidad porcentual de incremento en cobertura. En comparación con REF, la Corta de aclareo registró una reducción significativa de aproximadamente 12 %. El modelo mostró buen desempeño predictivo, con una r=0.98 entre la Pd observada y predicha y un Error absoluto medio de 1.16 mm. Los resultados indican que el manejo forestal modifica la estructura del dosel, y esta influye en la partición de la lluvia. Destaca la necesidad de estudios de largo plazo para orientar un manejo adaptativo.

Palabras clave: Balance hídrico, intercepción de lluvia, partición de la precipitación, precipitación total, servicios ecosistémicos hidrológicos, tratamientos silvícolas.

Forests under forest management provide a wide range of ecosystem services (Millennium Ecosystem Assessment, 2005). Beyond providing wood and pulp, they contribute to carbon sequestration (Nayak et al., 2022), habitat conservation (Pawson et al., 2008), as well as increased soil water storage capacity and the mitigation of surface runoff (Jeong et al., 2022). In addition, forests play a vital role in the water cycle by breaking down total or incident rainfall (TR) into three main components (Savenije, 2004; Yang et al., 2024): (1) Throughfall (TF), which is rainwater that passes through the canopy to the ground; (2) Stemflow (SF), in which water is conducted from the canopy down the stems to the ground; and (3) Canopy interception (CI), which refers to water retained by the vegetation that subsequently evaporates without reaching the ground (Crockford & Richardson, 2000; Van Stan et al., 2020).

TF is the primary source of water entering forest soils and depends on both the TR and the canopy structure (canopy architecture, leaf area index, height, and bark roughness); it regulates water storage and redistribution processes. In temperate forests, TF accounts for 60 % to 90 % of the TR, depending on the canopy structure and the rainfall characteristics (Barbier et al., 2009; Muzylo et al., 2009; Van Stan et al., 2020). However, information on the quantitative relationship between the TF and the TR in Mexico’s forest ecosystems is scarce, particularly for assessing the impact of forest management in temperate forests.

Measuring TF in forest ecosystems poses a methodological challenge due to the high spatial heterogeneity of the canopy, even when using a large number of rain gauges. Often, the first challenge is to obtain an accurate estimate of the TR. For practical reasons, this is most often done in an open area near the area of interest, rather than directly over the trees. The structural variability of the forest stand calls for multiple rain gauges to accurately measure the TF, as the error in estimating this parameter decreases as the number of rain gauges increases (Deguchi et al., 2006). Aussenac (1981) suggests using 12 to 16 rain gauges per hectare, while other authors believe the number depends on the measurement time scale. Therefore, a larger number of rain gauges is required for daily TF estimates than for monthly or annual estimations (Deguchi et al., 2006; Huber & Iroumé, 2001).

The role of the forest canopy in the redistribution of TF is recognized; however, there is limited empirical evidence on how different forest management practices under comparable climatic conditions alter the amount of TF reaching the ground (Gazol et al., 2025). This information is key to linking forestry with the provision of hydrological ecosystem services, particularly in areas with a forestry tradition that spans several decades, such as Chignahuapan, state of Puebla, Mexico. Within this context, the objective of this study was to evaluate the effect of silvicultural treatments on TF and to analyze its relationship with structural attributes (height, diameter at breast height, crown cover, basimetric area and density) in a temperate forest in central Mexico. The hypothesis posits that TF varies across silvicultural treatments due to the changes these cause in the canopy structure, which determines the amount of water that reaches the forest floor.

The research was conducted in Emiliano Zapata ejido, in Chignahuapan municipality, state of Puebla, Mexico, located between 19°39'42" and 19°58'48" North and -97°57'18" and -98°18'06" West. The ejido covers an area of 309.92 hectares, of which 305.3 are used for forestry. The climate is temperate subhumid, with rainfall in the summer (800-1 200 mm of annual precipitation), an average annual temperature of 13.4 °C, and an average relative humidity of 85 %. The vegetation consists of a temperate forest dominated by Pinus patula Schiede ex Schltdl. & Cham., with a smaller presence of other conifers such as Pinus ayacahuite Ehrenb. ex Schltdl. and Abies religiosa (Kunth) Schltdl. & Cham. (Lazcano-Hernández, 2006; Pérez-Miranda, 2014).

The study was conducted in managed forest stands during the third rotation cycle (2014-2024), using two silvicultural methods: the Silvicultural Development Method (SDM) and the Mexican Method for the Management of Irregular Forests (MMOBI) (Figure 1). Based on a review of the ejido’s Timber Forest Management Plan (PMFM by its Spanish acronym, 2013) and information provided by the ejido’s technical officer, a quasi-systematic sampling design was adopted, using circular sampling units (SU) of 1 000 m² spaced 100 m apart (Figure 1).

REF = Reference stand; T = Third thinning; LC = Liberation cut; RC = Regeneration cut from parent trees; SL = Selective logging. Puebla = State of Puebla; Tlaxcala = State of Tlaxcala.

Four silvicultural treatments were selected: Third thinning (T), Regeneration cutting using parent trees (RC), and Release cutting (LC) for the SDM, which were carried out in 2016, 2014, and 2014, respectively; as well as the Selective logging (SL) plot, which was managed in 2016 for the MMOBI (Figure 1), and a Reference stand (REF), which has not been managed in the last 30 years (Correa-Díaz et al., 2025). In each treatment, five sampling units (SU) were considered, with the exception of the RC, where only three SUs were established due to the stand’s spatial limitations and the need to maintain independence among sampling units. In the RC treatment, the removal of parent trees had already been completed at the time of the study; however, the stand still exhibited structural characteristics associated with the regeneration process.

In each sampling unit, the species and mensuration variables such as height (Suunto^® clinometer, % and 0°-90°), breast height diameter (BHD, 95 cm model Mantax Haglöf^® blue aluminum caliper), crown diameter (50 m Truper^® measuring tape) of trees with a breast height diameter greater than 7.5 cm, a threshold commonly used in forest inventories, were recorded. The variables’ measurements were based on the operational and mensuration criteria established by the National Forest and Soil Inventory (Comisión Nacional Forestal [Conafor], 2017).

The canopy cover was estimated from a digital analysis of hemispheric photographs taken above each rain gauge using a smartphone camera (Apple iPhone^® XS). The images were taken in portrait orientation, at a height of approximately 1.20 m above ground level, under uniform lighting conditions to minimize the effects of shadows and overexposure. In addition, an approximate radius of influence was estimated from the average crown height of the trees, which ranged between 4.5 and 6 m around each rain gauge. The original images taken under the canopy (Figure 2A) were binarized using the automatic global thresholding method, based on the Otsu algorithm implemented in Python using the OpenCV library (cv2.threshold with the THRESH_BINARY+THRESH_OTSU settings). In accordance with the resulting binary image, the percentage of canopy cover associated with each rain gauge was quantified as the proportion of pixels classified as vegetation relative to the total number of pixels in the image (Figure 2B) (Tichý, 2016). In addition, coverage measurements were taken using a spherical densitometer (model C concave Forest Densitometers^®) to validate the image-based method. The comparison showed an average difference of ~10 %, indicating good agreement between the two methods.

A = Original RGB image captured in the field; B = Binarized image resulting from digital thresholding, where dark areas represent vegetation cover and light areas represent the visible sky.

Figure 2. Binarization process of hemispherical photographs for crown cover estimation.

The TF was measured using True-Check^® (model 110800) rain gauges with a 150 mm diameter and a 0.1 mm accuracy, installed at a height of 1.20 m above ground level (Figure 3B), located 8.5 m from the center of each sampling unit in the north, east, south, and west directions; theoretical spatial pattern arranged in a cross shape was followed

(Ramos-Madrigal et al., 2024; Westfall & Edgar, 2022; Yim et al., 2015) in order to reduce operating costs and, where appropriate, explicitly account for the intra-cluster spatial variability of the TF within each sampling unit (Figure 3C) (Holwerda et al., 2006; Sadeghi et al., 2024). The practical implementation of this sampling design may result in a seemingly random pattern of some sample units (Ramos-Madrigal et al., 2024). The spatial variability of the TF across orientations within each sampling unit was assessed using the Coefficient of variation (CV) calculated on a weekly basis. Two-stage cluster sampling (the four rain gauges form a cluster within the management unit) works well when the elements within each group exhibit high variability and all groups are relatively similar (Scheaffer et al., 2011).

A= Open-air weather station; B = Second open-air measurement area; C = Rain gauge; D = Sampling design for measuring the TF under the canopy. The red circle represents a 1 000 m² circular sampling plot. Figure 3D: O = W.

Figure 3. Instrumentation and spatial arrangement of rain gauges for recording the TF at under-canopy sampling sites.

Rainfall at each rain gauge was measured on a weekly basis over a three-year period: from August 21^st to October 30^th, 2023 (9 weeks), from June 24^th to November 11^th, 2024 (20 weeks), and from July 1^st to October 27^th, 2025 (17 weeks). The sample size for the analysis consisted of 92 rain gauges (grouped in sets of four in each of the 23 study areas). TR and other climatic variables (temperature, relative humidity, radiation, and wind speed) were recorded using a model Pro2 Davis^® weather station equipped with a 0.2-mm resolution rain gauge, located in an area without tree cover adjacent to the study area (Figure 1 and Figure 3A). For the 2024 and 2025 periods, a second measurement area was added (Figure 1 and Figure 2C), where a model Vantage Vue Davis^® weather station was installed, along with two rain gauges designed to measure the total rainfall, with the aim of verifying the consistency of the TR measurements in the study area and providing an additional reference point to compare the records.

The comparison of TR measurements between the two measurement zones was performed using paired t-tests and Wilcoxon tests. TF measurements in which the rain gauges were blocked, as well as cases in which inconsistencies were detected (when the TF exceeded the TR during the week of measurement), were excluded from the statistical analysis. These observations accounted for 3 % of the total; therefore, they are not considered to have had a significant impact on the results (Anys & Weiler, 2024). A total of 4 232 spatially and temporally correlated TF observations were analyzed.

For the initial data analysis, descriptive statistics (mean, standard deviation, and Coefficient of variation) of the TF were calculated by aspect, treatment, and silvicultural method. The normality and homogeneity of the variances of the TF were assessed using the Shapiro-Wilk (Shapiro & Wilk, 1965) and Levene (Levene, 1960) tests, based on an initial model. However, because these assumptions were not met and the data were spatially and temporally correlated, the TF was modeled using a Generalized linear mixed model (GLMM) with a Gamma distribution and a logarithmic link function, which is well-suited to positive and asymmetric continuous variables (Zuur et al., 2009). The logarithmic link was used to ensure positive values in the TF predictions and stability in the process of estimating the model parameters. The TR was included as the primary explanatory covariate. In addition, various variables associated with the canopy structure were evaluated, including crown cover, basal area, stand density, and dominant species, as well as silvicultural treatment as a fixed categorical factor.

To assess the effects of management at different scales, models were tested at the management type and treatment levels. To adequately represent the spatial dependence of the observations (e. g., rain gauges nested within management and silvicultural treatment units), the management unit was treated as a random effect. In regard to the temporal correlation among observations from a single rain gauge, a first-order autoregressive model (AR(1)) was incorporated to model the correlation between weekly observations from a single rain gauge, defined within each sampling unit combination. In addition, the model dispersion was allowed to vary depending on the magnitude of the TR. Model fitting and comparison were performed in R version 4.4 (R Core Team, 2024) using the glmmTMB package, and model selection was evaluated using the Akaike information criterion (AIC). The model's performance was evaluated using the Coefficient of determination (pseudo-R²) and the Mean absolute error (MAE) as measures of fit and predictive accuracy, respectively.

Thus, according to Gbur et al. (2012), the proposed general model presents the following predictor (Equation 1):

= Assumption of distribution of the response variable (throughfall, mm), where

is the direct precipitation recorded in the rain gauge l in week k corresponding to the sampling unit i in treatment j; is the mean and

is the variance under the Gamma distribution assumption.

The total seasonal rainfall (rainy season) accumulated over the three-year assessment period was 1 599 mm, distributed as follows: 201 mm in 2023, 884 mm in 2024, and 514 mm in 2025, with differences primarily associated with the duration of the sampling period in 2023 compared to 2024 and 2025. The TR was assessed in two measurement zones without cover; however, no significant differences were detected between these (paired t-test, p=0.148; Wilcoxon test, p=0.115), the mean difference being of 3.04 mm. Each point in Figure 4 represents an individual TF record associated with a rain gauge and its location within the sampling unit, which made it possible to visualize the spatial variability of the throughfall across treatments and years of evaluation.

REF = Reference stand; T = Third thinning; LC = Liberation cut; RC = Regeneration cut from parent trees; SL = Selective logging.

Figure 4. Relationship between TF and TR on a weekly scale, by rain gauge orientation and type of silvicultural treatment.

Overall, a consistent positive correlation was observed between the TF and the TR, with a Correlation coefficient of r=0.92. By treatment, the correlation coefficients ranged between r=0.90 and 0.94. On the other hand, when evaluating this relationship by year, the highest correlations were registered in 2024 and 2025 (r=0.92 and 0.93, respectively), while a lower value was registered in 2023 (r=0.80), possibly due to the lower frequency and intensity of the rainfall events observed during that period, which was characterized by dry conditions. This pattern is consistent with the data recorded for that year, which show a decrease in the frequency and duration of the rainfall across large areas of central Mexico associated with meteorological drought conditions (Servicio Meteorológico Nacional [SMN], 2023) and heat waves (Cavazos, 2024). Figure 5 shows interannual differences in the magnitude of weekly TF, with peak events of nearly 50 mm in 2023, 150 mm in 2024, and 100 mm in 2025. In contrast, the medians and overall distribution of the weekly TF across treatments were relatively similar.

REF = Reference stand; T = Third thinning; RC = Regeneration cut by parent trees; LC = Liberation cut; SL = Selective logging.

Figure 5. Interannual variability in weekly rainfall under different silvicultural treatments.

The spatial variability of the throughfall across orientations, estimated using the Coefficient of variation (CV) on a weekly scale within each sampling unit, had an average value of 21.6 % (±13.1 %), with most events falling below 40 %. Extreme values were primarily associated with low-intensity rainfall events, in which small absolute differences between rain gauges result in disproportionate increases in relative variability (Bolaños-Sánchez et al., 2021).

In the silvicultural treatments evaluated using the SDM and MMOBI management methods, the dominant species was Pinus patula, which accounted for nearly 60 % of the total number of individuals recorded. By treatment, its relative abundance ranged between 41 % in REF and 72 % in LC. Abies religiosa was the second most common species in SL, while Pinus ayacahuite had a moderate presence, particularly in T (22 %) and RC (34 %). The structural characteristics of the treatments show distinct differences; in particular, T, REF, and SL had the highest canopy cover values, a variable directly related to the proportion of TF (Barbier et al., 2009; Van Stan et al., 2020). The RC treatment showed intermediate values for diameter and height (Table 1).

Management	Treatment	SU	Time since the last silvicultural intervention (years)	*BHD* (cm)	Total height (m)	Crown diameter (m)	Basal area (m² ha^-1)	Tree cover (%)	Density (trees SU^-1)
NA	REF	5	>30	20.90±12.50	16.24±6.56	3.34±1.90	18.70±3.51	75.55±4.01	40±11
SDM	T	5	8	26.51±15.25	20.44±9.08	4.03±1.79	33.32±16.66	73.21±3.68	45±18
	LC	5	9	21.38±9.47	16.90±4.70	3.28±1.59	18.19±7.53	69.21±3.85	42±15
	RC	3	9	20.14±8.92	17.39±5.45	3.44±1.20	25.26±12.71	69.48±4.57	66±26
MMOBI	SL	5	8	23.43±14.02	19.15±9.22	3.12±1.53	29.84±10.59	70.28±2.14	51±13

Note: values expressed as mean±standard deviation. NA = Not applicable; SDM = Silvicultural Development Method; MMOBI = Mexican Method for Managing Irregular Forests; REF = Reference stand; T = Thinning; LC = Release cutting; RC = Regeneration cutting; SL = Selective logging.

During the study period, the TF accounted for 76 % to 86 % of the total accumulated TR (1 599 mm) across all events recorded over the three years of monitoring. The highest values were recorded with the RC (86 %) and REF (84.5 %) treatments, which is consistent with the shorter tree height in these treatments, while T (76.0 %) and SL (77.4 %) showed lower proportions. The LC treatment had an intermediate value (81.3 %). In this regard, a general trend was observed whereby the TF/TR ratio decreased as the canopy cover increased. However, when analyzing the data on a weekly basis, the proportion of TF depended on precipitation amount (Figure 6). During weeks with low rainfall (0-10 mm), interception was higher and more variable; in contrast, during weeks with higher rainfall (>100 mm), interception decreased significantly —a pattern that held true across all treatments. This behavior is consistent with the canopy saturation effect, in which the canopy's storage capacity is exceeded during intense events, allowing more rain to reach the ground (Brasil et al., 2022; Cisneros-Vaca et al., 2018; Rutter et al., 1971; Van Stan et al., 2020).

The observed throughfall fraction falls within the range reported for temperate forests in Mexico and other ecologically similar regions. Previous studies have shown that between 14 % and 40 % of TF can be retained by the canopy, implying that the proportion of rainfall reaching the ground is comparable to the estimates in this study (Bolaños-Sánchez et al., 2021; Flores-Ayala et al., 2016; Huber & Iroumé, 2001; Llorens, 1997).

Based on the general predictor in Equation 1, candidate models with different combinations of structural variables and random structures were evaluated. The final model was selected based on the lowest AIC and the biological plausibility of the estimated effects (Table 2). The model with the best statistical performance used rainfall gauge-level observations and included weekly precipitation at the rainfall gauge level, canopy cover for each orientation, and the silvicultural treatment as explanatory variables; this enabled explicit capture of spatial variability within each sampling unit.

Model	Fixed variables	Random structure	*AIC*
M1	TR+Crown cover+Treatment	(1\|SU)+SU:Orientation+AR(1)	24 689.40
M2	TR+Basal area	SU:Orientation+AR(1)	24 698.75

AIC = Akaike information criterion; TR = Incident rainfall; SU = Sampling units; AR(1) = First-order autoregressive model.

The overall analysis of variance, which used a mixed-model approach based on Wald tests (Type II), revealed significant effects of TR (χ²=6 418.99, df=1, p<0.001), canopy cover (χ²=18.95, df=1, p<0.001), and silvicultural treatment (χ²=10.95, df=4, p=0.036) on throughfall. The model showed that weekly throughfall was closely associated with the TR, increasing by about 1.8 % on average for each additional millimeter of rain. For its part, canopy cover had a significant negative effect, such that increases in canopy cover reduced the amount of rainfall reaching the ground (Deguchi et al., 2006; Kaushal et al., 2017; Staelens et al., 2006).

When analyzing the effect of the treatments specifically, only T showed a significant reduction in TR of approximately 12 % compared with the reference treatment (Table 3). In contrast, the LC, RC, and SL treatments showed no statistically significant effects. In the SL, the most common species was Abies religiosa, whose crown structure has been linked to a greater ability to intercept (Flores-Ayala et al., 2016). Given this context, one might expect a lower throughfall in SL; however, the present study did not detect any statistically significant effects of the dominant species on throughfall (p=0.42).

Table 3. Model coefficients (M1) for the Gamma GLMM of the weekly throughfall.

Variable	Estimator (β)	Standard error	z	p
Intercept	3.383	0.184	18.39	<0.001
Total rainfall (mm)	0.0182	0.00023	80.12	<0.001
Canopy cover (%)	-0.0104	0.00240	-4.35	<0.001
T treatment	-0.1305	0.0410	-3.18	0.001
LC treatment	-0.0606	0.0434	-1.40	0.162
RC treatment	-0.0662	0.0491	-1.35	0.177
SL treatment	-0.0750	0.0426	-1.76	0.078

Treatments: T = Thinning; LC = Release cutting; RC = Regeneration cutting; SL = Selective logging. Note: the coefficients are plotted on a logarithmic scale; the interpretation of the percentage was obtained through exponential transformation

. The AR(1) autoregressive component showed a low temporal correlation (p=-0.02), with a variance of 0.0883 and a standard deviation of 0.2972.

The selected model proved statistically robust and was considered appropriate (AIC=24 689.4), as the observed values compared to the predicted values lay along a straight line passing through the origin, regardless of the treatments, without showing any systematic bias pattern. Furthermore, the proper application of the Gamma distribution is reflected in the positive values of the predicted throughfall. In fact, the correlation between observed and predicted throughfall was high across all treatments, with average values of r=0.98 and pseudo-R²=0.97 (Figure 7). Consistent with this, the mean absolute errors (MAE=1.16 mm) remained low and relatively homogenous, as did the Root mean square error (RMSE=2.28 mm), indicating that the model has an adequate predictive power for estimating the weekly TR regardless of the type of forest management.

Figure 7. Relationship between observed throughfall and the throughfall predicted by the Generalized linear mixed model (GLMM) at the rain gauge level for each silvicultural treatment.

The effect of the canopy structure—as measured by basal area—on throughfall was also assessed. In this case, the basal area had a significant negative effect on throughfall (β=-0.032, p<0.01), which implies that a higher structural density of the stand promotes interception processes. In both models, the TR remained the most influential variable (p<0.001), which highlights its dominant role in determining throughfall (Table 4).

Table 4. Model coefficients (M2) for the Gamma GLMM of the weekly throughfall.

Variable	Estimator (β)	Standard error	z	p
Intercept	2.648	0.0377	70.32	<0.001
Total rainfall (mm)	0.0182	0.00023	80.12	<0.001
Basal area (m² ha^-1)	-0.0319	0.0127	-2.50	0.012

Note: the coefficients are plotted on a logarithmic scale; the interpretation of the percentage was obtained through exponential transformation

. The autoregressive component AR(1) showed a low temporal correlation (p=-0.03), with a variance of 0.0884 and a standard deviation of 0.2973.

It should be noted that the stands evaluated were harvested between 2014 and 2016; therefore, at the time of the monitoring, they were in a structural recovery phase lasting approximately 8-10 years. Within this context, the effects of management on throughfall may be influenced not only by the type of treatment applied but also by the resulting tree density and canopy structure following the intervention (Kaushal et al., 2017; Monárrez-González et al., 2018; Muñoz-Villers et al., 2012; Staelens et al., 2006). In this regard, it is important to examine in greater detail the role of the time elapsed since the intervention on the interception dynamics under managed forest conditions.

The results show that the variability of the TF in managed forests depends largely on the TR; however, its variability at the stand level was determined primarily by the attributes of the canopy structure. In particular, the canopy cover and basal area showed significant negative effects, indicating that stands with greater structural development increase water retention and reduce the amount of water that reaches the ground. In contrast, the effect of the silvicultural treatment was secondary and inconsistent across categories, suggesting that its influence is primarily indirect, manifested through structural changes in the canopy, especially when considering recovery times following the intervention. The weekly descriptive analysis showed that interception depends on precipitation amount, with higher values during low-precipitation events and a gradual decrease as rainfall increases. Taken together, these results indicate that the hydrological response under the study conditions is better explained by the current canopy structure than by management type.

The authors would like to thank the Emiliano Zapata ejido in Chignahuapan, Puebla, for assistance with conducting fieldwork. Also, to the Instituto Nacional de Investigaciones Forestales, Agrícolas y Pecuarias, INIFAP (National Institute for Research on Forest, Agriculture and Livestock) for the funding provided through the project titled “Integrated management of forest resources for the sustainability of ecosystem services in the face of climate change”.

The authors declare that they have no conflicts of interest. Arian Correa-Díaz declares that he did not participate in the editorial process for this article.

Jesús Valentín Gutiérrez-García: analysis of results and drafting of the original manuscript; Efraín Velasco-Bautista: review of the methodology and contribution to the drafting of the manuscript; Arian Correa-Díaz: data analysis and assistance with the drafting of the manuscript; Antonio González-Hernández: information review and contribution to the drafting of the manuscript; Francisco Moreno-Sánchez: information review and contribution to the drafting of the manuscript; Bertha Patricia Zamora-Morales: revision and drafting of the manuscript; Vidal Guerra-De la Cruz: technical revision and editing of the manuscript; Eulogio Flores-Ayala: scientific consultation and final revision of the manuscript.

Anys, M., & Weiler, M. (2024). Rainfall interception by urban trees: Event characteristics and tree morphological traits. Hydrological Processes, 38(4), Article e15146. https://doi.org/10.1002/hyp.15146

Barbier, S., Balandier, P., & Gosselin, F. (2009). Influence of several tree traits on rainfall partitioning in temperate and boreal forests: a review. Annals of Forest Science, 66(6), 602. https://doi.org/10.1051/forest/2009041

Bolaños-Sánchez, C., Prado-Hernández, J. V., Silván-Cárdenas, J. L., Vázquez-Peña, M. A., Madrigal-Gómez, J. M., & Martínez-Ruíz, A. (2021). Estimating rainfall interception of Pinus hartwegii and Abies religiosa using analytical models and point cloud. Forests, 12(7), Article 866. https://doi.org/10.3390/f12070866

Brasil, J. B., de Andrade, E. M., de Queiroz-Palácio, H. A., Fernández-Raga, M., Ribeiro-Filho, J. C., Medeiros, P. H. A., & Guerreiro, M. S. (2022). Canopy effects on rainfall partition and throughfall drop size distribution in a tropical dry forest. Atmosphere, 13(7), Article 1126. https://doi.org/10.3390/atmos13071126

Cisneros-Vaca, C., van der Tol, C., & Ghimire, C. P. (2018). The influence of long-term changes in canopy structure on rainfall interception loss: a case study in Speulderbos, the Netherlands. Hydrology and Earth System Sciences, 22(7), 3701-3719. https://doi.org/10.5194/hess-22-3701-2018

Correa-Díaz, A., Villanueva-Díaz, J., Gutiérrez-García, J. V., Velasco-Bautista, E., Moreno-Sánchez, F., & Zamora-Morales, B. P. (2025). Efecto del clima y el manejo forestal en el crecimiento radial de un bosque de coníferas en Puebla, México. Madera y Bosques, 31, Article e312717. https://doi.org/10.21829/myb.2025.312717

Deguchi, A., Hattori, S., & Park, H.-T. (2006). The influence of seasonal changes in canopy structure on interception loss: Application of the revised Gash model. Journal of Hydrology, 318(1-4), 80-102. https://doi.org/10.1016/j.jhydrol.2005.06.005

Flores-Ayala, E., Guerra-De la Cruz, V., Terrazas-González, G. H., Carrillo-Anzures, F., Islas-Gutiérrez, F., Acosta-Mireles, M., & Buendía-Rodríguez, E. (2016). Intercepción de lluvia en bosques de montaña en la cuenca del río Texcoco, México. Revista Mexicana de Ciencias Forestales, 7(37), 65-76. https://doi.org/10.29298/rmcf.v7i37.52

Gazol, A., Pizarro, M., Hammond, W. M., Allen, C. D., & Camarero, J. J. (2025). Droughts preceding tree mortality events have increased in duration and intensity, especially in dry biomes. Nature Communications, 16, Article 5779. https://doi.org/10.1038/s41467-025-60856-5

Gbur, E. E., Stroup, W. W., McCarter, K. S., Durham, S., Young, L. J., Christman, M., West, M., & Kramer, M. (2012). Analysis of generalized linear mixed models in the agricultural and natural resources sciences. John Wiley & Sons. https://doi.org/10.2134/2012.generalized-linear-mixed-models

Holwerda, F., Scatena, F. N., & Bruijnzeel, L. A. (2006). Throughfall in a Puerto Rican lower montane rain forest: A comparison of sampling strategies. Journal of Hydrology, 327(3), 592-602. https://doi.org/10.1016/j.jhydrol.2005.12.014

Huber, A., & Iroumé, A. (2001). Variability of annual rainfall partitioning for different sites and forest covers in Chile. Journal of Hydrology, 248(1-4), 78-92. https://doi.org/10.1016/S0022-1694(01)00394-8

Jeong, S., Kume, T., Shinohara, Y., Farahnak, M., & Otsuki, K. (2022). Application of the reformulated gash analytical model for rainfall interception loss to unmanaged high-density coniferous plantations laden with dead branches. Forests, 13(5), Article 657. https://doi.org/10.3390/f13050657

Kaushal, R., Kumar, A., Alam, N. M., Mandal, D., Jayaparkash, J., Tomar, J. M. S., Patra, S., Gupta, A. K., Mehta, H., Panwar, P., Chaturvedi, O. P., & Mishra, P. K. (2017). Effect of different canopy management practices on rainfall partitioning in Morus alba. Ecological Engineering, 102, 374-380. https://doi.org/10.1016/j.ecoleng.2017.02.029

Llorens, P. (1997). Rainfall interception by a Pinus sylvestris forest patch overgrown in a Mediterranean mountainous abandoned area II. Assessment of the applicability of Gash’s analytical model. Journal of Hydrology, 199(3), 346-359. https://doi.org/10.1016/S0022-1694(96)03335-5

Millennium Ecosystem Assessment (Ed.). (2005). Ecosystems and human well-being. Opportunities and challenges for business and industry. World Resources Institute. http://pdf.wri.org/ma_synthesis_business.pdf

Monárrez-González, J. C., Pérez-Verdín, G., López-González, C., Márquez-Linares, M. A., & González-Elizondo, M. del S. (2018). Efecto del manejo forestal sobre algunos servicios ecosistémicos en los bosques templados de México. Madera y Bosques, 24(2), Artículo e2421569. https://doi.org/10.21829/myb.2018.2421569

Muñoz-Villers, L. E., Holwerda, F., Gómez-Cárdenas, M., Equihua, M., Asbjornsen, H., Bruijnzeel, L. A., Marín-Castro, B. E., & Tobón, C. (2012). Water balances of old-growth and regenerating montane cloud forests in central Veracruz, Mexico. Journal of Hydrology, 462-463, 53-66. https://doi.org/10.1016/j.jhydrol.2011.01.062

Muzylo, A., Llorens, P., Valente, F., Keizer, J. J., Domingo, F., & Gash, J. H. C. (2009). A review of rainfall interception modelling. Journal of Hydrology, 370(1-4), 191-206. https://doi.org/10.1016/j.jhydrol.2009.02.058

Nayak, N., Mehrotra, R., & Mehrotra, S. (2022). Carbon biosequestration strategies: a review. Carbon Capture Science & Technology, 4, Article 100065. https://doi.org/10.1016/j.ccst.2022.100065

Pawson, S. M., Brockerhoff, E. G., Meenken, E. D., & Didham, R. K. (2008). Non-native plantation forests as alternative habitat for native forest beetles in a heavily modified landscape. Biodiversity and Conservation, 17(5), 1127-1148. https://doi.org/10.1007/s10531-008-9363-y

R Core Team. (2024). R: A language and environment for statistical computing (Version 4.4) [Computer software]. R Foundation for Statistical Computing. https://www.r-project.org/

Ramos-Madrigal, R., de los Santos-Posadas, H. M., Valdez-Lazalde, J. R., Velasco-Bautista, E., Ángeles-Pérez, G., & Ortiz-Reyes, A. D. (2024). Evaluation of potential productivity in coniferous forests by integrating field data and aerial laser scanning in Hidalgo, México. Forest Systems, 33(3), Article 20886. https://doi.org/10.5424/fs/2024333-20886

Rutter, A. J., Kershaw, K. A., Robins, P. C., & Morton, A. J. (1971). A predictive model of rainfall interception in forests, 1. Derivation of the model from observations in a plantation of Corsican pine. Agricultural Meteorology, 9, 367-384. https://doi.org/10.1016/0002-1571(71)90034-3

Sadeghi, S. M. M., Epstein, J. M., Deljouei, A., Gorora, F. J., & Cohen, M. J. (2024). Stand age controls canopy and soil rainfall partitioning in slash pine forests. Forest Ecology and Management, 572, Article 122307. https://doi.org/10.1016/j.foreco.2024.122307

Savenije, H. H. G. (2004). The importance of interception and why we should delete the term evapotranspiration from our vocabulary. Hydrological Processes, 18(8), 1507-1511. https://doi.org/10.1002/hyp.5563

Shapiro, S. S., & Wilk, M. B. (1965). An analysis of variance test for normality (complete samples). Biometrika, 52(3-4), 591-611. https://doi.org/10.1093/biomet/52.3-4.591

Staelens, J., De Schrijver, A., Verheyen, K., & Verhoest, N. E. C. (2006). Spatial variability and temporal stability of throughfall water under a dominant beech (Fagus sylvatica L.) tree in relationship to canopy cover. Journal of Hydrology, 330(3), 651-662. https://doi.org/10.1016/j.jhydrol.2006.04.032

Tichý, L. (2016). Field test of canopy cover estimation by hemispherical photographs taken with a smartphone. Journal of Vegetation Science, 27(2), 427-435. https://doi.org/10.1111/jvs.12350

Van Stan, II, J. T., Gutmann, E., & Friesen, J. (Eds.). (2020). Precipitation partitioning by vegetation. A global synthesis. Springer Nature Switzerland AG. https://doi.org/10.1007/978-3-030-29702-2

Westfall, J. A., & Edgar, C. B. (2022). Addressing non-response bias in urban forest inventories: an estimation approach. Frontiers in Forests and Global Change, 5, Article 895969. https://doi.org/10.3389/ffgc.2022.895969

Yang, J., Wang, A., Shen, L., Dai, G., Liu, Y., Zhang, Y., Fei, W., & Wu, J. (2024). The impact of canopy on nutrient fluxes through rainfall partitioning in a mixed broadleaf and coniferous forest. Forests, 15(4), Article 623. https://doi.org/10.3390/f15040623

Yim, J.-S., Shin, M.-Y., Son, Y., & Kleinn, C. (2015). Cluster plot optimization for a large area forest resource inventory in Korea. Forest Science and Technology, 11(3), 139-146. https://doi.org/10.1080/21580103.2014.968222

Zuur, A. F., Ieno, E. N., Walker, N., Saveliev, A. A., & Smith, G. M. (2009). Mixed effects models and extensions in ecology with R. Springer-Verlag New York. https://doi.org/10.1007/978-0-387-87458-6

Todos los textos publicados por la Revista Mexicana de Ciencias Forestales –sin excepción– se distribuyen amparados bajo la licencia Creative Commons 4.0 Atribución-No Comercial (CC BY-NC 4.0 Internacional), que permite a terceros utilizar lo publicado siempre que mencionen la autoría del trabajo y a la primera publicación en esta revista.