Can plasticity make spatial structure irrelevant in individual-tree models?

Abstract

Background

Distance-dependent individual-tree models have commonly been found to add little predictive power to that of distance-independent ones. One possible reason is plasticity, the ability of trees to lean and to alter crown and root development to better occupy available growing space. Being able to redeploy foliage (and roots) into canopy gaps and less contested areas can diminish the importance of stem ground locations.

Methods

Plasticity was simulated for 3 intensively measured forest stands, to see to what extent and under what conditions the allocation of resources (e.g., light) to the individual trees depended on their ground coordinates. The data came from 50 × 60 m stem-mapped plots in natural monospecific stands of jack pine, trembling aspen and black spruce from central Canada. Explicit perfect-plasticity equations were derived for tessellation-type models.

Results

Qualitatively similar simulation results were obtained under a variety of modelling assumptions. The effects of plasticity varied somewhat with stand uniformity and with assumed plasticity limits and other factors. Stand-level implications for canopy depth, distribution modelling and total productivity were examined.

Conclusions

Generally, under what seem like conservative maximum plasticity constraints, spatial structure accounted for less than 10% of the variance in resource allocation. The perfect-plasticity equations approximated well the simulation results from tessellation models, but not those from models with less extreme competition asymmetry. Whole-stand perfect plasticity approximations seem an attractive alternative to individual-tree models.

Background

Distance-dependent individual-tree growth models, also known as spatially explicit individual-based models, have a long history in forestry ([Dudek and Ek 1980]; [Newnham and Smith 1964]; [Reventlow 1879]; [Staebler 1951]), and more recently have received considerable attention in plant ecology ([Grimm 1999]; [Grimm and Railsback 2005]; [Wyszomirski 1983]). In them, stem base or breast-height coordinates are used to compute indices that reflect the competitive status of each tree and predict growth and mortality. Although such models are valuable research tools, it has been generally found that tree locations contribute little to predictive power, and in practical forest management they have been almost entirely replaced by non-spatial approaches ([Burkhart and Tomé 2012]; [Weiskittel et al. 2011]).

One possible reason for the insensitivity to stem location is plasticity, the ability of trees to lean and/or to adjust crown development so as to occupy less contested spaces ([Longuetaud et al. 2013]; [Muth and Bazzaz 2003]; [Rouvinen and Kuuluvainen 1997]; [Schröter et al. 2012]; [Seidel et al. 2011]; [Stoll and Schmid 1998]; [Umeki 1995]). Strigul et al. ([2008]) simulated forest stand development combining ideas from SORTIE ([Pacala et al. 1993]) and from the canopy tessellation methods of Mitchell ([1969], [1975]), but allowing for crown displacements as in [Umeki (1995]). They proposed a perfect plasticity approximation (PPA) as a limit where crowns are free to move so as to equalize competition intensity along their periphery. It was found that the simulation results were close to the PPA predictions, which do not depend on tree coordinates.

The main objective of this study was to complement the findings of [Strigul et al. (2008]), at the same time simplifying and generalizing aspects of that work. Their results are influenced by a number of specific assumptions and design choices, the importance of which are difficult to assess. Those include a space tessellation based on the intersection of crowns of a certain shape, allometric relationships linking tree dimensions to dbh, and particular growth functions. Here a general spatial individual-plant modelling framework implemented in the siplab R package was used, testing several alternative assumptions about neighbouring tree interactions ([García 2014]). Simulations were run for three data sets with different species, tree sizes, and spatial structures. Unessential complications were avoided by using single-species even-aged stands, focusing instead on key mechanisms. Some extensions to mixed-species are discussed elsewhere (Lee and García, manuscript in preparation).

Spatial individual-tree models predict growth and mortality rates as functions of the target tree size, and of a competition or resource capture index that encapsulates neighbourhood effects. Observed correlations, however, do not imply causality, trees that grow faster because of genetic, microsite or other factors will be larger. Extrapolation of the individual variability in past growth rates, represented by current size, may be another reason for the prediction efficiency of aspatial models ([García 2014]). For the purposes of evaluating the effects of plasticity and spatial structure, we side-stepped the circular size-growth ambiguity issues by limiting the analysis to an index of effective resource capture (assimilation index, for short). Clearly, dependence or independence between the assimilation index and spatial structure implies the same for growth and survival. In addition, results will not depend on specific assumptions about growth relationships.

The next section describes the test data, the spatial individual-based models, perfect plasticity approximations, and the analysis of simulation output. Simulation results follow, focusing on how much of the assimilation variability is explained only by tree size (ignoring tree coordinates), with and without plasticity. The article continues with stand-level implications of perfect plasticity useful for whole-stand modelling and other applications, and ends with a Results and discussion and Conclusions. The Additional files 1 and 2 include computer code and additional details.

Methods

Data

Simulations were based on 3 plots from the Boreal Ecosystem-Atmosphere Study (BOREAS, Rich and Fourni [1999]). They were established in unmanaged natural stands in Manitoba and Saskatchewan, central Canada. Coordinates and diameter at breast height (dbh) were measured for all trees taller than 2 m on a 50 × 60 m area. Heights and crown dimensions were measured on a subsample, and estimated for all trees by regression on dbh. The stands were single species, approximately even aged, and situated on flat terrain. Plot characteristics are shown in Table 1; the crown base heights and crown widths are subsample averages.

Table 1 Data statistics

[García (2006]) analyzed the same data and includes additional details. From two similar jack pine plots, only the one in the southern research site was used for this study. The data sets are included and documented in the siplab package.

The jack pine attributes are intermediate between those of the other two plots. The aspen trees are larger, and their density somewhat lower. The spruce stand is much denser, with smaller trees and an irregular spatial pattern. Spatial distributions are shown in Figure 1 (top row).

Figure 1

Tree spatial distributions for the 3 sample plots. Top: actual stem-base locations. Bottom: simulated crown displacements (simulation runs with a=1, α→∞, and distance-decreasing use efficiency). Circle diameters are proportional to tree height and zone of influence. Colour or shading follows influence function values. Trees in the 5 m border were excluded from displacements and from the analysis.

Models and simulation

Following [Strigul et al. (2008]), the model is most easily visualized through the physical crown space interactions of Mitchell’s TASS model ([Mitchell 1975]); generalizations are introduced later. In TASS, trees have a radially symmetrical potential crown shape with lateral crown expansion stopping at the points of contact, tessellating the plane on a horizontal projection (Figure 2a). The shapes move upward with height growth, modifying the tessellation, and possibly over-topping and eventually causing the death of the smaller trees. Light extinction produces a constant depth of live foliage, so that light interception and growth are essentially proportional to the horizontal area occupied by each tree.

Figure 2

Plasticity mechanisms that tend to equalize the competition pressure from the tree neighbours. The curves represent potential crown profiles, or more generally, shading potential or competitive strength (influence function). Dashed lines join neighbour contact points. (a) No plasticity. (b) Leaning. (c) Differential branch growth. A combination of these mechanisms can be expected, in addition to a redistribution of foliage density. Produced with Asymptote (http://asymptote.org) using a=1.5 in Eq. (1).

Influence functions

It is often observed that crowns do not interlock, especially at higher latitudes, and a direct interpretation of the models of [Mitchell (1975]) and [Strigul et al. (2008]) is then unrealistic ([Fish et al. 2006]; [Goudie et al. 2009]). However, rather than as a crown surface, the shapes can be seen in a more abstract way as a shading potential or competition intensity function, possibly representing both above-ground and below-ground processes, which we call an influence function. Light arrives at various angles, especially the important diffuse light radiation, so that the influence function is likely to extend somewhat beyond the physical crown limits.

Regardless of interpretation, [Gates et al. (1979]) derived forms for a crown profile or influence function that ensure that the induced growing-space partition satisfies a number of reasonable properties. Their conditions, together with the assumption of shape preservation by upward movement through height growth (“gnomonic scaling”), imply that the surface height must follow the equation

$$z=H-b{R}^a,$$

(1)

where H is tree height, R is horizontal distance, and a and b are positive parameters ([García 2014]). Only positive values are used, otherwise the influence is taken as 0. The circle z=0 defines the tree zone of influence (ZOI).

The simulations used a=1, which gives a cone, and a=2 that corresponds to a paraboloid of revolution. The pointed convex crown shapes used by [Mitchell (1975]) and [Strigul et al. (2008]) are intermediate between these two. Figure 2 was drawn using Eq. (1) with a=1.5.

The parameter b determines the shape slenderness, the ZOI extent, and the height of the function intersections between competing trees. The choice of values for the simulations is discussed later in Section ‘Simulation parameters’.

Allotment

In TASS and in [Strigul et al. (2008]), the horizontal space is subdivided on an exclusive basis, with the tree having the largest influence function value taking all the resource (e.g., light) available at each point. Competition is completely asymmetric. A less extreme alternative is to assume that the resource is somehow shared among trees where their ZOIs overlap. Siplab implements a general allotment rule where at each point (or pixel) a tree with influence function value z _i captures a proportion

$$\frac{z_i^{\alpha }}{\sum \underset{j}{\overset{\alpha }{z}}}.$$

(2)

The sum is over all the trees (or over all trees with positive influence at that point), and the parameter α is a measure of local competition asymmetry. With α=1 resource capture is directly proportional to z. In the limit α→0 capture is fully symmetric, the same for all competing trees independently of z. For α→∞ one has the TASS tessellation where the largest z takes all. Simulations were run with α→∞ and with α=1.

Efficiency

As mentioned before, the area allocated to a tree can be used as a summary of the effect of the neighbours on its development. For a given α, the values of (2) are spatially integrated, assuming a uniform resource distribution with one unit per unit area. Siplab discretizes these calculations, a 10-cm square pixel was used. More generally, the integration can be weighted by an efficiency function, to produce an effective resource capture or assimilation index reflecting contributions that diminish with distance from the tree location. An efficiency function of the same form as the influence function was used, scaled by its value at the origin:

$$1-b{R}^a/H.$$

(3)

This efficiency is 1 at the tree location, and decreases to 0 at the edge of the ZOI. The tessellations of [Mitchell (1975]) and [Strigul et al. (2008]) correspond to the special case α→∞ with a flat efficiency function.

Plasticity

With plasticity, phototropism induces a displacement toward areas where more light is available. Leaning of the stem in the direction of canopy gaps is common (Figure 2b). Differential branch growth (Figure 2c) and redistribution of foliage density can also be important. In most cases a combination of all these mechanisms is likely. Either way, the result is a more even resource allocation and less dependence on the basal stem locations. Note that plasticity makes the height of the crown contact or influence function intersection points more uniform. Below ground, roots can follow similar asymmetric patterns ([Brisson and Reynolds 1994]).

To simulate plasticity, tree coordinates were iteratively displaced to the centroid of the tree efficiency-weighted pixel resource captures. This tends to equalize competitive pressure on opposite sides, in the spirit of [Umeki (1995]) and [Strigul et al. (2008]). Iterations terminated when all coordinates changed by less than 5 cm. To prevent unlimited drifting, limits on the maximum displacement from the original tree position were enforced (Section ‘Simulation parameters’). Distortions due to the absence of competitors beyond the plot boundaries were limited by excluding a 5 m border from coordinate changes and from the results (Figure 1).

This approach does not distort the influence profiles as in Figure 2c. However, it can be seen that only the vertical change in cross-sectional area near the contact height is important. Another simplification is that, as in [Strigul et al. (2008]), crowns are displaced but the cross-sections remain circular. In reality, in many tree species horizontal crown shape distortion can add significantly to the spatial regularization of the canopy, although crown displacements still seem to be the main factor ([Longuetaud et al. 2013]).

Perfect plasticity

In the perfect plasticity approximation (PPA) of [Strigul et al. (2008]), plasticity causes all the crown contact points to be at a common height z^∗ (Figure 2). Consider a general influence function, where the horizontal cross-sectional area for tree i is some function f _i of the distance from the top. Then, for α→∞, the area captured by the tree is

$$A_i=f_i(H_i-z^{\ast }).$$

(4)

Assuming full canopy closure (z^∗>0), these areas must add up to the total stand area, so that the mean area is

$$\overline{A}=\overline{f_i(H_i-z^{\ast })}=1/N,$$

(5)

where N is the number of trees per unit area. This equation determines the contact level z^∗; a numerical solution is generally necessary.

In the models of [Strigul et al. (2008]) all tree dimensions for a tree species are fixed functions of dbh, and f _i varies with species and dbh. With gnomonic scaling f _i is size-invariant, and for Eq. (1) $A_i=\pi R_i^2$ is

$$A_i=\frac{\pi }{b^{2/a}}{(H_i-z^{\ast })}^{2/a}$$

(6)

if H _i ≥z^∗, otherwise A _i =0. The parameters can be species-dependent.

From equations (5) and (6),

$$\overline{{(H_i-z^{\ast })}^{2/a}}=\frac{b^{2/a}}{\mathrm{\pi N}},$$

(7)

provided that all trees are taller than z^∗. Explicit expressions for z^∗ can be obtained for a=1 and a=2. In the paraboloid a=2, the tree area A _i is a linear function of H _i , and eq. (7) reduces to

$$z^{\ast }=\overline{H}-\frac{b}{\mathrm{\pi N}}.$$

(8)

In the cone a=1 the function is quadratic, and

$$\overline{{(H_i-z^{\ast })}^2}=\overline{{\left[\right(H_i-\overline{H})+(\overline{H}-z^{\ast }\left)\right]}^2}={\sigma }^2+{(\overline{H}-z^{\ast })}^2$$

(9)

where σ² is the height variance, giving

$$z^{\ast }=\overline{H}-\sqrt{\frac{b^2}{\mathrm{\pi N}}-{\sigma }^2}.$$

(10)

With the weighting of eq. (3), the effective resource capture or assimilation index A i′ can be obtained by integration over annular differentials of area 2π r dr:

$$A_i^{\prime }=2\pi \underset{0}{\overset{R_i}{\int }}(1-b{r}^a/H_i)r\mathrm{d}r,$$

(11)

giving

$$A_i^{\prime }=A_i\left[1-\frac{2}{a+2}\left(1-\frac{z^{\ast }}{H_i}\right)\right].$$

(12)

Some additional relationships are derived in Section ‘Stand-level implications’.

Finding an explicit assimilation PPA for non-tessellation models (α<∞) seems more complicated. In any case, resource capture under perfect plasticity, and the consequent predicted growth and mortality, depend on tree size and stand density but not on spatial coordinates.

Simulation parameters

It remains to choose values for the parameter b in eq. (1), and limits for the ZOI displacements.

The stands have closed canopies, signalled by a rising canopy base. Influence intersection heights should therefore lie mostly above the average green crown level. It seemed reasonable to choose b so that the PPA intersections are about 2 m above the crown base, for a foliage depth of approximately 2 m ([Mitchell 1975]). Values thus obtained from equations (8)–(9) and Table 1 are shown in Table 2. Other values of b gave qualitatively similar simulation results.

Table 2 Values of the parameter b used in the simulations

With regards to displacement limits, the plasticity literature usually reports means for the horizontal distances between crown centroid and stem base, often as relative displacements (displacement divided by mean crown radius), but maximum values are less common. [Muth and Bazzaz (2003]) show relative displacements less than 1 for mixed hardwoods. In old-growth European beech, [Schröter et al. (2012]) found a maximum displacement of 6.26 m. [Longuetaud et al. (2013]) gave relative displacements of up to 7.68 in mixed broadleaves. In Scots pine, [Vacchiano et al. (2011]) calculated stand means between 1.0 and 3.9 m across 4 sites in the Alps. Figure five of [Gatziolis et al. (2010]) shows horizontal deviations between surveyed stem base and LiDAR-assessed tree top exceeding 10 m in both conifers and hardwoods, although some of that may be due to measurement error.

Besides crown displacement, foliage distribution also affects the influence function, and the contribution of crown shape distortion ([Longuetaud et al. 2013]) is ignored in the model. Based on this information, for the main results the algorithm total ZOI displacement was limited by a seemingly conservative upper bound of 3 m in all cases. To assess the effects of this parameter, also partial results with more restrictive bounds of 1.5 m for pine and 1 m for spruce will be shown (the spruce narrower crowns and perhaps their closer spacing might justify the smaller bound; aspen exhibits more stem leaning than conifers).

Analysis

We are interested in to what extent allowing for plasticity in the simulations diminishes the effect of spatial structure. In other words, how good is a perfect plasticity approximation, which assumes that assimilation indices depend only on tree size. To that effect, assimilation indices were analyzed both in the absence of plasticity, i.e., with the original tree coordinates, and after convergence of the ZOI displacement algorithm.

As a direct visual evaluation, the spread in scattergrams of assimilation index over tree height was examined, for the various model variants and data sets. The simulated values were compared to the theoretical tessellation PPA relationships of Section ‘Perfect plasticity’. As a numerical summary, the r-squared from a quadratic polynomial regression indicated the proportion of variance accounted for by tree size alone.

Also of interest are the differences in total resource capture (sum of the assimilation indices). Values per unit area and tree averages were calculated.

In a tessellation, some smaller trees may be completely over-topped, receiving no resources. The centroid-driven ZOI displacement algorithm ignores these trees, affecting total resource utilization (this is not the case with α<∞). The proportion of such trees is shown in the results. The proportion of trees where the displacement was constrained by the bounds of Section ‘Simulation parameters’ was also computed. This indicates to what extent approaching the PPA was limited by those assumptions.

The trend toward a more regular spatial distribution under plasticity was assessed with the Clark and Evans aggregation index, calculated with spatstat ([Baddeley and Turner 2005]) using the default Donnelly edge correction ([Rouvinen and Kuuluvainen 1997]; [Schröter et al. 2012]). The index is 1 for a “random” (Poisson) pattern, values less than 1 indicate clustering, while more uniform spacings produce indices greater than 1.

Results and discussion

Simulation

Simulation results are shown graphically in Figures 3, 4, 5, and numerical summaries are given in Table 3. Note that the simple iterative algorithm generally progresses toward more balanced competition, but there is no guaranteed convergence to any sort of optimum; anomalies can occur.

Figure 3

Simulation results for jack pine. Without plasticity, much of the assimilation variability is due to the spatial structure of tree locations. Plasticity reduces the effect of tree coordinates, so that the variability is largely explained by tree height. Curves are the tessellation (α→∞) assimilation PPA. It is seen that the approximation is satisfactory in that case, but clearly it is not appropriate for α=1.

Figure 4

Simulation results for trembling aspen (see the legend of Figure3). Outliers correspond to trees close to a large gap.

Figure 5

Simulation results for black spruce (see the legend of Figure3). In this denser and more heterogeneous stand the effect of plasticity is weaker.

Table 3 Simulation results, maximum displacement 3m

Simulated plasticity resulted in good relationships between assimilation indices and tree size, in most instances explaining over 90% of the variance as indicated by the regression R². That is, less than 10% of the variability can be attributed to spatial structure. The assimilation PPA works well for the tessellation models on which it is based (α→∞), but the relationships for α=1 are substantially different.

The Clark-Evans aggregation index indicates clustering of the stem coordinates. Plasticity moves the ZOIs into less occupied areas, tending to a more regular pattern, See also Figure 1. For the same reason, the total resource capture (or mean assimilation) increases, although perhaps not as much as might have been expected. The outliers with high assimilation indices in the aspen graphs of Figure 4 correspond to trees near the large gap on the bottom-right corner of Figure 1.

In the tessellation models, trees that are underneath the influence functions of other trees have zero assimilation, and do not participate in the displacements unless movement of their competitors changes their circumstances. That can be seen in the graphs, and in the rows labelled “% with assimilation 0” in Table 3. These trees are relatively more numerous in the more heterogeneous spruce stand. With α<∞ all trees receive some resources.

About half of the pine and spruce trees, and 1/4 of the aspens, are constrained by the 3 m displacement bound at the last iteration. This and the visually tighter relationships for the aspen suggest that relaxing this bound would decrease further the importance of spatiality. To examine the effects of more severe constraints, simulations were run with displacement limits of 1.5 m for pine and 1 m for spruce. Results are shown in Table 4 for distance-weighted efficiency. As expected, the R² decreased, although size still accounts for much of the variation.

Table 4 Plasticity simulations with displacement bounds of 1.5 m for pine and 1 m for spruce

Stand-level implications

A number of global relationships can be derived from the perfect plasticity approximation. These can be useful for canopy depth prediction, distribution modelling, and to develop whole-stand models for complex stands.

Canopy depth

With perfect plasticity, a tree crown length equals the distance from the top to the contact level z^∗, plus the depth d of the foliage layer:

$$L_i=H_i-z^{\ast }+d,$$

(13)

and z^∗ is related to stand density N as described in Section ‘Perfect plasticity’. The mean canopy depth is $\overline{L}=\overline{H}-z^{\ast }+d$.

From eq. (7), for the crown profiles or influence functions of eq. (1)

$$\overline{{(L_i-d)}^{2/a}}=\frac{b^{2/a}}{\mathrm{\pi N}}.$$

(14)

As a first approximation, ignoring the height variability gives

$$\overline{L}\approx \frac{b}{{\left(\mathrm{\pi N}\right)}^{a/2}}+d.$$

(15)

A second order approximation can be obtained applying the delta method (a 2nd order Taylor expansion around $\overline{L}$) on the left-hand side of (11), leading to

$$\overline{L}\approx \frac{b}{{[1+(2-a)C^2/a^2]}^{a/2}}\frac{1}{{\left(\mathrm{\pi N}\right)}^{a/2}}+d,$$

(16)

where $C=\sigma /\overline{L}$ is the coefficient of variation of the crown lengths. Either way, assuming that C is relatively stable, the mean canopy depth is approximately linear in 1/N^a/2, or in terms of the average spacing $S=1/\sqrt{N}$, approximately linear in S^a. In the special case a=1 eq. (9) gives an exact relationship that is only slightly nonlinear, while for a=2 equations (8) and (12) coincide.

[Brown (1962]) and [Valentine et al. (1994]) used a similar argument with equal-sized conical crowns, but assuming that the crown base coincided with z^∗, to conclude that $\overline{L}$ should be proportional to S. [Beekhuis (1965]) and [Valentine et al. (2013]) found that including an intercept in a linear regression gave better results for pine plantation data. That agrees with eq. (12) for a=1 (Beekhuis’ canopy depth was based on stand top height instead of mean height, so that his intercept includes the difference between those height measures).

A nonlinear least-squares regression with Beekhuis’ data, in metric units, gives $\overline{L}=3.67+3.77S^{0.742}$ or $\overline{L}=6.26+1.97S$, with the exponent not significantly different from 1 (p=0.456). A similar analysis with the grouped summaries from Table two of [Amateis and Burkhart (2012]), excluding the age 5 data which might not have reached canopy closure, gives $\overline{L}=2.30+0.854S^{1.27}$ or $\overline{L}=1.54+1.41S$, again with the exponent not significantly different from 1 (p=0.429). The linear regression is similar to those of [Valentine et al. (2013]) for individual trees from the same experiment. These observations support a value of a≈1, at least for conifers.

Distributions

Perfect plasticity makes growth and mortality independent of tree locations, they only depend on tree size and stand density. The state of a stand is then fully characterized by a size distribution and N. In a discrete approximation, individual-based aspatial models can then simulate the development of a finite sample of trees, as in the traditional distance-independent individual-tree growth and yield models (Burkhart and Tomé [2012]; Weiskittel et al. [2011]).

[Strigul et al. (2008]) calculated the evolution of a continuous size distribution with a partial differential equation (PDE) known as the McKendrick or von Foerster equation. It extends the classical Liouville equation to include mortality and recruitment ([Picard and Franc 2004]). Introducing stochastic elements would produce the Fokker-Planck (or Kolmogorov forward) PDE. Stochastic differential equations are a sometimes advantageous alternative representation.

Another possibility is to avoid PDEs by describing the state through the distribution moments instead of a continuous function. In general, this would require an infinite sequence of equations, one for the rate of change of each moment. The problem is solved by moment closure methods, which ignore higher moments or approximate them as functions of lower-order moments (e.g., Milner et al. [2011]; Murrell et al. [2004]).

Most models have used a simple scalar measure of tree size, usually dbh or basal area. A biologically meaningful tree description, however, requires at least two variables such as height and volume or biomass, a two-dimensional size vector ([García 2014]). In principle, the methods above can be applied to vectors, although that has rarely been done.

It should be recognized that the size-dependent growth or mortality relationships apply only within a stand. Large dominant trees tend to grow faster than the stand average, but at the landscape level large trees may correspond to older stands with lower growth rates. Hierarchical statistical methods should give better results than the common practice of fitting simple regressions to data gathered from different stands.

A spatially uniform resource availability is also assumed. This is appropriate for light, but nutrients and moisture vary and tend to be spatially correlated. In fact, for these data sets neighbouring tree sizes were found to be positively correlated, instead of the negative correlation predicted by spatial competition models, or of the independence assumed when using distributions ([García 2006]). In aspen, clonal vegetative propagation adds to the positive size correlations through spatial clustering of genotypes and growth rates.

As mentioned in the Introduction, a growth-size correlation may not be due to large size causing faster growth, but to the fact that faster-growing trees are larger. Specially with diameter-driven models and under management or natural disturbances, this statistical confounding can be problematic ([García 2014]).

Whole-stand

Production can be expected to be approximately proportional to effective resource capture, so that there is interest in the total assimilation or assimilation per unit area. With a flat efficiency function, assimilation corresponds to the area allocated to the tree, and with full canopy closure the total is independent of stand density and spatial structure (although it can change with stand age or height, and with site productivity). If use efficiency decreases with distance, however, stand assimilation and biomass or volume growth should decrease with increasing average spacing ([García 1990]).

Under the PPA with α→∞ and ignoring size variability, equations (10) and (7) give the following approximation for the assimilation per unit area:

$$N\overline{A^{\prime }}\approx 1-\frac{2b}{(a+2){\pi }^{a/2}}\frac{S^a}{\overline{H}},$$

(17)

with $S=1/\sqrt{N}$ being the average spacing. As in Section ‘Canopy depth’, the delta method can be used to produce a more accurate approximation of the same form. This is similar to moment closure (Section ‘Distributions’), retaining only the first moment. [García (1990]) found in radiata pine that gross volume increment per hectare, adjusted for site quality, decreased linearly with S, pointing again to a≈1.

If this is still approximately true with α<∞ is an open question.

Conclusions

These simulations are static, representing one point in time, instead of simulating stand development over long periods of time as in [Strigul et al. (2008]). However, dynamically the plasticity adjustments correspond to fast variables that act on much shorter time scales than tree growth, leading approximately to a dynamic equilibrium. Dispensing with the details of full growth, mortality and regeneration models allows for more general inferences, appropriate to the time horizons typically considered in growth and yield prediction. This might not be sufficient for studying succession mechanisms, including natural regeneration and changes in species composition over several centuries ([Strigul et al. 2008]).

Ignoring crown distortion underestimates the effectiveness of plasticity in reducing the effects of spatial pattern. With circular cross-sections it is not possible to achieve in three dimensions the idealized equalizing of inter-tree competition pressures suggested by Figure 2. That would require the crown or influence function width to vary with azimuth, leading to irregular cross-sections similar to those commonly observed in practice ([Longuetaud et al. 2013]). Therefore, one would expect tighter relationships than those shown in Figures 3, 4, 5, unless the 3 m bound on displacements turns out to be far too high.

Intuition and experimental results suggest that the simple centroid-chasing algorithm tends to the PPA, although that is not mathematically proven. Even so, the algorithm does not necessarily always converge to a “best” solution, and improvements might be possible. At a more fundamental level, it is not entirely obvious why the uniform z^∗ of [Strigul et al. (2008]) might be biologically optimal or desirable. In fact, it can be shown that deviations from it can increase the total effective resource capture by the stand, so it may not be strategically optimal at the population level in the sense of [Parker and Smith (1990]). Balancing competitive pressure on all sides seems however a reasonable tree-level tactic.

The simulations supported the hypothesis that plasticity causes assimilation indices, and hence growth and mortality, to be affected much less by tree location than by tree size. The degree of spatial independence varied among the test stands, presumably mainly in relation to the regularity of the spatial pattern. The plastic capabilities of the trees, represented by the bounds on maximum displacement, can be important in the more irregular stands. These conclusions were robust, mostly insensitive to assumptions about influence function shape, competition asymmetry, or resource use efficiency. The PPA assimilation-size relationships derived from tessellation models were not appropriate with less extreme asymmetry (α=1); it would be interesting to find correct explicit equations for that case.

Plasticity can also help to explain the insensitivity to spacing rectangularity in forest plantations. It has been observed that planting pattern has little or no effect on tree sizes, mortality, or yield (e.g., Amateis and Burkhart [2012]). A simulation of equal-sized trees on a 1:3 rectangular spacing (Figure 6) showed that plasticity produced tree growing areas indistinguishable from those arising from square spacing (a small random coordinate perturbation was used to get the algorithm started). With plasticity the total assimilation was 15% higher than without plasticity.

Figure 6

Rectangular planting. Left: no plasticity, crowns centred on the rectangular planting pattern. Right: plasticity simulated with the centroidchasing algorithm, crowns are displaced into less contested spaces. Shading indicates effective resource capture, line segments show the displacement of the crown centroid from the stem base locations. Trees at 1.5×4.5 m spacing, H _i =14 m, a=1, b=3.5, α→∞, efficiency weighting.

Fertility and other spatial correlations, and growth-size statistical confounding, are additional reasons for the generally small contribution of tree coordinates to growth predictions (Section ‘Distributions’). More research on these topics is needed.

It appears that for spatial patterns that are not too irregular, and for sufficiently plastic tree species, there is little to be gained by including spatial structure in growth and yield forecasting models. However, spatial modelling is likely to be still important in relation to severe disturbances and gap dynamics, and their effects on natural regeneration and succession. Difficulties with distance-independent and other distribution-based models (Section ‘Distributions’) make whole-stand modelling approaches attractive. Extensions of the PPA and similar approximations to multiple species, and to cohorts in uneven-aged forests, can produce whole-stand level equations suitable for complex stands.