¹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.

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.

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.

Los bosques con manejo forestal proporcionan una amplia gama de servicios ecosistémicos (Millennium Ecosystem Assessment, 2005). Más allá del suministro de madera y pulpa, contribuyen al secuestro de carbono (Nayak et al., 2022), a la conservación de hábitat (Pawson et al., 2008), así como al aumento de la capacidad de almacenamiento de agua en el suelo y mitigación de la escorrentía superficial (Jeong et al., 2022). Además, los bosques desempeñan un papel esencial en el ciclo hidrológico al dividir la precipitación total o incidente (Pt) en tres componentes principales (Savenije, 2004; Yang et al., 2024): (1) precipitación directa (Pd), es el agua de lluvia que atraviesa el dosel hasta el suelo; (2) flujo por troncos o escurrimiento fustal (Ef), en el cual el agua es conducida desde la copa a lo largo de los tallos hasta el suelo; y (3) intercepción de lluvia del dosel (Id), que corresponde al agua retenida en la vegetación que posteriormente se evapora, sin llegar al suelo (Crockford & Richardson, 2000; Van Stan et al., 2020).

La Pd es la principal fuente de entrada de agua al suelo forestal y depende, tanto de la Pt como de la estructura del dosel (arquitectura de copa, índice de área foliar, altura y rugosidad de la corteza), y regula los procesos de almacenamiento y redistribución del agua. En bosques templados, la Pd representa entre 60 % y 90 % de la Pt, según la estructura del dosel y las características de la lluvia (Barbier et al., 2009; Muzylo et al., 2009; Van Stan et al., 2020). Sin embargo, en los ecosistemas forestales de México, particularmente en la evaluación del efecto del manejo forestal en bosques de clima templado existe un vacío de información sobre la relación cuantitativa entre Pd y Pt.

La medición de la Pd en ecosistemas forestales representa un reto metodológico debido a la alta heterogeneidad espacial del dosel, incluso si se utiliza una gran cantidad de pluviómetros. A menudo, el primer desafío consiste en obtener una estimación correcta de la Pt. Por razones prácticas, esto se hace con mayor frecuencia en un área abierta cercana al área de interés, en lugar de registrarla directamente sobre los árboles. La variabilidad estructural de la masa forestal hace necesario usar diversos pluviómetros para evaluar de manera correcta la Pd, ya que el error de estimación de este parámetro disminuye conforme el número de pluviómetros se incrementa (Deguchi et al., 2006). Aussenac (1981) sugiere utilizar de 12 a 16 pluviómetros por hectárea, mientras que otros autores consideran que el número depende de la escala temporal de medición. Entonces, se requiere una mayor cantidad de pluviómetros para estimaciones diarias de Pd que para las mensuales o anuales (Deguchi et al., 2006; Huber & Iroumé, 2001).

A pesar de reconocer el papel del dosel forestal en la redistribución de la Pd, persiste una limitada evidencia empírica que indique cómo distintos esquemas de manejo forestal, bajo condiciones climáticas comparables, modifican la Pd que alcanza el suelo (Gazol et al., 2025). Esa información resulta clave para vincular la silvicultura con la provisión de servicios ecosistémicos hidrológicos, particularmente, en zonas con tradición silvícola de varias décadas como ocurre en Chignahuapan, Puebla, México. En este contexto, el objetivo del presente estudio fue evaluar el efecto de los tratamientos silvícolas sobre la Pd y analizar su relación con los atributos estructurales (altura, diámetro a la altura del pecho, cobertura de copa, área basal y densidad) en un bosque templado del centro de México. La hipótesis planteada establece que la Pd difiere entre tratamientos silvícolas debido a las modificaciones que generan en la estructura del dosel, la cual condiciona la cantidad de agua que alcanza el suelo forestal.

La investigación se desarrolló en el ejido Emiliano Zapata, municipio Chignahuapan, Puebla, México, ubicado en los paralelos 19°39'42" y 19°58'48" de latitud norte y los meridianos -97°57'18" y -98°18'06" de longitud oeste. El ejido tiene una superficie de 309.92 ha, de las cuales 305.3 son de uso forestal. El clima es templado subhúmedo con lluvias en verano (800-1 200 mm de precipitación anual), temperatura promedio anual de 13.4 °C y humedad relativa promedio de 85 %. La vegetación corresponde a un bosque templado dominado por Pinus patula Schiede ex Schltdl. & Cham. y con menor presencia de otras coníferas como Pinus ayacahuite Ehrenb. ex Schltdl. y Abies religiosa (Kunth) Schltdl. & Cham. (Lazcano-Hernández, 2006; Pérez-Miranda, 2014).

El estudio se llevó a cabo en rodales con manejo forestal, en el tercer ciclo de corta 2014-2024, bajo dos métodos silvícolas: el Método de Desarrollo Silvícola (MDS) y el Método Mexicano de Ordenación de Bosques Irregulares (MMOBI) (Figura 1). A partir de la revisión del Programa de Manejo Forestal Maderable (PMFM, 2013) del ejido y de información proporcionada por la responsiva técnica de este, se consideró un diseño de muestro cuasi-sistemático con unidades de muestreo circulares (UM) de 1 000 m², equidistantes a 100 m (Figura 1).

CA = Corta de tercer aclareo; CL = Corta de liberación; CR = Corta de regeneración por árboles padres; CS = Corta de selección; REF = Rodal de referencia.

Se seleccionaron cuatro tratamientos silvícolas: corta de tercer aclareo (CA), corta de regeneración por árboles Padre (CR) y corta de liberación (CL) para el MDS, intervenidos en el 2016, 2014 y 2014, respectivamente; así como la corta de selección (CS) intervenido en el 2016 para el MMOBI (Figura 1) y un rodal de Referencia (REF), el cual no ha sido intervenido en los últimos 30 años (Correa-Díaz et al., 2025). En cada tratamiento se consideraron cinco unidades de muestreo (UM) a excepción de CR, donde se establecieron únicamente tres UM debido a las limitaciones espaciales del rodal y a la necesidad de mantener independencia entre unidades de muestreo. En el tratamiento CR, la remoción de árboles Padre se había realizada al momento del estudio; sin embargo, el rodal aún presentaba características estructurales asociadas al proceso de regeneración.

En cada unidad de muestreo se registró la especie y variables dasométricas como altura (clinómetro Suunto^® % y 0°-90°), diámetro a la altura del pecho (DAP; forcípula de aluminio 95 cm marca Haglöf^® Mantax blue), diámetro de copa (cinta métrica de 50 m marca Truper^®) de árboles mayores a 7.5 cm en diámetro a la altura del pecho, umbral comúnmente empleado en inventarios forestales. La medición de las variables se realizó a partir de los criterios operativos y dasométricos establecidos en el Inventario Nacional Forestal y de Suelos (Comisión Nacional Forestal [Conafor], 2017).

La cobertura de copa se estimó con base en el análisis digital de fotografías hemisféricas tomadas sobre cada pluviómetro, con una cámara de teléfono móvil (Apple iPhone^® XS). Las imágenes se hicieron en orientación vertical, a una altura aproximada de 1.20 m sobre el nivel del suelo, bajo condiciones de iluminación homogéneas para reducir efectos de sombra y sobreexposición. De manera complementaria, se calculó un radio de influencia aproximado a partir de la altura media de copa del arbolado, el cual se situó entre 4.5 y 6 m alrededor de cada pluviómetro. Las imágenes originales tomadas bajo dosel (Figura 2A) fueron binarizadas mediante el método de umbralización global automática, basado en el algoritmo de Otsu, implementado en el entorno Python utilizando la biblioteca OpenCV (cv2.threshold con la configuración THRESH_BINARY+THRESH_OTSU). A partir de la imagen binaria resultante, se cuantificó el porcentaje de cobertura de copa asociado a cada pluviómetro como la proporción de pixeles correspondientes a vegetación, respecto al total de pixeles de la imagen (Figura 2B) (Tichý, 2016). De manera complementaria, se realizaron mediciones de cobertura con un densitómetro esférico (Forest Densitometers^® modelo C cóncavo) para validar el método basado en imágenes. La comparación mostró una diferencia promedio de ~10 %, lo que indica una adecuada concordancia entre ambos métodos.

A = Imagen original en formato RGB obtenida en campo; B = Imagen binarizada resultante del proceso de umbralización digital, donde las áreas oscuras representan la cobertura vegetal y las áreas claras el cielo visible.

Figura 2. Proceso de binarización de fotografías hemisféricas para la estimación de cobertura de copa.

La Pd se midió con pluviómetros True-Check^® (modelo 110800) de 150 mm y 0.1 mm de precisión, instalados a una altura de 1.20 m sobre el nivel del suelo (Figura 3B), ubicados a 8.5 m del centro de cada UM en direcciones norte, este, sur y oeste, en una estructura espacial teórica conglomerada en forma de cruz (Ramos-Madrigal et al., 2024; Westfall & Edgar, 2022; Yim et al., 2015) a fin de reducir costos operativos y, en su caso, considerar explícitamente la variabilidad espacial intraconglomerado de la Pd dentro de cada unidad de muestreo (Figura 3C) (Holwerda et al., 2006; Sadeghi et al., 2024). La implementación práctica de este diseño de muestreo puede conducir a un patrón con apariencia aleatoria de algunas unidades de muestreo (Ramos-Madrigal et al., 2024). La variabilidad espacial de la Pd entre orientaciones en cada unidad de muestreo se evaluó mediante el coeficiente de variación (CV) calculado a escala semanal. El muestreo por conglomerados en dos etapas (los cuatro pluviómetros forman un conglomerado en la UM) funciona bien cuando los elementos de cada grupo presentan gran variabilidad y todos los grupos son relativamente similares (Scheaffer et al., 2011).

A = Estación meteorológica en condiciones abiertas; B = Segunda área de medición en condiciones abiertas; C = Pluviómetro; D = Diseño de muestreo para la medición de Pd bajo dosel. El círculo rojo representa una unidad de muestreo circular de 1 000 m².

Figura 3. Instrumentación y disposición espacial de los pluviómetros para el registro de la Pd en unidades de muestreo bajo dosel.

La precipitación en cada pluviómetro se midió de forma semanal durante tres años: del 21 de agosto al 30 de octubre de 2023 (9 semanas), del 24 de junio al 11 de noviembre de 2024 (20 semanas) y del 1 de julio al 27 de octubre de 2025 (17 semanas). El tamaño de muestra para el análisis consistió en 92 pluviómetros (anidados en grupos de cuatro en cada una de las 23 UM). La Pt y otras variables climáticas (temperatura, humedad relativa, radiación y velocidad del viento) se registraron con una estación meteorológica Davis^®Pro2 equipada con un pluviómetro de resolución 0.2 mm, ubicada en un lugar sin cobertura arbórea aledaña al área de estudio (Figura 1 y Figura 3A). Para los periodos 2024 y 2025 se incorporó una segunda área de medición (Figura 1 y Figura 2C), en la cual se instaló una estación meteorológica Davis^® Vantage Vue, junto con dos pluviómetros destinados a la medición de la precipitación total, con el objetivo de verificar la consistencia de la medición de Pt y tener un punto adicional de referencia para la comparación de los registros.

La comparación de la medición de Pt entre ambas zonas de medición, se realizó con pruebas t pareada y de Wilcoxon. Las mediciones de Pd registradas en los pluviómetros que estuvieron obstruidos, así como aquellos casos en los que se detectaron inconsistencias (cuando Pd superaba la Pt en la semana de registro) se excluyeron del análisis estadístico. Estas observaciones representaron 3 % del total, por lo que se considera que no influyeron de manera significativa en los resultados (Anys & Weiler, 2024). En total se consideraron 4 232 observaciones de Pd correlacionadas espacial y temporalmente.

Para la exploración inicial de los datos se calcularon estadísticos descriptivos (media, desviación estándar y coeficiente de variación) de la Pd, por orientación, tratamiento y método silvícola. La normalidad y homogeneidad de varianzas de Pd se evaluaron con las pruebas de Shapiro-Wilk (Shapiro & Wilk, 1965) y Levene (Levene, 1960), a partir de un modelo inicial. Sin embargo, al no cumplirse estos supuestos y por tratarse de datos correlacionados espacial y temporalmente, la Pd se modeló con un modelo lineal mixto generalizado (GLMM) con distribución Gamma y función de enlace logarítmica, adecuado para variables continuas positivas y asimétricas (Zuur et al., 2009). El enlace logarítmico se utilizó para garantizar positividad en las predicciones de la Pd y estabilidad en el proceso de estimación de los parámetros del modelo. La Pt se incluyó como covariable explicativa principal. Asimismo, se evaluaron diferentes variables asociadas a la estructura del dosel, que incluyeron la cobertura de copa, el área basal, la densidad del rodal y especie dominante, así como el tratamiento silvícola como factor categórico fijo.

Con el objetivo de evaluar los efectos del manejo a distintas escalas, se probaron modelos a nivel de tipo de manejo y de tratamiento. Para representar de manera adecuada la dependencia espacial de las observaciones (e. g., pluviómetros anidados en unidad de manejo y tratamiento silvícola) se consideró a la unidad de manejo como efecto aleatorio; respecto a la correlación temporal entre observaciones dentro de un mismo pluviómetro, se incorporó una estructura autorregresiva de primer orden (AR(1)) para modelar lo correlación entre observaciones semanales de un mismo pluviómetro, definida dentro de cada combinación de unidad de muestreo. Adicionalmente, se permitió que la dispersión del modelo variara en función de la magnitud de la Pt. El ajuste y la comparación de modelos se hicieron en el entorno R versión 4.4 (R Core Team, 2024), mediante el paquete glmmTMB, y la selección del modelo se evaluó utilizando el criterio de información de Akaike (AIC). El desempeño del modelo se evaluó con el coeficiente de determinación (pseudo-R²) y el error absoluto medio (MAE), como medidas de ajuste y precisión predictiva, respectivamente.

Entonces, conforme a Gbur et al. (2012), el modelo general postulado presenta el siguiente predictor (Ecuación 1):

= Supuesto de distribución de la variable respuesta (precipitación directa, mm), donde

es la precipitación directa registrada en el pluviómetro l en la semana k correspondiente a la unidad de muestreo i en el tratamiento j;

es la media y

la varianza bajo el suspuesto de distribución Gamma.

Comportamiento general de la precipitación total (Pt) y precipitación directa (Pd)

La Pt acumulada a nivel estacional (época de lluvias) en los tres años de evaluación fue de 1 599 mm, la cual se distribuyó de la siguiente manera: 201 mm en 2023, 884 mm en 2024 y 514 mm en 2025 con diferencias asociadas principalmente a la duración del periodo de muestreo de 2023, respecto a 2024 y 2025. La Pt se evaluó en dos zonas de medición sin cobertura; sin embargo, no se detectaron diferencias significativas entre ellas (t pareada, p=0.148; Wilcoxon, p=0.115), con una diferencia media de 3.04 mm. Cada punto de la Figura 4 representa un registro individual de Pd asociado a un pluviómetro y orientación dentro de la unidad de muestreo, ello permitió visualizar la variabilidad espacial de la precipitación directa entre tratamientos y años de evaluación.

REF = Rodal de referencia; CA = Corta de tercer aclareo; CL = Corta de liberación; CR = Corta de regeneración por árboles padres; CS = Corta de selección.

Figura 4. Relación entre la Pd y la Pt a escala semanal, por orientación de pluviómetro y tipo de tratamiento silvícola.

En general, se observó una correlación positiva consistente entre la Pd y la Ptcon un coeficiente de correlación de r=0.92. Por tratamiento, los coeficientes de correlación oscilaron entre r=0.90 y 0.94. Al evaluar esta relación por año, las correlaciones más altas se registraron en 2024 y 2025 (r=0.92 y 0.93, respectivamente), mientras que en 2023 el valor fue menor (r=0.80); posiblemente, asociado a la baja frecuencia e intensidad de eventos de precipitación observados durante ese periodo, caracterizado por condiciones secas. Este patrón concuerda con los datos documentados para ese año, con una reducción en la frecuencia y continuidad de las lluvias en amplias regiones del centro de México, asociada a condiciones de sequía meteorológica (Servicio Meteorológico Nacional [SMN], 2023) y olas de calor (Cavazos, 2024). En la Figura 5 se muestran diferencias interanuales en la magnitud de la Pd semanal, con eventos máximos cercanos a 50 mm en 2023, 150 mm en 2024 y 100 mm en 2025. En contraste, las medianas y la distribución general de la Pd semanal entre tratamientos fueron relativamente similares.

REF = Rodal de referencia; CA = Corta de tercer aclareo; CR = Corta de regeneración por arboles padres; CL = Corta de liberación; CS = Corta de selección.

Figura 5. Variabilidad interanual de la Pd semanal bajo distintos tratamientos silvícolas.

La variabilidad espacial de la precipitación directa entre orientaciones, estimada mediante el coeficiente de variación (CV) a escala semanal en cada unidad de muestreo, presentó un valor promedio de 21.6 % (±13.1 %), con la mayoría de los eventos concentrados por debajo de 40 %. Los valores extremos se asociaron, principalmente, a eventos de baja magnitud de precipitación, en los cuales pequeñas diferencias absolutas entre pluviómetros generan incrementos desproporcionados en la variabilidad relativa (Bolaños-Sánchez et al., 2021).

En los tratamientos silvícolas evaluados bajo los métodos de manejo MDS y MMOBI, la especie dominante fue Pinus patula, que concentró casi 60 % del total de individuos registrados. Por tratamiento, su abundancia relativa varió entre 41 % en REF y 72 % en CL. Abies religiosa fue la segunda especie más frecuente en CS, mientras que Pinus ayacahuite presentó una participación moderada, especialmente en CA (22 %) y CR (34 %). Las características estructurales de los tratamientos reflejan ciertas diferencias; en particular CA, REF y CS presentaron los mayores valores de cobertura de copa, variable directamente relacionada con la fracción de Pd (Barbier et al., 2009; Van Stan et al., 2020). El tratamiento de CR exhibió valores intermedios de diámetro y altura (Cuadro 1).

Manejo	Tratamiento	UM	Tiempo desde la última intervención silvícola (años)	*DAP* (cm)	Altura total (m)	Diámetro de copa (m)	Área basal (m² ha^-1)	Cobertura arbórea (%)	Densidad (árboles UM^-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
MDS	CA	5	8	26.51±15.25	20.44±9.08	4.03±1.79	33.32±16.66	73.21±3.68	45±18
	CL	5	9	21.38±9.47	16.90±4.70	3.28±1.59	18.19±7.53	69.21±3.85	42±15
	CR	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	CS	5	8	23.43±14.02	19.15±9.22	3.12±1.53	29.84±10.59	70.28±2.14	51±13

Nota: valores expresados como media±desviación estándar. NA = No Aplica; MDS = Método de Desarrollo Silvícola; MMOBI = Método Mexicano de Ordenación de Bosques Irregulares; REF = Rodal de referencia; CA = Corta de aclareo; CL = Corta de liberación; CR = Corta de regeneración; CS = Corta de selección.

Durante el periodo de estudio, la Pd representó entre 76 % y 86 % respecto a la Pt acumulada (1 599 mm) en el total de eventos registrados a lo largo de los tres años de monitoreo. Los valores más altos se registraron en los tratamientos CR (86 %) y REF (84.5 %), lo cual es congruente con la menor altura del arbolado en estos tratamientos, mientras que en CA (76.0 %) y CS (77.4 %) se observaron proporciones más bajas. El tratamiento CL tuvo un valor intermedio (81.3 %). En este sentido, hubo una tendencia general en la que a mayor cobertura de copa, la relación Pd/Pt disminuyó. No obstante, al analizar la información a escala semanal, la fracción de Pd dependió de la magnitud de la precipitación (Figura 6). En periodos semanales de baja intensidad (0-10 mm), la intercepción fue mayor y más variable; en cambio para las semanas de mayor magnitud (>100 mm), la intercepción disminuyó considerablemente, patrón que se mantuvo entre tratamientos. Dicho comportamiento es consistente con el efecto de saturación del dosel, en el cual la capacidad de almacenamiento se supera durante eventos intensos, lo cual favorece un mayor paso de la lluvia hacia el suelo (Brasil et al., 2022; Cisneros-Vaca et al., 2018; Rutter et al., 1971; Van Stan et al., 2020).

Figura 6. Variación de la intercepción en función de la magnitud de la precipitación semanal.

La fracción de precipitación directa observada está dentro del intervalo citado para bosques templados de México y otras regiones ecológicamente similares. Estudios previos han documentado que entre 14 % y 40 % de la Pd puede retenerse en el dosel, lo que implica proporciones de lluvia que alcanzan el suelo comparables a las estimadas en el presente estudio (Bolaños-Sánchez et al., 2021; Flores-Ayala et al., 2016; Huber & Iroumé, 2001; Llorens, 1997).

Efectos del manejo forestal y la estructura del dosel sobre la precipitación directa

Conforme al predictor general de la Ecuación 1, se evaluaron modelos candidatos con diferentes combinaciones de variables estructurales y estructuras aleatorias. El modelo final fue seleccionado con base en el menor AIC y en la coherencia biológica de los efectos estimados (Cuadro 2). En el modelo con mejor desempeño estadístico se consideraron las observaciones a nivel de pluviómetro, e incorporaron como variables explicativas, la cobertura del dosel en cada orientación, la Pt semanal a nivel de pluviómetro y el tratamiento silvícola; ello permitió capturar explícitamente la variabilidad espacial dentro de cada unidad de muestreo.

Cuadro 2. Comparación de modelos candidatos para la estimación de Pd mediante GLMM.

Modelo	Variables fijas	Estructura aleatoria	*AIC*
M1	Pt+Cobertura de copa+Tratamiento	(1\|UM)+UM:Orientación+AR(1)	24 689.40
M2	Pt+Área basal	UM:Orientación+AR(1)	24 698.75

AIC = Criterio de información de Akaike; Pt = Precipitación total; UM = Unidad de muestreo; AR(1) = Estructura autorregresiva de primer orden.

El análisis de varianza global, que consideró el modelo mixto, basado en pruebas de Wald (tipo II) indicó efectos significativos de la Pt (χ²=6 418.99, gl=1, p<0.001), de la cobertura del dosel (χ²=18.95, gl=1, p<0.001) y del tratamiento silvícola (χ²=10.95, gl=4, p=0.036) sobre la precipitación directa. El modelo evidenció que la precipitación directa semanal estuvo estrechamente asociada con la Pt, y aumentó en promedio alrededor de 1.8 % por cada milímetro adicional de lluvia. La cobertura del dosel ejerció un efecto negativo significativo, de modo que incrementos en la cobertura se tradujeron en una menor cantidad de lluvia que alcanza el suelo (Deguchi et al., 2006; Kaushal et al., 2017; Staelens et al., 2006).

Al analizar de manera específica el efecto de los tratamientos, solo la CA presentó una disminución significativa en la Pt de aproximadamente 12 % en comparación con el tratamiento de referencia (Cuadro 3). En contraste, los tratamientos CL, CR y CS no mostraron efectos estadísticamente significativos. En la CS la especie más común fue Abies religiosa, cuya arquitectura de copa se ha asociado con una mayor capacidad de intercepción (Flores-Ayala et al., 2016). Bajo este contexto, podría esperarse una menor precipitación directa en CS; sin embargo, en el presente estudio no se detectaron efectos estadísticamente significativos de la especie dominante sobre la precipitación directa (p=0.42).

Cuadro 3. Coeficientes del modelo (M1) GLMM Gamma para la precipitación directa semanal.

Variable	Estimador (β)	Error estándar	z	p
Intercepto	3.383	0.184	18.39	<0.001
Precipitación total (mm)	0.0182	0.00023	80.12	<0.001
Cobertura del dosel (%)	-0.0104	0.00240	-4.35	<0.001
Tratamiento CA	-0.1305	0.0410	-3.18	0.001
Tratamiento CL	-0.0606	0.0434	-1.40	0.162
Tratamiento CR	-0.0662	0.0491	-1.35	0.177
Tratamiento CS	-0.0750	0.0426	-1.76	0.078

CA = Corta de aclareo; CL = Corta de liberación; CR = Corta de regeneración; CS = Corta de selección. Nota: los coeficientes están estimados en escala logarítmica; la interpretación porcentual se realizó mediante la transformación exponencial

. El componente autorregresivo AR(1) presentó una correlación temporal baja (p=-0.02), con varianza de 0.0883 y desviación estándar de 0.2972.

El modelo seleccionado fue estadísticamente robusto y se consideró apropiado (AIC=24 689.4), ya que los valores observados frente a los predichos se ubicaron alrededor de una línea recta que pasa por el origen, independientemente de los tratamientos, sin evidenciar un patrón sistemático de sesgo. Asimismo, la adecuada aplicación de la distribución Gamma se refleja en la positividad de los valores predichos de precipitación directa. De hecho, la correlación entre la precipitación directa observada y la predicha fue alta en todos los tratamientos, con valores promedio de r=0.98 y pseudo-R²=0.97 (Figura 7). De manera congruente, los errores absolutos medios (MAE=1.16 mm) se mantuvieron bajos y relativamente homogéneos, al igual que la raíz del error cuadrático medio (RMSE=2.28 mm), lo que indica una adecuada capacidad predictiva del modelo para estimar Pt semanal independientemente del tipo de manejo forestal.

Figura 7. Relación entre la precipitación directa observada y la predicha por el Modelo lineal mixto generalizado (GLMM) a nivel de pluviómetro para cada tratamiento silvícola.

De manera complementaria, se evaluó el efecto de la estructura del dosel mediante el área basal, sobre la precipitación directa. En este caso, el área basal tuvo un efecto negativo significativo sobre la precipitación directa (β=-0.032, p<0.01), lo cual implica que una mayor densidad estructural del rodal favorece los procesos de intercepción. En ambos modelos, la Pt se mantuvo como la variable de mayor influencia (p<0.001), lo que evidencia su papel dominante en la determinación de la precipitación directa (Cuadro 4).

Cuadro 4. Coeficientes del modelo (M2) GLMM Gamma para la precipitación directa semanal.

Variable	Estimador (β)	Error estándar	z	p
Intercepto	2.648	0.0377	70.32	<0.001
Precipitación total (mm)	0.0182	0.00023	80.12	<0.001
Área basal (m² ha^-1)	-0.0319	0.0127	-2.50	0.012

Nota: los coeficientes están estimados en escala logarítmica; la interpretación porcentual se realizó mediante la transformación exponencial

. El componente autorregresivo AR(1) presentó una correlación temporal baja (p=-0.03), con varianza de 0.0884 y desviación estándar de 0.2973.

Cabe señalar que los rodales evaluados fueron intervenidos entre 2014 y 2016, por lo que al momento del monitoreo presentaban un periodo de recuperación estructural de 8-10 años. En ese contexto, es posible que los efectos del manejo sobre la precipitación directa estén condicionados no solo por el tipo de tratamiento aplicado, sino también por la densidad arbórea y la estructura del dosel resultante después de la intervención (Kaushal et al., 2017; Monárrez-González et al., 2018; Muñoz-Villers et al., 2012; Staelens et al., 2006). En este sentido, resulta relevante profundizar en el papel del tiempo transcurrido desde la intervención sobre la dinámica de intercepción en condiciones de bosques bajo manejo.

Los resultados evidencian que la variabilidad de la Pd en el bosque bajo manejo depende en gran medida de la Pt; sin embargo, su variabilidad a escala de rodal está determinada principalmente por los atributos estructurales del dosel. En particular, la cobertura de copa y el área basal muestran efectos negativos significativos, lo cual indica que rodales con mayor desarrollo estructural incrementan la retención de agua y reducen la fracción que alcanza el suelo. En contraste, el efecto del tratamiento silvícola es secundario y no consistente entre categorías, ello sugiere que su influencia se manifiesta principalmente de manera indirecta a través de los cambios estructurales inducidos en el dosel, en especial si se consideran los tiempos de recuperación posteriores a la intervención. El análisis descriptivo a escala semanal evidencia que la intercepción depende de la magnitud de la precipitación, con mayores valores en eventos de baja precipitación y una disminución progresiva conforme aumenta la lluvia. En conjunto, estos resultados indican que la respuesta hidrológica bajo las condiciones del estudio está mejor explicada por la estructura actual del dosel que por el tipo de manejo.

Los autores agradecen al ejido Emiliano Zapata, Chignahuapan, Puebla por las facilidades otorgadas para el desarrollo de las actividades de campo. Asimismo, al Instituto Nacional de Investigaciones Forestales, Agrícolas y Pecuarias (INIFAP) por el financiamiento otorgado a través del proyecto “Manejo integrado de recursos forestales para la sustentabilidad de los servicios ecosistémicos ante el cambio climático”.

Los autores declaran no tener conflicto de intereses. Arian Correa-Díaz declara no haber participado en el proceso editorial del artículo.

Jesús Valentín Gutiérrez-García: análisis de resultados y redacción del manuscrito original; Efraín Velasco-Bautista: revisión de la metodología y contribución a la redacción del manuscrito; Arian Correa-Díaz: análisis de la información y apoyo en la redacción del manuscrito; Antonio González-Hernández: revisión de la información y contribución a la redacción del manuscrito; Francisco Moreno-Sánchez: revisión de la información y contribución a la redacción del manuscrito; Bertha Patricia Zamora-Morales: revisión y redacción del manuscrito; Vidal Guerra-De la Cruz: revisión técnica y corrección del manuscrito; Eulogio Flores-Ayala: asesoría científica y revisión final del manuscrito.

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.