Revista Mexicana de Ciencias Forestales Vol. 9 (49)  

Septiembre-Octubre (2018)

Fecha de recepción/Reception date: 10 de diciembre de 2017

Fecha de aceptación/Acceptance date: 30 de julio de 2018

_______________________________

1Campo Experimental Valle del Guadiana, INIFAP. México.  correo-e: quinonez.geronimo@inifap.gob.mx

2Posgrado en Ciencias Forestales, Universidad Autónoma de Nuevo León. México.

3Facultad de Ciencias Forestales, Universidad Autónoma de Nuevo León. México.

Abstract

In modeling of taper functions and dominant height growth with time series data, the presence of heteroscedasticity and autocorrelation in residuals is common. Variance Functions (varFunc) and correlation structures (corStruct) were used to correct heteroscedasticity and autocorrelation; both were combined and evaluated through taper and height growth equations for Pinus teocote in Durango, Mexico. A dataset of 51 stems analysis with 768 taper observations and 634 height growth observations was used. The varFuncs applied were: 1) power function (varPower); 2) exponential function (varExp); 3) constant plus power function (varConstPower); and 4) a combination of power and exponential functions (varComb). The corStructs were: compound symmetry (corCompSymm), autoregressive of order 1 (corAR1), continuous-time autoregressive of order 1 (corCAR1), autoregressive-moving average (corARMA2-0), corARMA1-1, corARMA2-1, corARMA2-2, corARMA3-1 and corARMA3-2. To fit the equations, the generalized nonlinear least squares method was used and evaluated with a rating system through: RMSE, R2, AIC, BIC, LogLik, VC and average bias. According to the rating system, the best combinations for taper and height growth equations were 1-9, 2-5, 3-8 and 4-6 and 1-6, 2-9, 3-7 and 4-4, respectively. In the taper equation, only the combination 2-5 was homoscedastic with independent residuals, and the selected height growth equations were homoscedastic with independent residuals; the varFunc and corStruct had influence on the trajectories of site index curves.

Key words: Taper, dominant height, correlation structures, variance functions, Pinus teocote Schiede ex Schltdl. & Cham., residuals.

Resumen

En la modelación del ahusamiento y del crecimiento en altura dominante con datos de series de tiempo, es muy común la presencia de heterocedasticidad y autocorrelación de los errores. Funciones de varianza (varFunc) y estructuras de correlación (corStruct) para corregir la heterocedasticidad y modelar dependencia de los errores, respectivamente. Estas fueron combinadas y evaluadas en ecuaciones de ahusamiento y crecimiento en altura de Pinus teocote en Durango, México. La base de datos se obtuvo de 51 análisis troncales con 768 observaciones de ahusamiento y 634 de altura. Las varFunc utilizadas fueron: 1) función de potencia (varPower); 2) función exponencial (varExp); 3) función constante y de potencia (varConstPower); y 4) función combinada de potencia y exponencial (varComb). Las corStruct incluyeron: simetría compuesta (corCompSymm), autorregresiva de orden 1 (corAR1), autorregresiva continua (corCAR1), autorregresiva de media móvil (corARMA2-0), corARMA1-1, corARMA2-1, corARMA2-2, corARMA3-1 y corARMA3-2. Las ecuaciones se ajustaron por mínimos cuadrados generalizados no lineales; y se evaluaron con un sistema de calificación con los estadísticos de ajuste: RMSE, R2, AIC, BIC, LogLik, CV y sesgo promedio. Con base en la calificación, las mejores combinaciones para el ahusamiento y crecimiento en altura fueron 1-9, 2-5, 3-8 y 4-6 y 1-6, 2-9, 3-7 y 4-4, respectivamente. En el ahusamiento solo la combinación 2-5 fue homocedástica con residuales independientes al igual que las ecuaciones de altura seleccionadas y las varFunc y corStruct presentaron influencia en la trayectoria de las curvas de índice de sitio construidas.

Palabras clave: Ahusamiento, altura dominante, estructuras de correlación, funciones de varianza, Pinus teocote Schiede ex Schltdl. & Cham., residuales.

Introducción

Para la planeación, ejecución y monitoreo del manejo forestal sustentable, se requiere la generación de investigación que apoye la toma de decisiones y la evaluación de los objetivos establecidos. La estimación de las existencias maderables y productividad de los rodales es un objetivo principal en los sistemas de manejo forestal; por tanto, resulta imprescindible conocer el crecimiento de especies comerciales (Aguirre-Calderón, 2015; Salas et al., 2016).

Las investigaciones sobre la estimación del volumen, crecimiento e incremento son herramientas clave para comprender la dinámica de los ecosistemas con manejo; por lo anterior, estos enfoques siguen siendo necesarios para la planeación y ejecución de las actividades forestales. Las ecuaciones de ahusamiento y de crecimiento en altura dominante han sido temas muy explorados (Castillo et al., 2013; Paulo et al., 2015; Rodríguez-Carrillo et al., 2015; Krisnawati, 2016; Corral-Rivas et al., 2017; Fierros-Mateo et al., 2017; Tamarit et al., 2017), porque representan herramientas biométricas que permiten caracterizar el perfil de los árboles, estimar el volumen total y el volumen comercial, así como la productividad de terrenos forestales a través del índice de sitio, respectivamente (Crecente-Campo et al., 2013; Quiñonez-Barraza et al., 2015; Özçelik y Crecente-Campo, 2016).

Dado a que las bases de datos que se utilizan para ajustar las ecuaciones de ahusamiento y crecimiento en altura son series de tiempo que se obtienen a partir de la observación de variables medidas en el mismo árbol, derivadas de análisis de ahusamiento y análisis troncales, es razonable asumir que las observaciones en cada árbol estén correlacionadas y por consecuencia, las perturbaciones o residuales de las ecuaciones ajustadas (Arias-Rodil et al., 2015; Quiñonez-Barraza et al., 2015; Corral-Rivas et al., 2017); además las predicciones podrían mantener variación en los niveles de las variables independientes, lo que se conoce como heterocedasticidad (Gujarati y Porter, 2011).

El término autocorrelación se refiere a la correlación de los residuales de un modelo de regresión, cuando se trabaja con series de observaciones ordenadas en el tiempo, como en datos de series de tiempo, o en el espacio y en datos de corte transversal. Los modelos de regresión lineal y no lineal tienen el supuesto teórico que los residuos tienen la misma varianza y por tanto son homocedásticos. La presencia de autocorrelación y heterocedasticidad da lugar a estimaciones de parámetros de varianza no mínima e intervalos de predicción poco fiables, especialmente en la construcción de ecuaciones de ahusamiento y volumen (Fortin et al., 2007; Xu et al., 2014; Tang et al., 2016). En consecuencia, las pruebas usuales de  o  no son válidas; por lo tanto, el uso de estimadores por mínimos cuadrados generalizados (gnls), con funciones de varianza y estructuras de correlación son una alternativa para generar los mejores estimadores lineales insesgados (Gujarati y Porter, 2011).

En estudios de modelos de ahusamiento y crecimiento en altura dominante, se han utilizado estructuras de correlación y funciones de potencia para corregir la autocorrelación y heterocedasticidad de los residuales (Quiñonez-Barraza et al., 2014; Quiñonez-Barraza et al., 2015; Sharma et al., 2015; Özçelik y Crecente-Campo, 2016; Corral-Rivas et al., 2017; Tamarit et al., 2017), en estos estudios se destacan estructuras autorregresivas continuas de los errores de orden uno, dos o tres y funciones de potencia con varianza conocida.

Por la importancia de mejorar la capacidad predictiva y la interpretación de las propiedades estadísticas en el ajuste de ecuaciones, el objetivo del presente trabajo fue evaluar la combinación de funciones de varianza con estructuras de correlación, para modelar la heterocedasticidad y la dependencia de los errores en ecuaciones de ahusamiento y de crecimiento en altura dominante de Pinus teocote Schiede ex Schltdl. & Cham. en Durango, México.

Materiales y Métodos

Área de estudio y descripción de datos experimentales

Se utilizó información de análisis troncales de 51 árboles de Pinus teocote recolectados en rodales mezclados de la Unidad de Manejo Forestal (Umafor) 1005 Santiago Papasquiaro y Anexos, en la región noroeste de Durango México. El área de estudio se ubica en la Sierra Madre Occidental, entre 24°48’16.98’’ y 25°13’47.25’’ LN y 105°53’9.81’’ y 106°12’52.58’’ LO. El predio forestal fue el ejido San Diego de Tezains, con una superficie de 61 098.25 ha, de ellas 26 636.09 ha son de producción maderable con manejo forestal (Quiñonez-Barraza et al., 2014). Los tipos de climas predominantes son templado cálido húmedo y templado subhúmedo, con precipitación media anual de 1 375 mm. Las temperaturas medias varían de 8 °C en las zonas más altas a 24 °C en las zonas bajas, en las cuales la altitud media llega a 600 m (García, 1981; Quiñonez-Barraza et al., 2015).

Los datos se obtuvieron con un diseño completamente al azar para los rodales de la superficie de producción maderable, y se consideró una distribución normal para las categorías de diámetro. Los árboles se derribaron y seccionaron para registrar la información de crecimiento en altura dominante y ahusamiento por alturas relativas. La primera medición correspondió a la altura del tocón, posteriormente a longitudes de 0.30 m, 0.60 m y 1.3 m, además a cada 2 m se tomaron rodajas de crecimiento (Quiñonez-Barraza et al., 2015), de las que se anotaron los diámetros, alturas y número de anillos. Las combinaciones de diámetro-altura (ahusamiento) y altura-edad fueron de 768 y 685, respectivamente. En el Cuadro 1 se presentan los estadísticos descriptivos de las variables analizadas.

Cuadro 1. Estadísticos descriptivos de las variables analizadas para ajustar el ahusamiento y crecimiento en altura dominante de Pinus teocote Schiede ex Schltdl. & Cham.

D = Diámetro normal (cm); d = Diámetro a la altura comercial h (cm); H = Altura total (m); h = Altura comercial (m); hb = Altura del tocón (m); dt = Diámetro a la altura del tocón (cm); Ec = Edad a la altura comercial h (años); Et = Edad total (años); DE = Desviación estándar de la media.

Modelos utilizados

El ahusamiento se modeló con la ecuación segmentada desarrollada por Fang et al. (2000) que ha sido utilizada en diferentes estudios para generar sistemas compatibles de ahusamiento y volumen comercial en especies de interés (Quiñonez-Barraza et al., 2014; Uranga-Valencia et al., 2015; Özçelik y Crecente-Campo, 2016; Corral-Rivas et al., 2017; Tamarit et al., 2017). La ecuación de ahusamiento segmentada se anota a continuación:

               (1)

Donde:

;

 = Diámetro  del árbol  a la altura comercial  (cm)

 = Altura total del árbol  (m)

 = Altura comercial  del árbol  (m)

 = Diámetro normal del árbol  (cm)

 = Altura del tocón del árbol  (m)

 = Parámetros de volumen total (i = 1, 2, 3)

 = Parámetros del ahusamiento (i = 1, 2, 3)

 = Error del diámetro  en el árbol

Para modelar la altura dominante como una ecuación intrínseca de índice de sitio, se utilizó la ecuación dinámica en diferencia algebraica generalizada (GADA), derivada por Quiñonez-Barraza et al. (2015), la cual es dada por la ecuación 2.

                       (2)

Donde:

 = Altura  del árbol  a la edad  del árbol  en el estado 2 (m)

 = Altura  del árbol  a la edad  del árbol  en el estado 1 (m)

 = Edad  del árbol  en el estado 1 (años)

 = Edad  del árbol  en el estado 2 (años)

 = Parámetros a estimar (i = 1, 2, 3)

 = Función exponencial

 = Logaritmo natural

 = Error de la altura  en el árbol

Autocorrelación y heterocedasticidad

Combinaciones de funciones de varianza con estructuras de correlación fueron usadas en la ecuación de ahusamiento (Ec. 1) y la de altura dominante (Ec. 2). Las varFunc y corStruct fueron definidas de acuerdo con Pinheiro y Bates (2000).

Funciones de varianza. Las funciones de varianza utilizadas fueron: 1) función de potencia (varPower); 2) función exponencial (varExp); 3) función constante y de potencia (varConstPower); y 4) función combinada de potencia y exponencial (varComb). Las funciones de varianza fueron usadas para modelar la variabilidad entre las mediciones de cada árbol  con las covariables de altura comercial () para la ecuación de ahusamiento y la altura dominante () en el estado  para la ecuación de crecimiento en altura. La estructura general de las funciones de potencia para modelar la heterocedasticidad considera dos argumentos para la mayoría de las varFunc, el valor y la forma. El primero especifica el valor del parámetro de la varianza () y el segundo especifica la covariable de la varianza (), la cual en este caso puede ser considerada una variable estratificada para cada árbol, es decir se puede obtener un parámetro de varianza para cada conjunto de datos de los diferentes árboles (Pinheiro y Bates, 2000).

El modelo de varianza ( o ) para varPower es representado por la ecuación 3 y corresponde a la función () de varianza de la ecuación 4. Para la ecuación de ahusamiento (Ec. 1) se utilizó la covariable , mientras que para la ecuación de crecimiento en altura dominante (Ec. 2) .

        (3)

                    (4)

El modelo de varianza de varExp es representado por la ecuación 5 y la función correspondiente por la ecuación 6, para las mismas covariables de la ecuación de ahusamiento y crecimiento en altura dominante, definidas previamente. El modelo de varianza de varConstPower está definido en la ecuación 7 y la función de varianza en la ecuación 8. El modelo de varianza de varComb (varExp y varPrower), está definido en la ecuación 9, con la función respectiva dada en la ecuación 10 (Pinheiro y Bates, 2000). En todos los casos se utilizaron las mismas covariables definidas anteriormente. En el de varConstPower,  representa el parámetro constante y  el parámetro de potencia correspondiente y para el de varComb,  corresponde el parámetro de varExp, mientras que  a varPower.

                 (5)

                     (6)

                         (7)

                         (8)

                         (9)

                             (10)

Estructuras de correlación. Se utilizaron: 1) simetría compuesta (corCompSymm); 2) autorregresiva de orden 1 (corAR1); 3) autorregresiva continua de orden 1 (corCAR1); 4) autorregresiva de media móvil 2, 0 (corARMA2-0); 5) autorregresiva de media móvil 1, 1 (corARMA1-1); 6) autorregresiva de media móvil 2, 1 (corARMA2-1); 7) autorregresiva de media móvil 2, 2 (corARMA2-2); 8) autorregresiva de media móvil 3, 1 (corARMA3-1); y 9) autorregresiva de media móvil 3, 2 (corARMA3-2). Las estructuras de correlación se usaron para modelar la dependencia entre los residuales de cada árbol, con datos de series de tiempo (Pinheiro y Bates, 2000). En este estudio se modeló la dependencia entre las mediciones de diámetro y alturas en el mismo árbol para lograr independencia en los residuales de la ecuación de ahusamiento (Ec. 1), y la de crecimiento en altura dominante (Ec. 2). La estructura general de correlación entre grupos para un solo nivel de agrupación está expresado como la ecuación 11 (Pinheiro y Bates, 2000).

                                 (11)

Donde:

 = Residual de diámetro o altura  del árbol

 = Residual de diámetro o altura con un retardo en el árbol

 = Función de correlación , que toma valores entre  y -1

 = Distancia entre los vectores de posición  y , que para el caso de la ecuación ahusamiento fueron las diferencias entre la altura comercial  y ; mientras que, para la ecuación de crecimiento en altura dominante fueron las alturas definidas en dos posiciones  y

 = Vector de parámetros de correlación

En la estructura de correlación corCompSymm, se asume una correlación igual para todos los errores dentro del árbol pertenecientes al mismo grupo, el modelo de correlación es dado por la ecuación 12. El modelo de corAR1 está representado en la ecuación 13, mientras que el modelo carCAR1 por la ecuación 14. La estructura autorregresiva de media móvil (corARMA2-0) está definida por la ecuación 15 y la corARMA1-1 por la ecuación 16, la cual es la generalización para las estructuras corARMA p, q; es decir, corARMA2-0; corARMA1-1; corARMA2-1; corARMA2-2; corARMA3-1; y corARMA3-2; están representadas por los cambios en los valores de p para la estructura regresiva (AR) y q para la estructura de media móvil (MA). El término nuevo de error () define los residuales dependientes en el modelo de regresión; en la estructura de media móvil (MA), el residual está dado por .

                            (12)

                               (13)

                       (14)

                       (15)

                                 (16)

Ajuste de las ecuaciones y estadísticos de bondad de ajuste

Las ecuaciones se ajustaron por mínimos cuadrados no lineales generalizados (GNLS) del paquete estadístico NLME (modelos de efectos mixtos no lineales) (Pinheiro et al., 2015), en el ambiente R (R Core Team, 2017). Para evaluar el ajuste de las ecuaciones de ahusamiento y crecimiento en altura dominante, se utilizaron combinaciones de las funciones de varianza con estructuras de correlación, como las estudiadas por Pinheiro y Bates (2000). La bondad de ajuste se evaluó con seis estadísticos: raíz del cuadrado medido del error (RMSE), coeficiente de determinación ajustado (R2a), criterio de Información de Akaike (AIC), criterio de Información Bayesiano (BIC), logaritmo de la verosimilitud (logLik), coeficiente de variación (CV) y sesgo promedio absoluto (MB). Con ellos, se generó un sistema de calificación para seleccionar las mejores combinaciones varFunc y corStruct. A cada estadístico se le asignó un valor de 1 a 9; el 1 correspondió a la combinación con el mejor estadístico y el 9 para el peor (Sakici et al., 2008; Tamarit et al., 2014). Para evaluar la corrección de la autocorrelación, se utilizó el estadístico Durbin Watson (Dw) (Durbin y Watson, 1971), con una modificación robusta (DwM), como el promedio del Dw obtenido entre grupos, ya que los errores se consideran dependientes en las mediciones de cada árbol, y no en la base de datos general.

El estadístico modificado se presenta en la ecuación 17. La homogeneidad de varianzas se verificó con la prueba de Bartlett (Bartlett, 1954; Layard, 1973; Hidding et al., 2013), se consideró la variación entre los grupos y un valor de significancia de 1 % (), para la ecuación de ahusamiento y de crecimiento en altura dominante. El estadístico de prueba de homogeneidad de varianzas (hipótesis nula, H0) es dado por la ecuación 18; por tanto, valores mayores a 0.01 del valor de probabilidad de la prueba asumen homogeneidad de varianzas y de manera contraria (hipótesis alterna, Ha), los menores implican que al menos dos varianzas son diferentes.

                                              (17)

Donde:

 = Estadístico de Durbin Watson modificado

 = Residual de diámetro o altura  para el árbol

 = Número de observaciones en cada árbol

 = Número de árboles i en la base de datos

                   (18)

Donde:

 = Muestras con tamaño

 = Varianza de las muestras

 (estimación de la varianza combinada)

Resultados y Discusión

Las combinaciones de funciones de varianza con estructuras de correlación generaron 36 opciones de ajuste de ahusamiento y 36 para crecimiento en altura, basadas en las ecuaciones 1 y 2, respectivamente. Los estadísticos de ajuste y la calificación total (Ct) mostraron las bondades de ajuste de las ecuaciones (Cuadro 2 para el ahusamiento y Cuadro 3 para crecimiento en altura). También, se presentan los estadísticos DwM y la probabilidad de la prueba de homogeneidad de varianzas de Bartlett (valor-P).

El sistema de calificación mostró las combinaciones de funciones de varianza con estructuras de correlación por orden jerárquico (Sakici et al., 2008). El valor de Ct menor, fue la combinación estadísticamente mejor y a la Ct mayor le correspondió la peor, a partir de la suma de las calificaciones para cada estadístico de ajuste (Tamarit et al., 2014).

En el caso del ahusamiento, la función varPower combinada con corARMA3-2 fue la mejor y la combinación peor resultó con corCAR1. Para todos los casos de las estructuras de correlación combinadas con varPower, se obtuvieron valores del estadístico DwM entre 1.634 y 1.862; sin embargo, para la prueba de homogeneidad de varianzas, todas las combinaciones fueron heterocedásticas. En la función de varianza varExp, todas las combinaciones con las estructuras de correlación, excepto corCompSymm, fueron homocedásticas (valores de P>0.01) y la estructura corARMA1-1 resultó la mejor de acuerdo con la Ct y un valor del DwM de 1.553.

Las funciones de varConstPower y varComb combinadas con las estructuras de correlación, aunque con estadísticos confiables y residuales independientes (valores de DwM entre 1.722 y 1.973), presentan heterocedastidad en la prueba de Bartlett. Las estructuras corARMA3-1 y corARMA2-1 fueron las mejores para varConstPower y varComb, respectivamente.

El crecimiento en altura, en las combinaciones de las funciones de varianza y estructuras de correlación generaron estadísticos del DwM alrededor de 1.12; sin embargo, para la mayoría de las ecuaciones registraron varianzas homogéneas. La función varPower se combina estadísticamente de manera similar con las estructuras corARMA2-1, corARMA2-2 y corARMA3-2, ya que la Ct fue de 22 para los tres casos. La estructura corCAR1 tuvo el ajuste más bajo con varianzas desiguales. En el Cuadro 3 se puede apreciar el comportamiento de las funciones de varianza varExp, varConstPower y varComb con las combinaciones de las estructuras de correlación.

Cuadro 2. Estadísticos de ajuste de las ecuaciones de ahusamiento para las combinaciones de funciones de varianza con estructuras de correlación.

H1 = varPower; H2 = varExp; H3 = varConstPower; H4 = varComb; A1 = corCompSymm; A2 = corAR1; A3 = corCAR1; A4 = corARMA2-0; A5 = corARMA1-1; A6 = corARMA2-1; A7 = corARMA2-2; A8 = corARMA3-1; A9 = corARMA3-2; RMSE = Raíz del cuadrado medio del error; R2a = Coeficiente de determinación ajustado; AIC = Criterio de Información de Akaike; BIC = Criterio de Información Bayesiano; logLik = Logaritmo de la verosimilitud; CV= Coeficiente de variación; BM = Sesgo promedio absoluto; Ct = Calificación total del sistema; Valor-P = Probabilidad del estadístico de Bartlett.

Cuadro 3. Estadísticos de ajuste de las ecuaciones de crecimiento en altura con las combinaciones de funciones de varianza con estructuras de correlación.

H1 = varPower; H2 = varExp; H3 = varConstPower; H4 = varComb; A1 = corCompSymm; A2 = corAR1; A3 = corCAR1; A4 = corARMA2-0; A5 = corARMA1-1; A6 = corARMA2-1; A7 = corARMA2-2; A8 = corARMA3-1; A9 = corARMA3-2; RMSE = Raíz del cuadrado medio del error; R2a = Coeficiente de determinación ajustado; AIC = Criterio de Información de Akaike; BIC = Criterio de Información Bayesiano; logLik = Logaritmo de la verosimilitud; CV= Coeficiente de variación; BM = Sesgo promedio absoluto; Ct = Calificación total del sistema; Valor-P = Probabilidad del estadístico de Bartlett.

En todas las ecuaciones ajustadas, los parámetros estimados fueron estadísticamente diferentes de cero a un nivel de significancia de 5 % (P < 0.05). Los parámetros estimados para las mejores combinaciones de cada función de varianza con las estructuras de correlación se resumen en el Cuadro 4, para el ahusamiento (A) y crecimiento en altura (CA). En el ahusamiento las combinaciones de varExp&corARMA3-2, varExp&corARMA2-2 y varConstPower&corARMA3-1 representaron 14 parámetros estimados y la prueba de homogeneidad de varianzas reveló varianzas desiguales; mientras que, la combinación varExp&corARMA1-1, con 12 parámetros estimados registró varianzas constantes.

Los parámetros estimados que definen los cambios de las formas dendrométricas del modelo segmentado de ahusamiento se determinaron para el cambio del neiloide a paraboloide entre 4.0 % y 7.2 %; mientras que, el cambio de paraboloide a cono ocurre en promedio a 70 % de la altura total de los perfiles de los fustes de los árboles; resultados similares han sido documentados por investigadores para diferentes especies del género Pinus (Uranga-Valencia et al., 2015; Özçelik y Crecente-Campo, 2016; Corral-Rivas et al., 2017; Tamarit et al., 2017) con valores de entre 4 % y 7 % para el primer punto y entre 55 % y 75 % para el segundo punto de inflexión de la ecuación. Además, se utilizaron estructuras autorregresivas continuas de primer y segundo orden y funciones de potencia y se asumió una varianza conocida. No obstante, los estudios no incluyen la prueba correspondiente de homogeneidad de varianzas.

Las tendencias de los residuales por altura relativa (h/H; %) de las cuatro combinaciones mejores de funciones de varianza, con estructuras de correlación para el ahusamiento se muestran en forma de gráficos de cajas y alambres (Figura 1); aunque las tendencias son muy similares, solo la combinación varExp&corARMA1-1 (H2-A5) presentó varianza constante y se aprecia que la dispersión de los valores extremos son más cercanos a los valores de los cuantiles correspondientes.

Las combinaciones en las ecuaciones de crecimiento en altura, con más parámetros estimados (9) fueron varExp&corARMA3-2 y varConstPower&corARMA2-2; las combinaciones varExp&corARMA2-1 y varComb&corARMA2-0 tuvieron siete parámetros estimados.

El uso de ecuaciones con muchos parámetros puede ocasionar un riesgo de sobreparametrización y las predicciones no son las más eficientes (Gregoire y Schabenberger, 1996). En los cuatro casos las varianzas fueron constantes, lo que garantiza que los estimadores de los parámetros son insesgados y eficientes y con varianza mínima (Gujarati y Porter, 2011; Tang et al., 2016).

Las combinaciones seleccionadas de las ecuaciones de crecimiento en altura generaron una dispersión de residuales cercana a la línea del cero, las cuales se muestran en forma de cajas y alambres por altura relativa (h/H; %) en la Figura 2. Aunque, existe una sobrestimación entre 1 % y 10 % de la altura relativa y una subestimación 10 % a 25 %, las varianzas fueron constantes para los cuatro casos, según la prueba de homogeneidad de varianzas de Bartlett (Bartlett, 1954; Hidding et al., 2013).

Respecto a los residuales, a pesar de que la prueba DwM representó valores alrededor de 1.2, se consideraron independientes dentro de las mediciones de altura y edad en el mismo árbol (Crecente-Campo et al., 2013; Quiñonez-Barraza et al., 2015), por lo que los estimadores son insesgados, consistentes y eficientes (Williams et al., 2013).

H1 = varPower; H2 = varExp; H3 = varConstPower; H4 = varComb; A5 = corARMA1-1; A6 = corARMA2-2; A8 = corARMA3-1; A9 = corARMA3-2.

Figura 1. Gráficas de cajas y alambres para la distribución de los residuales de ahusamiento por altura relativa para las combinaciones de varFunc y corStruct.

En ecuaciones de crecimiento de altura dominante e índice de sitio con enfoques de ecuaciones en diferencia algebraica (ADA) o generalizada (GADA), se han utilizado funciones de potencia o exponenciales para la corrección de la heterocedasticidad, en las que se asumen varianzas desiguales en el proceso de ajuste por mínimos cuadrados generalizados (Castillo et al., 2013; Quiñonez-Barraza et al., 2015; Rodríguez-Carrillo et al., 2015; González et al., 2016); además de estructuras autorregresivas continuas de los errores para corregir la autocorrelación con datos de análisis troncales de diferentes especies en bosques naturales y plantaciones forestales comerciales, como Pinus durangensis Martínez, P. arizonica Engelm., Juniperus deppeana Steud, P. teocote Schiede ex Schltdl. & Cham. y P. ayacahuite Ehrenb. ex Schltdl. El enfoque desarrollado en este artículo considera la corrección de la heterocedásticidad y de la autocorrelación como combinaciones de funciones (Rodríguez et al., 2013; Williams et al., 2013).

H1 = varPower; H2 = varExp; H3 = varConstPower; H4 = varComb; A4 = corARMA2-0; A6 = corARMA2-2; A7 = corARMA2-2; A9 = corARMA3-2.

Figura 2. Gráficas de cajas y alambres para la distribución de los residuales de crecimiento en altura por altura relativa para las combinaciones de varFunc y corStruct.

Cuadro 4. Parámetros estimados de las combinaciones mejores de las funciones de varianza y estructuras de autocorrelación en las ecuaciones de ahusamiento (A) y crecimiento en altura (CA).

Est = Estadístico;  = Parámetro;  = Estimador del parámetro;  = Parámetros de las estructuras de autocorrelación continuas (corCompSymm, AR y CAR) (i=1, 2, 3);  = Parámetros de las estructuras de autocorrelación de media móvil (MA) (i=1, 2);  = Parámetros de las funciones de potencia (i= 1, 2); * = Combinaciones con varianzas constantes de acuerdo con la prueba de homogeneidad de varianzas de Bartlett. H1 = varPower; H2 = varExp; H3 = varConstPower; H4 = varComb; A4 = corARMA2-0; A5 = corARMA1-1; A6 = corARMA2-1; A7 = corARMA2-2; A8 = corARMA3-1; A9 = corARMA3-2.

En la Figura 3 se contrastan las predicciones de las ecuaciones de ahusamiento y crecimiento en altura con las combinaciones de funciones de potencia y estructuras de correlación seleccionadas, se presenta la tendencia observada del perfil de un árbol, la ecuación sin corregir heterocedasticidad y autocorrelación (SC) y las combinaciones H1-A9, H2-A5, H3-A8 y H4A6; solo la combinación H2-A5 presentó varianzas constantes. Con respecto a la ecuación de crecimiento en altura, gráficamente se presenta la comparación de tres curvas de IS, para la ecuación sin corregir (IS10, IS14 e IS18), y las ecuaciones con las combinaciones H1-A6, H2-A9, H3-A7 y H4-A4; además, se ilustran las tendencias de crecimiento de los árboles de la base de datos, en este caso todas las combinaciones evidenciaron mostraron varianzas constantes, lo cual preserva las propiedades deseables en los parámetros estimados (Beale et al., 2010; Vogelsang, 2012).

H1 = varPower; H2 = varExp; H3 = varConstPower; H4 = varComb; A4 = corARMA2-0; A5 = corARMA1-1; A6 = corARMA2-1; A7 = corARMA2-2; A8 = corARMA3-1; A9 = corARMA3-2.

Figura 3. Gráficas de predicción del ahusamiento para el perfil de un árbol, y curvas de crecimiento en altura por índice de sitio (IS 10, 14 y 18 m a la edad base de 60 años) con las combinaciones de funciones de varianza y estructuras de autocorrelación, y las ecuaciones sin estructuras de varianza ni autocorrelación (SC).

Conclusiones

Las funciones de varianza combinadas con las estructuras de correlación corrigieron los supuestos de heterocedasticidad de varianzas y autocorrelación de los errores en las ecuaciones de ahusamiento y crecimiento en altura dominante, a través del ajuste de mínimos cuadrados generalizados no lineales; con ello, se obtuvieron parámetros insesgados con varianza mínima.

Las predicciones de las ecuaciones seleccionadas para el ahusamiento son más eficientes, con estadísticos de ajuste confiables, la combinación de la función exponencial de varianza con una estructura de correlación autorregresiva de media móvil (corARMA1-1) produce varianzas constantes por categorías de alturas relativas de los perfiles del fuste. El modelo de altura dominante e índice de sitio, a una edad base de 60 años, muestra predicciones realistas. Las combinaciones seleccionadas (H1-A6, H2-A9, H3-A7 y H4-A4) presentan varianzas constantes por categorías de altura relativa a edades específicas. En ambos, los residuales son independientes y se garantizan las propiedades de las pruebas de hipótesis de los parámetros estimados.

El uso de las ecuaciones compatibles de ahusamiento y volumen comercial, así como las de crecimiento en altura e índice de sitio, está definido por la utilización de los parámetros intrínsecos de cada ecuación; por lo que, los parámetros de las funciones de varianza y estructuras de correlación solo constituyen los indicadores estadísticos para que los ajustes sean más eficientes.

Agradecimientos

Al ejido San Diego de Tezains, Santiago Papasquiaro, Durango, México por la disponibilidad para el uso de la información de ahusamiento y crecimiento en altura utilizada en este estudio.

Conflicto de intereses

Los autores declaran no tener conflicto de intereses.

Contribución por autor

Gerónimo Quiñonez-Barraza: análisis de datos, ajuste de ecuaciones y elaboración y revisión del manuscrito; Guadalupe Geraldine García-Espinoza: análisis de datos, ajuste de ecuaciones y elaboración y revisión del manuscrito; Oscar Alberto Aguirre-Calderón: análisis de datos, ajuste de ecuaciones y elaboración y revisión del manuscrito.

Referencias

Aguirre-Calderón, O. A. 2015. Manejo Forestal en el Siglo XXI. Madera y Bosques 21(21):17-28.

Arias-Rodil, M., U. Diéguez-Aranda, F. Rodríguez P., C. A. López-Sánchez, E. Canga L., A. Cámara O. and F. Castedo-Dorado. 2015. Modelling and localizing a stem taper function for Pinus radiata in Spain. Canadian Journal of Forest Research 45(6):647-658.

Bartlett, M. S. 1954. A Note on the Multiplying Factors for Various χ2 Approximations. Journal of the Royal Statistical Society. Series B (Methodological) 16(2):296-298.

Beale, C. M., J. J. Lennon, J. M. Yearsley, M. J. Brewer and D. A. Elston. 2010. Regression analysis of spatial data. Ecology Letters 13(2):246-264.

Castillo L., A., B. Vargas-Larreta, J. J. Corral R., J. A. Nájera L., F. Cruz C. y F. J. Hernández. 2013. Modelo compatible altura-índice de sitio para cuatro especies de pino en Santiago Papasquiaro, Durango. Revista Mexicana de Ciencias Forestales 4(18):89-103.

Corral-Rivas, J. J., D. J. Vega-Nieva, R. Rodríguez-Soalleiro, C. A. López-Sánchez, C. Wehenkel, B. Vargas-Larreta, J. G. Álvarez-González and A. Ruiz-González. 2017. Compatible system for predicting total and merchantable stem volume over and under bark, branch volume and whole-tree volume of Pine species. Forests 8(11):417.

Crecente-Campo, F., J. G. Álvarez-González, F. Castedo-Dorado, E. Gómez-García and U. Diéguez-Aranda. 2013. Development of crown profile models for Pinus pinaster Ait. and Pinus sylvestris L. in northwestern Spain. Forestry: An International Journal of Forest Research 86(4):481-491.

Durbin, J. and G. S. Watson. 1971. Testing for serial correlation in least squares regression. III. Biometrika 58(1):1-19.

Fang, Z., B. E. Borders and R. L. Bailey. 2000. Compatible volume-taper models for Loblolly and Slash pine based on a system with segmented-stem form factors. Forest Science 46(1):1-12.

Fierros-Mateo, R., H. M. De los Santos-Posadas, M. A. Fierros-González y F. Cruz-Cobos. 2017. Crecimiento y rendimiento maderable en plantaciones de Pinus chiapensis (Martínez) Andresen. Agrociencia 51(2):201-214.

Fortin, M., G. Daigle, C.-H. Ung, J. Bégin and L. Archambault. 2007. A variance-covariance structure to take into account repeated measurements and heteroscedasticity in growth modeling. European Journal of Forest Research 126(4):573-585.

García, E. M. 1981. Modificaciones al Sistema de Clasificación Climática de Köppen. México, D.F., México. 500 p.

González M., M., F. Cruz C., G. Quiñonez B., B. Vargas L. y J. A. Nájera L. 2016. Modelo de crecimiento en altura dominante para Pinus pseudostrobus Lindl. en el estado de Guerrero. Revista Mexicana de Ciencias Forestales 7(37):7-20.

Gregoire, T. G. and O. Schabenberger. 1996. A non-linear mixed-effects model to predict cumulative bole volume of standing trees. Journal of Applied Statistics 23(2-3):257-272.

Gujarati, D. N. y D. C. Porter. 2011. Econometria Básica. AMGH Editora. México, D.F., México. pp. 365, 371, 388 y 413.

Hidding, B., J. P. Tremblay and S. D. Côté. 2013. A large herbivore triggers alternative successional trajectories in the boreal forest. Ecology 94(12):2852-2860.

Krisnawati, H. 2016. A compatible estimation model of stem volume and taper for Acacia mangium Willd. plantations. Indonesian Journal of Forestry Research 3(1):49-64.

Layard, M. 1973. Robust large-sample tests for homogeneity of variances. Journal of the American Statistical Association 68(341):195-198.

Özçelik, R. and F. Crecente-Campo. 2016. Stem taper equations for estimating merchantable volume of Lebanon Cedar trees in the Taurus Mountains, Southern Turkey. Forest Science 62(1):78-91.

Paulo, J. A., J. H. N. Palma, A. A. Gomes, S. P. Faias, J. Tomé and M. Tomé. 2015. Predicting site index from climate and soil variables for cork oak (Quercus suber L.) stands in Portugal. New Forests 46(2):293-307.

Pinheiro, J. and D. Bates. 2000. Mixed-effects models in S and S-PLUS. Statistics and computing. Springer-Verlag. New York, NY USA. 528 p.

Pinheiro, J., D. Bates, S. DebRoy, D. Sarkar and R. C. Team. 2015. nlme: Linear and Nonlinear Mixed Effects Models. R package version 3.1–120. https://CRAN.R-project.org/package=nlme (4 de diciembre de 2017).

Quiñonez-Barraza, G., H. M. De los Santos-Posadas, J. G. Álvarez-González y A. Velázquez-Martínez. 2014. Sistema compatible de ahusamiento y volumen comercial para las principales especies de Pinus en Durango, México. Agrociencia 48(5):553-567.

Quiñonez-Barraza, G., H. M. De los Santos-Posadas, F. Cruz-Cobos, A. Velázquez-Martínez, G. Ángeles-Pérez y G. Ramírez-Valverde. 2015. Índice de sitio con polimorfismo complejo para masas forestales de Durango, México. Agrociencia 49(4):439-454.

R Core Team. 2017. R: A language and environment for statistical computing. Editorial. Vienna, Austria. 3475 p.

Rodríguez-Carrillo, A., F. Cruz-Cobos, B. Vargas-Larreta and F. J. Hernández. 2015. Compatible dominant height-site index model for juniper (Juniperus deppeana Steud.). Revista Chapingo. Serie Ciencias Forestales y del Ambiente 21(1):97-108.

Rodríguez, F., I. Lizarralde and F. Bravo. 2013. Additivity on nonlinear stem taper functions: A case for corsican Pine in Northern Spain. Forest Science 59(4):464-471.

Sakici, O. E., N. Misir, H. Yavuz and M. Misir. 2008. Stem taper functions for Abies nordmanniana subsp. bornmulleriana in Turkey. Scandinavian Journal of Forest Research 23(6):522-533.

Salas, C., T. G. Gregoire, D. J. Craven y H. Gilabert. 2016. Modelación del crecimiento de bosques: estado del arte. Bosque 37(1):3-12.

Sharma, M., N. Subedi, M. Ter-Mikaelian and J. Parton. 2015. Modeling climatic effects on stand height/site index of Plantation-Grown Jack Pine and Black Spruce trees. Forest Science 61(1):25-34.

Tamarit U., J. C., E. Rojas D., G. Quiñonez B., C. Ordoñez P. y J. C. Monárrez G. 2017. Sistema de cubicación para árboles individuales de Quercus sp. en bosques bajo manejo de Puebla, México. Revista Mexicana de Ciencias Forestales 8(40):69-88.

Tamarit U., J. C., H. M. De los Santos P., A. Aldrete, J. R. Valdez L., H. Ramírez M. y V. Guerra D. 2014. Sistema de cubicación para árboles individuales de Tectona grandis L. f. mediante funciones compatibles de ahusamiento-volumen. Revista Mexicana de Ciencias Forestales 5(21):58-74.

Tang, X., C. Pérez-Cruzado, L. Fehrmann, J. G. Álvarez-González, Y. Lu and C. Kleinn. 2016. Development of a compatible taper function and stand-level merchantable volume model for Chinese fir plantations. PloS one 11(1):e0147610.

Uranga-Valencia, L. P., H. M. De los Santos-Posadas, J. R. Valdez-Lazalde, J. López-Upton y H. Navarro-Garza. 2015. Volumen total y ahusamiento para Pinus patula Schiede ex Schltdl. et Cham. en tres condiciones de bosque. Agrociencia 49(7):787-801.

Vogelsang, T. J. 2012. Heteroskedasticity, autocorrelation, and spatial correlation robust inference in linear panel models with fixed-effects. Journal of Econometrics 166(2):303-319.

Williams, M. N., C. A. Gómez G. and D. Kurkiewicz. 2013. Assumptions of multiple regression: correcting two misconceptions. Practical Research & Evaluation 18(11):1-14.

Xu, H., Y. Sun, X. Wang, Y. Fu, Y. Dong and Y. Li. 2014. Nonlinear mixed-effects (NLME) diameter growth models for individual China-fir (Cunninghamia lanceolata) trees in Southeast China. PloS one 9(8):e104012.

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.