Original article
Allometric equations to predict the total aboveground biomass of radiata pine trees
Centre for Timber Engineering, Edinburgh Napier University ,
Merchiston Campus, 10 Colinton
Road, Edinburgh,
EH10 5DT,
UK
^{*} Corresponding author:
j.moore@napier.ac.uk
Received:
14
October
2009
Accepted:
3
June
2010
• Radiata pine (Pinus radiata D. Don) is the main exotic plantation tree species grown in New Zealand for wood production and as such represents a significant component of the terrestrial carbon cycle.
• Using data for 637 trees collected in 13 different studies, a series of equations was developed that enable the total aboveground biomass of individual radiata pine trees to be estimated from information about height and diameter. A mixedeffects modelling approach was used when fitting these equations in order to account for random fluctuations in model parameters between studies due to site and methodological differences. Linear models were fitted to logarithmically transformed data, while weighted linear and nonlinear models were fitted to data on the original arithmetic scale.
• Based on a modified likelihood statistic (Furnival's Index of Fit), models fitted to transformed data were found to perform slightly better than weighted models fitted to data on the original arithmetic scale; however, the latter do not require a means for correcting for the bias that occurs when estimates of biomass obtained from transformed models are back transformed to the original scale.
• Recommendations for further development of these models including additional data collection priorities are given.
Key words: radiata pine / allometric equation / biomass / carbon sequestration
1. Introduction
Radiata pine (Pinus radiata D. Don) is a common plantation species in a number of Southern Hemisphere countries, including New Zealand, Chile, Australia, and South Africa, with more than 4 million hectares of plantations having been established (Lavery and Mead, 1998). In New Zealand, the total area of radiata pine plantations is approximately 1.6 million ha, which is 89 per cent of the total area of exotic plantations (Ministry of Agriculture and Forestry, 2007). Since the 1960s a considerable number of biomass studies have been conducted in this species in New Zealand, Australia and South Africa. These were generally undertaken to estimate forest productivity (as measured by dry matter production), nutrient budgets (e.g. Baker and Attiwill, 1985; van Laar and van Lill, 1978; Will, 1964) and how these are affected by a range of site and management factors (e.g. Barton, 1984; Cromer et al., 1985; Mead et al., 1985), as well as the partitioning of biomass among foliage, branches, stems and roots (e.g. Beets and Pollock, 1987). Much of this earlier research into radiata pine biomass has been summarised by Madgwick (1994).
Recently, there has been a renewed interest globally in biomass research in most tree species, including radiata pine, due to the need to predict forest carbon stocks and the potential amount of biomass available as a source of energy. One of the main drivers for this has been the Kyoto Protocol of the United Nations Framework Convention on Climate Change (UNFCCC), which stipulates mechanisms whereby post1990 carbon storage in forests can be used as an allowable carbon sink to offset some greenhouse gas emissions (IPCC, 2003). Under Article 3.3 of the Kyoto Protocol, net carbon stock changes on land subject to afforestation, reforestation, and deforestation must be estimated each year over the defined commitment period from information on land area and its corresponding carbon density. To help meet its obligations under the Kyoto Protocol, the New Zealand Government established the Emissions Trading Scheme (ETS), which introduced a price for greenhouse gases in order to provide an incentive to reduce emissions and increase forest sinks (New Zealand Government, 2008). Forestry is the first sector to enter the ETS and therefore a robust method is needed for determining the carbon stocks in forests registered under this scheme.
A number of approaches exist for estimating tree biomass. These include physiologicallybased models, such as 3PG (Landsberg and Waring, 1997) that use information on temperature, solar radiation and rainfall to predict total biomass; biomass expansion factors (e.g. Lehtonen et al., 2004) that convert stem volume into wholetree biomass; and allometric equations that predict total biomass from easily measured variables, such as diameter at breast height and height (e.g. Parresol, 1999; TerMikaelian and Korzukhin, 1997; Zianus et al., 2005). Allometric equations have been developed for predicting the biomass of individual radiata pine trees and stands within a number of previous biomass studies (see Keith et al. 2000 for a catalogue of equations). However, many of these equations are stand or sitespecific (e.g. Baker et al., 1984; Dargavel, 1970) and so potentially introduce high uncertainty when applied at a national level where a wide range of site types, stand structures and age classes exist. A number of individualtree equations were developed by Snowdon et al. (2000) using a dataset consisting of 351 trees compiled from various Australian and New Zealand biomass studies. However, additional data are known to exist from studies performed in New Zealand, which could be used to develop nationallyapplicable functions for predicting biomass. Many of these data have already been used to develop functions for different biomass components (e.g. Madgwick, 1983; 1985), but not total aboveground biomass.
Therefore, the objective of this paper was to develop functions that can be used to predict the total aboveground biomass of individual radiata pine trees growing in New Zealand. Mixed effects models are used to account for any peculiarities of the stands investigated or different measurement methodologies used in the different studies, following a similar approach taken by Wirth et al. (2004).
2. Materials and methods
2.1. Radiata pine biomass data
A database was compiled consisting of biomass measurements for 637 individual radiata pine trees from 13 different studies (Tab. I). Only those trees that had both diameter at breast height (D) and height (H) measurements were considered for subsequent analysis. These trees ranged in age from 2 to 42 y, in D from 7 to 806 mm and in height from 1.47 to 42.1 m. However, there was a relatively large number of smaller trees in the dataset (Fig. 1). Of the 637 trees in the dataset, 620 had information available on the stem, bark, branch and foliage components in addition to total aboveground biomass. The methods used to measure biomass differed between studies. In general, the biomass of smaller trees was determined by subdividing the whole tree into its components (i.e., stem wood, bark, branches and foliage) and weighing each of these components. In larger trees, the basal area of every branch in the tree was typically measured and a subsample of these branches were weighed. The relationship between branch component mass and branch basal area was determined and used to estimate these biomass components for the whole tree.
Summary of the 13 datasets used to develop the radiata pine biomass functions. Missing values indicate that data were either not recorded or were not available.
Figure
1 Distribution of (a) diameter at breast height and (b) total height of the 637 radiata pine trees used to develop the biomass equations. 
The majority of records in the database come from studies that have been carried out in the North Island of New Zealand where rainfall is in excess of 1000 mm per annum. Approximately half of the observations are from the studies conducted by Beets and Pollock (1987) and Smith et al. (1994). Most of the studies in eastern areas of the South Island, where rainfall is often less than 700 mm, have measured the biomass of trees less than 10 years old (e.g. Bandara, 1997; Clinton, 1990; Mason, unpublished data). In addition, a number of other studies have reported stand level estimates of biomass (see Madgwick, 1994 for details) or refer to individual tree measurements (e.g. Will, 1965); however, it was not possible to obtain these individual tree measurements in the course of assembling the database.
2.2. Development of biomass equations
A preliminary visual inspection of the data was undertaken to check for outliers and to investigate potential predictive relationships. Regression equations were developed to relate total aboveground tree mass (AGB) to different measures of tree size, including diameter at breast height (D [cm]), height (H [m]) and D^{2}H. The first of these equations was based on the traditional nonlinear allometric relationship which has the following form: where: AGB = total biomass or one of its components, X_{l}= tree dimension variable (l = 1,..., p), and β_{l} are the model parameters (l = 0,1,..., p). In equation (1) either an additive error term or a multiplicative error term can be assumed (Parresol, 1999). Because biomass data generally exhibit nonconstant variance, a common approach taken when developing predictive equations is to convert equation (1) into a linear equation with an additive error term by taking the natural logarithm of both sides, i.e., where ε is the residual error term. This transformation linearises the allometric relationships and stabilises the error variance (Baskerville, 1972). In this paper, the first set of equations fitted to the radiata pine data followed this approach, which also enabled parameter estimates and goodness of fit statistics to be compared with those from previous biomass equations developed for radiata pine (e.g. Snowdon et al., 2000) and other species (e.g. TerMikaelian and Korzukhin, 1997; Parresol, 1999; Zianus et al., 2005). The following equations with this form were fitted to the data: To account for methodological and site differences between studies, random fluctuations in the model parameters were included, such that or in vector notation: where AGB_{ij} is the above ground biomass at the ith measurement in the jth study, b_{lj} is the random deviation from the mean regression coefficient β_{l} due to study j, and ε_{ij} is the residual error ( ∼ N(0, σ^{2}). Models were fitted using the linear mixed effects model function within the R software program (R Project Core Team, 2009). Following the approach taken by Wirth et al. (2004), random effects representing differences between studies were only included for the intercept and coefficients of the first order terms in Equations (3)–(9). Random effects were assumed to be independent of the residual error and normally distributed with covariance matrix Ψ ( b _{j} ∼ N( 0 , Ψ )). While random effects for the different studies were assumed to be independent, random effects for an individual study were not assumed to be independent, and therefore the variancecovariance matrix for the random effects was positive definite symmetric (Pinheiro and Bates, 2000). Equations (3), (5) and (8) were also refitted to data from trees older than five years to allow direct comparison with the equations developed by Snowdon et al. (2000). Predictions from the two sets of equations were compared using Lin’s concordance correlation coefficient (Lin, 1989; 2000), which considers both the correlation and agreement between the sets of predictions.
To correct for the bias that occurs when estimates from logarithmic equations are back transformed to arithmetic units, the modified smearing estimate developed for linear mixed models by Wirth et al. (2004) was used. In matrix notation, the modified smearing estimate (i.e. the estimated expected value) of the total above ground biomass for a new tree (AGB_{new}) is given by: where X _{ new } = (1,X_{1new},X_{2new},.....,X_{pnew}), β̂ is the vector of estimated fixed effects from the model and cf_{1} and cf_{2} are correction factors that are given by: where b̂ _{ j } is the vector of estimated random effects associated with study j, Z _{ new } is the vector of predictor variables for which random fluctuations were calculated, and ε̂_{ij} is the vector of estimated residuals from the model. From equation (12b) it can be seen that the value of cf_{2} depends on the vector Z _{ new }. For this study the mean value of cf_{2} was calculated for the underlying dataset used to develop the models.
As an alternative to correcting for the bias that occurs when converting estimates from logarithmic models back to arithmetic units, a series of linear and nonlinear models were fitted to data on the original arithmetic scale. Specifically, the following models were fitted to the data: Weighted analyses were used when fitting the models given by equations (13)–(20) due to the nonconstant error variance. Weightings for each observation were determined by assuming that in each model the variance (Var(ε_{ij})) was a power of the absolute value of the variance covariate (v_{ij}), i.e. . When fitting equations (13) and (18), D was chosen as the variance covariate, while in the remaining equations it was D^{2}H. The value of the variance function coefficient δ was determined using varPower constructor within the R software program (R Development Core Team, 2009) and a common value was assumed across all studies. The linear models (Eqs. (18)–(20)) were fitted following the same approach taken for the models fitted to the transformed data, with random effects included for both parameters in each model. For the nonlinear models (Eqs. (13)–(17)), random effects were included for each parameter following the same approach taken for the linear models; however, in order for the algorithm to converge, the variancecovariance matrices for the random effects had to be simplified by assuming a block diagonal structure, rather than a symmetric structure (Pinheiro and Bates, 2000).
Figure 2 Relationship between total above ground biomass (AGB) for an individual tree and tree height (H), diameter at breast height squared (D^{2}) and D^{2}H. 
2.3. Model evaluation
Comparisons between models were made using several goodnessoffit statistics. Models with the same response variable (i.e., those based on the transformed allometric equation or those fitted on the arithmetic scale using weighted analysis) were compared using the Akaike Information Criterion (AIC). The AIC was calculated from the results of model fitting as − 2 loglikelihood + 2 P, where P is the number of parameters in the model (including elements of the variancecovariance matrix). In order to calculate AIC, models were refitted and the parameter estimates obtained using the maximum likelihood method (see Pinheiro and Bates, 2000). Using the AIC statistic, the ’best’ model was taken to be the one with the smallest value of AIC. Because AIC depends on the scale of the response variable, it cannot be used to compare models with different response variables. Therefore, Furnival’s index (FI) (Furnival, 1961) was used to compare transformed models with weighted linear and nonlinear models. This index was calculated as: where f^{′}(Y) is the derivative of the dependent variable with respect to total aboveground biomass (the square brackets signify the geometric mean) and RMSE is the root mean square error of the fitted equation. For untransformed models, FI is equivalent to the RMSE and a rescaled RMSE for transformed or weighted models. A fit index (R^{2}) was also calculated for each model using the following equation (Parresol, 1999): where TSS is the total sum of squares and RSS is the residual sum of squares. These were calculated as: where is the mean value of above ground biomass across all trees in the dataset and AĜB_{i} is the predicted value for the ith tree. In calculating these fit indices, model predictions were based on the fixed effects and predictions from the logarithmic models were back transformed and corrected for bias.
Models were also evaluated by examination of residuals, with an emphasis on investigating the existence of any possible trends with age, stand density and wood density as these factors could contribute to site level differences in total above ground carbon for a given tree diameter and height. Mean wood density for the whole stem was calculated by dividing total ovendry stem mass by the total stem volume under bark as determined from sectional measurements.
3. Results
Visual examination of the data (Figs. 2 and 3) revealed the existence of strong relationships between AGB and both D^{2} and D^{2}H on the original scale of measurement. The relationship between AGB and tree height was clearly nonlinear and there was also evidence of increasing variance in AGB with increasing values of all three tree size variables. On the logarithmic scale strong relationships were observed between log(AGB) and log(H), log(D) and log(D^{2}H) with a reduced heterogeneity in variance compared with the arithmetic scale (Fig. 3).
Figure 3
Relationship between total above ground biomass (AGB) for an individual tree and tree height (H), diameter at breast height squared (D^{2}) and D^{2}H on a logarithmic scale for x and y axes. 
3.1. Models fitted to transformed data
On the transformed scale, the models predicting AGB as a linear function of D (Eq. (3)) or the combined variable D^{2}H (Eq. (8)) were the poorest performed as measured by AIC (Tab. II). Comparison of the residuals from Equation (3) against predicted values (Fig. 4) (and similarly for Eq. (8)–not shown) revealed evidence of a systematic trend which indicated that the relationships were possibly not linear and that an additional higher order term could be required. The corresponding models containing additional quadratic terms (i.e., Eqs. (4) and (9)) had lower values of AIC than their linear counterparts. Further improvement to the model containing linear and quadratic diameter terms was obtained by the inclusion of tree height terms (Eqs. (6) and (7)). Overall, the ‘best’ model on the logarithmic scale in terms of AIC (i.e., the one with the lowest value) was that given by Equation (7) (Tab. II). Values for the first correction factor for the smearing estimate (cf_{1}) ranged from 1.0254 to 1.0474, while the mean value of cf_{2} for the underlying dataset ranged from 1.0527 to 1.0970. Values of the fit indices (R^{2}) based on the corrected back transformed predictions from the models ranged from 0.92 for the model based on a linear function of D (Eq. (3)) up to 0.95 for the models given by Equations (6) and (7), indicating that the latter are able to explain more of the variation in total above ground biomass.
Parameter estimates, their standard errors and goodness of fit statistics for the models developed to predict total above ground biomass. Values of root mean square error (RMSE) and the Akaike information criterion (AIC) can only be used for comparisons with the same class of model (i.e. logarithmic transformed models or weighted arithmetic models).
Figure 4
Comparison of residuals against predicted values of above ground biomass (on the logarithmic scale) for Equation (3) which contains diameter at breast height as the sole predictor. The solid line represents a locallyweighted smoothing function (lowess) fitted to the data. 
3.2. Weighted arithmetic models
On the arithmetic scale, the linear model (Eq. (20)) containing the combined variable (D^{2}H^{0.5}) had the lowest value of AIC (5237). The linear models containing D^{2} and D^{2}H (Eqs. (18) and (19)) had the highest values of AIC (5721 and 5571, respectively). However, the higher weighting (δ = 1.00) given to results from smaller trees in the model containing D^{2} (Eq. (18)), coupled with some apparent nonlinearity in the relationship with AGB, meant that this model considerably underpredicted AGB in larger trees (Fig. 5a). It was possible to rectify this problem by adding a higher order term to the model (either D^{3} or D^{4}). The resulting models exhibited more satisfactory behaviour for larger trees (Fig. 5b). The three and four parameter nonlinear models containing D and H (Eqs. (16) and (17)) also had relatively low values of AIC (5347 and 5349, respectively). For these models the value of the variance function coefficient δ was approximately 0.5 indicating the that variance was linearly proportional to D^{2}H. However, with this weighting there was some indication that Equations (16) and (17) underpredicted when AGB was greater than 1500 kg (data not shown), although there were fewer than 30 observations for trees with a total biomass greater than this value.
Figure 5 Comparison of measured aboveground biomass with that predicted by (a) Equation (18) and (b) the same Equation but with a higher order diameter term (D^{3}) added. 
3.3. Model evaluation
Comparison of the values of Furnival’s index showed that overall the transformed models performed better than the weighted arithmetic models (Tab. II), with Equations (6) and (7) having the lowest value of FI (10.67 and 10.66, respectively) of all the models developed. These models were able to explain approximately 95% of the variation in AGB, while the best weighted model (Eq. (20)) was able to explain approximately 97% of the variation. All three of these models had both tree height and diameter as inputs. Of the three models (Eqs. (4), (13) and (18)) that only had diameter as an input, Equation (4) was the best performing and was able to account for 95% of the variation in AGB.
Figure 6 Comparison of residuals from the model given by Equation (7) against tree age, stand density and wood density. 
Correction factors for the modified smearing estimates used when back transforming mixed effects models fitted on the logarithmic scale.
For the best overall model in terms of Furnival’s index (Eq. (7)), there was some evidence that it overpredicted biomass in stands younger than five years old (Fig. 6). Many of these trees were sampled by Bandara (1997) and this study had the largest random deviations from the model parameters of all the studies considered. An attempt was made to include a stand age term in this model, but this term was not significant (p > 0.05). This was possibly due to the high correlation between stand age and tree height (r = 0.92) which meant that stand age was already implicitly included in the model. There was no clear evidence of a trend in the residuals with stand density (Fig. 6). Mean wood density was only able to be calculated for 258 of the 637 trees in the biomass dataset (Tab. IV), with most of these trees from the studies conducted by Beets and Pollock (1987) and Smith et al. (1984). For these 258 trees, mean stem density ranged from 297 kg m^{3} up to 534 kg m^{3} with an average value of 344 kg m^{3}. There was some indication that the model underpredicted AGB for the nine trees with the highest values of wood density, but there were few observations at higher values of wood density upon which to assess the presence of any systematic trend.
Values of mean wholetree wood density from the seven studies where this parameter could be calculated.
3.4. Comparison with previous models
Refitting Equations (3), (5) and (8) to data from trees greater than five years of age resulted in the following parameter estimates: These were very similar to the parameter estimates for the equations developed by Snowdon et al. (2000). Comparison of the predictions from the two sets of equations for the range of tree sizes in the database assembled for this study showed that there was almost perfect agreement between them, as values of the concordance correlation coefficient all exceeded 0.99.
4. Discussion and conclusions
This metaanalysis built on the work of Magwick (1983) and Snowdon et al. (2000) and has produced a series of models that can be readily applied to predict the total aboveground biomass of individual radiata pine trees. These equations are intended to be applied at a national level, primarily for carbon accounting purposes. At a national level, unbiased estimates of carbon stocks are required, but biomass estimates obtained by applying these national level equations to individual stands may be biased as the relationship between tree dimensions and biomass may differ between stands. In developing these equations a mixed modelling approach was used to account for random perturbations caused by methodological and site differences, which had not been done in the studies by Magwick (1983) and Snowdon et al. (2000) that used ordinary least squares regression. However, comparisons showed that the functions developed in this study were very similar to those developed by Snowdon et al. (2000) using a dataset containing measurements from 351 radiata pine trees greater than five years in age.
The mixed modelling approach has been used by Wirth et al. (2004) to develop functions for Norway spruce (Picea abies (L.) Karst.) using data obtained from studies conducted by different authors. Wirth et al. (2004) assumed that differences between studies reflected methodological and site differences. In the radiata pine dataset, methodological differences were often related to tree size, with the biomass of small trees determined from weighing the entire tree, whereas subsampling was normally undertaken to determine the biomass of older larger trees. However, such methodological differences may not always be captured by the random effects structure assumed in the models as several of the studies sampled trees from across an age series (e.g. Beets and Pollock, 1987; Madgwick et al., 1977). Site differences may reflect different silvicultural treatments such as initial spacing and thinning intensity as well as fertiliser treatments; many of the data available were from studies investigating the effects of such treatments on biomass accumulation and nutrient uptake (e.g. Barton, 1984; Beets and Pollock, 1987). In some instances different treatments were applied within a study, but this level of intrastudy variation was not included as a random effect as records linking individual trees to a particular experimental treatment were incomplete.
Site differences could also be due to broad scale differences in stem form and wood density. In New Zealand, wood density of radiata pine can vary by up to 25–30% across the latitudinal range over which it is grown (Cown et al., 1991), and is related to variation in air temperature and the ratio of carbon to nitrogen in the soil (Beets et al., 2007a). While there was no strong evidence of bias across the range of wood density values for the best model developed in this study, there was some indication of an underestimate of AGB at high wood density and biomass data were only available for a subset of the latitudinal range on which radiata pine is grown in New Zealand. In particular, data were lacking from older trees grown in the southern parts of the South Island, i.e., Canterbury, Otago and Southland, where wood density is lower. Therefore, if the equations developed in this paper were applied to these stands, they could potentially overpredict aboveground biomass. For example, when the equations developed by Madgwick (1983) were applied to an independent dataset from Victoria, Australia they were found to underestimate stem biomass by 19% (Baker et al., 1984). One reason given for this underprediction was the high wood density of the Australian trees (Baker et al., 1983). It may be possible to capture some of the variation in the factors affecting wood density through the inclusion of appropriate sitelevel variables such as latitude, elevation and site index in biomass equations. Elevation and site index were included in the equations developed for Norway spruce branch and needle biomass by Wirth et al. (2004) and a similar approach could be taken for radiata pine; however, more information is required about the characteristics of those sites for which data are available and additional biomass data are required from a broader range of sites. In addition, explicitly representing these factors in national level biomass functions could reduce the potential bias that may occur when applying these functions to individual stands or regions.
The functions developed in this paper do not explicitly account for the effect of live crown pruning on total aboveground biomass. Approximately, 62 percent of the radiata pine forest estate in New Zealand is, or is expected to be, pruned to a height of at least four metres (Ministry of Agriculture and Forestry, 2007). Therefore, ignoring the effect of green crown pruning when estimating the total aboveground biomass of young radiata pine trees could result in aboveground biomass being significantly overestimated as the ratio of crown biomass to total aboveground biomass is relatively high (Cromer et al., 1985; Madgwick et al., 1977). In order to explicitly account for pruning, functions need to be developed to predict foliage and branch mass that include some measure of crown size, e.g. crown length. Using this approach, separate functions would be developed for each component of total aboveground biomass and the predicted values from these functions combined to give the predicted total for the entire tree. Ideally, these equations would be additive so that the sum of predicted biomass for tree components would equal the predicted total for the whole tree (Parresol, 1999).
Although estimates of below ground biomass are often required for estimating total live biomass, only aboveground biomass was considered in the analysis reported here. One of the main reasons for this was that data on root biomass were only available for 31 of the 637 trees in the dataset. However, it is possible to estimate below ground biomass using previously developed allometric functions (Jackson and Chittenden, 1981) or through the ratio of aboveground to belowground biomass (Beets et al., 2007b).
The models developed in this paper require validation against independent datasets. As is often the case in studies in which biomass functions have been developed for different species, full use was made of all available data. While it is possible to evaluate the accuracy of biomass equations by various crossvalidation techniques, this was not done here. The priority for further work is to incorporate data from other studies that have been conducted into the biomass database. Many of the papers that present standlevel biomass estimates describe the measurement of biomass for individual trees that was carried out in order to develop these standlevel estimates. It would also be worthwhile to combine the dataset of measurements from Australian studies that was used in the analysis performed by Snowdon et al. (2000). This would enable an even broader metaanalysis to be performed that would allow the possible significance of site and stand level factors in biomass equations to be tested.
The predictions from the biomass equations presented in this paper should also be compared with those obtained from the C_change model developed by Beets et al. (1999). This model predicts the amount of carbon in the stem, crown, roots and forest floor by using functions for partitioning tree growth into different components (i.e., stem, bark, branches and foliage) coupled with functions that model the mortality and decay of these components. This approach is being used to predict the carbon stocks and stock changes in New Zealand’s planted forests as part of the Land Use and Carbon Analysis System (LUCAS), which is used to fulfil New Zealand’s reporting requirements under the Kyoto Protocol (Stephens et al., 2008). It is therefore important that the two approaches give consistent results when they are used to calculate carbon storage in forests.
Acknowledgments
I would like to thank those researchers who willingly contributed their data which made it possible to develop the functions contained in this report. Drs Herb Madgwick and Craig Trotter, Richard Woollons and Mark Kimberley provided numerous useful advice in the preparation of this report. Funding for this project was provided by the New Zealand Ministry for the Environment under contract 03/040271L and the New Zealand Ministry of Agriculture and Forestry under contract MAF POL 080919.
References
 Baker T.G. and Attiwill P.M., 1985. Aboveground nutrient distribution and cycling in Pinus radiata D. Don and Eucalyptus obliqua L’H’erit. forests in southeastern Australia. For. Ecol. Manage. 13: 41–52. [CrossRef] (In the text)
 Baker T.G.,Attiwill P.M., and Stewart H.T.L., 1984. Biomass equations for Pinus radiata in Gippsland, Victoria. N. Z. J. For. Sci. 14: 89–96. (In the text)
 Bandara G.D., 1997. Ecophysiology of clonal and seedling trees of Pinus radiata D. Don in an agroforestry system. Ph.D. thesis, Lincoln University. (In the text)
 Barton P.G., 1984. Preliminary results from spray irrigating septic tank effluent into Pinus radiata forests at Warkworth. In: Paterson A. (Ed.), 16th New Zealand Biotechnology Conference, Forest Industries and Biotechnology, Palmerston North, pp. 258–275. (In the text)
 Beets P.N. and Madgwick H.A.I., 1988. Aboveground dry matter and nutrient content of Pinus radiata as affected by lupin, fertiliser, thinning, and stand age. N. Z. J. For. Sci. 18: 43–64.
 Beets P.N. and Pollock D.S., 1987. Accumulation and partitioning of dry matter in Pinus radiata as related to stand age and thinning. N. Z. J. For. Sci. 17: 246–271. (In the text)
 Beets P.N.,Kimberley M.O., and McKinley R.B., 2007a. Predicting wood density of Pinus radiata annual growth increments. N. Z. J. For. Sci. 37: 241–266. (In the text)
 Beets P.N.,Pearce S.H.,Oliver G.R., and Clinton P.W., 2007b. Root/shoot ratios for deriving belowground biomass of Pinus radiata stands. N. Z. J. For. Sci. 37: 267–288. (In the text)
 Beets P.N.,Robertson K.A.,FordRobertson J.B.,Gordon J., and Maclaren J.P., 1999. Description and validation of C_change: a model for simulating carbon content in managed Pinus radiata stands. N. Z. J. For. Sci. 29: 409–427. (In the text)
 Clinton P.W., 1990. Competition for nitrogen and moisture in a Pinus radiata  pasture agroforestry system. Ph.D. thesis, University of Canterbury, Christchurch. (In the text)
 Cown D.J., McConchie D.L., and Young G.D., 1991: Radiata pine wood properties survey. New Zealand Ministry of Forestry, FRI Bulletin No. 50 (revised edition). (In the text)
 Cromer R.N.,Barr N.J.,Williams E.R., and McNaught A.M., 1985. Response to fertiliser in a Pinus radiata plantation. 1. Aboveground biomass and wood density. N. Z. J. For. Sci. 15: 59–70. (In the text)
 Dargavel J.B., 1970. Provisional tree weight tables for radiata pine. Aust. For. 34: 131–140. (In the text)
 Dyck W.J., Hodgkiss P.D., Oliver G.R., and Mees C.A., 1991. Harvesting sanddune forests: impacts on secondrotation productivity. In: Dyck W.J. and Mees C.A. (Eds.), Longterm Field Trials to Assess Environmental Impacts of Harvesting, FRI Bulletin 161, Forest Research Institute, Rotorua. (In the text)
 Furnival G.M., 1961. An index for comparing equations used in constructing volume tables. For. Sci. 7: 337–341. (In the text)
 Intergovernmental Panel on Climate Change, 2003. Good Practice Guidance for Land Use, LandUse Change and Forestry. Institute for Global Environmental Strategies, Kanagawa.
 Jackson D.S. and Chittenden J., 1981. Estimation of dry matter in Pinus radiata root systems. 1. Individual trees. N. Z. J. For. Sci. 11: 164–182. (In the text)
 Keith H., Barrett D., and Keenan R., 2000. Review of allometric relationships for estimating woody biomass for New South Wales, the Australian Capital Territory, Victoria, Tasmania, and South Australia. National Carbon Accounting System Technical Report 5B. Australian Greenhouse Office, Canberra. 114 p. (In the text)
 Landsberg J.J. and Waring R.H., 1997. A generalised model of forest productivity using simplified concepts of radiationuse efficiency, carbon balance and partitioning. For. Ecol. Manage. 95: 209–228. [CrossRef] (In the text)
 Lavery P.B. and Mead D.J., 1998. Pinus radiata: a narrow endemic from North America takes on the world. In: Richardson D.M. (Ed.), Ecology and Biogeography of Pinus, Cambridge University Press, Cambridge, pp. 432–449. (In the text)
 Lehtonen A., Makipaa Heikkinen J.,Sievanen R., and Liski J., 2004. Biomass expansion factors (BEF) for Scots pine, Norway spruce and birch according to stand age for boreal forests. For. Ecol. Manage. 188: 211–224. (In the text)
 Lin L.I., 1989. A concordance correlation coefficient to evaluate reproducibility. Biometrics, 45: 255–268. [CrossRef] [PubMed] (In the text)
 Lin L.I., 2000. A note on the concordance correlation coefficient. Biometrics 56: 324–325. [CrossRef] (In the text)
 Madgwick H.A.I., 1983. Estimation of the ovendry weight of stems, needles and branches of individual Pinus radiata trees. N. Z. J. For. Sci. 13: 108–109. (In the text)
 Madgwick H.A.I., 1985. Dry matter and nutrient relationships in stands of Pinus radiata. N. Z. J. For. Sci. 15: 324–336. (In the text)
 Madgwick H.A.I., 1994. Pinus radiata – Biomass, Form and Growth. H.A.I. Madgwick, 36 Selwyn Road, Rotorua, 428 p. (In the text)
 Madgwick H.A.I.,Jackson D.S., and Knight P.J., 1977. Aboveground dry matter, energy and nutrient contents of trees in an age series of Pinus radiata plantations. N. Z. J. For. Sci. 7: 445–468. (In the text)
 Mead D.J.,Draper D., and Madgwick H.A.I., 1984. Dry matter production of a young stand of Pinus radiata: Some effects of nitrogen fertiliser and thinning. N. Z. J. For. Sci. 14: 97–108.
 Ministry of Agriculture and Forestry, 2007. A national exotic forest description as at 1 April 2006. Ministry of Agriculture and Forestry, Wellington, 62 p. (In the text)
 New Zealand Government, 2008. Climate Change Response (Emissions Trading) Amendment Act 2008. Public Act 2008 No. 85. (In the text)
 Parresol B.R., 1999. Assessing tree and stand biomass: a review with examples and critical comparisons. For. Sci. 45: 573–593. (In the text)
 Pinheiro J.C. and Bates D.M., 2000. MixedEffects Models in S and SPlus. SpringerVerlag, New York, 528 p. (In the text)
 R Development Core Team, 2009. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. ISBN 3900051070, URL http://www.Rproject.org. (In the text)
 Smith C.T.,Dyck W.J.,Beets P.N.,Hodgkiss P.D., and Lowe A.T., 1994. Nutrition and productivity of Pinus radiata following harvesting disturbance and fertilization of coastal sand dunes. For. Ecol. Manage. 66: 5–38. [CrossRef] (In the text)
 Snowdon P., Eamus D., Gibbons P., Khanna P.K., Keith H., Raison R.J., and Kirschbaum M.U.F., 2000. Synthesis of allometrics, review of root biomass and design of future woody biomass sampling strategies. National Carbon Accounting System Technical Report 17, Australian Greenhouse Office, Canberra. (In the text)
 Stephens P., Smith B., and Lincoln R., 2008. New Zealand: Forest carbon reporting and the role of forestry in climate change policy. In: Hendrick E. and Black K.G. (Eds.), Forests, Carbon and Climate Change – Local and International Perspectives, Proceedings of the COFORD conference held at the Glenview Hotel, Co Wicklow on 19 September 2007. COFORD, Dublin. (In the text)
 TerMikaelian M.E. and Korzukhin M.D., 1997. Biomass equations for sixtyfive NorthAmerican tree species. For. Ecol. Manage. 97: 1–24. [CrossRef] (In the text)
 van Laar A. and van Lill W.S., 1978. A biomass study in Pinus radiata D. Don. S. Afr. For. J. 107: 71–76. (In the text)
 Webber B. and Madgwick H.A.I., 1983. Biomass and nutrient content of a 29yearold Pinus radiata stand. N. Z. J. For. Sci. 13: 222–228. (In the text)
 Will G.M., 1964. Dry matter production and nutrient uptake by Pinus radiata in New Zealand. Commonwealth For. Rev. 43: 57–70. (In the text)
 Wirth C.,Schumacher J., and Schulze E.D., 2004. Generic biomass functions for Norway spruce in Central Europe – a metaanalysis approach toward prediction and uncertainty estimation. Tree Physiol. 24: 121–139. [PubMed] (In the text)
 Zianis D., Muukkonen P., Mäkipää R., and Mencuccini M., 2005. Biomass and stem volume equations for tree species in Europe. Silva Fenn. Monogr. 4: 63. (In the text)
All Tables
Summary of the 13 datasets used to develop the radiata pine biomass functions. Missing values indicate that data were either not recorded or were not available.
Parameter estimates, their standard errors and goodness of fit statistics for the models developed to predict total above ground biomass. Values of root mean square error (RMSE) and the Akaike information criterion (AIC) can only be used for comparisons with the same class of model (i.e. logarithmic transformed models or weighted arithmetic models).
Correction factors for the modified smearing estimates used when back transforming mixed effects models fitted on the logarithmic scale.
Values of mean wholetree wood density from the seven studies where this parameter could be calculated.
All Figures
Figure
1 Distribution of (a) diameter at breast height and (b) total height of the 637 radiata pine trees used to develop the biomass equations. 

In the text 
Figure 2 Relationship between total above ground biomass (AGB) for an individual tree and tree height (H), diameter at breast height squared (D^{2}) and D^{2}H. 

In the text 
Figure 3
Relationship between total above ground biomass (AGB) for an individual tree and tree height (H), diameter at breast height squared (D^{2}) and D^{2}H on a logarithmic scale for x and y axes. 

In the text 
Figure 4
Comparison of residuals against predicted values of above ground biomass (on the logarithmic scale) for Equation (3) which contains diameter at breast height as the sole predictor. The solid line represents a locallyweighted smoothing function (lowess) fitted to the data. 

In the text 
Figure 5 Comparison of measured aboveground biomass with that predicted by (a) Equation (18) and (b) the same Equation but with a higher order diameter term (D^{3}) added. 

In the text 
Figure 6 Comparison of residuals from the model given by Equation (7) against tree age, stand density and wood density. 

In the text 