Allometric Bark Biomass Model for Daphne bholua in the Mid-Hills of Nepal

Bark of Daphne bholua is an important non-timber forest product and makes a substantial contribution to the Nepalese economy. A precise estimate of the amount of D. bholua bark in mountain forests is possible using a biomass model. We developed an allometric bark biomass model for naturally grown D. bholua in Baglung District in the mid-hills of Nepal. The model was based on data from 101 destructively sampled D. bholua on 20 sample plots representing different growth stages (regeneration, established, and matured), site qualities, and stand densities, and we used diameter and height–diameter ratio as predictors. Among 9 functions evaluated, a simple power function showed the best fit to the data. This model described most of the variations in bark biomass with no substantial trends in the residuals. Leave-one-out cross-validation also confirmed the high precision of this model, because it described most of the variations in bark biomass with no substantial trends in the prediction errors. The model can be applied for a precise prediction of bark biomass for individuals of D. bholua with diameters and height–diameter ratios similar to those used in this study. It is site-specific, and its application should therefore be limited to sites with growth stage, site quality, stand density, and species distribution similar to those that formed the basis of this study. Further validation and verification of this model, with a larger dataset collected from sites with a wider range of these characteristics, is recommended.


Introduction
Plant biomass can be estimated by using direct or indirect methods (Vogt et al 1998). Direct methods, traditionally preferred, require destructive sampling-that is, felling an individual and separating and weighing its parts. This method is the most accurate, but it is costly and time consuming. Generally, it is more suitable for small plants and sample sizes, which require fewer resources than similar testing for big trees and large sample sizes. Since forest area is often large, estimation of tree biomass for an entire forest using destructive sampling is almost impossible. An indirect approach using allometric models is the most appropriate option in such a situation (Chaturvedi et al 2010(Chaturvedi et al , 2012Kuyah et al 2012;Subedi and Sharma 2012;Chaturvedi and Raghubanshi 2015). Allometric models serve as tools to relate variables such as diameter, height, height-diameter ratio, crown dimension, and wood density to biomass of the whole individual or its part. Felling and measurement of various components, including bark, of a number of sampled individuals of a plant species in any population is needed to establish a precise allometric biomass model, which can then be applied to predict biomass for the rest of the individuals in the population.
Because of diverse chemical and physical properties, tree bark could emerge as a potential raw material for various biomass-based industries. Bark offers a substantial added value to the biomass assortments and utilities. Because of the higher concentrations of several elements in the bark, it has higher metabolic activity than that of wood (Kenney et al 1990;Adler et al 2005). This is the reason that bark can be used for various purposes. For example, bark fibers from coppiced branches of Tilia cordata are used in Europe to make cordage for fishing and trap construction (Myking et al 2005). Bark of Betula species can be used as a covering in the making of canoes, as a drainage layer in roofs, and in the manufacture of shoes and backpacks (Adney and Howard 1964). Many commercial products are made from bark, including cork, Mountain Research and Development (MRD) An international, peer-reviewed open access journal published by the International Mountain Society (IMS) www.mrd-journal.org

MountainResearch
Systems knowledge cinnamon, and quinine; it is also used in tanning. Bark chips are used as mulch to protect soil moisture and nutrients (Adney and Howard 1964;Aronson et al 2009). The inner bark (phloem) of some species is also edible; for example, the Sami people of northern Europe use the bark of Scots pine (Pinus sylvestris) as a staple food (Zackrisson et al 2000). Similarly, bark of Cinnamomum tamala is used as a food flavoring (Subedi and Sharma 2012).
Non-timber forest products including bark play an important role in efforts to reduce the poverty of forestdependent people, contributing to livelihoods and improvements in food security, income, and health (FAO 1995;Falconer 1997;Ahenkan and Boon 2008;Uprety et al 2010). Most of the people residing in rural parts of developing countries like Nepal depend heavily on nontimber forest products for food, medicine, and income (Olsen 1998;Olsen and Larsen 2003;Chhetri and Gupta 2006;Kunwar and Bussmann 2008;Piya et al 2011). These products also play an important role in the development of entrepreneurship and contribute substantially to local and national economies, especially in mountain areas (Shrestha et al 2003;Gauli and Hauser 2009;Uprety et al 2010).
Daphne bholua (family Thymelaeaceae), locally known as lokta, is an important non-timber forest product that grows naturally in the mid-hills of Nepal. It is an erect 1 to 3 m tall shrub found at elevations between 1600 and 4000 m (Khadgi et al 2013). Two species of Daphne are common in Nepal: D. bholua and D. papyracea (Khadgi et al 2013). Both species have strong fibers in their bark, which is used to make handmade Nepali paper; because of its higher fiber density, D. bholua is preferred (Paudyal 2004;Kharal et al 2011;Khadgi et al 2013). The handmade paper is used to prepare various value-added products, which have strong markets in Europe, the United States, and Japan (Banjara 2007). D. bholua has also been identified in Nepal's Master Plan for the Forestry Sector as 1 of the 7 non-timber forest products available in Nepal (MPFS 1989). D. bholua covers 2.9 million ha in 55 districts in Nepal, and total bark stock has been estimated at 110,481 metric tons, which could support paper production of more than 950 metric tons every year (FSRO 1984). Several enterprises that produce high-quality Nepali paper from the bark of D. bholua provide employment to local people (Paudyal 2004).
Realizing its economic importance, the Nepal government and some private enterprises have attempted to promote sustainable management of D. bholua (Kharal et al 2011). However, the total amount of available D. bholua bark remains unknown, because of lack of bark biomass models and scientific forest inventory methods suited to local conditions. Even though several studies have been conducted on Daphne species (Paudyal 2004;Biggs and Messerschmidt 2005;Banjara 2007;Kharal et al 2011;Khadgi et al 2013), none of them developed a bark biomass model. This study developed such a model using data from 101 destructively sampled D. bholua in Baglung District in Nepal's mid-hills. The proposed model can be useful for predicting the total amount of D. bholua bark available in a forest and can thus support effective management of Daphne species.

Study area
This study was conducted in the Dekhothulojanti Community Forest of Lekhani village development committee (VDC), Baglung District, Nepal ( Figure 1). This district experiences subtropical to alpine climates, with mean annual rainfall ranging from 1500 to 3000 mm and mean annual temperature from 19 to 278C. Most of the forests in this district, including Dekhothulojanti, are managed by local user groups. Dekhothulojanti covers 129.12 ha and is a natural, uneven-aged, mixed-species forest. The main tree species are Abies spp, Rhododendron spp, Quercus lamellosa, Quercus semecarpifolia, Aesculus indica, Alnus nepalensis, Saurauia nepalensis, Semecarpus anacardium, and Daphniphyllum himalense (DFO 2012). Daphne bholua or Nepalese paper plant, the species on which this study focused, grows throughout this forest.

Sampling and measurement
The forest was divided into 6 strata based on stand conditions, growth stages, and management modalities. Stand condition could be even-aged or uneven aged, monospecies or mixed species, and sparse or dense. We took 3 growth stages into account (defined in the next paragraph). Management modality was defined as management practices performed by user groups-ie, priority given either to timber or fodder production or biodiversity conservation. We applied stratified systematic line plot sampling (sampling intensity about 0.1%), but only in the forest area where D. bholua grows (effective area). The sampling design applied in this study is presented in Table 1. D. bholua was absent in the first stratum, therefore only 5 strata were used for allocation of sample plots. In total, 20 sample plots (5 3 5 m) were used, allocated to each stratum using the proportional allocation method (Chaturvedi and Khanna 2011). We avoided sampling damaged individuals. We identified 3 growth stages based on the diameter (at 30 cm above ground) and height of D. bholua individuals: 1. Regeneration was defined as a seedling with diameter ,1 cm and height up to 1 m. 2. Established individuals were defined as those with diameter 1-3 cm and height .2 m or diameter .3 cm and height 1-2 m. 3. Matured individuals were defined as those with diameter . 3 cm and height . 2 m.
For each sample plot, regeneration and established individuals were counted and 2 from each category were felled; matured individuals were counted and all were felled. The total height of each felled individual was measured. We felled 101 individuals. Because of the greater abundance of smaller (regeneration and established) individuals in each sample plot, the number of felled individuals in those categories was higher than that for matured plants, even though all matured individuals were felled ( Table 2).
The stems of all felled individuals were debarked. Where applicable, large branches (.1 cm diameter) from established and matured individuals were also debarked and the bark added to the stem-bark samples. The total green weight of the bark of each felled individual was measured using a digital weighing machine with a precision of 10 mg. Of the 101 felled individuals, we randomly selected 75 and took 75% of the total bark sample of each to the laboratory, where it was dried in an oven at 808C until a constant weight was achieved. The weight of the sample was recorded after 24 hours and again at 6-hour intervals. It took 30 hours to get constant weights for most of the samples, and only a few samples needed to be kept in the oven longer (30-36 hours) in order to get constant weights. From this, we derived an average dry-to-fresh-weight ratio, which we used to estimate the dry weight of the remaining 26 individuals. Summary statistics for these data are presented in Table 2, and the relation of bark biomass to diameter is displayed in Figure 2. We plotted height-diameter ratio (HDR) against diameter ( Figure 3) and used the results to inform our modeling. The HDR is also known as a tree slenderness quotient and is defined as total height divided either by diameter at breast height or diameter at the base of tree or other smaller individual. As HDR is a measure of slenderness, it is an appropriate characteristic to describe the form of an individual plant (Opio et al 2000;Sharma et al 2016). We developed and compared models based on the first 2 approaches, with the second approach adjusted to use HDR instead of height.
In preliminary analysis, we fitted 9 different functions (Table 3) applying the first approach (modeling bark biomass as a function of diameter alone). However, the resulting models performed poorly, and this approach was dropped from further evaluation. Then we applied the second approach using the parameter prediction method (Clutter 1983;Chapagain et al 2014), in which parameter b 1 of each function in Table 3 was modeled as a nonlinear function of HDR. A simple power function (F1 in Table 3) with b 1 modeled as a nonlinear function of HDR best fitted the data and therefore was selected (Equation 1) for further evaluation.
where W i ¼ bark biomass of individual i (g); D i ¼ diameter of individual i (cm); HDR i ¼ height-diameter ratio of individual i (m/cm); b 1 , b 2 , a 1 , a 2 ¼ parameters; and E i ¼ residual error for individual i, which was assumed to be independent with zero mean and variance r 2 .

Model estimation and evaluation
We estimated models using PROC MODEL in SAS (SAS Institute 2012). The models were evaluated using root mean squared error (RMSE) and adjusted coefficient of determination (R 2 adj ) (Montgomery et al 2001). We examined graphs of residuals against the fitted biomass and predictors. We used a 1% level of significance in our analyses.
Cross-validation increases the reliability and credibility of models; a common method of doing this is use of split data (Vanclay 1994;Kozak and Kozak 2003;Yang et al 2004). However, our dataset was too small for this to be feasible. Instead, we evaluated our model by applying the leave-one-out cross-validation (LOOCV)   . This involved exclusion of 1 observation from the full dataset, fitting the model to the remaining data, and predicting the biomass of the left-out observation. This was replicated for all 101 observations. The RMSE and R 2 were calculated using the difference between predicted and observed biomasses; graphs of the mean prediction errors were also examined.

Results
Among 9 functions evaluated, a simple power function (Equation 1) emerged as the best, describing a larger part of the variations in bark biomass of D. bholua (R 2 adj ¼ 0.8873, RMSE ¼ 20.98) than the other functions (Table 3). Inclusion of HDR significantly improved the fit of statistics over that of the model with diameter as a single predictor. Without inclusion of HDR, the model fitted relatively poorly (R 2 adj ¼ 0.8516, RMSE ¼ 25.39). All other functions fitted with diameter as single predictor also showed poorer fit statistics than those of a simple power function (F1 in Table 3). All parameter estimates of each function in Table 3 were significant (p , 0.01) and biologically plausible. The allometric bark biomass model with estimated parameters is All symbols and abbreviations are the same as in Equation 1. No heteroscedasticity problem in the residuals was observed. Also, no serious systematic trend in the residuals was observed within the ranges of observed diameters and heights of the individuals. The LOOCV did   Huxley and Teissier (1936) 20.9858 0.8873 Rizvi et al (2008) 23.9807 0.8545 This study 22.5321 0.8700  not show substantial trends in the prediction errors for the observed ranges of diameter and height (Figure 4). The model described a large part of the variations in the observed biomass of the left-out individuals in the LOOCV (R 2 ¼ 0.8653, RMSE ¼ 23.03), and all parameter estimates from each iteration in the LOOCV were significant (p , 0.01). The model predicted the biomass data for larger individuals less precisely than that for smaller individuals. The bark biomass curves simulated using Equation 2 adequately covered the observed data with a few exceptions ( Figure 5).

Discussion
Even though the model described a large part of the bark biomass variations, some residuals seemed larger, but were distributed randomly. The absence of serious trends in the residuals and prediction errors within the observed range confirms the high precision of the model ( Figure 4). As in other studies (Subedi and Sharma 2012;Chapagain et al 2014), the model we developed using diameter as a single predictor performed poorly compared to the model that combined diameter and HDR and therefore could not be applied for precise prediction of all individuals in a stand. This is because D. bholua individuals with similar diameters within the same stand may have different heights, HDRs, and bark densities and consequently different bark biomasses. Variation in HDR was much greater for small diameters than for large diameters (Figure 3). This might be due to the sampling design, whereby individuals were selected from stands with large variations in competition. Greater competition results in thinner and taller individuals-that is, larger HDRs (Nyk€ anen et al 1997; Opio et al 2000;Sharma et al 2016). To describe more variation in the bark biomass, inclusion of HDR as an additional predictor is thus justifiable. A clear differentiation of the curves within the observed data range, even for individuals with similar diameters, is due to the large effect of HDR on the model. Because of differing individual heights (Table 2), their predicted biomasses are expected to be different, even for individuals with similar diameters, as shown in Figure 5. Adequate covering of the observed data by simulated curves also suggests that our model is precise for all possible sizes of individuals. Smaller variations of the residuals and prediction errors for smaller individuals (Figure 4) suggest that our model could be more precise for smaller individuals. Model precision depends largely on sample size, the predictors used, and application of the model to growth conditions of individuals. It may not be relevant to compare biomass models developed for different species, components, and locations, because several factors determine the accuracy of the model. Our model is as precise as a recent bark biomass model for Cinnamomum tamala (Subedi and Sharma 2012), which included both diameter and height as predictors. Our sample size (101 individuals) is larger than those used by other recent biomass studies (Ajit et al 2011;Sharma 2011;Subedi and Sharma 2012;Chaturvedi and Raghubanshi 2013;Chapagain et al 2014), which ranged from 27 to 46 individuals.
Total bark biomass of D. bholua in the forest can be estimated using our model. This information can be useful for sustainable management and effective commercial utilization of this economically important species. Several enterprises are producing high-quality handmade paper from the bark of D. bholua and providing employment to people residing in the mid-hills of Nepal (Paudyal 2004). Our model can also serve as a tool for predicting the total quantity of the raw material (bark) in the forest, which will directly contribute to feasibility studies for new enterprises to be established to meet growing demand for handmade paper, which is used in the preparation of various value-added products that have strong markets in Europe, the United States, and Japan (Banjara 2007).
Our model is site-specific and therefore may not be applicable to a wider geographical range, where growth conditions for D. bholua may differ. To make the model more comprehensive, accurate, and broadly applicable, it needs to incorporate measurements of site quality (eg site index), stand density (eg competition index), topographic features such as aspect, slope, and altitude; and climate and soil features such as temperature, precipitation, and soil properties-elements for which we did not have data.
Furthermore, we did not validate our model with data from other sites and regions, which would be necessary to confirm whether the model can be used for precise prediction under different growing conditions. Several studies have found that validating a model by applying LOOCV can yield less reliable results than validating it against external data (Molinaro et al 2005;Nord-Larsen et al 2009;Timilsina and Staudhammer 2013;Rushing et al 2015). Collecting a large amount of biomass data or data external to the main study site is often difficult due to lack of resources and time (Vanclay 1994). With a few exceptions (eg Brown et al 1993;Chave et al 2005), biomass modeling studies requiring destructive sampling often use small sample sizes from a small geographic area.

Conclusions
An allometric bark biomass model was developed using HDR and diameter of D. bholua individuals. The model described most of the variations in bark biomass (R 2 adj ¼ 0.8873, RMSE ¼ 20.98) with no substantial trends in the residuals. LOOCV confirmed the high precision of the model, because it described most of the variations in bark biomass (R 2 adj ¼ 0.8653, RMSE ¼ 23.03) with no substantial trends in the prediction errors.
Since site quality measures (eg site index), stand density measures (eg competition index), and measures of other factors that influence growth and development of D. bholua were not included in the model, it is site-specific. For high precision in prediction, application of the model should therefore be limited to site, growth stage, and stand conditions similar to those that formed the basis of this study.
The advantage of this model is that its application only requires measurement of outside bark diameter (0.3 m above ground) and total height of individuals to predict the total bark biomass of D. bholua in a forest, information which can then be used in formulation of a management plan. Further research is recommended to validate and verify our model using a larger dataset with a wider range of values for site quality, stand density, growth stage, and species distribution.  (Equation 2), overlaid on the observed biomass of D. bholua. The curves were produced using HDRs at 0.2 intervals, with the highest and lowest curves belonging to biomass for individuals with HDR of 0.2 and 1.6, respectively. The dots represent the same observed values shown in Figure 2.