Using a stand-level model to predict light absorption in stands with vertically and horizontally heterogeneous canopies

David I Forrester,Rubén Guisasola,Xiaolu Tang,Axel T Albrecht,Tran Lam Dongand Guerric le Maire

Using a stand-level model to predict light absorption in stands with vertically and horizontally heterogeneous canopies

David I Forrester1*,Rubén Guisasola1,Xiaolu Tang2,Axel T Albrecht3,Tran Lam Dong4,5and Guerric le Maire6

Background:Forest ecosystem functioning is strongly influenced by the absorption of photosynthetically active radiation(APAR),and therefore,accurate predictions of APAR are critical for many process-based forest growth models.The Lambert-Beer law can be applied to estimate APAR for simple homogeneous canopies composed of one layer,one species,and no canopy gaps.However,the vertical and horizontal structure of forest canopies is rarely homogeneous.Detailed tree-level models can account for this heterogeneity but these often have high input and computational demands and work on finer temporal and spatial resolutions than required by stand-level growth models.The aim of this study was to test a stand-level light absorption model that can estimate APAR by individual species in mixed-species and multi-layered stands with any degree of canopy openness including open-grown trees to closed canopies.

Complex forests;Mixed-species;Stand structure;Extinction coefficient;Lambert-Beer law;Light absorption

Background

The absorption of photosynthetically active radiation (APAR)by trees is an important determinant of their growth and accurate estimates of APAR are often critical for process-based growth models.The fraction of PAR(f)that is absorbed by homogeneous canopies can be estimated using the Lambert-Beer equation(Monsi and Saeki 1953;or in English as Monsi and Saeki 2005)

where k is the light extinction coefficient for the considered growth period,and L is the leaf area index (m2m-2).However,most forests and plantations do not have homogeneous canopies;they may consist of multiple species or the canopies may contain gaps,such as in young stands or resulting from thinning and following natural disturbances.Tree-level light absorption models have been developed to deal with this canopy heterogeneityand some have been shown to give comparable predictions to field measurements of APAR(Norman and Welles 1983;Oker-Blom et al.1989;Wang and Jarvis 1990; Bartelink 1998;Brunner 1998;Canham et al.1999;Bartelink 1998;Courbaud et al.2003;Gersonde et al.2004;Abraha and Savage 2010;Ligot et al.2014b).Inputs for these models may be the leaf area of each tree,vertical and horizontal leaf area distributions,leaf angle distribution,leaf and soil optical properties,and x and y coordinates to indicate the tree positions.Generally,the accuracy of these tree-level models increases with the level of detail used to describe the tree crowns and canopy structure(Brunner 1998;Sinoquet et al. 2000;Parveaud et al.2008).

Some of these models,such as the Maestra(or Maestro) model(Grace et al.1987;Wang and Jarvis 1990;Medlyn 2004),have been shown to provide realistic predictions of absorbed light in mixed-species stands(le Maire et al. 2013;Charbonnier et al.2013).However,these tree-level models often require extensive input data and can have high computational demands.In contrast,several equations based on Equation 1 have been developed to predict APAR at the stand level using less information for parameterisation,and models using these can run faster than tree-level models(Duursma and M?kel? 2007; Charbonnier et al.2013;Forrester 2014).Such models are therefore useful for stand-level process-based models that need to be easy to use and parameterise and quick to run at large temporal and spatial scales(Ligot et al.2014a).

Stand-level light absorption models,such as Equation 1 or modifications of it,are often used in process-based forest growth models but they are rarely tested for the heterogeneous canopies being modelled,even when the same growth models are thoroughly tested for their predictions of growth,transpiration,carbon partitioning and nutrient availability.This project had three aims. The first was to test the APAR predictions of a standlevel light absorption model,similar to that described by Forrester(2014),for a wide range of canopy structures. The model used by Forrester(2014)has been tested against a detailed tree-level model but not using data from real forests.The second aim was to compare this stand-level model to several other approaches that have been used to predict APAR by mixed-species canopies. The third aim was to determine which parameters this model is most sensitive to.

To achieve these aims,stand-level APAR predictions were compared with those of a detailed tree-level model, Maestra(Grace et al.1987;Wang and Jarvis 1990;Medlyn 2004),which was run using data collected from five experiments:(i)thinned and unthinned Eucalyptus nitens plantations in southeastern Australia,(ii)monospecific and mixed-species plantations of Eucalyptus grandis and Acacia mangium in the tropics of Brazil,(iii)mixedspecies plantations of Hopea odorata and Acacia hybrid (A.mangium×A.auriculiformis)in Vietnam,(iv)mixedspecies forests and Cunninghamia lanceolata plantations in the subtropics of China,and(v)temperate mixedspecies forests composed of Abies alba,Picea abies and Fagus sylvatica at six sites in southwestern Germany,each of which was managed according to a range of thinning treatments.

Methods

Description of study sites

The stands in this study were composed of monocultures or mixtures of up to four species and all but one of the five data sets contained different tree spacing treatments,resulting in different vertical and horizontal complexity of the canopies(Table 1).

In most of the mixed-species stands the mixing of species was relatively uniform;i.e.tree-by-tree mixtures rather than a spatial clustering of species.However,the mixtures of Acacia hybrid and H.odorata in Vietnam (Stand 3,Table 1)included two different nurse crop designs.One of these consisted of rows of H.odorata (about 1.5 m tall)planted between rows of much taller Acacia hybrid(about 14 m tall).The second design was a 22-m diameter circular gap inside an Acacia hybrid plantation that had been planted with H.odorata trees. The stand-level model assumes that species mixtures are intimate.Therefore,this planting design allows for testing of possible errors that can result when this assumption is not met.Details about each stand are provided in Table 1.

Description of the stand-level model

The stand-level light model is based on the model described by Forrester(2014)and calculates APAR using monthly PAR,and species(or age class,or dominance class)-specific information including the number of trees per ha and mean crown characteristics such as livecrown length(LCL,m),crown diameter(m),tree height (m),tree leaf area(LA,m2),and crown shapes,which are approximated by simple geometric forms like ellipsoids, half-ellipsoids,spheres or cones,depending on the respective species.Table 2 contains a list of parameters.

When a stand contains more than one species the canopy is divided into horizontal layers that do not overlap vertically(Figure 1).A given layer may consist of a single species,or several species whose crowns do overlap vertically.A given species is only present in a single layer and,as shown in Figure 1,the crowns are not split between different layers.However,an exception to this is when there are two age classes of a given species between which there is no vertical overlap.In that case, each age class is treated as a different species.

Equation 2 is used to calculate the APAR by each layer, starting with the highest layer,and then reducing theavailable PAR for each layer by the APAR of all higher layers.

Table 2 A list of the abbreviations and symbols used in the text and in equations

The fraction of PAR intercepted by the whole canopy layer is the term inside the square brackets.kHis defined as an extinction coefficient for a homogeneous canopy. Typical long-term(monthly,annual)extinction coefficients,k(Equation 1)vary for a given species as the crown architecture and canopy structures(e.g.the space between neighbouring crowns)change with factors such as age and resource availability(Binkley et al.2013; Bryars et al.2013).Duursma and M?kel?(2007)showed that this variability in k could be accounted for by first replacing k in Equation 1 with an extinction coefficient for a homogeneous canopy,kH.This canopy is characterized as being composed of trees of the same height,with the same live-crown length,having box-shaped crowns that fit together perfectly(no space between crowns), and of the same leaf-area density(LAD,leaf area per crown volume,m2m-3),leaf angle distribution,leaf reflectance and leaf transmittance.In such stands the kHis independent of trees per ha.

Figure 1 Example canopy structures that can be modelled using the stand-level light absorption model.(a)and(b)show simple structures with only one layer and one species but varying in stand density.(c)shows a single-layered canopy containing three species whose crowns overlap vertically,while(d)and(e)show canopies with two layers,each containing either one or two species.Modified from Forrester(2014).

The λviis an empirical parameter that partitions the absorbed light between each of the species(i)within a given layer based on the vertical structure of the layer as well as the extinction coefficients and leaf area.The λhjis used to account for the effects of horizontal heterogeneity within the jthcanopy layer.This replaces an empirical parameter(?,Appendix A)that was proposed by Duursma and M?kel?(2007)and used in the stand-level mixedspecies model by Forrester(2014).However,preliminary analyses indicated that ? can be difficult to calculate for some crown architectures and it is occasionally not calculable for others due to its asymptotic form.

λvican be calculated for any species and stand structure using Equation 3(Forrester 2014).

where hm,iis the mid-crown-height of species i,hmis the mid-crown-height of the layer in which species i belongs.For Equation 3 R2=0.75 and P＜0.0001.The mid-crown-height is the height of a point half way between tree height and the height to crown base(Hb), where the Hbis the height where the lowest live branch joins the stem.Thus,the mid-crown-height of a given species is calculated as(height-Hb)/2+Hb, and the mid-crown-height of a given layer is calculated using the height of the tallest species in that layer,and the minimum Hbof all species in that layer.The sum of λvfor all species in a given layer is then normalized to unity such that

Adjusted

The horizontal heterogeneity λhis calculated as a function of mean midday solar zenith angle,kH,the ratio of mean tree leaf area(LA)to mean tree crown surface area (SA)and vfrac,which is the sum of the crown volume (m3)of all crowns within a given layer within one hectare divided by the total volume of that layer(m3)over one hectare.The total volume of a layer is(height of the tallest species-minimum Hb)×100 m×100 m.

The vfracparameter is used to quantify canopy openness and therefore the potential inter-tree shading. Inter-tree shading increases as the mean midday zenith angle increases due to season and latitude.The zenith angle used to calculate λhis an adjusted value,zadj,because while mean midday solar zenith angle is a sineshaped function of Julian day,at latitudes of about＜23°, the sine shape is distorted when solar zenith angles decline to 0 and increase again,instead of continuing to decline below zero,which would maintain the sine shape (see Appendix B).The zadjallows that part of the curve to be negative(×-1).The λhis calculated using two empirical equations(5a and 5b)because zadjdoes not influence λhuntil it is greater than 30°.These equations are general and can be applied to any species and stand structure.

Table 3 Parameters used for the Maestra model in each of the five stands

The standard errors of parameter values in Equations 5a and 5b are in parentheses below the estimate and P＜0.0001.These equations were fitted using the monospecific data set described by Forrester(2014).Monocultures were used so that the effects of horizontal heterogeneity could be quantified in the absence of any vertical heterogeneity that occurs in mixtures due to inter-specific differences in the vertical distribution of LAD,even if the height and height to crown base does not vary between species.When Equations 5a and 5b are used for mixed-species layers the kH×LA/SAis a weighted mean of all species within the given layer.The mean is weighted by the contribution that each species makes to the sum of kH×L.

The monospecific data contained simulated monospecific stands with crown characteristics that varied in terms of LAD(0.27 to 2.62 m2m-3),monthly kH(0.04 to 1.48),live-crown lengths(2 to 18 m),LA(10 to 390 m2),crown diameters(2.7 to 7.6 m),LA/SA(0.16 to 1.94),mean leaf angles(20°to 70°from horizontal),crown shapes(cones,ellipses,and half-ellipses),and densities ranged from 10 to 5300 trees per ha.About 548 stands were replicated at five latitudes between 0 and 65°and f was calculated for each month using Maestra.The kHof these stands was calculated using simulations of homogeneous canopies consisting of box-shaped crowns.

Description of detailed tree-level light model

Maestra is a three dimensional tree-level model that calculates APAR by individual trees based on their crown architectures and any shading from neighbouring trees.To account for shading from neighbours the canopy is represented as an array of tree crowns(shaped as spheres, cones,ellipses,half-ellipses)whose positions are defined by x and y coordinates.Each crown is divided into horizontal layers and each layer is divided into several points, and for each,LAD,leaf angle distributions and leaf optical properties are used to calculate APAR.The penetration of radiation through the canopy is calculated using the radiative transfer model of Norman and Welles(1983).Transmission and absorption of diffuse radiation is modelled using the method of Norman(1979),while the nonintercepted radiation reaching a canopy point is calculated in the sun direction for direct beam according to the hourly zenith and azimuth angles of the sun and of various azimuth and zenith directions for the diffuse sky radiation.

Table 4 Homogeneous extinction coefficients(kH)for each species and typical extinction coefficients(k)for all monospecific stands

Maestra has been tested against field measurements in several stands,and performed well.For example,Maestra calculations of diffuse radiation were compared with field measurements in agroforestry plantations of widely-spaced Erythrina poeppigiana trees and a dense understorey of Coffea arabica plants(Charbonnier et al.2013).Transmittance through the overstorey layer as well as through both layers,was compared and the goodness of fit(R2)was＞0.75 in all cases.Maestra slightly underestimated values at the high end of the range of diffuse transmittance and Charbonnier et al.(2013)suggested that this may have resulted from errors relating to 1)field data collection,2)estimates of Maestra parameters,and 3)assumptions used by Maestra,such as the simplification of crown architecture descriptions.

Maestra was also tested in an experiment containing monocultures and mixtures of Eucalyptus grandis and Acacia mangium(Stand 2 in Table 1;le Maire et al.2013). Gap fractions were slightly overestimated by Maestra in the monocultures and underestimated in the mixtures, which resulted in underestimates of APAR of 3.4%for A.mangium monocultures,4.5%for E.grandis monocultures and overestimates of 4.6%for 1:1 mixtures.

Model runs and parameterization

Previous studies have already used Maestra to calculate APAR for four of the data sets used in this study,and in the other stand the crown architecture has been examined in detail(Table 1).Therefore,the same Maestra parameters were used in this study(Table 3).Input data to the Maestra model included individual tree x and y coordinates,and individual tree characteristics including crown diameter,LCL and LA.These were estimated for each tree or were predicted using site-and species-specific allometric equations described in the studies cited in Table 1.One exception was Stand 5,where the LAallometric equations were not specific to the sites. As APAR of trees within the plot could still be influenced by trees beyond the buffers,which were up to 10 m wide (Table 1),a 25-m wide buffer with mean stand characteristics was created around the outside of the plot to avoid potential estimation bias.This buffer had the mean characteristics of the plot.Daily PAR data was estimated from solar radiation measured at the specific sites or calculated from other climatic data collected from the sites or nearby.Maestra parameters that differed between stands and species are provided in Table 3.Other parameters were constant for all stands.These included the number of individual tree crown layers to integrate over:nolay(6), pplay(12),nzen(5),naz(11);the distribution of diffuse radiation incident from the sky(difsky;0);the number of time steps per day(khrsperday)was 24;crown shape was a half-ellipse(cshape=ELIP);there was only one age considered for the beta distribution of LAD;there was noclumping of foliage into shoots.Further details about the meaning and significance of these parameters/settings can be found in the Maestra manual(see Medlyn 2004).

Figure 2 The vertical canopy structure within a mixed-species forest containing four different species(Sp.)as considered by(a)the stand-level model based on Equation 2 in this study,and(b)the assumption made by some forest gap models that all of the leaf area is distributed at the top of the crown.

The PAR data used for the stand-level model were monthly totals of the daily data used for Maestra.The stand-level model was parameterised using the kH,mean height,mean LCL,mean crown diameter and mean LAfor all trees of a given species in each plot.The actual number of trees per ha was also used.That is,plot specific means were used.For each species,kHwas estimated using Maestra and Equation 7,where fHis the fraction of PAR that is absorbed by a homogeneous canopy

Figure 3 Comparison between the predicted APAR from Maestra and predictions by the stand-level model(Equation 2)for each species within a stand or the total APAR of all species.The data correspond to(a)Stand 1,(b)Stand 2,(c,d)Stand 3,(e,f)Stand 4,and (g,h)Stand 5 in Table 1.(f,h)shows the predictions for the same stands as(e,g)but where the leaf area of a given species is positioned at the top of the crown and therefore each species is considered to be in a separate canopy layer to the others(as in Figure 2b).The solid lines are 1:1 lines and the dashed lines are lines fitted to the data that pass through the origin.

The kHused for each species is shown in Table 4 and were calculated using Maestra.As explained above,kHiscalculated for a given species by simulating a stand with a homogeneous canopy such that all trees have the same height,live-crown length,LAD,leaf angle distribution, and leaf optical properties,and have box-shaped crowns that fit perfectly together(no space between crowns). The varying kHvalues for the Acacia hybrid(Stand 3)result from the different LAD and leaf angles that were found in each of the experiments.In the same stand the high kHfor H.odorata is related to the relatively high LAD and low to medium leaf inclination angle(Tables 1 and 3).

Some of the stands contained deciduous species.The months during which leaf fall and leaf development occurred were excluded from the analyses.These were March and October for L.formosana in Stand 4 and May and October for F.sylvatica in Stand 5.During the leafless period,the deciduous species were given LAvalues of zero,while all other species continued to absorb light.

Model comparisons and sensitivity analysis

The predictions of f from the stand-level model were compared with those from the Maestra model using criteria such as the relative average error(average bias, e%,Equation 8),the relative mean absolute error (MAE%,Equation 9),and the mean square error(MSE, Equation 10)(Janssen and Heuberger 1995;Vanclay and Skovsgaard 1997).

where O is the observed calculations from Maestra and P is the predicted values from the stand-level model, and~O and~P are the means.

The accuracy test described by Freese(1960)was used to test whether the calculations from Maestra and those from the stand-level model differed by＞10%with a 95% confidence limit(α=0.05).All statistical analyses were performed using R 3.0.2(R Core Team 2013).

The APAR predictions of the stand-level model based on Equation 2,were also compared with three other approaches.The first was to simply use Equations 7 and 11.

where ζ is the fraction of ground covered by the canopy.Both of these equations ignore any vertical heterogeneity within the canopy and inter-specific differencesin shading abilities.Equation 7 also ignores horizontal heterogeneity in terms of gaps between trees,while Equation 11 takes this into account by considering the fraction of the ground that is actually covered by the canopy of that species.The second approach was to use Equation 12,which ignores vertical or horizontal heterogeneity(in terms of canopy gaps)but considers the inter-specific differences in shading abilities.

Table 5 Statistical information that describes the relationship between the predicted f from Maestra and that predicted by the stand-level models using Equation 2 or Equation 12

In Equation 12,instead of using the empirical parameter λvto partition the PAR absorbed by a given canopy layer to each species within that layer,the contribution made by each of n species is weighted based on their L and k(Rimmington 1984;Sinoquet and Bonhomme 1991;Sinoquet et al.2000).The third approach was taken from some forest gap models,which assume that all of the leaf area of a species(or age-class)is at the top of its crown(Figure 2b;Bugmann 2001),as opposed to the stand-level model of this study where the leaf area is distributed between the top and bottom of the crowns(Figure 2a;Equation 2).A potential problem with this gap model approach is shown in Figure 2.For example,species 1 is taller than species 2,but species 2 probably shades species 1 a lot more than vice versa.In the five stands in this study,all species had different heights and so for this third approach,each species was assumed to occupy its owncanopy layer,even if it was only a few centimetres taller than another species,and Equation 2 was used without the λvparameter.

Figure 4 A comparison between the predicted f from the stand-level model based on Equation 2 and the earlier version of this stand-level model described by Forrester(2014),which required an empirical parameter that is difficult to calculate for some crown architectures.The data correspond to(a)Stand 1,(b)Stand 2,(c)Stand 3,(d)Stand 4,and(e)Stand 5 in Table 1.The solid lines are 1:1 lines and the dashed lines are lines fitted to the data that pass through the origin.The relationships did not differ depending on whether the mixtures or monocultures were used,so for stands containing mixtures(b-e)the figure only contains data from the mixtures.

In all cases kHwas used instead of k because k could not be determined for all species using this data set (monocultures are required)and it was assumed that kHwas a good approximation of k.It is also problematic to determine a k for each species because it varies with stand density and L,unlike kH.To increase the range of stand structures and inter-specific differences that were examined the comparisons for the three different approaches were made using data from Stands 4 and 5,or using the simulated dataset developed by Forrester(2014) that contained 548 monocultures and 495 mixtures replicated at 3-5 latitudes between 0 and 65°.In these stands crown characteristics varied in terms of LAD(0.27 to 2.62 m2m-3),monthly kH(0.04 to 1.60),live-crown lengths(2 to 18 m),LA(10 to 390 m2),crown diameters(2.7 to 8.0 m),LA/SA(0.16 to 2.05),mean leaf angles(20°to 70°from horizontal),crown shapes(cones, ellipses,and half-ellipses),and densities ranged from 10 to 5300 trees per ha.

Finally,a sensitivity analysis was used to examine how much the APAR predictions using Equation 2 change in response to a 10%change in kH,LA,crown diameter and live-crown length.

Results and discussion

The stand-level light absorption model provided similar predictions of f or APAR to the much more detailed tree-level model(Maestra)in four of the five stands (Figure 3,Table 5).In contrast,f or APAR was not predicted well in Stand 2 where the crown architectures (e.g.kH)were more variable and it appears unlikely that ficould be predicted in such stands using tree-or stand-level models that do not allow different trees of the same species to have different kHor other architectural variables such as LAD.

The stand-level model based on Equation 2 and the model from which it was developed described in Forrester (2014)both consider the horizontal and vertical heterogeneity within the canopy,and gave very similar predictions of fi(Figure 4).In terms of horizontal heterogeneity,as L declines and canopy openness increases,the total stand APAR declines but the APAR per tree increases because there is less inter-tree shading.The effect of openness or canopy gaps was taken into account by considering the zenith angle of the sun,canopy openness and crown architecture(Equation 5),which are all summarised using the λhparameter(Equation 2).In contrast, the openness could not be sufficiently accounted for by only considering the reduction in L(e.g.Equation 1 or 7; Figures 5a and 6a)or the fraction of ground area covered by the canopy by using Equation 11(Figure 5b). Equation 11 was also applied to mixtures and resulted in very biased predictions of APAR that were similar to those for monocultures in Figure 5b(data not shown).

Figure 5 Comparison between the predicted f from Maestra and that predicted by(a)Equation 7,(b)a variation of that equation,Equation 11,and(c)the predicted f from Equation 2. The solid lines are 1:1 lines and the dashed lines are lines fitted to the data that pass through the origin.The data used were the monoculture simulations used to fit Equations 5a,b.

Figure 6 Comparison between the predicted APAR by Maestra and predictions using(a)Equation 7(or 1)in Stand 5,(b)Equation 12 in Stand 4,(c)Equation 12 at the total stand level in Stand 5 and(d)Equation 12 at the species level in Stand 5.The solid lines are 1:1 lines and the dashed lines are lines fitted to the data that pass through the origin.

For some species,the use of Equation 12 to predict APAR was as good as,or even better,than Equation 2 (Table 5),even though it does not consider the vertical or horizontal heterogeneity within the canopy layers. Nevertheless,Equation 2 outperformed Equation 12 more often than not,and APAR predictions from Equation 12 were nearly always more variable,as indicated by the higher MSE in Table 5,Figure 7 and also Figure 3e compared with Figure 6b.

Vertical heterogeneity is initially considered in the model by dividing the canopy into layers,each of which contains species with crowns that overlap vertically(Figure 1).Then the vertical heterogeneity within a layer is considered by using the constant λv(Equations 2 &3),which is determined by the height of the midpoint of the crowns compared with the midpoint of the given layer.The same constant is used to partition the total layer f to each of the species within that layer based on their kH×L.

Partitioning and consideration of the vertical heterogeneity is particularly important in mixtures because if this is ignored the model will give one species an unrealistic competitive advantage for light absorption and hence growth.Biased estimates can be obtained if this vertical heterogeneity is quantified by only using the mean height of each species as shown in Figure 2b(Sinoquet et al. 2000).This will exaggerate the asymmetry of competition for light such that if one tree or species is able to overtop another only by a few metres or even centimetres it gains a complete competitive advantage in terms of light because shorter trees can only intercept the light that is transmitted through the crowns of the taller trees.The bias that can result from assuming the canopy structure of Figure 2b is shown in Figure 3.Figure 3e and 3g show the APAR estimates when considering the relative positions of the whole crown length(i.e.Figure 2a),not just the height of each species.In contrast,Figure 3f& 3h assume that all of the leaf area of a given species is positioned at the top of the crown(i.e.Figure 2b).This results in an overestimate of APAR for the taller species,such as P.abies,and a corresponding underestimate for the shorter species,A.alba or F.sylvatica (Figure 3h).Even though there is a difference in height, the majority of the length of the crowns in Stands 4 and 5 overlap in a similar way to species 1 and species 3 in Figure 2(Table 1).If one species is only slightly shorter than another,but has a similar or shorter LCL, or a higher LAor extinction coefficient then it may even have a greater shading effect on the other species than vice versa.This may not be a problem when all species have very low kHvalues,but it could lead to biased estimates where the taller species within a given layer have higher leaf areas or higher extinction coefficients.Total layer APAR should be unaffected.

Figure 7 Comparison between the predicted fiof each species in synthesised stands with a single canopy layer that contains four species using Maestra and(a)predictions by the stand-level model based on Equation 2,and(b)the stand-level model described using Equation 12.The solid lines are 1:1 lines and the dashed lines are lines fitted to the data that pass through the origin.

The model is more sensitive to changes in kHand LAthan to changes in crown length and width(Figure 8). This sensitivity depends on the mean midday solar zenith angle and therefore varies between months,especially at higher latitudes(data not shown).For a given species and month,the sensitivity tends to increase with declining L because at higher stand densities(or L)there is already a lot of inter-tree shading and changes in crown parameters have a smaller influence on APAR.

In this study it was assumed that there is only minor variability in kHfor a given species,and a single kHvalue was used for each species.That is,it is assumed that leaf angle distributions or LAD or LA/SAdo not change enough with age,resource availability or climatic conditions to significantly influence kH,and instead these variables influence tree APAR via changes in crown length,crown diameter or LA.This is based on the results of other studies that have found that LAD is less sensitive to spacing and resource availability than LA,crown length or crown width (Forrester et al.2013;Dong 2014;Guisasola 2014)and that inter-specific variability in LAD tends to be greater than intra-specific variability(Ligot et al.2014a).However, there was considerable variability in kHin Stands 1 and 2(Table 4).

In Stand 1 the variability in kHwas probably an artifact of assuming a half-ellipsoidal crown shape,when the crowns were not quite half-ellipsoidal.By making this assumption,the LAD estimates for trees in unthinned plots were 0.67 m2m-3while those in thinned plots were 0.57 m2m-3,however more detailed calculations that did not assume any crown shape resulted in no significant thinning effect on LAD(Forrester et al.2013).It is not clear what caused the large variability in kHand LAD in Stand 2.This level of variability is also large in comparison to each of the other stands in this study.Stand 2 differs from the others in that it was examined for a whole rotation period,which spanned the age of 1 to 6 years.During this time the mean tree heights were between 2.5 to 24.4 m and the stand L increased to its peak at about age 3 years before stabilizing and finally declining towards the end of the rotation(le Maire et al.2013).During this rotation period there may be significant changes in leaf display,however,there was no clear relationship between age and kHfor either species.

It may be hard to predict APAR in stands containing species with such variable crowns using any tree-or standlevel model that does not allow individuals(or cohorts)of a given species to vary in terms of kH,LAD or variables that describe the within-crown architecture.It is worth noting that while there was clearly some variability in kHin the five stands examined in this study,this variability,in terms of standard deviations,was only about half that for typical(non-uniform canopy)extinction coefficients,k (Table 4).This is because the crowns of trees within a canopy do not fit together perfectly,which results in some empty space between individual tree crowns and this space will influence the estimates of k.The variability in kHis lower because in a homogeneous canopy there are no spaces and the kHis only influenced by the crown architecture and leaf characteristics.

It is important to note,that when calculating kHvalues using tree-level models,it is critical to use models that do not require all trees of a given species to have the same within-canopy k or LAD.Maestra is an example of a suitable model for this purpose because it does not assume any particular k and even though it allows the vertical and horizontal shapes of the LAD distribution to be defined,the absolute values of the LAD at any givenheight can vary between individual trees in response to their leaf area and crown sizes.

Figure 8 The relative change in APAR(or f),estimated using Equation 2,in response to 10%increases in(a)kHor LA(m2), (b)crown diameter(d,m)and(c)LCL(m)for all species in Stands 1,4 and 5.X-axis names indicate the Stand(S)number and the species name(first letter of genus followed by first two letters of species name).Includes data from each month and all plots in the relevant stands.

In most of the stands the trees were roughly evenly spaced.However,in Stand 3 there was a very uneven distribution such that the H.odorata trees were planted within a 22-m wide circular gap within an Acacia hybrid plantation(Dong 2014).This spatial distribution is taken into account by the detailed tree-level model but not by the stand-level model.As a result the stand-level model assumes that the Acacia hybrid trees are evenly spaced,which would increase APAR per tree because there is less shading from other Acacia hybrid trees.However,in reality Acacia hybrid trees are grouped (and are all outside the gap)so the stand-level model overestimates their APAR(Figure 3)because it underestimates the shading effect.Despite this,the overestimate was minor and if the model is to be used for these designs the bias could be determined and accounted for.

Another simplification in terms of stand structure is that the stand-level model ignores slopes and aspects, whereas the detailed tree-level model takes these into account.At Stands 4 and 5 there were some plots with very steep slopes of＞50%and sometimes＞100% (Table 1).However,the slope and aspect do not appear to have any influence on the f estimates,which is consistent with other studies that have suggested that crown architecture and canopy structure are more important determinants of f than slope(Courbaud et al.2003; Duursma et al.2003).This may be because if a tree gains an advantage by being upslope of another tree,it is disadvantaged to a similar degree by being downslope of other trees.

Conclusions

The stand-level model provided adequate predictions of APAR for a wide range of crown architectures,species compositions,species proportions and stand densities.It uses two empirical parameters to account for the vertical and horizontal heterogeneity within the canopy and the model could be improved further by improving this part of the model.It can be parameterised with speciesspecific information about mean crown length,diameter, height and leaf area as well as extinction coefficients for a simulated homogeneous canopy.These extinction coefficients can be estimated using the mean crown length,diameter,height and leaf area information to run a treelevel light model.It should be noted that the stand-level model did not perform well in one of the five stands examined because of the high variability in kHwithin that stand.The prediction of APAR in such circumstances may require(tree-or stand-level)models that allow different individuals or cohorts of the same species to have different kH,LAD or analogous parameters.The stand-level model could be used to examine light dynamics in complex canopies and in stand-level growth models.

Appendix A

The stand-level model of this study is based on a model described by Forrester(2014)where the fraction of PAR (f)that is intercepted by species i within a given layer is calculated using Equation A1.

where keff,iis an effective extinction coefficient derived by Duursma and M?kel?(2007)and is given as

where ? is an empirical parameter that depends on the mean solar zenith angle and therefore considers latitude and season.Duursma and M?kel?(2007)and Forrester (2014)calculated ? using empirical equations.However, the equation used by Duursma and M?kel?(2007)was developed for a single kHand while the equation used by Forrester(2014)included a range of kH,it approached an asymptote and preliminary analyses showed that it could not provide adequate estimates of ? for some crown architectures and was occasionally not calculable for others due to its asymptotic form.

Appendix B

Figure 9 The relationship between the day of the year(Julian day)and mean midday solar zenith angle(a,c)or the adjusted mean midday solar zenith angle(b,d)and at four different latitudes(lat.).Note that(a)and(b)are for the southern hemisphere and(c)and(d) are for the northern hemisphere.

The relationship between mean midday solar zenith angle and the day of the year for four different latitudes is shown in Figure 9a(southern hemisphere)and Figure 9c(northern hemisphere).At latitudes betweenabout-23°and 23°the sine-shaped function is distorted because the solar zenith angle is never negative. However,a smooth sine-shaped function of solar zenith angles is easier to use when predicting the effect of solar zenith angle on the horizontal heterogeneity parameter λhin Equation 5.Therefore an adjusted zenith angle zadjis used,which allows the curve to be negative(×-1)to retain the sine shape.To determine when the curve is to be multiplied by-1 for a given latitude,it is necessary to calculate the two days of the year where the Zenith angle=0.The first time during the year when Zenith angle=0 is given by Equation B1 and the second time during the year when Zenith angle=0 is given by Equation B2.Note that when developing Equation 5,where zadjis used,Forrester (2014)incorrectly calculated the solar zenith angles such that those for the northern hemisphere were used for the southern hemisphere and vice versa and this reversal was retained in this study in order to retain Equation 5.

Competing interests

The authors declare that they have no competing interests.

Authors’contributions

All authors carried out the field work and processed the data.DIF developed the stand-level model,ran the light absorption models,analysed the data and drafted the manuscript.All authors jointly revised the manuscript.All authors read and approved the final manuscript.

Acknowledgements

This study was part of the Lin2Value project(project number 033 L049) supported by the Federal Ministry of Education and Research(BMBF, Bundesministerium für Bildung und Forschung).Thanks to B.Medlyn and R. Duursma for maintaining the Maestra model.We are also grateful to those who established and maintain the experiments.Hancock Victorian Plantations Pty.Ltd.provided the site for Stand 1 and established the plantation.Stand 2 measurements and model parameterisation were carried out within the ANR(Agence Nationale de la Recherche)SYSTERRA program, ANR-2010-STRA-004(Intens&fix).N.V.Phan and N.V.Minh provided the site for Stand 3 and established the plantation.We also thank Director An’guo Fan,Mr.Bailing Ding and Miss Yue’e Chu from the Shitai Forest Bureau for useful advice and assistance when organising the fieldwork at Stand 4,and also Mr.Xiaozhu Wang and Mr.Hongxing Ruan for fieldwork support.Plots in Stand 5 are maintained and monitored by the Forest Research Institute of Baden-Württemberg and the Forest Districts where the field sites are located. Thanks to Dr.Lutz Fehrmann,who provided useful comments that improved the manuscript and also helped with data collection and data processing for Stand 4.Lastly,we would like to thank two anonymous reviewers who provided useful comments that improved manuscript.

Author details

1Chair of Silviculture,Faculty of Environment and Natural Resources,Freiburg University,Tennenbacherstr.4,Freiburg 79108,Germany.2Chair of Forest Inventory and Remote Sensing,Georg-August-Universit?t G?ttingen, Büsgenweg 5,G?ttingen 37077,Germany.3Forest Research Institute of Baden-Württemberg,Wonnhaldestr.4,Freiburg 79100,Germany.4Tasmanian Institute of Agriculture(TIA),Private Bag 98,Hobart,Tasmania 7001,Australia.5Vietnamese Academy of Forest Sciences,Silviculture Research Institute,Duc Thang,Bac Tu Liem,Hanoi,Vietnam.6CIRAD,UMR Eco&Sols,2 Place Viala, Montpellier 34060,France.

Received:13 May 2014 Accepted:22 August 2014

Abraha MG,Savage MJ(2010)Validation of a three-dimensional solar radiation interception model for tree crops.Agr Forest Meteorol 139:636-652

Bartelink HH(1998)Radiation interception by forest trees:a simulation study on effects of stand density and foliage clustering on absorption and transmission.Ecol Model 105:213-225

Binkley D,Campoe OC,Gspaltl M,Forrester DI(2013)Light absorption and use efficiency in forests:Why patterns differ for trees and forests.For Ecol Manage 288:5-13,doi:10.1016/j.foreco.2011.11.002

Brunner A(1998)A light model for spatially explicit forest stand models.For Ecol Manage 107:19-46

Bryars C,Maier C,Zhao D,Kane M,Borders B,Will R,Teskey R(2013)Fixed physiological parameters in the 3-PG model produced accurate estimates of loblolly pine growth on sites in different geographic regions.For Ecol Manage 289:501-514

Bugmann H(2001)A review of forest gap models.Clim Change 51:259-305

Canham C,Coates KD,Bartemucci P,Quaglia S(1999)Measurement and modeling of spatially explicit variation in light transmission through interior cedar-hemlock forests of British Columbia.Can J Forest Res 29:1775-1783

Charbonnier F,Gl M,Dreyer E,Casanoves F,Christina M,Dauzat J,Eitel JUH,Vaast P,Vierling LA,Roupsard O(2013)Competition for light in heterogeneous canopies:Application of MAESTRA to a coffee(Coffea arabica L.)agroforestry system.Agr Forest Meteorol 181:152-169

R Core Team(2013)R:A Language and Environment for Statistical Computing.R Foundation for Statistical Computing,Vienna,Austria,http://www.R-project.org/

Courbaud B,de Coligny F,Cordonnier T(2003)Simulating radiation distribution in a heterogeneous Norway spruce forest on a slope.Agr Forest Meteorol 116:1-18

Dong TL(2014)Using Acacia as a Nurse Crop for Re-Establishing Native-Tree Species Plantation on Degraded Lands in Vietnam.PhD Thesis.School of Land and Food.University of Tasmania,Hobart,p 216

Duursma RA,M?kel? A(2007)Summary models for light interception and lightuse efficiency of non-homogeneous canopies.Tree Physiol 27:859-870

Duursma RA,Marshall JD,Robinson AP(2003)Leaf area index inferred from solar beam transmission in mixed conifer forests on complex terrain.Agr Forest Meteorol 118:221-236

Forrester DI(2014)A stand-level light interception model for horizontally and vertically heterogeneous canopies.Ecol Model 276:14-22

Forrester DI,Albrecht AT(2014)Light absorption and light-use efficiency in mixtures of Abies alba and Picea abies along a productivity gradient.For Ecol Manage 328:94-102

Forrester DI,Collopy JJ,Beadle CL,Baker TG(2013)Effect of thinning,pruning and nitrogen fertiliser application on light interception and light-use efficiency in a young Eucalyptus nitens plantation.For Ecol Manage 288:21-30

Freese F(1960)Testing accuracy.For Sci 6(2):139-145

Gersonde R,Battles JJ,O’Hara KL(2004)Characterizing the light environment in Sierra Nevada mixed-conifer forests using a spatially explicit light model.Can J Forest Res 34:1332-1342

Grace JC,Jarvis PG,Norman JM(1987)Modelling the interception of solar radiant energy in intensively managed stands.New Zealand J For Sci 17:193-209

Guisasola R(2014)Allometric biomass equations and crown architecture in mixed-species forests of subtropical China.In:Masters Thesis.Chair of Silviculture. Albert-Ludwigs University of Freiburg,Freiburg,p 48

Janssen PHM,Heuberger PSC(1995)Calibration of process-oriented models.Ecol Model 83:55-66

le Maire G,Nouvellon Y,Christina M,Ponzoni FJ,Gon?alves JLM,Bouillet J-P, Laclau J-P(2013)Tree and stand light use efficiencies over a full rotation ofsingle-and mixed-species Eucalyptus grandis and Acacia mangium plantations. For Ecol Manage 288:31-42

Ligot G,Balandier P,Courbaud B,Claessens H(2014a)Forest radiative transfer models:which approach for which application?Can J Forest Res 44:385-397

Ligot G,Balandier P,Courbaud B,Jonard M,Kneeshaw D,Claessens H(2014b) Managing understory light to maintain a mixture of species with different shade tolerance.For Ecol Manage 327:189-200

Medlyn BE(2004)A MAESTRO Retrospective.In:Mencuccini M,Moncrieff J, McNaughton K,Grace J(eds)Forests at the Land-Atmosphere Interface.CABI Publishing,Wallingford UK,pp 105-122,http://maespa.github.io/

Monsi M,Saeki T(1953)über den Lichtfaktor in den Pflanzengesellschaften und seine Bedeutung für die Stoffproduktion.Jpn J Botany 14:22-52

Monsi M,Saeki T(2005)On the factor light in plant communities and its importance for matter production.Ann Bot 95:549-567

Norman JM(1979)Modeling the complete crop canopy.In:Barfield BJ,Gerber JF (eds)Modification of the Aerial Environment of Plants.American Society of Agricultural Engineers,Michigan,pp 249-277

Norman JM,Welles JM(1983)Radiative transfer in an array of canopies.Agro J 75:481-488

Oker-Blom P,Pukkala T,Kuuluvainen T(1989)Relationship between radiation interception and photosynthesis in forest canopies:effect of stand structure and latitude.Ecol Model 49:73-87

Parveaud C-E,Chopard JRM,Dauzat J,Courbaud BT,Auclair D(2008)Modelling foliage characteristics in 3D tree crowns:influence on light interception and leaf irradiance.Trees Struct Funct 22:87-104

Rimmington GM(1984)A model of the effect of interspecies competition for light on dry-matter production.Aust J Plant Phys 11:277-286

Sinoquet H,Bonhomme R(1991)A theoretical analysis of radiation interception in a two-species plant canopy.Mathematical Biosci 105:23-45

Sinoquet H,Rakocevic M,Varlet-Grancher C(2000)Comparison of models for daily light partitioning in multispecies canopies.Agr Forest Meteorol 101:251-263

Vanclay JK,Skovsgaard JP(1997)Evaluating forest growth models.Ecol Model 98:1-12

Wang YP,Jarvis PG(1990)Description and validation of an array model-MAESTRO. Agr Forest Meteorol 51:257-280

Cite this article as:Forrester et al.:Using a stand-level model to predict light absorption in stands with vertically and horizontally heterogeneous canopies.Forest Ecosystems 2014 1:17.

Submit your manuscript to a SpringerOpen journal and benef i t from:

? Convenient online submission

? Rigorous peer review

? Immediate publication on acceptance

? Open access: articles freely available online

? High visibility within the fi eld

? Retaining the copyright to your article

Submit your next manuscript at ? springeropen.com

10.1186/s40663-014-0017-0

*Correspondence:david.forrester@waldbau.uni-freiburg.de

1Chair of Silviculture,Faculty of Environment and Natural Resources,Freiburg University,Tennenbacherstr.4,Freiburg 79108,Germany

Full list of author information is available at the end of the article

?2014 Forrester et al.;licensee Springer.This is an Open Access article distributed under the terms of the Creative Commons Attribution License(http://creativecommons.org/licenses/by/4.0),which permits unrestricted use,distribution,and reproduction in any medium,provided the original work is properly credited.

Methods:The stand-level model was compared with a detailed tree-level model that has already been tested in mixed-species stands using empirical data.Both models were parameterised for five different forests,including a wide range of species compositions,species proportions,stand densities,crown architectures and canopy structures.

Results:The stand-level model performed well in all stands except in the stand where extinction coefficients were unusually variable and it appears unlikely that APAR could be predicted in such stands using(tree-or stand-level)models that do not allow individuals of a given species to have different extinction coefficients, leaf-area density or analogous parameters.

Conclusion:This model is parameterised with species-specific information about extinction coefficients and mean crown length,diameter,height and leaf area.It could be used to examine light dynamics in complex canopies and in stand-level growth models.