La investigación se desarrolló en los municipios de Cuatro Ciénegas, Parras de la Fuente, Ramos Arizpe, San Pedro de las Colonias, Sierra Mojada y Viesca, correspondientes a la parte centro y sur del Estado de Coahuila (Fig. 1); abarcando 10 poblaciones naturales localizadas en predios con programa de manejo forestal vigente de candelilla. Estas localidades se caracterizan por presentar climas del tipo árido y semiárido (Bs) y muy árido (Bw), con temperaturas medias anuales de 18 °C a 22 °C y precipitación anual de 100 mm a 400 mm (Instituto Nacional de Estadística y Geografía [Inegi], 2007). La altitud está comprendida entre 950 m y 1500 m. La vegetación correspondiente al área es matorral desértico micrófilo y rosetófilo (Inegi, 2016).

Figura 1
Distribución geográfica de las poblaciones de Euphorbia antisyphilitica Zucc., bajo estudio en Coahuila, México.

Los datos se colectaron mediante un muestreo selectivo de forma destructiva, donde a cada individuo se le midió altura media (A, cm), diámetro menor y mayor de la base (cm) y diámetro menor y mayor de la cobertura de copa (cm). Posteriormente la planta fue cosechada para obtener la biomasa total verde (Btv) (en kilogramos por individuo “kg individuo^-1”), con una báscula Torrey de 0.005 kg de resolución (Winck et al., 2015; Hernández-Ramos et al., 2017; Flores et al., 2018; Villavicencio-Gutiérrez, Hernández-Ramos, Aguilar-González y García-Cuevas, 2018). La muestra se conformó de 1083 individuos distribuidos de forma que abarcaran la mayor variabilidad de condiciones de crecimiento de la especie y las categorías diamétricas existentes en las diez poblaciones muestreadas (Tabla 1).

Tabla 1
Distribución de la muestra de las poblaciones de Euphorbia antisyphilitica Zucc., bajo estudio en Coahuila, México.

Agrupación de datos. La agrupación de los municipios con atributos similares se realizó mediante un análisis de varianza (Andeva), donde el factor primario fue el municipio y la variable de respuesta la Btv. Posterior al análisis de significancia entre los seis municipios, se realizó el contraste de medias mediante la técnica de la diferencia mínima significativa o método LSD (least significant difference), que consiste en la prueba de hipótesis por parejas de la forma H₀: µ_i = µ_j; H₁: µ_i ≠ µ_j. Por lo tanto, las medias fueron estadísticamente diferentes si${(\stackrel{-}{\mathrm{y}}}_i-{\stackrel{-}{\mathrm{y}}}_j)>LSD$, donde la mínima diferencia significativa (LSD) está dada por:$LSD={\tau }_{\alpha /2;}N-1\sqrt{{{\widehat{s}}_R}^2(\frac{1}{n_i}+\frac{1}{n_j})}$, siendo n_i y n_j el número de observaciones correspondientes a cada medida, N - 1 el número de grados de libertad y t_α/2;N-1 el valor crítico de la distribución t con N - 1 grados de libertad, que deja una probabilidad a su derecha igual a α/2 (Williams y Abdi, 2010). Finalmente, se realizó la agrupación de medias por grupo de municipios (Statistical Analysis Systems Institute Inc., [SAS], 2017).

Selección de variables para el ajuste de modelos. El diámetro promedio de la cobertura de copa (DPC, cm) y el diámetro promedio de la base (DPB, cm) resultaron de la media de las medidas de la planta tomadas en campo para cada variable; la cobertura arbustiva (CA, m²) fue calculada con la expresión$CA=\frac{\pi }{40000}*{\left(\raisebox{1ex}{$d_1+d_2$}\!\left/ \!\raisebox{-1ex}{$2$}\right.\right)}^2$ , donde d₁ y d₂ son los diámetros de cobertura medidos a la planta en centímetros. Las variables antes generadas y la A fueron relacionadas con la Btv mediante un análisis de correlación de Pearson, seleccionando para el ajuste de los modelos únicamente aquellas que presentaron los valores de correlación (r) más altos y las que aportaron mayor variabilidad en el modelo lineal ajustado mediante la técnica Stepwise (Volke, 2008). Además, se graficaron las variables independientes contra la Btv para observar su tendencia.

Depuración de bases de datos. Para la detección de los valores atípicos de la muestra, previo al ajuste de los modelos, se empleó el método de eliminación directa mediante el análisis gráfico de la dispersión de los datos y con los valores de los residuales studentizados; estos últimos se calcularon en SAS 9.4 ® con PROC GLM (Modelos lineales generalizado por sus siglas en ingles) y la opción r-Influence (SAS, 2017). Se consideraron valores atípicos los residuales mayores al valor critico de ± 2, con el supuesto de que al menos 5% de la muestra no pertenecen a la población (Coenders y Saez, 2000) y se eliminaron estos valores de la base de datos para evitar tener resultados sesgados en los análisis de regresión (Chen, 1978).

Ajuste estadístico de los modelos. Se emplearon tres modelos seleccionados de la literatura, que presentaron resultados sobresalientes al estimar la biomasa de diversas especies arbóreas y arbustivas (Gómez-García, Crecente-Campo y Diéguez-Aranda, 2013; Soriano-Luna et al., 2015; Hernández-Ramos et al., 2017; Villavicencio-Gutiérrez et al., 2018) (Tabla 2). Una vez seleccionadas las variables más relacionadas con la Btv, se realizó el ajuste de los modelos, mediante el método de mínimos cuadrados ordinarios con PROC MODEL en SAS 9.3 ® (SAS, 2017) para cada uno de los grupos de municipios resultantes de la agrupación.

Tabla 2
Modelos seleccionados para el ajuste de la biomasa verde de Euphorbia antisyphilitica Zucc., en seis municipios de Coahuila.

La selección de los modelos de mejor ajuste al conjunto de datos de Btv se efectuó de acuerdo con los mayores valores del coeficiente de determinación ajustado por el número de parámetros (R²_Ajus), el cual mide la variabilidad explicada por el modelo, y los menores de la raíz del cuadrado medio del error (RCME), que proporciona la media de las diferencias entre los valores estimados y los observados (Tlaxcala-Méndez, Santos-Posadas, Hernández-de la Rosa y López-Ayala, 2016), además de la significancia de sus parámetros y del coeficiente de determinación entre los valores predichos y observados (Tlaxcala-Méndez et al., 2016; Flores et al., 2018). Los supuestos de regresión del modelo se analizaron con la prueba de Shapiro-Wilk para la normalidad de los datos, la prueba de White para detectar la homocedastisidad de los residuales y la prueba de Durbin-Watson para detectar autocorrelación de los errores. Por último, a los modelos se les evaluó la capacidad predictiva mediante el sesgo promedio absoluto (Ē) en kilogramos de Btv y la diferencia agregada en porcentaje (DA%) (Trincado y Leal, 2006; Tlaxcala-Méndez et al., 2016; Hernández-Ramos et al., 2017). El análisis estadístico se realizó en el programa SAS 9.3 ®.

Donde:$Y_i$ es el valor observado,${\widehat{Y}}_i$ es el valor predicho con el modelo,${\stackrel{-}{Y}}_i$ es la media de la variable dependiente y$n$ es el número de observaciones.

El análisis de varianza indica diferencias significativas entre los municipios y entre las poblaciones, el error Tipo III en ambos es significativo (α = 0.0001). Mientras que los resultados obtenidos con LSD muestran una agrupación de los municipios en tres grupos de acuerdo con la Btv por individuo. El municipio de Parras de la Fuente (El Amparo) formó el primer grupo (GI) con la mayor media registrada igual a 1.722 kg individuo^-1, el segundo grupo (GII) lo conformaron las poblaciones de los municipios de San Pedro (Mala Noche), Cuatro Ciénegas (Estanque de León, Las Morenas y Tanque Nuevo), Ramos Arizpe (Las Coloradas, Plan de Guadalupe y Cosme) y Viesca (La Fe), con medias similares (1.046, 0.924, 0.917 y 0.809) kg individuo^-1, respectivamente; por último, el grupo tres (GIII) se conformó por la población colectada en el municipio de Sierra Mojada con 0.657 kg individuo^-1, siendo este valor de Btv el más bajo de las poblaciones evaluadas (Fig. 2).

Figura 2
Biomasa verde promedio por municipio (arriba) y por población (abajo) de Euphorbia antisyphilitica Zucc.

El municipio de Parras de la Fuente presentó valores altos de biomasa (GI) y Sierra Mojada con la población de San José de Madero presentó los valores más bajos (GIII), comparados con el resto del listado, por ello, se clasificó en un grupo específico a cada una de estas poblaciones con el fin de disminuir la variabilidad de estimación al momento de ajustar los modelos de regresión. En cuanto al GII, se integró con cuatro municipios con condiciones de biomasa similares. Las poblaciones que presentaron los mayores valores de biomasa fueron El Amparo y Las Coloradas con 1.722 kg individuo^-1 y 1.507 kg individuo^-1, respectivamente. Las poblaciones menos productivas fueron San José de Madero, Estanque de León y Plan de Guadalupe con 0.657 kg individuo^-1, 0.526 kg individuo^-1 y 0.420 kg individuo^-1, para cada una; el resto de las poblaciones oscilaron entre 1.178 kg individuo^-1 a 0.809 kg individuo^-1.

Las variables correlacionadas con la Btv son la CA (0.90), la A por DPC (0.88), el DPC (0.87) y el DPB (0.78), considerando una relación de tipo lineal. Sin embargo, la dispersión de los datos es de tipo sigmoidal por lo que se analizaron con una regresión Stepwise, indicando que las variables A, A por DPB (ADPB), DPC, ADPC y cobertura arbustiva son significativas (0.05) dentro del modelo. Las dos pruebas coincidieron en que las variables A, DPC, DPB y la CA son adecuadas para estimar la Btv (Fig. 3), con ellas se ajustaron los modelos seleccionados, además se integraron las iteraciones entre las variables.

Figura 3
Relación entre la biomasa verde (kg) y las variables explicativas con mayor influencia en el modelo.

La representación gráfica de los valores observados de Btv contra las variables independientes A, DPC y DPB indican una tendencia exponencial que es representada con el tipo de modelos no lineales propuestos, caso contrario para la cobertura, que tiene una relación lineal con la Btv, razón por la que no se consideró en el ajuste de los modelos, además de que es una variable calculada con el DPC, con el cual esta correlacionado.

El ajuste de los modelos se realizó para cada grupo (GI, GII y GIII). Los modelos seleccionados presentaron una R²_Ajus superior a 0.86, en los primeros dos grupos y, para GIII, fue de 0.64. La RCME fue inferior a 0.425 y los parámetros de los modelos fueron todos significativos (p<0.05). Los modelos seleccionados para estimar la Btv de cada grupo quedaron de la siguiente forma: para GI el modelo 2 con las variables DPB y DPC (Tabla 3); para GII los modelos con mejor ajuste y factibles para estimar la Btv de candelilla fueron 1 y 2 con las variables DPC, DPB y A (Tabla 4); finalmente, para el GIII fue el modelo 7, que emplea el DPC como variable independiente; aunque en este grupo los modelos restantes, a excepción del 6, presentaron valores de los estadísticos superiores al modelo 7, estos se descartaron por no mostrar significancia en uno de sus parámetros (Tabla 5).

Tabla 3
Estadísticos de ajuste de los modelos analizados para los datos del conjunto GI, para estimar la biomasa verde de Euphorbia antisyphilitica Zucc.

Tabla 4
Estadísticos de ajuste de los modelos analizados para los datos del conjunto GII, para estimar la biomasa verde de Euphorbia antisyphilitica Zucc.

Tabla 5
Estadísticos de ajuste de los modelos analizados para los datos del conjunto GIII, para estimar la biomasa verde de Euphorbia antisyphilitica Zucc.

En los cuatro modelos seleccionados se verificaron los supuestos de regresión, con la prueba de Shapiro-Wilk se confirmó que los datos son tendentes a la distribución normal, ya que el valor del estadístico fue W = 0.87 (Pr < W = 0.0001, para los modelos seleccionados), y la prueba de White comprobó la homogeneidad de varianzas. La inexistencia de autocorrelación de los errores se observó con la prueba de Durbin-Watson, al obtener valores de 2.17, 1.64 y 1.81 para GI, GII y GIII, respectivamente. Además, el análisis gráfico de los valores predichos con los cuatro modelos y de los observados muestran una tendencia recta, con una correlación lineal entre estos de 0.923 para GI, 0.862 y 0.881 de los modelos seleccionados 1 y 2 de GII y 0.653 para GIII, esto confirma un ajuste adecuado de los modelos seleccionados para estimar la Btv de candelilla en la región muestreada (Fig. 4).

Figura 4
Valores de biomasa verde (Btv) estimada contra la observada para las regiones GI (a), GII (b) y GIII (c), empleando los modelos seleccionados.

La capacidad predictiva de los modelos seleccionados mostró un Ē y DA% de 0.0197 kg (1.144%) para GI, 0.0134 kg (1.437%) y 0.018 kg (1.927%) en GII y para GIII fue de -0.01684 kg (-2.566%) individuo^-1. En los primeros dos grupos los modelos presentan una sobreestimación de la Btv no mayor a 2% y para el grupo restante fue ligeramente mayor por la variabilidad que existe en las poblaciones evaluadas.

^*Autor de correspondencia. cele64@gmail.com