Introduction
Juniperus L. is one of the nine genera of Mexican conifers, among which Juniperus coahuilensis (Martínez) Gaussen ex R. P. Adams, commonly known in Mexico as táscate, is one of the most important species in the arid regions of northern Mexico. In addition to this region, the species is also distributed in parts of the southwestern United States, including central and southeastern Arizona, southwestern New Mexico, and western Texas (Farjon & Filer, 2013). In Durango, J. coahuilensis is found near the lower distribution limit of the pinyon pine (Pinus cembroides Zucc.) forest and appears as an ecotone with grasslands and, in some cases, with montane chaparral (Encina-Domínguez et al., 2019).
Juniperus coahuilensis wood varies in color from yellow to reddish-brown and tends to darken with age. It is considered moderately heavy, with a basic density of approximately 0.46 g∙cm-3, and shows good dimensional stability (Rendón et al., 2021). Although the wood is well-suited for pencil manufacturing, its limited availability means it is primarily used locally for constructing rural homes, crafting handmade furniture, and producing charcoal (Niembro et al., 2010). In the semi-desert region of Durango, táscate is also used for poles and firewood, and has been reported as an ecological niche and food source for some wildlife species. Due to its ecological, economic, and social importance in Durango's semi-arid zones, J. coahuilensis is classified as a species of Least Concern on the IUCN Red List (Rodríguez-Trejo & Vázquez-Soto, 2021).
While volume equation systems have been developed for a variety of temperate forest species (Vargas-Larreta et al., 2017) and tropical species (López-Martínez et al., 2020)-including Juniperus deppeana Steud (Vargas et al., 2012)-there is still a lack of volumetric modeling information for timber species in Mexico's arid and semi-arid regions (Rodríguez-Carrillo et al., 2015; Silva-García et al., 2018). Specifically, no studies have been conducted to estimate the volume of J. coahuilensis in the state of Durango. The only available reference comes from the Regional Forest Management Unit 1013 (FMU 1013), which reports an authorized annual harvest volume of 25 000 m3∙yr-1 (5 864 ha) until 2030, even higher than the authorized volume for Pinus cembroides in the same region (1 387 ha, 13 527 m3∙yr-1) (Prestadores de Servicios Profesionales, Agroforestales y Empresariales de Durango [PSPAED], 2009). However, that report does not detail the methodology or equations used to derive these volume estimates.
The volume of interest may include the entire stem or only a portion between the stump and a specific point along the stem. Equations that estimate volume up to a specific diameter or height along the stem are known as merchantable volume equations. Traditionally, these have been developed using two main approaches: (i) taper functions and (ii) ratio equations, the latter estimating merchantable volume as a proportion of the tree's total volume (Flores et al., 2021). Although ratio equations are relatively simple to develop and apply, taper functions are currently preferred.
Most studies involving taper functions have focused on temperate forest species, with a limited number addressing tropical species (Cruz-Cobos et al., 2023; Hernández et al., 2023; Li et al., 2024; López-Martínez et al., 2020; Lumbres et al., 2016). However, no specific equations have yet been published for native forest species in the semi-desert region of Durango. Therefore, the objective of this study was to develop the first compatible taper-volume system to describe stem profile and provide accurate estimates of stem, branch, and total tree volume for J. coahuilensis in the state of Durango.
Materials and Methods
Study area
The study was conducted in the municipalities of Hidalgo, El Oro, Indé, and Mapimí, which are part of the Regional Forest Management Unit (FMU) 1013, known as “Semidesierto Duranguense”, in the state of Durango (Figure 1). This FMU covers a total area of 4 736 722 ha, of which 2 753 704 ha (58.14 %) have forest cover. The forested area is primarily composed of 2 486 288 ha of arid zone vegetation (90.29 %), 147 843 ha of hydrophilic and halophilic vegetation (5.37 %), and 119 573 ha of forest (4.34 %) (PSPAED, 2009).

Figure 1 Location of the study area at the Regional Forest Management Unit (FMU) 1013 “Semidesierto Duranguense” in the state of Durango.
According to the Köppen climate classification system, as modified by García (1981), the region falls under the dry climate group (B), with subgroups S (semi-arid) and W (arid). Temperatures range from -3 °C to 28 °C, and average annual precipitation varies between 350 and 400 mm. The main commercially utilized timber species in the area include P. cembroides (pinyon pine) and Prosopis laevigata (Humb. et Bonpl. ex Willd) M. C. Johnst. (mesquite), both primarily used for charcoal production. Juniperus coahuilensis also plays a significant role in the region, covering approximately 20 000 ha of pure stands (Secretaría de Medio Ambiente y Recursos Naturales [SEMARNAT], 2013).
Database
Data were collected from 153 J. coahuilensis trees, selected subjectively to ensure a representative distribution per diameter and height classes. The following variables were measured: diameter at breast height with bark (DBH, cm), total height (Ht, m), crown ratio (cr, m; defined as the ratio of crown length to total height), diameter with bark (d i , cm) at each section based on its height above the ground (h i , cm), and diameter with bark of all branches with a basal diameter greater than 5 cm. Two sections were obtained: one at 0.30 m above the stump, and the next corresponding to the diameter at breast height (1.30 m). Additional sections were taken at 1 m intervals (or another commercial measurement) up to the tree top. Table 1 presents a summary of the number of observations, along with the mean, minimum, maximum, and standard deviation values for the measured variables.
Table 1 Descriptive statistics of the dataset (n = 153) for Juniperus coahuilensis.
| Variable | No. observations | Minimum | Maximum | Mean ± SD |
|---|---|---|---|---|
| Logs | 1 073 | 5.000 | 11.000 | 8.100 ± 1.160 |
| Stem volume (m3) | 153 | 0.003 | 0.184 | 0.022 ± 0.023 |
| Branch volume (m3) | 942 | 0.001 | 0.091 | 0.014 ± 0.015 |
| Total volume (m3) | 153 | 0.008 | 0.275 | 0.039 ± 0.034 |
| Diameter at breast height (cm) | 153 | 4.000 | 31.500 | 8.500 ± 4.212 |
| Total height (m) | 153 | 2.910 | 7.140 | 4.210 ± 0.632 |
| Crown ratio | 153 | 0.481 | 0.998 | 0.753 ± 0.122 |
| Stump height (m) | 153 | 0.050 | 0.400 | 0.086 ± 0.640 |
Total tree volume = stem volume + volume of thick branches. SD: standard deviation.
The trees were measured in sections using Smalian’s formula (Avery & Burkhart, 2002), with the tree top modeled as a cone. The total stem volume with bark was obtained by summing the volumes of the individual sections and the top. The volume of the branches was calculated using the same method, while the total tree volume was determined by adding the total volume of the stem and branches. Merchantable volume (v i , m3) was also calculated at a specific point where the diameter was equal to the merchantability threshold (d i ), along with the relative diameter (d i /DBH) and relative height (h i /Ht) of the stem. To detect potential outliers, a local quadratic non-parametric fit (Figure 2) was performed using local regression (LOESS) (Bi, 2000) with the SAS/ETS® statistical package (SAS Institute, 2021).
Fitted equation system
The function proposed by Fang et al. (2000), part of a group of segmented models, was fitted to describe the three possible stem shapes of a tree: neiloid, paraboloid, and cone (Husch et al., 1982). This function has demonstrated good performance in describing the stem profile of various forest species (Cruz-Cobos et al., 2023; Guzmán-Santiago et al., 2022; López-Martínez et al., 2020):
where: p 1 = h 1/Ht and p 2 = h 2/Ht are the relative heights at which the two inflection points assumed by the model occur; q = h i /Ht, k = π/40000, h st is the stump height (m) and a 0-a 2, b 1-b 3, p 1 and p 2 are parameters to be estimated.
The model developed by Fang et al. (2000) also includes equations for merchantable volume (vᵢ) and total stem volume (v f ) derived through direct integration of the stem profile function. These expressions are as follows:
As a first step, linear and nonlinear models were fitted to estimate branch volume (vᵣ), using DBH, Ht and cr as independent variables. The best results were obtained using DBH and cr; therefore, two new equations were included, where rᵢ are parameters to be estimated:
Fitting procedure
The fitting of compatible equation systems requires meeting certain statistical assumptions, such as ensuring compatibility and additivity. Additivity means that the sum of the estimated volumes for the stem and branches must equal the volume obtained directly from the total volume equation. To achieve this, the equations v f , v r and v t were fitted simultaneously using the Iterated Seemingly Unrelated Regression (ITSUR) method, implemented through the MODEL procedure in SAS/ETS® (SAS Institute, 2021).
Correction of heteroscedasticity and autocorrelation
Heteroscedasticity was addressed using weighted regression, with weights equal to the inverse of the variance of each observation (Parresol, 1999). The weighting factors used were 1/DBH 2 cr for vr and 1/DBH 2 Ht for vf and vt.
Autocorrelation in the longitudinal data was corrected using generalized nonlinear least squares, by expanding the error term through a second-order continuous autoregressive model [CAR(2)]:
where,
e ij = j-th ordinary residual of the i-th tree
e ij-k = j-th ordinary residual of the i-kth tree
l k = 1 for j > k and 0 for j ≤ k
ρ k = autoregressive parameter of order k to be estimated
h ij-h ij-k = distance between the jth-kth observations within each tree h ij > h ij-k
𝜀 ij = new error term under the assumption of independence.
The error structure expressed in the equation e ij was simultaneously fitted along with the mean structure of the stem function using the MODEL procedure of the SAS/ETS® statistical software (SAS Institute, 2021).
Model fit evaluation
The goodness of fit of the models was assessed through both numerical and graphical analysis of the residuals. Bias, root mean square error (RMSE), and the coefficient of determination (R2) estimated for nonlinear regression were examined (Ryan, 1997). In addition, model behavior was evaluated by analyzing d i and v i residuals at different relative heights, as well as by observing the evolution of bias and RMSE per diameter classes.
Results
Stem, branch, and total tree volume equations
Initially, the system of equations was fitted without expanding the error term (ρ), in order to account for autocorrelation. A clear trend in d i residuals was observed (Figure 3a). After correcting for autocorrelation, the trend disappeared (Figure 3c), and highly significant parameters (p < 0.0001) were obtained.

Figure 3 Diameter residuals (d i , cm) of Juniperus coahuilensis vs. lagged residuals (LAG) for the Fang et al. (2000) model fitted without accounting for autocorrelation (a) and using a continuous first- and second-order autoregressive model (b and c).
The parameter estimates and their standard errors, as well as the goodness-of-fit statistics for the simultaneously fitted equations, are presented in Table 2. The final form of the stem, branch, and total tree volume equations were as follows:
Table 2 Parameter estimates and goodness-of-fit statistics for the simultaneously fitted volume equations for Juniperus coahuilensis.
| Parameters | Equation | ||||
|---|---|---|---|---|---|
| di (cm) | vi (m3) | vf (m3) | vr (m3) | vt (m3) | |
| b1 | 0.000013 (0.00000) | 0.000013 (0.00000) | |||
| b2 | 0.000020 (0.00000) | 0.000020 (0.00000) | |||
| b3 | 0.000025 (0.00001) | 0.000025 (0.00001) | |||
| p1 | 0.095110 (0.00489) | ||||
| p2 | 0.525960 (0.03210) | ||||
| a0 | 0.000260 (0.00000) | 0.000260 (0.00000) | |||
| a1 | 1.552646 (0.00705) | 1.552646 (0.00705) | |||
| a2 | 0.658212 (0.01610) | 0.658212 (0.01610) | |||
| r0 | 0.001465 (0.00010) | 0.001465 (0.00010) | |||
| r1 | 1.26646 (0.03360) | 1.26646 (0.03360) | |||
| r2 | 0.610324 (0.10880) | 0.610324 (0.10880) | |||
| Bias | -0.08 | -0.00005 | 0.00017 | -0.00024 | 0.00026 |
| RMSE | 1.6994 | 0.0019 | 0.00013 | 0.00023 | 0.00015 |
| R2 | 0.921 | 0.989 | 0.857 | 0.453 | 0.826 |
Standard errors of the estimates are shown in parentheses. di = diameter, vi = merchantable volume, vf = stem volume, vr = branch volume, vt = total tree volume (stem volume + volume of large branches). RMSE: Root Mean Square Error.
All parameters were significant at the 1 % level. The percentage of variance explained by the models was 85.7 %, 45.3 %, and 82.6 % for stem volume, branch volume, and total tree volume, respectively. The RMSE values were 0.00013 m3 (v f ), 0.00023 m3 (v r ) and 0.00015 m3 (v t ).
Figure 4 shows the bias and RMSE of the volume estimates for the stem, branches, and total volume per diameter class. A similar pattern was observed for all three volume types, that is, a shift in bias from positive to negative as the diameter class increased, with a positive bias for the thickest trees. The RMSE values for branch volume estimates increased in the 5 to 15 cm classes and then decreased as tree diameter increased. In contrast, the RMSE for stem volume and total volume increased steadily up to the 25 cm class and then decreased in the 30 cm class.

Figure 4 Bias and root mean square error (RMSE) of the volume estimates for stem (vf), branches (vr), and total volume (vt) per diameter class for Juniperus coahuilensis.
Table 2 also presents the parameters and goodness-of-fit statistics (RMSE and R²) for d i and v i . The equations explained 92 % and 98 % of the variation in diameter and volume along the stem, respectively. The system's accuracy in estimating d i was 1.69 cm, while for v i it was 0.002 m3. Figure 5 shows the residuals of di and v i per relative height category, where no trend indicating heteroscedasticity issues was observed.

Figure 5 Residuals of diameter (d i ) and volume (v i ) against relative height classes (%) for Juniperus coahuilensis. The red point represents the mean prediction error.
The bias showed a similar trend in both equations, with values generally remaining around the zero line, except at relative heights above 80 %, where a positive bias was observed. Additionally, the accuracy of the d i equation improved as the height above the stump increased, with the mean error near the base of the tree close to 3 cm, while at the top of the tree it was less than 0.5 cm. The accuracy of the merchantable volume equation remained consistent in all relative height categories, with an RMSE below 0.005 m3 in all cases (Figure 6).
Discussion
During the development of this study, no volume equation systems were found for similar species in arid or semi-arid zones; however, models have been developed for other species. For example, in pine-oak forests in Durango, Cruz et al. (2012) and Vargas et al. (2012) generated stem volume equations for J. deppeana, with R² values of 0.94 and 0.96 and RMSEs of 0.101 and 0.052 m3, respectively. The quality of their model fits was superior to that observed in the present study (Table 2). Later, Vargas-Larreta et al. (2017) developed compatible taper-volume equations for more than 80 pine and oak species, as well as for 11 tropical species, based on the model by Fang et al. (2000), and reported goodness-of-fit values similar to those in this study. Similar results were obtained by Corral-Rivas et al. (2017) and Simental-Cano et al. (2017). On the other hand, in tropical forests of southeastern Mexico, López-Martínez et al. (2020) developed a compatible volume system for the 11 most important species using the same modeling approach, and their model fits yielded results similar to those obtained for J. coahuilensis.
The model developed by Fang et al. (2000) includes two inflection points (p 1 and p 2) and describes both the diameter along the stem and the height to a specific diameter. The first inflection point is located near breast height diameter, and the second occurs higher up the stem. In J. coahuilensis, the parameter p 1 was located at 10.3 % of the tree's total height, a value higher than those reported for Pinus durangensis Martínez (4.3 %; Cruz-Cobos et al., 2023), Pinus cooperi var. ornelasii (Martínez) Blanco (4.6 %; Corral-Rivas et al., 2017), and other pine species (4.7 %; Quiñonez-Barraza et al., 2014), but similar to that obtained for P. sylvestris L. (10.3 %; Diéguez-Aranda et al., 2006). The second inflection point (p 2) was found at 64 % of the tree's total height, which is consistent with the values observed for tropical species in southeastern Mexico (55-92 %; López-Martínez et al., 2020), P. sylvestris (60.7 %; Diéguez-Aranda et al., 2006), P. leiophylla Schltdl. & Cham. (60 %; Quiñonez-Barraza et al., 2014), P. cooperi (69.9 %; Cruz-Cobos et al., 2023), and P. durangensis (75 %; Cruz-Cobos et al., 2023). Although no similar studies were found for J. coahuilensis, the two inflection points appear to be suitable for describing the profile of its stem.
The taper function estimated the bark diameter in the lower part of J. coahuilensis trees with less accuracy, likely because this species sometimes has an irregular trunk or multiple stems, which may explain the biases associated with this portion of the tree. However, this lack of precision in the lower part of the tree has also been observed in other species. For example, Rodríguez and Broto (2003) reported lower accuracy and higher bias in the lower sections of Quercus pirenaica Willd., Populus euramericana (Dode) Guinier and Fagus sylvatica L. in Castilla and León, Spain, which they attributed to the relatively larger diameters in this part of the stem and the high variability in tree shape. Modeling this portion of the stem has also proven challenging for conifer species; for instance, Diéguez-Aranda et al. (2006) found reduced accuracy and increased bias at relative heights below 15 % in taper models for P. sylvestris. As in the present study, the accuracy of the taper equation improved in the portion of the stem between 10 and 60 % of total tree height (Figure 5), which corresponds to the most cylindrical part of the trees, showing a mean bias < 0.05 cm and average RMSE values < 1 cm.
The RMSE associated with estimating volume up to a specified diameter or height (commercial volume) remained consistent for all diameter classes (maximum RMSE < 0.005 m3). Although the goodness-of-fit values were lower than those reported in other studies-primarily for conifer species (Guzmán-Santiago et al., 2022; Tamarit-Urías et al., 2014; Tang et al., 2016)-the taper function by Fang et al. (2000) generally performed well for estimating both commercial and total stem volume for all tree size classes included in the sample.
On the other hand, modeling branch volume is challenging (González-Benecke et al., 2022), mainly due to the wide variety of crown structures, sizes, and shapes. As reported in other studies (García-Espinoza et al., 2018; Vega-Nieva et al., 2015), the accuracy of the branch volume equation was low, with RMSE = 0.00023 m3 and R² = 0.45. These goodness-of-fit values are similar to those found by Corral-Rivas et al. (2017), Gómez-García et al. (2015), Jiménez et al. (2013), and López-Martínez et al. (2020). Unlike previous studies that used only DBH as the independent variable in the branch volume equation, the equation developed in this study also includes crown ratio (cr). This additional variable more accurately explains the variability in branch volume, as crown architecture is directly influenced by stand density, as well as by slope and aspect (Hernández-Ramos et al., 2022).
Conclusions
This study presents the first compatible system of equations that simultaneously estimates stem volume, branch volume, commercial volume, and total tree volume. It also describes the stem taper of individual Juniperus coahuilensis trees in arid regions. The choice of equations will depend on their intended application. The total volume equation simplifies calculations and is recommended when product classification by commercial use is not required. The new commercial volume equation will be useful for species management, as the authorized volume allocation for táscate is typically divided into 70 % for secondary products and 30 % as waste. This equation will enable more accurate validation of that proportion and will provide technical support to the authority (SEMARNAT) in issuing legal documentation. Additionally, the equation system can enhance carbon monitoring in arid ecosystems by reducing errors in volume estimates.










texto en 





