Introducción
La proyección del rendimiento maderable en el tiempo permite planificar el manejo forestal, en especial con silvicultura intensiva. Los modelos de proyección son útiles en los proyectos de plantaciones forestales comerciales, para establecer la valoración en el tiempo y prever el abastecimiento de cada mercado (Clutter et al., 1983; Burkhart y Tomé, 2012).
La empresa Agroproductores Forestales de Zacualpan, S. P. R. de R. L. ha establecido 3500 ha en 15 años. Al inicio el objetivo fue restaurar áreas degradadas, y fue convertirlas en plantaciones forestales comerciales. Ahora en la región de Zacualpan no hay información cuantitativa para un manejo forestal eficiente y se realizan actividades silvícolas basadas en la experiencia cotidiana. Además, el mercado local de productos maderables está basado en una selección diamétrica, por lo cual es importante tener un sistema de proyección robusto que permita tomar decisiones de forma oportuna.
Para modelar el AB se utilizan estructuras matemáticas influenciadas por la altura dominante (A), el número de árboles por unidad de superficie, la edad (E) y en ocasiones el índice de sitio (IS) (Gonzalez-Benecke et al., 2012; Santiago-García et al., 2013). La variable número de árboles residuales (NA) se explica en función de la edad (E) (Rose et al., 2004), además de indicadores de estado del sitio. Asimismo, el volumen (V) se puede explicar al considerar como variables independientes el área basal, la altura dominante o índice de sitio y la edad (Magaña et al. 2008).
Los sistemas de crecimiento y rendimiento maderable (SCRM) son de tipo explícito o implícito. Un sistema explícito puede ajustarse en forma simultánea o individual, y el sistema implícito es un complemento del sistema explícito y permite estimar las variables de estado por categoría diamétrica en forma dinámica (Torres-Rojo et al., 1992). Ambos procedimientos son buenos y ello dependerá de la especie, lugar y tipo de datos utilizados para el ajuste.
La función usada para el sistema implícito es el modelo Weibull de tres parámetros (Weibull 3P). Éste se implementó para proyectar condiciones del rodal en Pinus rudis (Magaña et al., 2008), P. patula (Santiago-García et al., 2014) y P. arizonica, P. duranguensis, P. teocote, P. leiophylla, P. lumholtzii, P. ayacahuite (Quiñonez et al., 2015a). En otros estudios se usan modelos de distribución como el de Burr, distribución logística, y percentiles libres o proyección tabulada del rodal (Stand-Table-Projection) (Borders et al., 1987; Knowe et al., 1997; Lindsay et al., 1996; Wang y Rennolls, 2005).
El método de percentiles se usa para predecir los parámetros del modelo Weibull 3P (localización, escala y forma) y describir la distribución de las categorías diamétricas. El percentil utiliza como variable dependiente el diámetro medio cuadrático (Dq, cm) (Torres-Rojo et al., 1992; Santiago-García et al., 2014; Quiñonez et al., 2015a).
Los SCRM son herramientas para programar el plan de cortas, el abastecimiento de materia prima a la industria, implementar aclareos, estimar ingresos económicos futuros y sus usos se han expandido para cuantificar biomasa y carbono (Clutter et al., 1983; Pienaar y Rheney, 1993, Galán et al., 2008).
El objetivo de este estudio fue construir un SCRM dinámico, que incorpore modelos de rodal y de distribución diamétrica. La hipótesis fue que parcelas de observación permanente remedidas en el tiempo describen el patrón de crecimiento y rendimiento lo que, integrado a una función de distribución de probabilidades, permite proyectar el volumen maderable por categoría diamétrica (Kiviste y Kiviste, 2009; Gonzalez-Benecke et al., 2012).
Materiales y Métodos
Área de estudio y datos dasométricos
El área de estudio pertenece al municipio de Zacualpan, Veracruz, México (20° 25’ y 20° 31’ N, 98° 24’ y 98° 25’ O, con 1670 m de altitud). El clima es templado húmedo con una temperatura promedio de 18 °C y una precipitación media anual de 1900 mm; el suelo es de tipo Feozem.
En el 2011 se establecieron 50 parcelas de observación permanente, las cuales se remidieron en 2012. De estas, en el 2015 se encontraron 14 parcelas aclareadas, 26 con mortalidad parcial (competencia intra-específica), siete con aprovechamiento total y tres con muerte total de las plantas por dominancia de maleza. Las parcelas son cuadradas, de 20x20 m (400 m2). En cada período se midió el diámetro normal de todos los árboles (DN en cm) y la altura de los cinco árboles dominantes (A en m) y se registró la edad de la plantación (E en años).
Las variables de estado por hectárea estimadas fueron: altura dominante promedio (A, m), volumen total (V, m3), diámetro cuadrático (Dq, cm), área basal (AB, m2) y número de árboles residual (NA), todas referidas al índice de sitio (IS) como criterio para clasificar la calidad de estación (Clutter et al., 1983).
Sistema de crecimiento explícito
El SCRM explícito se estructuró integrando modelos de predicción y proyección para las variables A, AB, V y NA. El resultado fueron modelos de estimar el rendimiento futuro por hectárea y por calidad de estación determinado por un índice de sitio (Clutter et al., 1983; Burkhart y Tomé, 2012). La construcción del sistema inició con la selección de un modelo genérico para altura dominante, y después por el método de diferencia algebraica (ADA) se derivó una condición futura en función de una situación inicial (Pienaar y Rheney, 1993; Jordan et al., 2006). La altura dominante fue ajustada con regresión no lineal en forma independiente y el AB, V y NA como sistema compatible (Santiago-García et al., 2013).
Altura dominante
En el ajuste de los modelos en altura dominante (A) e índice de sitio (IS) se utilizaron 76 pares de datos no traslapados que corresponden al promedio de las alturas dominantes con sus edades que oscilaron de 3 a 21 años. A partir del modelo de predicción se derivaron los modelos de proyección e índice de sitio. El método de diferencia algebraica (ADA) genera familias de curvas anamórficas y polimórficas de la forma A 2=f (A 1, E 1, E 2, β) donde: A 2 es la altura dominante proyectada o índice de sitio a la edad E 2 (edad de proyección), A 1 es la altura dominante (a una edad base de 15 años) medida a la edad inicial (E 1) y β es el vector de parámetros de regresión (De los Santos-Posadas et al., 2006; Quiñonez-Barraza et al., 2015). La expresión de Schumacher-Korf se utilizó en su forma de curva guía (1) y polimórfica-2 (2) con hipótesis de crecimiento en el parámetro a2 del modelo genérico (1) (Schumacher, 1939).
donde A: altura dominante de predicción, A1: altura dominante a la edad inicial, A 2: altura dominante a la edad proyectada, α0, α1, α2: parámetros de la regresión, E1: Edad inicial o edad base (15 años), E2: Edad de proyección, e: función exponencial y ln: logaritmo neperiano.
Área basal
Veinticuatro parcelas de observación permanente fueron segregadas de la muestra del 2015 por la realización de aclareos, cortas definitivas o muerte total. Para el ajuste de área basal, volumen y número de árboles por hectárea se usaron 73 pares de datos sin traslape. El modelo de área basal propuesto (3) es una función dependiente de la altura dominante, la edad y el número de árboles por hectárea. Con esta expresión se generó el modelo (4) para estimar el área basal a la edad proyectada (AB 2) (Gonzalez-Benecke et al., 2012). Las expresiones son:
donde AB1, A1, E1, NA1: área basal, altura dominante, edad y número de árboles por ha a la edad inicial, AB2, A2, E2, NA2: área basal, altura dominante, edad y número de árboles por ha a la edad proyectada y β0, β1, β2, β3: parámetros a ser estimados.
Volumen
Para el volumen se usó el par de ecuaciones (5) y (6) que describen el volumen de predicción y de proyección respectivamente:
donde V1, AB1, A1: son volumen, área basal y altura dominante en su condición inicial, V2, AB2, A 2: volumen, área basal y altura dominante en la condición proyectada y γ0, γ 1, γ 2: parámetros a ser estimados.
Número de árboles residual
En el modelo de número de árboles residual destaca un parámetro que representa la tasa de mortalidad (ω1); este parámetro multiplica directamente a la edad de referencia. La sobrevivencia puede expresarse en términos de predicción (7) y proyección (8) de la forma siguiente:
donde NA1, E1, NA2, E2: Número de árboles por hectárea y edad en la condición inicial y proyectada respectivamente y ω0, ω1: parámetros de la regresión.
Sistema de crecimiento implícito mediante la fdp Weibull
Weibull (1951) publicó una función de distribución de probabilidad (fdp) y su potencia es tal que ha incursionado incluso en la dasonomía (Villaseñor-Alba y Villanueva-Morales, 2000). Para ajustar la función Weibull y modelar la distribución de clases diamétricas existen métodos como el de regresión lineal y no lineal, percentiles y máxima verosimilitud (Torres-Rojo et al., 1992; Torres-Rojo et al., 2000). El método de estimación de parámetros a través de percentiles ha mostrado buenos resultados (Santiago-García et al., 2014).
Las variables de estado propuestas fueron: área basal por hectárea (AB, m2 ha-1) y diámetro medio cuadrático (Dq, cm) dado por:
La altura de los árboles no medidos se obtuvo en forma indirecta con la ecuación H i =0.880151xA 1.026997x(1-e (-1.075746xDN) ) que relaciona la altura dominante (A) y el diámetro normal (DN) para estimar la altura total individual. El volumen total árbol con corteza se obtuvo del modelo Schumacher: Vix0.000114xDN 1.718158xH 0.939592 (Uranga-Valencia et al., 2015) dónde H: altura total del árbol y DN: diámetro normal. Los percentiles propuestos fueron:
donde p0, p63, p95: corresponden al percentil cero, 63 y 95, Dq: diámetro medio cuadrático, AB: área basal y A: altura dominante, φ0, φ1, φ2, ω0, ω1, ξ0, ξ 1, ξ 2 son los parámetros a ser estimados en el sistema de percentiles. La función de distribución de probabilidades Weibull 3P está dada por la función probabilística:
(a≤ x≤ ∞); =0, de otra forma, donde a es el parámetro de localización, b el de escala y c el de forma. La función define el valor de densidad de probabilidad asociada a una variable aleatoria, que para este estudio fue el diámetro normal (x=DN) (Santiago-García et al., 2014). La distribución acumulada en su forma cerrada es:
Los valores de los parámetros b y c son positivos, en tanto que a puede ser negativo cero o positivo. Sin embargo en la modelación de distribución de clases diamétricas a no debe ser negativo (Clutter et al., 1983). Conociendo los parámetros puede obtenerse la función de densidad para un intervalo dado por:
donde P es la proporción de árboles de la clase diamétrica asignada, L y U son el límite inferior y superior de la clase respectivamente, e es la función exponencial y X es el centro de la categoría diamétrica.
Los parámetros fueron estimados por el método de momentos (Pienaar y Rheney, 1993), el cual ha sido probado con resultados satisfactorios (Santiago-García et al., 2014; Quiñonez-Barraza et al., 2015b), utilizando las siguientes ecuaciones:
El parámetro a de localización se determinó con:
El parámetro c está dado por c=ln [((-ln(1-0.95))/(-ln(1-0.63)))/ln(p95-a)/p63-a))]. Para el parámetro b se utilizó la siguiente expresión del segundo momento de la distribución Weibull:
donde Γ1= Γ(1+1/c), Γ2= Γ(1+2/c), Γ(.) es la función Gamma, las demás variables ya fueron definidas previamente.
Variables de estado por categoría diamétrica
El número de árboles por clase diamétrica se obtuvo de la forma siguiente: NA (CD) =P(L<X<U)xNAHA, donde NA (CD) es la distribución de la clase diamétrica por hectárea, P(L<X<U) es la función de densidad de probabilidad para el intervalo de esa misma categoría (CD) y NAHA: número de árboles por hectárea residuales a una edad determinada. El volumen por hectárea por categoría diamétrica se obtuvo con la expresión V (CD) =NA (CD) xVi, donde NA (CD) es el número de árboles de la distribución de la clase diamétrica.
Técnica y prueba de bondad de ajuste
El ajuste del modelo de altura dominante e índice de sitio se realizó mediante regresión no lineal con el procedimiento MODEL del software SAS/ETS® (SAS, 2014). El coeficiente de determinación ajustado por el número de parámetros (R2 adj), la suma de cuadrados del error (SCE), el cuadrado medio del error (CME), el nivel de significancia de los parámetros (Pr>|t|), y la modelación lógica de las curvas con respecto al gráfico de datos observados en campo fueron los criterios de bondad para validar el modelo (De los Santos-Posadas et al., 2006).
Los modelos compatibles de área basal, volumen y mortalidad fueron ajustados simultáneamente con regresión aparentemente no relacionada (SUR, siglas en inglés), esta técnica asume que las variables son independientes entre ellas. La prueba de Kolmogorov-Smirnov (KS) fue utilizada para calificar la bondad de ajuste de los percentiles. Para descartar parcelas se usó un nivel de significancia de α=0.05 y 0.1 (Santiago-García et al., 2014; Quiñonez et al., 2015a).
Resultados y Discusión
Altura dominante y rendimiento explícito
Para altura dominante los modelos de Hossfeld IV, Schumacher-Korf y Chapman- Richards tienen los mejores ajustes y de ellos derivan familias de índice de sitio (IS) (De los Santos-Posadas et al., 2006; Quiñonez-Barraza et al., 2015). El modelo (2) ajustado para predecir la altura dominante describe coherentemente el patrón de crecimiento de las parcelas de observación permanente (Figura 1). Este permitió derivar cuatro calidades de estación etiquetadas por el índice de sitio (IS, m) a una edad base de 15 años, siendo IS=10 para rodales de baja productividad e IS=25 para rodales de mayor potencial productivo. El IS=20 fue considerado como la condición más representativa. El estimador del parámetro a0 relacionado con la asíntota tiene un valor de rechazo para el modelo (1) esto implicaría que la altura de los árboles alcanzara 125 m. Por el contrario el modelo (2) tiene un valor proyectado de 30.8 m de altura (Cuadro 1).
Área basal, volumen y número de árboles residual
Pienaar y Rheney (1993), Galán et al. (2008) y Magaña et al. (2008) generaron sistemas de crecimiento y rendimiento maderable (SCRM) relacionando modelos de altura dominante (A, m), área basal (AB, m2 ha-1), volumen (V, m3 ha-1) y mortalidad o número de árboles (NA), con diferentes técnicas estadísticas y los resultados fueron buenos. El sistema de ecuaciones del presente estudio generó resultados adecuados. Los criterios de bondad de ajuste (Cuadro 2) son aceptables y con altos niveles de significancia en sus parámetros (Cuadro 3). Con el sistema de ecuaciones y los parámetros estimados en combinación con la familia de curvas polimórficas se construyó la tabla de rendimiento maderable por calidad de estación y calidad promedio (Cuadro 4). El valor de R2 del modelo de volumen es bajo lo que implica una variación muy alta entre los valores de densidad inicial. Al año uno se fijó el parámetro ω0 con 1100 árboles por hectárea, que es la densidad inicial a la que se realiza la mayoría de las plantaciones de P patula en Zacualpan, Veracruz. La proyección de las curvas de crecimiento describe congruentemente el patrón observado en las parcelas de observación permanentes (Figura 2A).
Pr>|t|: prueba de significancia de parámetros; <0.001 es altamente significativo. β0, β1, β2, γ0, γ1, ω0, ω1. Parámetros de la regresión.
NA †: número de árboles por hectárea, IS: índice de sitio en metros a la edad base de 15 años y Edad en años.
El turno técnico para la calidad de estación con IS 10, 15, 20 y 25 ocurre a los 21, 18, 14 y 11 años, respectivamente. El incremento medio anual para estas calidades es de 7, 12, 18 y 26 m3 ha-1 año-1 (Figura 2B).
El tener tres mediciones sobre la misma unidad muestral y parcelas con edades representativas en las plantaciones, resultó en un buen ajuste de los modelos. En consecuencia permitió describir de manera adecuada el patrón de crecimiento de las variables de interés.
En bosques del municipio de Zacualtipán, Hidalgo, regenerados de manera natural con P patula, se ha registrado a la edad de 20 años, un rendimiento en volumen de 10.2 m3 ha-1 año-1 en calidades de estación con IS de 26 m a una edad base de 40 años (Santiago-García et al. 2013). Lo anterior significa que las plantaciones muestran mejor rendimiento comparado con bosques naturales.
La tasa de cambio (parámetro) que se traduce en mortalidad del número de árboles por hectárea es resultado de la competencia intra-específica del rodal. En Zacualtipán se reportó un parámetro de -0.0337 (Santiago-García et al. 2013) y en nuestra investigación fue un poco mayor (-0.04). El crecimiento en área basal en combinación con una secuela de aclareos es útil para simular modelos de respuesta a aclareos.
Sistema de crecimiento y rendimiento implícito con la fdp Weibull
Los estadísticos de bondad de ajuste y los valores estimados en cada parámetro permiten aplicar los modelos propuestos. La prueba KS con un α=0.05 [C(α)=1.36] indicó el rechazo de cuatro de las 50 parcelas (8 %) y con un α=0.1 [C(α)=1.22] un rechazo de cinco (10 %). Esros resultados son similares a los de Quiñonez et al. (2015a) y Santiago-García et al. (2013), estos últimos con 1.58 y 3.96 % de rechazo para α=0.05 y α=0.1, respectivamente. La prueba KS permitió confirmar que los datos de parcelas tienen una distribución Weibull. Al respecto, a mayor nivel de significancia el contraste de las desviaciones se reduce (Santiago-García et al., 2014).
El procedimiento SUR utiliza las correlaciones de los errores de los modelos para mejorar el ajuste. Para modelar Weibull 3P se requieren al menos dos percentiles, y la elección de estos dependen de la muestra (Santiago-García et al., 2014; Borders et al., 1987). La eficiencia del percentil 63 y 95 se muestra en los Cuadros 5 y 6. Al estar centrados se reduce el sesgo de predicción. En estudios similares han utilizado percentiles 50 y 90, 40 y 82 y 24 y 93 (Santiago-García et al., 2014).
La distribución diamétrica a la edad de 20 años indica que rodales con calidad de estación promedio (IS 20) tendrán 23 % más árboles para aserrío (DN>25 cm) comparado con sitios de calidad pobre (IS 10), lo cual repercute en el valor de la producción, a mayor diámetro el precio de venta del producto es más alto (Figuras 3 y 4). En IS 25 a la edad de 10 años se proyecta obtener 51 árboles con DN>30 cm; para clase de sitio pobre se proyecta cero árboles con esa misma categoría diamétrica.
Resultados promisorios se esperan para IS de 20 m; porque a la edad de 15 años se pronostican 429 árboles ha-1 con DN<30 cm y 151>30 cm ha-1, que corresponde a un volumen de 155 y 103 m3 ha-1, respectivamente. Santiago-García et al. (2014) presentan para la misma especie (sitio promedio) pero en bosque natural bajo manejo de Zacualtipán, Hidalgo, una proyección a 20 años de 70 árboles con DN>30 cm y 750 individuos con DN<30 cm.
Los resultados muestran una amplia variación en las variables de estado, lo anterior motivado por la calidad de estación. En IS 20, a la edad del turno técnico (14 años) el 85 % del volumen será comercial, mientras que en IS 10 el turno técnico se proyecta a los 21 años con similar producción maderable.
Conclusiones
A partir de las proyecciones por calidad de sitio se podrán programar las cortas finales a una edad conocida. El SCRM permite generar zonificación de las áreas ya plantadas en la región para establecer prioridades en los planes de expansión de esta actividad productiva, y también permite contar con elementos técnicos para realizar análisis financieros de los proyectos de plantaciones en la región de Zacualpan, Veracruz. La modelación explícita proyecta el volumen sin considerar categorías diamétricas. La forma implícita agrupa por clases diamétricas, por tanto el volumen será mayor en el primer sistema que en el segundo. El volumen por hectárea a 30 años en IS 25 producirá 522 y 740 m3 ha-1 con rendimiento explícito e implícito, respectivamente. El crecimiento libre (sin cortas intermedias) expresado en los sistemas de crecimiento aporta herramientas para el manejo como el turno, el volumen a futuro y la distribución diamétrica para contemplar ingresos por venta de madera.