Introducción
Los recursos naturales en México además de admirados son conocidos por mostrar una gran riqueza vegetal que puede ser explicada por la combinación de los factores clima, elevación, orografía y edafología, entre otros. Los científicos en el pasado mencionaron que la sinergia dada por la fisonomía y composición durante millones de años provocaron los tipos de vegetación presentes en el país, de ellos, el de mayor proporción es el matorral xerófilo que ocupa aproximadamente 40% de la superficie total. Este representa las zonas áridas y semiáridas nacionales donde se desarrollan actividades silvoagropecuarias como la recolección de subsistencia realizada en condiciones de gran variación climática por campesinos que habitan en esos lugares (Velazco 1991; Rzedowski, 2006; Villarreal, 2006; Martínez, 2013).
La actividad campesina de cosecha de subsistencia de productos vegetales de clima árido y semiárido se distribuye en el Altiplano Mexicano, en Querétaro, Guanajuato, Aguascalientes, Zacatecas, San Luis Potosí, Durango, Chihuahua, Coahuila, Nuevo León, Sonora y la península de Baja California. En menor escala en los estados de Oaxaca, Puebla, Hidalgo, Estado de México y Tamaulipas, sin que sea menos importante para el beneficio de las personas que la realizan (Beltrán et al., 1964; Rzedowski, 2006; Villarreal, 2006; Martínez, 2013). Dicho aprovechamiento se concentra en vegetación silvestre como Euphorbia antisyphilitica Zucc. (candelilla), Agave lechuguilla Torr. (lechuguilla), Lippia spp. (orégano), Opuntia spp. (nopal), Nolina spp. (pamilla), Yucca carnerosana Trel. (yucca) y Dasylirion spp. (sotol), entre otras (Martínez, 2013).
Las actividades que realiza el hombre en estas superficies han producido cambios en el paisaje, lo cual ha impactado en mayor medida a la vegetación y el suelo, un ejemplo se ha dado en E. antisyphilitica, la cual ocupa el tercer lugar en aprovechamiento a nivel nacional con 14 millones de hectáreas aproximadamente (Beltrán et al., 1964; Norma Oficial Mexicana, 2003; Rzedowski, 2006; Canales et al., 2006; CITES, 2009). La cosecha de E. antisyphilitica se hace en forma manual, se toma al individuo de la parte reproductiva llamada macollo y se le arranca del suelo con todo y raíz, si hay remanentes, se recuperan lentamente en períodos mayores a cinco años. Posteriormente, por un proceso de extracción en campo, se obtiene una cera altamente valorada por la industria de los cosméticos (jabones, lápices labiales, cremas corporales, preparaciones para el cabello y protecciones odontológicas, entre otros).
Este hecho permite que hay una metodología alternativa que evita la destrucción de individuos, en ella se utilizan modelos estadísticos que estiman el peso verde o biomasa de la planta a partir de la cuantificación de variables morfométricas fáciles de identificar y rápidas de medir en campo, en conjunto la literatura designa a dicho procedimiento, muestreo, medición y aplicación de métodos estadísticos como alometría. Su utilidad puede aplicarse para la obtención de tablas de rendimiento de la especie, acción deseable para inventarios, caracterización y conservación de E. antisyphilitica y en general de otras especies recolectadas por campesinos (Beltrán et al., 1964; Norma Oficial Mexicana, 2003; Cano et al., 2005; Rzedowski, 2006; Tapia y Reyes, 2008; CITES, 2009).
Por tal motivo, el objetivo del presente estudio fue seleccionar un modelo lineal, alométrico o geométrico que estime el peso verde o biomasa de la planta de E. antisyphilitica para elaborar una tabla de rendimiento en base a variables morfométricas como diámetro aéreo mayor y menor, altura de la planta y cobertura aérea de la planta en poblaciones silvestres al norte de Zacatecas.
Materiales y métodos
Área de estudio
Las poblaciones silvestres de E. antisyphilitica objeto de esta investigación se encuentran en el Ejido El Rodeo, municipio de Mazapil, al norte de Zacatecas. La localización geográfica se ubicó en las coordenadas 102° 11’ 00” longitud oeste y 24° 54’ 00” latitud norte a una elevación de 1 990 m (CONABIO, 2014), la superficie total del predio fue de 10 317.553 ha (Figura 1), esta área forma parte del Desierto Chihuahuense y gran proporción del predio está cubierto por matorral desértico rosetófilo. En el lugar se observaron especies como A. lechuguilla (lechuguilla), D. cedrosanum Trel. (sotol), Hechtia glomerata Zucc. (guapilla) y Agave stricta Salm-Dyck (espadín), entre otras, en ocasiones se mezclaron con matorral inerme parvifolio asociado con Larrea tridentata Cov. (gobernadora) y Y. carnerosana Trel. (Villa, 1981).
Muestreo de plantas con características de madurez de cosecha de E. antisyphilitica
Se localizaron puntos al azar en el área de estudio distanciados a 100 m uno de otro, así se ubicaron 48 sitios en el terreno. En cada uno se trazaron cuatro cuadrantes tomando la dirección norte-sur como eje principal y la este-oeste como secundaria. A partir del centro, se tomó la planta más cercana con características de madurez, tales como, diámetro aéreo mayor del individuo mayor a 0.25 m y altura mínima de 0.3 m (Norma Oficial Mexicana, 2003; Velasco et al., 2009).
Se midieron a cuatro individuos por sitio para tener representatividad en el muestreo y poder evaluar la variación en conjunto, se incluyeron plantas con medidas menores y mayores a las anteriores y que en forma práctica se sabe son cosechadas. Para ser congruentes con la normatividad vigente se contó con el aviso de aprovechamiento y el estudio técnico respectivo mencionando las existencias reales de la planta de candelilla en el lugar de muestreo (Zar, 1999; Norma Oficial Mexicana, 2003). En total se midieron n= 192 plantas para lograr representatividad estadística y cumplir con el objetivo del presente estudio.
Variables morfométricas estudiadas
Se midió en forma individual la h= altura de la planta (m); dma= diámetro aéreo mayor (m); y dme= diámetro aéreo menor (m), especificadas en la NOM-018-SEMARNAT-1999 (NOM-018, 2003). Para definir los diámetros aéreos de cada planta, se consideró una figura geométrica en forma de elipse, de la cual el dma (m) se obtuvo sobreponiendo un flexómetro en la parte más ancha de la cobertura aérea en cuestión, lo que representó el eje mayor de la elipse. El dme (m) se midió sobreponiendo el flexómetro en la parte más angosta de la cobertura aérea de la planta, lo que representó el eje menor.
Posteriormente, se obtuvo el semieje mayor (dma/2) y el semieje menor (dme/2); la coba= cobertura aérea (m2) de cada planta fue el resultado de veces el producto de los dos semiejes. La variable pv= peso verde (kg) individual, se obtuvo posterior a la medición de las variables anteriores. Cada planta se extrajo en su totalidad (incluyendo la raíz), se sacudió cada muestra para disminuir las impurezas (en su mayoría tierra), se ataron con cintas y se identificaron con una etiqueta cada una, posteriormente, se trasladaron a un centro de almacenamiento donde se pesaron con una báscula digital (Canales et al., 2006; CITES, 2009; Martínez, 2013).
Análisis estadístico de los datos
La información numérica del presente trabajo fue ordenada en una base de datos de doble entrada en Excel conformada por las variables, h (m); dma (m); dme (m); coba (m2); y pv (kg), esta última representó la biomasa de cada planta incluyendo raíz. Una vez conformada la información digital en una matriz de filas y columnas, se exportó al software R CRAN, se realizó un análisis de exploración de datos a las variables obteniendo los estimadores de la media muestral, desviación muestral, valores mínimos y máximos, intervalo ͞x ± s (donde: ͞x= media muestral y s= desviación muestral). Se aplicó la prueba de Shapiro-Wilk (W) para normalidad en cada una de las variables (Everitt y Hothorn, 2014).
Para conocer la relación lineal entre los pares de variables posibles de la matriz, se estimó el coeficiente de correlación simple de Pearson= r, se obtuvo la significación estadística focalizando en valores de p< 0.01. Debido a que el interés principal fue la estimación del pv (kg) de E. antisyphilitica, se atendieron las asociaciones de esta con respecto a las demás. De esta forma se promediaron todos los datos en diez rangos (0.1-1 m) de la variable con mayor correlación con pv (kg) y menor valor de P (Walpole et al., 2007) en las variables de la base de datos.
Se ajustaron modelos estadísticos (Cuadro 1) para estimar pv (kg) con las variables morfométricas más significativas (p< 0.01) de la base de datos. Los criterios de selección fueron, mayor valor de coeficiente de determinación ajustado (R2 aj), valor mínimo de cuadrado medio del error (CME), valor mínimo porcentual de coeficiente de variación (CV), menor valor de probabilidad de cometer el error tipo I (P) en el modelo de regresión y menor valor de la información de Akaike (AIC) (Crawley, 2007; Walpole et al., 2007; Anderson, 2008; The R Core Team, 2017). Después se construyeron bandas de confianza con el modelo (p= 0.05) de tal forma que los datos inmiscuidos de pv (kg) estuvieran incluidos a un 95% de confiabilidad de la estimación.
No. | Ecuación |
1 | pv= a +b (coba) |
2 | pv= a +b (dma) |
3 | pv= a +b (dme) |
4 | pv= abcoba |
5 | pv= abdma |
6 | pv= abdme |
7 | pv= a(coba)b |
8 | pv= a(dma)b |
9 | pv= a(dme)b |
pv = peso (kg); a y b = estimadores de los parámetros de la regresión; h = altura de la planta (m); dma = diámetro aéreo mayor (m); dme= diámetro aéreo menor (m); coba= cobertura aérea (m2).
Tabla de rendimiento de pv (kg) para E. antisyphilitica
Una vez seleccionado el modelo, se elaboró la tabla de rendimiento a partir de rangos (0.1-1 m) de la variable morfométrica con mayor significación estadística de la correlación lineal entre pares de variables. Esto apoyó la separación de individuos en madurez de cosecha a través de un criterio cuantitativo respaldado en la normatividad vigente (Norma Oficial Mexicana, 2003). La estimación del rendimiento tabular se basó en los valores incluidos del espacio de las bandas de confianza (95%) del modelo seleccionado (Walpole et al., 2007; Anderson, 2008; The R Core Team, 2017).
Resultados y discusión
Análisis estadístico y correlación entre variables morfométricas de E. antisyphilitica
El análisis de exploración de datos de las n= 192 plantas, mostraron los estimadores reportados en el Cuadro 2. La correlación entre las variables de estudio (Cuadro 3) exhibió la relación lineal y la significación estadística (p< 0.01). Todos los pares de variables fueron significativos (p< 0.01), los valores más altos de los coeficientes de correlación fueron, entre coba (m2) y dme (m) (r = 0.90718) y coba (m2) y dma (m) (r= 0.88426). Sin embargo, el pv (kg) también los tuvo (p< 0.01) con el resto de las variables morfométricas evaluadas, aunque menores que las primeras mencionadas. Debe destacarse que las correlaciones de interés están en la última fila del (Cuadro 3), las cuales consideran el pv (kg) con respecto a la h (m); dma (m); dme (m); y coba (m2), respectivamente.
h (m) | dma (m) | dme (m) | pv (kg) | |
Intervalo | 0.578 ±0.156 | 0.394 ±0.122 | 0.335 ±0.124 | 0.923 ±0.594 |
Mínimo | 0.2 | 0.2 | 0.15 | 0.15 |
Máximo | 0.98 | 1 | 0.9 | 3.89 |
W< P | 0.1086 | 1.079 E-08 | 2.013 E-09 | 8.071 E-15 |
h= altura de la candelilla; dma= diámetro aéreo mayor; dme= diámetro aéreo menor; W= Estadístico Shapiro-Wilk.
h | dma | dme | coba | pv | |
h | 1 | ||||
dma | 0.36805** | 1 | |||
dme | 0.26291** | 0.71582** | 1 | ||
coba | 0.32912** | 0.88426** | 0.90718** | 1 | |
pv | 0.36573** | 0.63625** | 0.54529** | 0.63599** | 1 |
Así que, al tomar r= 0.63625 de la relación lineal entre pv (kg) y dma (m) como valor de correlación más alto, fue posible la formación de intervalos con una asignación de 0.1-1 m. Posteriormente, se obtuvieron los gráficos de dispersión de la Figura 2 para dos variables inmiscuidas. Puede notarse que al elegir dma (m) para la reorganización de la información (Figura 2), mejoró la asociación de h (m); r= 0.8626 (antes r= 0.36573); dme (m); r= 0.9861 (antes r= 0.54529); la de dma (m); r= 0.9761 (antes r= 0.63625) y coba (m2); r= 0.9650 (antes r= 0.63599) con respecto a pv (kg).
Lo anterior permitió encontrar mejores ajustes para morfometrías similares en la especie estudiada, razón soportada y recomendada por la teoría básica de la regresión (Zar, 1999; Walpole et al., 2007) la que coincidió con trabajos en otras especies de interés forestal; por ejemplo, leñosas y leguminosas (Ares et al., 2002; Iglesias y Barchuk, 2010). Además, incluir un análisis de correlación, proporcionó los estimadores cuantitativos y gráficos de las asociaciones lineales (Cuadro 3 y Figura 2) que mostró contrastes con otros estudios similares en la carencia de una metodología reproducible con el evidente interés de realizar estimaciones con modelos designados como alométricos, al revisar la literatura se refiere a curvas exponenciales o geométricas que se linealizan con la aplicación de leyes de logaritmos y se basan en la aplicación del método de mínimos cuadrados donde los coeficientes son lineales, tienen el respaldo estadístico y no siempre es sencillo aplicarlos (Ares et al., 2002; Walpole et al., 2007; Fonseca et al., 2013).
En los estudios de Ares et al. (2002); Segura y Andrade (2008) se aplicó una metodología no destructiva en los individuos vegetales para estimar biomasa con modelos alométricos, el presente estudio también lo hizo, además de mostrar un orden a través de una metodología reproducible y dirigida a la selección del mejor modelo bajo criterios estadísticos definidos y que coincidieron con otras investigaciones (Segura y Andrade 2008; Iglesias y Barchuk, 2010; Fonseca et al., 2013).
Modelos ajustados y estimación de la tabla de rendimiento de pv (kg) de E. antisyphilitica
La estimación a través de modelos para las variables morfométricas organizados en intervalos de 0.1-1 m y los criterios de selección se concentran en el Cuadro 4. El ajuste de la regresión de cada uno de los modelos mostró que los tres primeros que representan líneas son los que tuvieron mejores estimadores estadísticos, los siguieron los geométricos y finalmente los alométricos, exhibieron valores menores de R2 aj y mayores en el AIC.
Modelo | R2 aj | CME | CV (%) | P | AIC |
1 | 0.9226 | 0.0926 | 24.16 | 6.28E-06 | 8.3657 |
2 | 0.9468 | 0.0636 | 20.03 | 1.39E-06 | 4.6063 |
3 | 0.9689 | 0.0371 | 15.29 | 1.59E-07 | 0.7909 |
4 | 0.7515 | 0.0873 | 136.39 | 0.000714 | 57.2297 |
5 | 0.6236 | 0.1324 | 167.95 | 0.004011 | 51.4759 |
6 | 0.5973 | 0.1416 | 173.7 | 0.005322 | 51.1438 |
7 | 0.2908 | 0.2495 | 230.52 | 0.062194 | 21.1029 |
8 | 0.4198 | 0.2041 | 208.51 | 0.025409 | 4.2571 |
9 | 0.3019 | 0.2456 | 228.72 | 0.057898 | 6.9131 |
pv= peso (kg); a y b= estimadores de los parámetros de la regresión; h= altura de la planta (m); dma= diámetro aéreo mayor (m); dme= diámetro aéreo menor (m); coba= cobertura aérea (m2).
Se eligió al modelo tres, pv= -0.1882 + 3.7831 (dme) (Cuadro 4) para la estimación de la tabla de rendimiento de candelilla, en el Cuadro 5 se muestra la estimación a través de intervalos de 0.1-1 m de dme (m). Al utilizar la ecuación lineal tres, todas las observaciones quedaron incluidas en el espacio de estimación de las bandas de confianza al 95% (Figura 3). Si se toma como base los criterios de selección de modelos para las n= 192 plantas evaluadas, la estimación de la tabla con el modelo lineal tres; pv= -0.1882 + 3.7831 (dme), en forma práctica solo utiliza una variable morfométrica medible en campo; dme (m), de esta forma se obtendrá pv (kg) incluyendo la raíz, lo cual representó una ventaja, debido a que ya no sería necesario extraer individuos y en consecuencia se evitaría su destrucción, resultado similar al de Iglesias y Barchuk (2010) en leguminosas.
Los criterios estadísticos, indicaron que el modelo en promedio estimó el valor de la biomasa (kg) para valores de la variable independiente dme (m) que sean sustituidos en la ecuación, dicha predicción cae dentro de las bandas de confianza al 95% del pv (kg). Resultados divergentes fueron obtenidos por Fonseca et al. (2013), quizás por la inclusión de modelos para estimar biomasa y carbono con criterios para coeficientes lineales en curvas exponenciales sin mostrar exploración de datos y análisis de correlación, dando mayor importancia a la distribución normal de los modelos utilizados, cuando teóricamente es incorrecto. De acuerdo con Zar (1999); Crawley (2007); Walpole et al. (2007); Everitt y Hothorn (2014); The R Core Team (2017) la normalidad es posible para las variables continúas evaluadas o el error aleatorio de un modelo lineal y, en su caso también para uno linealizado con ayuda de logaritmos, esta fue la razón por la cual se aseguró un tamaño de n= 192 plantas para cumplir con una confiabilidad de 95%.
dme (m) | pv (kg) |
0.1 | 0.19006 |
0.2 | 0.56838 |
0.3 | 0.9467 |
0.4 | 1.32501 |
0.5 | 1.70333 |
0.6 | 2.08165 |
0.7 | 2.45997 |
0.8 | 2.83828 |
0.9 | 3.2166 |
1 | 3.59492 |
Si se pone atención en el modelo lineal uno, pv= a +b (coba), tuvo estimadores estadísticos adecuados, sin embargo, la variable coba (m2) es el resultado de la medición de dos diámetros aéreos (el mayor y el menor) y posteriormente, la obtención de los semiejes mayor y menor y el producto de ellos con. Medir a coba (m2) no proporciona facilidad en los cálculos para la obtención del modelo respectivo de una tabla de rendimiento, razón por la cual no se eligió dicha ecuación, además, se alejó de lo encontrado por Iglesias y Barchuk (2010) sobre la sencillez de la estimación, abonando a que son dos mediciones y un cálculo intermedio.
Con respecto al modelo lineal dos, pv = a + b (dma), sus estimadores son adecuados; sin embargo, el ajuste que tuvo el modelo elegido fue mejor en forma operativa y por lo citado en la NOM-018-SEMARNAT-1999 (NOM-018, 2003), deberán medirse en campo tanto dma (m) como dme (m); los análisis estadísticos del presente estudio mostraron que la primera da orden en rangos de 0.1-1 m a todas las variables, siendo la segunda variable la que mejoró el ajuste del modelo lineal para la estimación de la tabla de rendimiento. Los resultados mostraron simpleza en la elección de un modelo, similar a lo realizado por Iglesias y Barchuk (2010) quienes hicieron la misma referencia en leguminosas.
Los modelos alométricos cuatro, cinco y seis y los geométricos siete, ocho y nueve, tuvieron estimadores que se alejaron de los criterios estadísticos, además, para su estimación y ajuste se requiere de mayor conocimiento en análisis de regresión, si a esto se adiciona obtener la variable coba (m2), dificultará la obtención de este tipo de ecuaciones. En el presente estudio dichos modelos fueron superados por los lineales y a pesar de que los alométricos son preferidos (Iglesias y Barchuk, 2010), no fue el caso para los datos morfométricos de E. antisyphilitica en el área de estudio. Por otro lado, diferimos de Fonseca et al. (2013) ya que sus resultados van encaminados a la misma variable con una clara tendencia al uso de modelos exponenciales del tipo geométrico.
Los llamados ajustes de modelos para la estimación de biomasa se han centrado en curvas del tipo alométricas y geométricas, no en todos los casos incluyen a los lineales (Fonseca et al., 2013, Iglesias y Barchuk 2010); sin embargo, los criterios estadísticos utilizados recaen en estimadores lineales (Zar, 1999; Walpole et al., 2007; Everitt y Hothorn, 2014; The R Core Team, 2017). Esa preferencia va acompañada de la llamada linealización que es la aplicación de las leyes de los logaritmos a un modelo que así lo requiera, distintas investigaciones evitan la correlación entre variables inmiscuidas como parte de un análisis anterior y exploratorio al fenómeno estudiado (Zar, 1999; Walpole et al., 2007; Fonseca et al., 2013).
En esta investigación se encontró que una sola variable estimó el peso pv (kg), en este caso dme (m) con el modelo tres de regresión lineal simple, debido a que entre más sencillo resulte en estructura estadística y en la conformación de la misma ecuación, la sustitución de valores será directa y rápida, tal y como lo corroboró Iglesias y Barchuk (2010) con respaldo estadístico, Zar (1999), Crawley (2007); Walpole et al. (2007); Anderson (2008); Everitt y Hothorn (2014); The R Core Team (2017), entre otros. Otro beneficio fue la construcción de una tabla con dos columnas, la primera compuesta con valores de la variable dme en rangos de 0.1-1 m; la segunda la formó la variable estimada pv (kg), así el modelo, pv= -0.1882 + 3.7831 (dme), representó la biomasa de la planta, lo cual cumplió con el objetivo planteado en este trabajo.
La coincidencia de la normatividad Norma Oficial Mexicana (2003) con la teoría de correlación mostrada en Zar (1999); Everitt y Hothorn (2014); Crawley (2007); Walpole et al. (2007); The R Core Team (2017), se reflejó en campo ya que dio la posibilidad de dejar distribuidas al 20% de las plantas en el área de explotación en madurez de cosecha (dma> 0.25 m, h≥ 0.3 m), representando no intervención en esos individuos, asegurando la recuperación de E. antisyphilitica para el futuro y por otro lado, en el sentido administrativo y en cumplimiento de la normatividad vigente, se cuidó que el estudio técnico tuviera el aviso de aprovechamiento reportando las existencias reales de la planta (Norma Oficial Mexicana, 2003).
Conclusiones
Fue posible estimar pv (kg) o biomasa de la planta de E. antisyphilitica a partir del dme (m), utilizando el modelo lineal, pv= - 0.1882 + 3.7831 (dme), tiene la característica de incluir la raíz de las plantas; los criterios utilizados para elegir al modelo fueron, mayor valor en el coeficiente de determinación ajustado (R2 aj= 0.9689), valor mínimo del cuadrado medio del error (CME= 0.0371); valor mínimo porcentual de coeficiente de variación (CV= 15.29%); menor valor de probabilidad de cometer el error tipo I (p< 0.01) en la regresión y menor valor de la Información de Akaike (AIC= 0.7909).
El análisis de correlación permitió identificar la variable de mayor significación estadística con la biomasa pv (kg) la cuál fue dma (m), permitiendo asignar rangos de 0.1-1 m, con la posibilidad de diferenciar plantas en madurez de cosecha (dma> 0.25 m, h ≥ 0.3 m) y de esta manera no intervenir en 20% de ellas, la tabla de rendimiento solo se recomienda para hacer estimaciones de pv (kg) en el área de estudio ubicada al norte del estado de Zacatecas.