Área de estudio

El estudio se realizó en rodales de P. patula localizados en bosques de la comunidad de Ixtlán de Juárez, Oaxaca, entre las coordenadas geográficas 17° 23’ 0.50’’ y 17° 23’ 0.58’’ latitud N y 96° 28’ 45’’ y 96° 28’ 53’’ longitud O, a una altitud promedio de 2780 m. Pinus patula es la especie que tiene mayor distribución en el predio, con aproximadamente 5000 ha como especie dominante, es de crecimiento rápido y posee un interés comercial alto para la región (^{Castellanos-Bolaños et al., 2008}). El predio, gestionado de forma comunitaria, abarca 19,310.14 ha de superficie.

Datos dasométricos

La base de datos que se utilizó fue generada a partir de tres mediciones en 66 sitios permanentes de muestreo de 400 m², las cuales se realizaron en 2015, 2016 y 2017; además, se utilizó una base de datos de 80 sitios de intervalo establecidos en 2007 y remedidos en 2008, con la finalidad de tener mayor representatividad del área de estudio y ampliar el intervalo de edades, que varió de cinco a 76 años.

Las variables registradas en campo fueron diámetro normal de todos los árboles dentro de la parcela con una forcípula (Mantrax Blue^®, Haglöf, Sweden) y la altura total de al menos ocho individuos por sitio con clinómetro electrónico (HCC Clinometer Compass^®, Haglöf, Sweden), de los cuales cuatro se identificaron como dominantes, de acuerdo con el concepto de altura dominante que emplea los 100 árboles más gruesos por hectárea (^{Alder, 1980}; ^{Assmann, 1970}).

Los 146 sitios permitieron obtener un intervalo amplio de distribución de la muestra (5-76 años) con heterogeneidad alta, lo cual permitió mejor apalancamiento de los modelos. El conjunto de datos se caracterizó por tener: una edad promedio de 28.1 ± 16.5 años, con valores mínimo de cinco y máximo de 76 años, la altura dominante promedio fue de 21.5 ± 8.9 m, con un mínimo de 4.2 y de 39.7 m como altura máxima.

Ecuaciones dinámicas

En este estudio se utilizaron ecuaciones polimórficas con asíntotas comunes y variables (tipo ADA y GADA). Los modelos base que se ensayaron han sido documentados en investigaciones sobre biometría forestal; son flexibles, tienen un punto de inflexión y tendencia a alcanzar una asíntota horizontal (^{Hernández-Cuevas et al., 2018}; ^{Kiviste et al., 2002}; ^{Krumland y Eng, 2005}; ^{Rojo-Alboreca et al., 2017}). En los Cuadros 1 y 2 se presentan las ecuaciones ajustadas a los datos.

Ajuste de los parámetros

La técnica de ajuste empleada en la estimación de parámetros de los modelos fue un método invariante con respecto a la edad de referencia (base-age invariant, BAI por sus siglas en inglés), denominado método iterativo anidado (^{Tait et al., 1988}), el cual se utilizó como criterio para lograr la estabilización de los parámetros globales, de tal forma que el error medio cuadrático del modelo entre dos iteraciones fuera menor de 0.0001 (^{Vargas-Larreta et al., 2010}). El método estima los efectos específicos del sitio y asume que siempre existen errores aleatorios y de medición en los datos, los cuales deben ser modelados (^{Cieszewski, 2003}). El ajuste de las ecuaciones se realizó con el procedimiento MODEL de SAS/ETS^® 9.3 (^{SAS Institute, 2011}).

Para reducir los efectos por heterocedasticidad, la varianza del error fue asumida como una función de potencia de la altura predicha (^{Diéguez-Aranda et al., 2005}; ^{Huang et al., 2000}). El factor de ponderación utilizado fue w _i = Ĥ^k ₁, donde k es una constante que tomó el valor de 0.5, para lograr la homogeneidad de los residuales y la consistencia de los estimadores de los parámetros (^{Quiñonez-Barraza et al., 2015}).

Indicadores de ajuste

La bondad de ajuste de las ecuaciones se analizó mediante el análisis numérico y gráfico de los residuales. Los estadísticos utilizados fueron coeficiente de determinación ajustado por el número de parámetros (R ² _adj ), raíz del error cuadrático medio (RECM), sesgo promedio absoluto (Ē) y criterio de información de Akaike (AIC) (^{Basuki et al., 2009}; ^{Zhao et al., 2007}); la autocorrelación de los residuales de los modelos fue medida con el estadístico ^{Durbin-Watson (dw) (Durbin y Watson, 1951}).

La significancia de los parámetros se evaluó mediante la prueba t-Student, bajo H ₀ : β_i = 0, b_i = 0 y H ₁ : β_i ≠ 0, b_i ≠ 0 con un valor de significancia de α = 0.05 (^{Infante y Zárate, 2012}) y se realizó un análisis gráfico de las ecuaciones ajustadas con la finalidad de verificar su capacidad predictiva (^{Goelz y Burk, 1992}; ^{Sharma et al., 2011}).

Ajuste de modelos

Los estadísticos de bondad de ajuste para las seis ecuaciones dinámicas mostraron resultados satisfactorios, las ecuaciones explican 99.8 % de la varianza total observada en la altura dominante, con un error aproximado de 0.33 y 0.34 m para las ecuaciones ADA y GADA, respectivamente; asimismo, todos los estimadores de los parámetros resultaron diferentes de cero a un nivel de significancia del 5 % (P < 0.0001) (Cuadro 3).

Los estadísticos de bondad de ajuste (RECM, R ² _adj , AIC) mostraron que el modelo de Korf (M1) polimórfico (β ₁ ) es el más apropiado para modelar la altura dominante bajo el enfoque ADA; sin embargo, las curvas de crecimiento generadas (Figura 1a) no cubren en su totalidad los datos en las calidades de estación altas.

Figura 1 Curvas de índice de sitio para Pinus patula a la edad base de 40 años: modelos ADA (M1-M3) y Modelos GADA (M4-M6). 

A pesar del ajuste preciso, los modelos Korf (M1) polimórfico (β ₁ ) y Hossfeld IV (M2) polimórfico (β ₂ ) presentaron algunas inconsistencias en la estimación de las alturas a edades jóvenes y, particularmente, en calidades de sitio altas (M2 = IS: 39 m). Adicionalmente, las curvas trazadas no se apegaron a la tendencia de los datos observados; por tal razón, se descartaron estas ecuaciones para modelar el crecimiento en altura dominante en los rodales de P. patula (Figuras 1a y 1b). Resultados similares fueron obtenidos por ^{Corral-Rivas et al. (2004)} y ^{Diéguez-Aranda et al. (2005)} al emplear estas ecuaciones en especies como P. cooperi C. E. Blanco, P. durangensis Martínez, P. engelmannii Carrière y P. sylvestris L.

Las curvas obtenidas con la ecuación de Levakovic II (M3) polimórfica (β ₁ ) fueron mejores al tener una tendencia similar a la observada y características deseables en una curva de crecimiento (Figura 1c); por esta razón, en el enfoque ADA se eligió dicha curva como la más apropiada para representar el patrón de crecimiento en altura dominante de P. patula. Estos resultados son similares a los reportados por ^{Santiago-García et al. (2017)} y ^{Elfving y Kiviste (1997)}, quienes destacaron que el modelo es flexible y predice de manera adecuada los datos observados en un sistema de crecimiento y rendimiento maderable. El estadístico de Durbin-Watson mostró que el modelo M3 (dw = 2.3) bajo el enfoque ADA no presenta problemas de autocorrelación; es decir, los estimadores de los parámetros son eficientes e insesgados; además, las curvas de crecimiento presentaron una tendencia biológica más realista al cubrir adecuadamente la variabilidad de la muestra.

Las curvas obtenidas con los modelos bajo el enfoque GADA mostraron un patrón de crecimiento más apegado a la trayectoria de los datos observados que las ecuaciones ADA. Por otro lado, la distribución de los residuales estandarizados presentó un comportamiento adecuado para ambos enfoques (Figura 2); sin embargo, los indicadores de bondad de ajuste indicaron que los modelos de Bertalanffy-Richards (M4) y de Hossfeld IV (M5) fueron los mejores al presentar características razonables y una tendencia biológica deseable en las curvas, además de cubrir completamente el intervalo de alturas de la muestra (Figuras 1d y 1e); no obstante, el modelo de BertalanffyRichards (M4) mostró mayor capacidad para modelar el patrón de crecimiento en altura dominante, y la corrección de la heterocedasticidad generó como resultado que los residuales estandarizados fueran homogéneos (Figura 2d); además, la prueba de Durbin-Watson reveló que los residuales no se encuentran correlacionados (dw = 2.3) (Cuadro 3). Estos resultados son similares a los reportados en los trabajos de ^{Corral-Rivas et al. (2004)}, ^{Vargas-Larreta et al. (2013)} y ^{Rojo-Alboreca et al. (2017)}, quienes obtuvieron curvas ADA y GADA con tendencia biológica deseable para estimar la productividad forestal en especies como P. cooperi C. E. Blanco, P. pseudostrobus Lindl., P. halepensis Mill., entre otras.

Figura 2 Gráficos de residuales estandarizados de las ecuaciones dinámicas ajustadas: modelos ADA (M1-M3) y modelos GADA (M4-M6). 

En los últimos años las metodologías ADA y GADA se han empleado para la modelación del crecimiento en altura dominante en el área de estudio; ^{Hernández-Cuevas et al. (2018)} ajustaron modelos de crecimiento ADA y GADA a partir de datos de análisis troncal de P. ayacahuite Ehrenb. ex Schltdl., y concluyeron que los modelos de Hossfeld IV bajo el enfoque GADA y Chapman-Richards polimórfico tipo ADA (β ₂ ) fueron los más apropiados para clasificar la productividad de esta especie. Por otro lado, ^{Castillo-López et al. (2018)} generaron curvas dinámicas de índice de sitio tipo GADA con el modelo de Von ^{Bertalanffy (1949}; ¹⁹⁵⁷⁾ para P. oaxacana Mirov, P. douglasiana Martínez, P. patula Schiede ex Schltdl. & Cham. y P. pseudostrobus Lindl. a nivel de la Unidad de Manejo Forestal (UMAFOR) 2001 que abarca la región de la Sierra Norte del estado de Oaxaca, México. Estos estudios demuestran que el enfoque GADA es una de las metodologías más utilizadas en la actualidad para modelar el crecimiento en altura dominante; asimismo, el enfoque ADA ha brindado resultados satisfactorios porque las ecuaciones dinámicas construidas con esta metodología para P. patula han permitido construir tablas de rendimiento de densidad variable y simular infinidad de escenarios silvícolas (^{Santiago-García et al., 2017}).

Las curvas generadas con los modelos M3 (ADA) y M4 (GADA) describen el crecimiento en altura dominante de manera apropiada porque cubren adecuadamente la trayectoria de los datos observados; además, las curvas generadas bajo el enfoque GADA son curvas polimórficas más flexibles, con tasa de crecimiento y potencialidades variables (Figura 3).

Figura 3 Curvas polimórficas obtenidas bajo el enfoque ADA (M3) y de polimorfismo complejo con GADA (M4) a una edad 

Con base en los indicadores de bondad de ajuste, las curvas de crecimiento de IS generadas, los gráficos de residuales, y por la estructura compacta de la ecuación, el modelo M4 es el más apropiado para modelar el patrón de crecimiento en altura dominante de los rodales de P. patula en los bosques de Ixtlán de Juárez, Oaxaca (Cuadro 3; Figuras 1d y 2d).

De esta manera, las curvas de índice de sitio generadas con el procedimiento GADA (Figura 1d) permitirán clasificar los bosques naturales de P. patula de acuerdo con su potencial productivo, y planificar las actividades silvícolas (^{McKenney y Pedlar, 2003}) al hacer posible la predicción del crecimiento y rendimiento maderable a nivel de rodal, por categoría diamétrica y árbol individual (^{Santiago-García et al., 2014}; ^{Santiago-García et al., 2017}).

Comparación de ecuaciones de altura dominante para Pinus patula

La obtención de datos recientes permitió desarrollar una ecuación dinámica de crecimiento de altura dominante para clasificar el potencial productivo de los rodales de P. patula en Ixtlán de Juárez, Oaxaca. La ecuación propuesta por ^{Pérez-López et al. (2017)} (Modelo previo), que expresa el índice de sitio, utiliza una edad base de 40 años, por lo que para efectos comparativos en este estudio se utilizó la misma edad de referencia; la comparación gráfica entre el modelo previo y el modelo M4 (actual) se muestra en la Figura 4. Al respecto, el modelo actual generó un patrón de crecimiento más coherente con los datos observados que el modelo previo (Figura 4).

Figura 4 Comparación de las curvas de crecimiento en altura dominante con ecuaciones obtenidas en investigaciones previas y la actual a una edad base de 40 años. 

Las curvas polimórficas y con polimorfismo asintótico generadas por ^{Pérez-López et al. (2017)} y ^{Castillo-López et al. (2018)}, respectivamente, donde las primeras subestiman y las segundas sobreestiman la altura de los árboles antes y después de la edad base (Figura 4). Esto era de esperarse para el modelo previo, porque tales curvas se construyeron con menor información bajo la metodología ADA, en tanto que en el modelo de ^{Castillo-López et al. (2018)} se empleó el enfoque GADA con información de análisis troncales, pero un hecho importante a destacar es que se desarrolló a nivel de UMAFOR, por lo que la escala espacial es mucho mayor que la usada en este estudio; por tanto, se esperaría que una ecuación a nivel local permita modelar mejor las condiciones específicas del sitio.

Es importante definir la ecuación que describa mejor el comportamiento del crecimiento en altura dominante, puesto que, al interrelacionarse con un sistema de crecimiento y rendimiento permite predecir la edad óptima de cosecha; así mismo, cuando se combina con herramientas de manejo de la densidad permite definir posibilidades maderables en aclareos y optimizar el turno de los rodales (^{Santiago-García et al., 2013}).

Con el modelo de ^{Castillo-López et al. (2018)} hay una sub- y sobreestimación a edades jóvenes; sin embargo, el modelo propuesto por ^{Pérez-López et al. (2017)} es ligeramente mejor en calidades bajas y altas (Figura 4). Este tipo de situaciones es común en la modelación forestal, debido a que los extremos climáticos anormales pueden influir más en la altura total de árboles jóvenes que en árboles adultos, y con frecuencia es determinado por factores distintos a la calidad de sitio (^{Huang, 1999}); al respecto, ^{Vargas-Larreta et al. (2013)} señalaron que no es deseable seleccionar una edad de referencia demasiado joven para clasificar la calidad del sitio de los rodales.

La ecuación actual (M4) genera predicciones mejores en la amplitud de la base de datos utilizada, en comparación con las ecuaciones previas contrastadas (^{Castillo-López et al., 2018}; ^{Pérez-López et al., 2017}) (Figura 4). El propósito mayor de modelar el crecimiento en altura dominante no es la predicción de la altura, sino la selección de un patrón de crecimiento de la altura dominante que el rodal puede seguir durante todo su desarrollo (^{Clutter et al., 1983}).