#### Abstract

Basal area has shortcomings as a measure of stand density, but it is often preferred for operational assessments because it is easy to measure. Previous work has demonstrated that an additive version of Reineke's stand density index can be estimated by a simple tree count using a modified horizontal point sampling technique. We show that this technique can be extended further to estimate a mixed-species density measure that has been developed for complex stands in the northeastern United States, using wood specific gravity to harmonize the density contributions of different species. The sampling technique provides design-unbiased estimates of stand density from a weighted tree count, where the weights depend on specific gravity but not on diameter. Rounding the specific gravity values for different species in the calculation of estimates introduces a trivial amount of bias but streamlines the procedure for rapid use in the field.

#### 1. Introduction

The ability to measure forest stand density is fundamental for both description and prescription in applied ecology including silviculture, forest health, and wildlife habitat management [1, 2]. Basal area is the most commonly used density metric in practice, primarily because it is very fast and simple to measure in the field [3–5]. However, the basal area at which key ecological behaviors occur (such as crown closure, the onset of self-thinning, and asymptotic density) depends on the size of the trees in the stand, the species of those trees, and also potentially on the age and site quality of the stand. Other metrics, such as Reineke’s [6] Stand Density Index (SDI) and related metrics such as those of Curtis [7] and Long and Daniel [8], offer improvements over basal area for monospecific stands. Reineke’s [6] SDI is

where is number of trees per hectare, and is stand quadratic mean diameter. SDI has units of trees/ha. Zeide ([2, page 5]) writes that SDI is “still the best [metric] we have” for assessing stand density.

Recent research efforts have attempted to extend Reineke’s [6] basic work to more complex stands than the even-aged, monospecific subjects of the original SDI. Long and Daniel [7] developed an additive version of SDI closely related to an index explored by Curtis [8]:

where and are the number of trees per hectare and diameter of size classes of trees, and the summation is taken over an exhaustive, mutually-exclusive set of classes. (For continuous diameter distributions, the summation in (2) can be replaced by an integral and can be replaced by the product of and an appropriate probability density function; see [9, 10], e.g.) Several authors [1, 7, 11] have suggested using ASDI for multicohort stands. More recently, Woodall et al. [12, 13], building on earlier mechanistic work [14], suggested that the maximum ASDI of mixed species stands could be estimated as a function of the mean specific gravity of the wood of trees on the plots and proposed a modified density measure. Ducey and Knapp [15], building on [12], suggested a relative density metric of the form

where is the number of trees per hectare in the th diameter class of the th species, and are parameters to be estimated from plot data, and is the specific gravity of wood for the th species.

A challenge for any stand density metric more complex than basal area is that of obtaining a fast and simple estimate in the field. That challenge seems more and more difficult as stand density metrics such as that in (3) attempt to encompass stands that are more complex in diameter distribution and species mixture. Rapid field estimates are particularly valuable for diagnosis and description of stand conditions when a density metric will be used in situ in conjunction with expert opinion, and for controlling the marking of thinnings, shelterwoods, single-tree selection, and other partial harvests. Basal area would seem to have a distinct advantage in these situations, because a simple count of trees using a prism or other angle gauge, multiplied by the basal area factor of that gauge, provides a direct estimate of basal area per hectare [16, 17]. However, a modified prism technique can provide estimates of ASDI through a simple count as well [18].

The goal of this study is to show how the modified prism technique [18] can be modified to provide direct estimates of the Ducey and Knapp [15] relative density for mixed species stands. We provide tables and instructions for implementing the approach and illustrate its application with data from a mixed-species stand in southeastern New Hampshire.

#### 2. Theory

##### 2.1. Mixed-Species Density Measure

Ducey and Knapp [15] proposed RD as a stand density metric for the mixed-species forests of New England and New York in the USA. They used data from the Forest Inventory and Analysis program of the U.S. Forest Service to estimate and , using specific gravity at 12% moisture content for each species. RD is a unitless fraction, scaled so that represents a completely empty stand, and represents a stocking level consistent with “normal” or A-line stocking using conventional stocking guides. For example, if for a particular stand, that stand has 50% of the stocking level implied by “normal” or A-line stocking.

Now, consider (3). It can be rewritten as

which makes clear that RD is simply a weighted form of ASDI, where the weights for different species are determined by their specific gravity.

Ducey and Valentine [18] developed a technique analogous to sampling for basal area with a prism, in which tree count times a factor gives ASDI. We will show that it is possible to modify that technique to give RD as a sum of weights for the tallied trees. First, we review the modified prism technique for ASDI.

##### 2.2. Direct Sampling Technique

Recall that in ordinary horizontal point sampling (HPS) with a prism, relascope, or other angle gauge, each tallied tree represents the same amount of basal area per hectare [16]. That quantity is called the basal area factor or BAF, and that quantity is constant for a given gauge because the area of the inclusion zone (the area within which a specific tree will be tallied) is directly proportional to the basal area of the individual tree. This result follows from the geometry of the gauge: the radius of the inclusion zone (also called the limiting distance) equals the gauge constant (or ratio of length to width) times the diameter of the tree. Mathematical proofs are provided by [19], while [20] provides a thorough discussion of HPS in its many applications.

If the inclusion zone radius could be made proportional to , rather than proportional to , then inclusion zone area would be proportional to and each tallied tree would represent a constant amount of ASDI [18]. This quantity is called the stand density factor (SDF), by analogy to the BAF of ordinary HPS. Recall that ASDI has units of trees/ha; so SDF also has units of trees/ha. If the method of [18] is used, and if trees are tallied at points, then the estimate of ASDI is simply

If the SDF is chosen with reasonable care, then the limiting distances for a given SDF are very close to those for the BAF that would ordinarily be used in similar stand conditions [18] (Figure 1). The slight departure, depending only on tree diameter, could easily be adjusted in the field by “pushing the point” or “pulling the point.” Pushing the point [21] consists of holding the prism not directly over the sample point but extended toward the tree of interest. In ordinary HPS pushing the point is a common field error, and it creates biased estimates. Pulling the point consists of holding the prism back and sighting across the sample point to determine whether a tree should be tallied. A simple table of the distance to push or pull the point, in relation to for individual trees, is adequate to provide a sample with a constant SDF [18]. In practice, it is not usually necessary to know the exact of an individual tree to determine whether it should be included in the sample. Because the degree of push or pull changes slowly with , an initial ocular estimate is usually sufficient. should be measured on borderline trees, and if an estimate of ASDI is desired for different diameter classes, then can be measured on the tallied trees as well.

**(a)**

**(b)**

Now suppose that we wish to modify the technique to estimate RD, building on (4). Let be the number of trees of the th species tallied, using sample points as before. Then

provides an estimate of the ASDI of species , and

provides a natural estimator of the summand in (4). Rearranging slightly, and substituting, we have

as a natural estimator for RD. In other words, RD can be estimated as a weighted tree count, in which each tree at a single sample point represents units of RD. Recall that RD is a unitless fraction; therefore the weight or contribution of an individual tree is also a unitless fraction. Proofs of unbiasedness of (5), (6), and (8) are given in the appendix.

In practice, if “exact” (or at least, unrounded) values of are used, the quantity will likely be a cumbersome fraction. This is not a problem if, for example, field computers are used to provide nearly instantaneous feedback on sample results. However, for the cruiser or timber marker who desires a simple, rapid field procedure, it may be more convenient to group species and use an approximate value so that is always a nearly round percentage. This requires advance tabulation of species into groups, which we undertake hereinafter. Furthermore, it does introduce a small amount of bias due to rounding error. As we show, this is likely to be negligible for samples taken at typical operational intensities.

#### 3. Application

One of the first steps in direct sampling for ASDI is to choose a BAF and corresponding SDF [18]. In the northeastern United States, a BAF of 2.3 m^{2}/ha is a common choice but can lead to nondetection errors; so a larger BAF is desirable [22, 23]. Prisms with a BAF of 4 m^{2}/ha are widely available and fall within the range recommended for the region [24]. For a given BAF, the SDF should be chosen such that the limiting distances are similar. This allows the push or pull to be done quickly and accurately. A good match to BAF 4 m^{2}/ha is an SDF of 50 trees/ha: the limiting distances differ by less than 2 m over a broad range of DBH, and “pulling the point” is only necessary for rare trees with cm. The limiting distances and push/pull for an SDF 50 using a BAF 4 are shown in Table 1.

Ducey and Knapp [15] tabulate SG for all the species found in New York and New England in the USFS FIA database. However, for fast application in the field it is useful to tabulate rounded values of as weights. This is done in Table 2, with weights multiplied by 100 so that they reflect percents rather than fractions.

Armed with Tables 1 and 2, it is straightforward to conduct field sampling and to calculate estimates. For example, suppose that we have used our prism and Table 1 to determine which trees to sample at a single point. We tallied 7 trees: 3 *Betula papyrifera*, 1 *Picea rubens*, 2 *Abies balsamea*, and 1 *Fagus grandifolia*. We can immediately determine that the ASDI at this point is , and RD is , or 0.44.

#### 4. Field Demonstration

To illustrate the calculations associated with the method, we conducted sampling in a 4.3 ha mixed-species stand located in Durham, NH. 12 sample points were located systematically on a 60 m by 60 m grid. We used a BAF 4 m^{2}/ha prism, and the push-pull values from Table 1, to obtain a sample of trees with an SDF of 50 trees/ha. Both the exact weights as calculated from the specific gravities tabulated in [15], and rounded weights from Table 2, were used to estimate RD for comparison.

A total of 133 trees were tallied on the 12 plots. Of these, 49 were *Tsuga canadensis*, 33 were *Pinus strobus*, 31 were *Acer rubrum*, 10 were *Quercus rubra*, 5 were *Fagus grandifolia*, 4 were *Betula allegheniensis*, and 1 was *Quercus alba.* RD as calculated using the exact weights was 65.6%, with a standard error of 7.6%; the coefficient of variation of the sample points was 40%. RD as calculated using the rounded weights was 63.3%, with a standard error of 7.7%; the coefficient of variation of the sample points was 42%. In this stand at least, the difference caused by rounding the weights is trivial in comparison with the sampling variability of RD.

#### 5. Discussion and Conclusions

Metrics of stand density can be viewed in two ways: as hypotheses about or models of stand dynamics and as operational aids to silviculture. From the first perspective, stand basal area has long been rejected. Even in single-species, single-cohort stands, the inadequacy of basal area for comparing density in stands of different tree size and species was known long before Reineke’s [6] work and provided much of the impetus for his efforts. Even in those stands, Reineke’s work represents a simplification, obscuring potentially important interspecific differences in allometry and growth [25]. Extension of Reineke’s ideas to multicohort stands and other stands with complex structure has been an important area of recent work [7, 26–28]. While empirical tests of those efforts in comparison to previous measures have been relatively rare [9, 29], some progress has been made at least in understanding how one might approach the design of experiments [10].

A variety of studies have explored adapations of SDI to mixed-species stands [12, 13, 30–36]. However, the Ducey and Knapp [15] formulation has two properties that are both required for the straightforward use of the sampling method developed here. The first property is additivity [9]. An index is said to be additive if, for two stand components and , the density of and considered together equals the density of considered alone, plus the density of considered alone. Additivity is required for unbiased estimation of any stand density metric using traditional methods (such as fixed-area plot sampling or HPS, or any method that can be identified as using a Horvitz-Thompson [37] estimator) [38, 39]. Among previous mixed-species indices, only that of Stout and colleagues [31, 32] is additive; that of Woodall et al. [12, 13] meets the criterion of separability but not additivity [9]. By analogy with site index curves, we may also say that a density measure is anamorphic if, for two different species, the density contribution of an individual of species divided by that of an individual of the same but of species is a constant that does not depend on . The density measure of Stout and colleagues [12, 13] is not anamorphic. If a density measure is not anamorphic, then a different “push” or “pull” would be required for each species or species group. Thus, the Ducey and Knapp [15] formulation has certain practical advantages from a measurement standpoint. We caution, however, that simplicity and practicality may be advantageous but do not imply biological correctness [9]. Both additivity and anamorphism are likely incorrect from a strictly biological and mechanistic standpoint [9, 25, 40]. On the other hand, detecting departures from additivity and anamorphism implies the availability of an additive, anamorphic measurement from which to construct a baseline.

At its heart, a stand density index is a model, hypothesis, or conjecture about the relationship between measurable quantities and ecological behaviors. Statistician G.E.P. Box [41, page 2] wrote “Models, of course, are never true, but fortunately it is only necessary that they be useful.” The durability of basal area as a practical measure of stand density points to the importance of rapid, simple measurement as a component of utility. Box [42, page 792] also wrote “Since all models are wrong the scientist must be alert to what is importantly wrong. It is inappropriate to be concerned about mice when there are tigers abroad.” One tiger abroad in stand density measurement is sampling variability. Fast, repeatable measurement with relatively low sampling variability has been a reality for basal area since the seminal work of Bitterlich [16] on HPS. However, the recent generalization of HPS for ASDI [18], and the approach for RD as explored here, may close the gap somewhat.

In theory, the modified HPS approach [18] and this method could be used as a general sampling approach, if tree diameters are measured. However, the need to push (or pull) the point does open the door to field errors. The results from ordinary HPS should be more precise and less prone to operator error, even when a final variable of interest is RD. Rather, the strength of the approach we present here is for near-instantaneous estimation of stand density in three primary situations.

The first situation is to allow quick assessment and description of stand conditions. Fast calculation of RD would allow more effective coupling between quantitative measurement and subjective, expert opinion [18]. Density indices and related diagrams are best viewed as complementary to qualitative silvicultural expertise [43]; a rapid measurement technique can make quantitative information available to the forester or ecologist on-site.

The second one is to enable control of residual density while marking partial harvests such as thinnings and shelterwoods.

The third one is to provide a simple point density estimate for use in ecological studies that would be available immediately for stratification or subsampling purposes.

Equation [7] provides an estimating equation for RD when trees are sampled using the modified point sampling technique. In this appendix, we present a formal proof of design-unbiasedness for this estimating equation.

Suppose that we are interested in the population of trees occurring in a tract of area . Let the trees be indexed by , where is the total number of trees. is fixed but unknown, as are the diameters of the trees at the outset of sampling. Without loss of generality, let us locate a single point at random and with uniform density within the tract. This point will be our sample point, and from it we will include all trees located at a distance no greater than , where is a constant chosen in advance. Denote the inclusion area of the th tree by ; is the area where the sample point may fall and the th tree will be included. From basic geometry we have

and the probability that the th tree will be included in the sample (ignoring boundary slopover, which can easily be corrected using the walkthrough method [44]) is

We can rewrite (2) as

Now let be an indicator variable taking the value 1 if the th tree is included in the sample, and 0 otherwise; the expectation is . The total number of trees included in the sample is

and its expectation is

from which it follows directly that if , then mSDF is a design-unbiased estimator of ASDI, since

and if we collect data at n such sample points, the sample mean will also provide a design-unbiased estimate. This proves the unbiasedness of the estimator in (5). Furthermore, if we restrict the population to consist only of trees of species and replace with and ASDI with in the above equations, the equations still hold; this proves the unbiasedness of the estimator in (6).

In general, if there is an attribute associated with the th tree, an estimator of the form

is a Horvitz-Thompson [37] estimator and provides a design-unbiased estimate of the population total . Now, consider (4), rewritten as a population total:

where is the specific gravity of the th tree, so that

Now, if we substitute this into the basic Horvitz-Thompson estimator, we obtain

Now, recognizing that all trees of the th species have the same specific gravity for our purposes, and letting SG_{j} represent the specific gravity of the th species as before, we see that

is the design-unbiased Horvitz-Thompson estimator associated with a single sample point, and the mean of several sample points is identical to (8), which must also be design-unbiased.

#### Acknowledgments

Colleen Didas assisted with the field work. Development of the relative density index was supported under a grant from the National Science Foundation, “Collaborative Research: The northeastern carbon sink: enhanced growth, regrowth, or both?” to M. J. Ducey. This paper is a scientific contribution of the New Hampshire Agricultural Experiment Station.