Revista Mexicana Ciencias Agrícolas   volume 13   number 3   April 01 - May 15, 2022

DOI: https://doi.org/10.29312/remexca.v13i3.2908

Article

Biomass yield table for Euphorbia antisyphilitica in the north of Zacatecas

Héctor Darío González López

Salvador Valencia Manz

Dino Ulises González Uribe§

Autonomous Agrarian University Antonio Narro. Road Antonio Narro 1923, Buenavista, Saltillo, Coahuila, Mexico. CP. 25315. Tel. 844 4110299, ext. 2702. (hectordarioua@gmail.com; svalenciam2016@gmail.com).

§Corresponding author: digon-mx@yahoo.com.

Abstract

The candelilla plant is a natural resource of importance in the arid and semiarid areas of Mexico, it provides marginal economic benefits to the inhabitants of the places where it grows in the wild. The overexploitation and lack of management plans for its conservation motivated to conduct this research in 2019 and as an objective, to select a model to estimate the pv= green weight (kg) or plant yield in the north of Zacatecas. The following morphometric variables were measured; h= average plant height (m); dma= largest aerial diameter (m); dme= smallest aerial diameter (m); coba= aerial cover (m2). The model selected was pv= -0.1882 + 3.7831 (dme); with an adjusted coefficient of determination; R2adj= 0.9689; mean square error; MSE= 0.0371; coefficient of variation; CV= 15.29% and Akaike information; AIC= 0.7909, in this way it was possible to estimate the yield table for the plant.

Keywords: Akaike information, morphometric variables, regression.

Reception date: November 2021

Acceptance date: February 2022

Introduction

The natural resources in Mexico in addition to admired are known for showing a great plant wealth that can be explained by the combination of the factor’s climate, elevation, orography and soil, among others. Scientists in the past mentioned that the synergy given by the physiognomy and composition for millions of years caused the types of vegetation present in the country, of them, the one with the highest proportion is the xerophilous shrubland, which occupies approximately 40% of the total area. This represents the arid and semiarid national areas where forestry and agricultural activities are carried out, such as subsistence gathering carried out in conditions of great climatic variation by peasants who live in these places (Velazco, 1991; Rzedowski, 2006; Villarreal, 2006; Martínez, 2013).

The peasant activity of subsistence harvest of vegetable products of arid and semiarid climate is distributed in the Mexican Plateau, in the states of Querétaro, Guanajuato, Aguascalientes, Zacatecas, San Luis Potosí, Durango, Chihuahua, Coahuila, Nuevo León, Sonora and the peninsula of Baja California. On a smaller scale in the states of Oaxaca, Puebla, Hidalgo, State of Mexico and Tamaulipas, without it being less important for the benefit of the people who carry it out (Beltrán et al., 1964; Rzedowski, 2006; Villarreal, 2006; Martínez, 2013). This use focuses on wild vegetation such as Euphorbia antisyphilitica Zucc. (candelilla), Agave lechuguilla Torr. (lechuguilla), Lippia spp. (oregano), Opuntia spp. (nopal), Nolina spp. (palmilla), Yucca carnerosana Trel. (yucca) and Dasylirion spp. (sotol), among others (Martínez, 2013).

The activities carried out by man in these areas have produced changes in the landscape, which has impacted to a greater extent on the vegetation and the soil, an example has occurred in E. antisyphilitica, which ranks third in use nationwide with approximately 14 million hectares (Beltrán et al., 1964; Norma Oficial Mexicana, 2003; Rzedowski, 2006; Canales et al., 2006; CITES, 2009). The harvest of E. antisyphilitica is done manually, the individual is taken from the reproductive part called tiller and pulled from the ground with everything and root, if there are remnants, they recover slowly in periods of more than five years. Subsequently, by a field extraction process, a wax highly valued by the cosmetic industry is obtained (soaps, lipsticks, body creams, hair preparations and dental protections, among others).

This fact allows seeing that there is an alternative methodology that avoids the destruction of individuals, in which statistical models are used, which estimate the green weight or biomass of the plant from the quantification of morphometric variables easy to identify and quick to measure in the field, as a whole, the literature designates this procedure, sampling, measurement and application of statistical methods as allometry. Its usefulness lies in the fact that it can be applied to obtain yield tables of the species, a desirable action for inventories, characterization and conservation of E. antisyphilitica and in general of other species collected by peasants (Beltrán et al., 1964; Norma Oficial Mexicana, 2003; Cano et al., 2005; Rzedowski, 2006; Tapia and Reyes, 2008; CITES, 2009).

For this reason, the objective of this study was to select a linear, allometric or geometric model that estimates the green weight or biomass of the E. antisyphilitica plant to construct a yield table based on morphometric variables such as plant height, largest aerial diameter, smallest aerial diameter and aerial cover of the plant in wild populations in the north of Zacatecas.

Materials and methods

Study area

The wild populations of E. antisyphilitica object of this research are in the Ejido El Rodeo, municipality of Mazapil, north of Zacatecas. The geographical location was located at the coordinates 102° 11’ 00” west longitude and 24° 54’ 00” north latitude at an altitude of 1 990 m (CONABIO, 2014), the total area of the farm was 10 317.553 ha (Figure 1), this area is part of the Chihuahuan Desert, and a large proportion of the farm is covered by desert rosette-like shrubland. Species such as A. lechuguilla (lechuguilla), D. cedrosanum Trel. (sotol), Hechtia glomerata Zucc. (guapilla) and Agave stricta Salm-Dyck (espadín), among others, were observed in the place, they were sometimes mixed with non-thorny small-leaved shrubland associated with Larrea tridentata Cov. (Creosote bush) and Y. carnerosana Trel. (Villa, 1981).

Figure 1. Location of the study area.

Sampling of E. antisyphilitica plants with harvest maturity characteristics

Random points were located 100 m apart from each other in the study area, thus representing 48 sites in the field. In each of them, four quadrants were drawn, taking the north-south direction as the main axis and the east-west direction as the secondary. From the center, the nearest plant with characteristics of harvest maturity, such as largest aerial diameter of the individual greater than 0.25 m and minimum height of 0.3 m, was taken (Norma Oficial Mexicana, 2003; Velasco et al., 2009).

Four individuals per site were measured to have representativeness in the sampling and to be able to evaluate the variation, plants with smaller and greater measures than the previous ones and that are known to be harvested in a practical way were included. To be consistent with the regulations in force, there was a notice of use and the respective technical study, mentioning the actual stocks of the candelilla plant at the sampling site (Zar, 1999; Norma Oficial Mexicana, 2003). In total, n= 192 plants were measured to achieve statistical representativeness and meet the objective of this study.

Morphometric variables studied

The following was measured: h= plant height (m); dma= largest aerial diameter (m); and dme= smallest aerial diameter (m), specified in NOM-018-SEMARNAT-1999 (NOM-018, 2003). To determine the aerial diameters of each plant, an ellipse-shaped geometric figure was considered, from which the dma (m) was obtained by placing a tape measure on the widest part of the aerial cover in question, which represented the major axis of the ellipse. The dme (m) was measured by placing the tape measure on the narrowest part of the plant’s aerial cover, which represented the minor axis.

Subsequently, the semi-major axis (dma/2) and the semi-minor axis (dme/2) were obtained; the coba= aerial cover (m2) of each plant was the result of times the product of the two semi-axes. The variable pv= individual green weight (kg) was obtained after the measurement of the previous variables. Each plant was extracted in its entirety (including the root), each sample was shaken to reduce impurities (mostly soil), tied with tapes and each one was identified with a label, subsequently they were moved to a storage center where they were weighed with a digital scale (Canales et al., 2006; CITES, 2009; Martínez, 2013).

Statistical analysis of data

The numerical information of the present work was ordered in a double-entry database in Excel made up of the variables, h (m); dma (m); dme (m); coba (m2); and pv (kg), the latter represented the biomass of each plant including root. Once the digital information was formed in a matrix of rows and columns, it was exported to the R CRAN software, a data exploration analysis was performed on the variables, obtaining the estimators of the sample mean, sample deviation, minimum and maximum values, interval ͞x ± s (Where: ͞x= sample mean and s= sample deviation). The Shapiro-Wilk test (W) for normality was applied in each of the variables (Everitt and Hothorn, 2014).

To know the linear relationship between the possible pairs of variables of the matrix, Pearson’s simple correlation coefficient= r was estimated, statistical significance was obtained, focusing on values of p< 0.01. Because the main interest was the estimation of the pv (kg) of E. antisyphilitica, the associations of this with respect to the others were addressed. In this way, all the data were averaged in ten ranges (0.1-1 m) of the variable with the highest correlation with pv (kg) and the lowest P-value (Walpole et al., 2007) in the variables of the database.

Statistical models (Table 1) were adjusted to estimate pv (kg), considering the most significant morphometric variables (p< 0.01) of the database. The selection criteria were the highest value of the adjusted coefficient of determination (R2adj), minimum value of the mean square error (MSE), minimum percentage value of the coefficient of variation (CV), lowest probability value of making the type I error (P) in the regression model and the lowest value of the Akaike information (AIC) (Crawley, 2007; Walpole et al., 2007; Anderson, 2008; The R Core Team, 2017). Subsequently, confidence bands were constructed with the model (p= 0.05) in such a way that the data of pv (kg) involved were included at 95% reliability of the estimate.

Table 1. Models used to estimate pv (kg) of E. antisyphilitica.

pv = weight (kg); a and b = estimators of the regression parameters; h = plant height (m); dma = largest aerial diameter (m); dme= smallest aerial diameter (m); coba= aerial cover (m2).

Table of yield of pv (kg) for E. antisyphilitica

Once the model was selected, the yield table was created from ranges (0.1-1 m) of the morphometric variable with the highest statistical significance of the linear correlation between pairs of variables. This supported the separation of individuals at harvest maturity through a quantitative criterion supported by current regulations (Norma Oficial Mexicana, 2003). The tabular yield estimate was based on the included values of the space of the confidence bands (95%) of the selected model (Walpole et al., 2007; Anderson, 2008; The R Core Team, 2017).

Results and discussion

Statistical analysis and correlation between morphometric variables of E. antisyphilitica

The data exploration analysis of the n= 192 plants showed the estimators reported in Table 2. The correlation between the study variables (Table 3) showed the linear relationship and statistical significance (p< 0.01). All pairs of variables were significant (p< 0.01), the highest values of the correlation coefficients were, between coba (m2) and dme (m) (r= 0.90718), and coba (m2) and dma (m) (r= 0.88426). However, the pv (kg) also had them (p< 0.01) with the rest of the morphometric variables evaluated, although lower than the first mentioned. It should be noted that the correlations of interest are in the last row of Table 3, which consider the pv (kg) with respect to h (m); dma (m); dme (m); and coba (m2), respectively.

Table 2. Descriptive statistics for the morphometric variables of E. antisyphilitica, n= 192 plants.

h= height of the candelilla plant; dma= largest aerial diameter; dme= smallest aerial diameter; W= Shapiro-Wilk statistic.

Table 3. Correlation coefficients for the morphometric variables of E. antisyphilitica, n= 192 plants.

So, by taking r= 0.63625 of the linear relationship between pv (kg) and dma (m) as the highest correlation value, it was possible to form intervals with a distribution of 0.1-1 m. Subsequently, the scatter plots of Figure 2 were obtained for two involved variables. It can be noted that by choosing dma (m) for the reorganization of the information (Figure 2), the association of h (m); r= 0.8626 (before r= 0.36573); dme (m); r= 0.9861 (before r= 0.54529); that of dma (m); r= 0.9761 (before r= 0.63625) and coba (m2); r= 0.965 (before r= 0.63599) with respect to pv (kg) improved.

Figure 2. Scatter plots showing the correlation (r) between pairs of morphometric variables in E. antisyphilitica plants.

This allowed finding better adjustments for similar morphometries in the species studied, a reason supported and recommended by the basic theory of regression (Zar, 1999; Walpole et al., 2007), which coincided with works on other species of forest interest, for example, woody and legumes (Ares et al., 2002; Iglesias and Barchuk, 2010). In addition, including a correlation analysis provided the quantitative estimators and graphs of the linear associations (Table 3, Figure 2) that showed contrasts with other similar studies in the lack of a reproducible methodology with the evident interest of making estimates with models designated as allometric, when reviewing the literature, it refers to exponential or geometric curves that are linearized with the application of laws of logarithms and are based on the application of the least squares method where the coefficients are linear, have statistical support and it is not always easy to apply them (Ares et al., 2002; Walpole et al., 2007; Fonseca et al., 2013).

In the studies by Ares et al. (2002); Segura and Andrade (2008), a non-destructive methodology was applied to plant individuals to estimate biomass with allometric models, the present study also did so, in addition to showing an order through a reproducible methodology and aimed at selecting the best model under defined statistical criteria and that coincided with other studies (Segura and Andrade 2008; Iglesias and Barchuk, 2010; Fonseca et al., 2013).

Adjusted models and estimation of the table of yield of the pv (kg) of E. antisyphilitica

The estimation through models for the morphometric variables organized in intervals of 0.1-1 m and the selection criteria are concentrated in Table 4. The adjustment of the regression of each of the models showed that the first three that represent lines are the ones that had the best statistical estimators, followed by the geometric ones and finally the allometric ones, they exhibited lower values of R2adj and higher values in the AIC.

Table 4. Adjusted regression models for the morphometric variables of E. antisyphilitica.

pv= weight (kg); a and b= estimators of the regression parameters; h= plant height (m); dma= largest aerial diameter (m); dme= smallest aerial diameter (m); coba= aerial cover (m2).

Model three, pv= -0.1882 + 3.7831 (dme) (Table 4), was chosen for the estimation of the candelilla yield table, Table 5 shows the estimate through intervals of 0.1-1 m dme (m). Using linear equation three, all observations were included in the estimation space of the 95% confidence bands (Figure 3). If the model selection criteria for the n= 192 plants evaluated are taken as a basis, the estimation of the table with linear model three; pv= -0.1882 + 3.7831 (dme) practically only uses a field-measurable morphometric variable; dme (m), in this way pv (kg) will be obtained including the root, which represented an advantage, because it would no longer be necessary to extract individuals and consequently their destruction would be avoided, a result similar to that of Iglesias and Barchuk (2010) in legumes.

The statistical criteria indicated that the model on average estimated the value of biomass (kg) for values of the independent variable dme (m) that are substituted in the equation, this prediction falls within the 95% confidence bands of pv (kg). Divergent results were obtained by Fonseca et al. (2013), perhaps due to the inclusion of models to estimate biomass and carbon with criteria for linear coefficients in exponential curves without showing data exploration and correlation analysis, giving greater importance to the normal distribution of the models used, when theoretically it is incorrect. According to Zar (1999); Crawley (2007); Walpole et al. (2007); Everitt and Hothorn (2014); The R Core Team (2017), normality is possible for the continuous variables evaluated or the random error of a linear model and, in its case, also for one linearized with the help of logarithms, this was the reason why a size of n= 192 plants was ensured to meet a reliability of 95%.

Table 5. Yield table for the estimation of pv (kg) in E. antisyphilitica.

Figure 3. Confidence bands for the regression model pv= -0.1882 + 3.7831 (dme). Where: BCI= lower confidence band; BCS= upper confidence band, at 95%.

If attention is paid to linear model one, pv= a +b (coba), it had adequate statistical estimators, however, the variable coba (m2) is the result of the measurement of two aerial diameters (the largest and the smallest) and subsequently, the obtaining of the major and minor semi-axes and the product of them. Measuring coba (m2) does not provide ease in the calculations for obtaining the respective model of a yield table, which is why this equation was not chosen, in addition, it was far from what was found by Iglesias and Barchuk (2010) about the simplicity of the estimate, adding that they are two measurements and an intermediate calculation.

With respect to linear model two, pv= a + b (dma), its estimators are adequate, however, the adjustment of the chosen model was better; operationally and as cited in NOM-018-SEMARNAT-1999 (NOM-018, 2003), both dma (m) and dme (m) should be measured in the field; the statistical analyses of the present study showed that the first gives an order in ranges of 0.1-  m to all variables, the second was the one that improved the adjustment of the linear model for the estimation of the yield table. The results showed simplicity in the choice of a model, similar to what was made by Iglesias and Barchuk (2010), who made the same reference in legumes.

The allometric models four, five and six and the geometric models seven, eight and nine had estimators that were far from the statistical criteria, in addition, for their estimation and adjustment, greater knowledge of regression analysis is required, if the obtaining the variable coba (m2) is added to this, it will make it difficult to obtain this type of equations. In the present study, these models were surpassed by linear models and although allometric models are preferred (Iglesias and Barchuk, 2010), this was not the case for morphometric data of E. antisyphilitica in the study area. On the other hand, we differ from Fonseca et al. (2013) since their results are aimed at the same variable with a clear tendency to use exponential models of the geometric type.

The so-called model adjustments for biomass estimation have focused on allometric and geometric curves, linear ones are not included in all cases (Fonseca et al., 2013, Iglesias and Barchuk 2010), however, the statistical criteria used rely on linear estimators (Zar, 1999; Walpole et al., 2007; Everitt and Hothorn, 2014; The R Core Team, 2017). This preference is accompanied by the so-called linearization, which is the application of the laws of logarithms to a model that requires it, different studies avoid the correlation between involved variables as part of a previous and exploratory analysis of the phenomenon studied (Zar, 1999; Walpole et al., 2007; Fonseca et al., 2013).

In this research, it was found that a single variable estimated the weight pv (kg), in this case dme (m) with the simple linear regression model three, because the simpler it is in statistical structure and in the conformation of the same equation, the substitution of values will be direct and fast, as corroborated by Iglesias and Barchuk (2010) with statistical support, Zar (1999); Crawley (2007); Walpole et al. (2007); Anderson (2008); Everitt and Hothorn (2014); The R Core Team (2017), among others. Another benefit was the construction of a table with two columns, the first composed of values of the variable dme in ranges of 0.1-1 m; the second was formed by the estimated variable pv (kg), so the model, pv = -0.1882 + 3.7831 (dme), represented the biomass of the plant, which met the objective proposed in this work.

The coincidence of the regulations of the Norma Oficial Mexicana (2003) with the correlation theory shown in Zar (1999); Everitt and Hothorn (2014); Crawley (2007); Walpole et al. (2007); The R Core Team (2017), was reflected in the field since it gave the possibility of leaving 20% of the plants distributed in the area of exploitation at harvest maturity (dma > 0.25 m, h ≥ 0.3 m), representing no intervention in these individuals, ensuring the recovery of E. antisyphilitica for the future and on the other hand, in the administrative sense and in compliance with current regulations, attention was paid to the fact that the technical study had the notice of use, reporting the actual stocks of the plant (Norma Oficial Mexicana, 2003).

Conclusions

It was possible to estimate pv (kg) or biomass of the E. antisyphilitica plant from the dme (m), using the linear model, pv = - 0.1882 + 3.7831 (dme), it has the characteristic of including the root of the plants; the criteria used to choose the model were the highest value in the adjusted coefficient of determination (R2adj= 0.9689), minimum value of the mean square error (MSE= 0.0371); minimum percentage value of the coefficient of variation (CV= 15.29%); lowest probability value of making the type I error (p< 0.01) in the regression and the lowest value of the Akaike Information (AIC= 0.7909).

The correlation analysis allowed identifying the variable with the highest statistical significance with the biomass pv (kg), which was dma (m), allowing assigning ranges of 0.1-1 m, with the possibility of differentiating plants at harvest maturity (dma > 0.25 m, h ≥ 0.3 m) and in this way not intervening in 20% of them, the yield table is only recommended to make estimates of pv (kg) in the study area located in the north of the state of Zacatecas.

Cited literature

Anderson, D. R. 2008. Model based inference in the life sciences: a primer on evidence. Springer. Colorado, USA. 184 p. Doi: 10.1007 / 978-0-387-74075-1.

Ares, A.; Boniche, J.; Quesada, P.; Yost, R.; Molina, E. y Smith, T. J. 2002. Estimación de biomasa por métodos alométricos, nutrimentos y carbono en plantaciones de palmito en Costa Rica. Agron. Costarric. 26(2):19-30. https://www.redalyc.org/articulo.oa?id=43626202.

Beltrán, E.; Miranda, F.; Hernández, E.; González, C. A.; Scheffey, W.; Hernández, X. E. y Rzedowski, J. 1964. Las Zonas Áridas del Centro y Noreste de México y el aprovechamiento de sus recursos. Instituto Mexicano de Recursos Naturales Renovables AC. Ciudad de México, México. 186 p.

Canales, E.; Canales, M. V. y Zamarrón, E. M. 2006. Candelilla, del desierto mexicano hacia el mundo. CONABIO. Biodiversitas. 69(1):1-5.

Cano, P. A.; Berlanga, R.; Castillo, C. A.; Martínez, D. y Zárate, A. 2005. Análisis dimensional y tablas de producción de sotol (Dasylirion cedrosanum Trel.) para el estado de Coahuila. INIFAP CIRNE. Campo Experimental Saltillo. Folleto técnico núm. 18. Coahuila, México. 24 p. http://www.inifapcirne.gob.mx/Biblioteca/Publicaciones/182.pdf.

CITES. 2009. Convención sobre el comercio internacional de especies amenazadas de fauna y flora silvestres. Evaluación del estatus de Euphorbia antisyphilitica en México dentro de los apéndices de la CITES. Decimoctava reunión del comité de flora Buenos Aires, Argentina. 13 p. https://cites.org/sites/default/files/common/com/pc/18/X-PC18-Inf10.pdf.

CONABIO. 2014. Comisión Nacional para el Conocimiento y Uso de la Biodiversidad. Localidades rurales y urbanas 1 2010, escala: 1:1. Comisión Nacional para el Conocimiento y Uso de la Biodiversidad. Datos estadísticos de 2010 del Instituto Nacional de Estadística y Geografía (INEGI). Ciudad de México. http://www.conabio.gob.mx/informacion/ metadata/gis/loc2010gw.xml?-xsl=/db/metadata/xsl/fgdc-html.xsl&-indent=no.

Crawley, M. J. 2007. The R book. Chichester: John Wiley & Sons. West Sussex, Inglaterra.  323-487 pp. Doi: 10.1002 / 9780470515075.

Everitt, B. T. and Hothorn, T. 2014. A handbook of statistical analyses using R. CRC Press. 3ra Edición. Florida, USA. 456 p. Doi: 10.1201/b17081.

Fonseca, W.; Ruíz, L.; Rojas, M. y Alice, F. 2013. Modelos alométricos para la estimación de biomasa y carbono en Alnus acuminata. Rev. Cienc. Amb. 46(2):37-50. Doi: 10.15359/rca.46-2.4.

Iglesias, M. y Barchuck, A. H. 2010. Estimación de la biomasa aérea de seis leguminosas leñosas del chaco árido (Argentina). Ecología Austral. 20(1):71-79. https://ri.conicet.gov.ar/ handle/11336/100072.

Martínez, S. M. 2013. Ecología y usos de especies forestales de interés comercial de las zonas áridas de México. Instituto Nacional de Investigaciones Forestales, Agrícolas y Pecuarias (INIFAP)-Centro de Investigación Regional Norte Centro. Sitio Experimental La Campana, Cd. Aldama. Libro técnico núm. 5. Chihuahua, México. 216 p. https://www.researchgate.net/profile/Martin-Martinez-Salvador/publication/261323618-Libro-zonas-aridas/links/00463533dceec7f707000000.pdf.

NOM-018. 2003. Norma Oficial Mexicana NOM-018-SEMARNAT-1999. Que establece los procedimientos, criterios y especificaciones técnicas y administrativas para realizar el aprovechamiento sostenible de la hierba de candelilla, transporte y almacenamiento del cerote. http://www.semarnat.gob.mx/leyes-y-normas/normas-oficiales-mexicanas.

Rzedowski, J. 2006. Vegetación de México. Comisión Nacional para el Conocimiento y Uso de la Biodiversidad (CONABIO). 1ra. (Ed.). Digital. Ciudad de México. 504 p. https://www.biodiversidad.gob.mx/publicaciones/librosDig/pdf/VegetacionMx-Cont.pdf.

Segura, M. y Andrade, H. J. 2008. ¿Cómo hacerlo? ¿Cómo construir modelos alométricos de volumen, biomasa o carbono de especies leñosas perennes? Agroforestería de las Américas. 46:89-96. https://www.worldcocoafoundation.org/wp-content/uploads/files-mf/segura 2008a.pdf.

Tapia, T. E. y Reyes, C. R. 2008. Productos forestales no maderables en México: aspectos económicos para el desarrollo sustentable. Maderas y Bosques. 14(3):95-112. Doi: 10.21829/myb.2008.1431208.

The R Core Team. 2017. R: a language and environment for statistical computing reference index version 3.4.2 (2017-09-28). Vienna, Austria. 1206-1564 pp. https://cran.r-project.org/.

Velasco, B. E.; Arredondo G. A.; Zamora M. M.; y Moreno, F. 2009. Modelos predictivos para la producción de productos forestales no maderables: Lechuguilla. Instituto Nacional de Investigación Forestales Agropecuarias y Pecuarias. Manual técnico núm. 2. Ciudad de México, México. 56 p.

Velazco, H. 1991. Las Zonas Áridas y Semiáridas. Sus características y manejo. (Ed.). Limusa. Ciudad de México, México. 27-70 pp.

Villa, A. 1981. Investigación Sobre Candelilla (Euphorbia antisyphilitica) en el área de CIFNE. Rev. Cienc. Forest. 6(32):3-18.

Villarreal, J. A. 2006. Introducción a la botánica forestal. Editorial Trillas. 3ra. (Ed.). Ciudad de México, México. 151 p. 

Walpole, E. R.; Myers, R. H.; Myers, L. S. y Ye, K. 2007. Probabilidad y estadística para ingeniería y ciencias. Pearson educación. 8a. (Ed.). Ciudad de México, México. 816 p. https://vereniciafunez94hotmail.files.wordpress.com/2014/08/8va-probabilidad-y-estadistica-para-ingenier-walpole-8.pdf.

Zar, J. H. 1999. Biostatistical analysis. Prentice Hall. 4th (Ed.). New Jersey, USA. 663 p.