Recovering the Metabolic, Self-Thinning, and Constant Final Yield Rules in Mono-Specific Stands

Competition among plants of the same species often results in power-law relations between measures of crowding, such as plant density, and average size, such as individual biomass. Yoda's self-thinning rule, the constant final yield rule, and metabolic scaling, all link individual plant biomass to plant density and are widely applied in crop, forest, and ecosystem management. These dictate how plant biomass increases with decreasing plant density following a given power-law exponent and a constant of proportionality. While the exponent has been proposed to be universal and thus independent of species, age, environmental, and edaphic conditions, different theoretical mechanisms yield absolute values ranging from less than 1 to nearly 2. Here, eight hypothetical mechanisms linking the exponent to constraints imposed on plant competition are featured and contrasted. Using dimensional considerations applied to plants growing isometrically, the predicted exponent is −3/2 (Yoda's rule). Other theories based on metabolic arguments and network transport predict an exponent of −4/3. These rules, which describe stand dynamics over time, differ from the “rule of constant final yield” that predicts an exponent of −1 between the initial planting density and the final yield attained across stands. The latter can be recovered from statistical arguments applied at the time scale in which the site carrying capacity is approached. Numerical models of plant competition produce plant biomass-density scaling relations with an exponent between −0.9 and −1.8 depending on the mechanism and strength of plant-plant interaction. These different mechanisms are framed here as a generic dynamical system describing the scaled-up carbon economy of all plants in an ecosystem subject to differing constraints. The implications of these mechanisms for forest management under a changing climate are discussed and recent research on the effects of changing aridity and site “quality” on self-thinning are highlighted.


INTRODUCTION
Power-law relations in ecology remain a subject of fascination and research interest given their simultaneous ubiquity and practical significance (Thompson, 1942;Vogel, 1988;Niklas, 1994;Brown and West, 2000;Farrior et al., 2016;West, 2017). That a complex phenomenon such as competition among plants may be succinctly summarized by a power-law expression between measures of plant size (e.g., biomass) and crowding (e.g., density) is arguably one of the most important examples prominently featured in ecological textbooks and research articles alike (Perry et al., 2008). In terrestrial ecology, two powerlaw relations have emerged between biomass and density, both developed for dense mono-specific stands (Shinozaki and Kira, 1956;Yoda, 1963): the self-thinning or Yoda's rule and the constant final yield rule. The usage of the term "rule" reflects the extensive experimental evidence supporting the universal character of the exponents of the size-density relations. The significance of these power-law relations to crop production, forestry and ecosystem management is rarely in dispute and has been reviewed elsewhere (Willey and Heath, 1969;Flewelling, 1977, 1979;White, 1981;Westoby, 1984;Peet and Christensen, 1987;Friedman, 2016). However, the ecological mechanisms responsible for their apparent universal character remains a subject of inquiry and debate since their inception in 1864 (Spencer, 1864). This debate frames the scope of this review.

The Self-Thinning Rule
Self-thinning, depicted in Figure 1A, describes a natural process in a single stand whereby the number of plants per unit area (p) decreases as average plant (or mean individual) aboveground weight (w) increases as time t progresses. That is, the relation between w(t) and p(t) is associated with transient dynamics initiated when p(t) begins to decline from its initial value with increasing time due to overcrowding. Self-thinning is, by definition, a process arising from space-filling where vegetation has covered the whole area under consideration. Self-thinning is presumed to be a process intrinsic to many managed and unmanaged terrestrial plant communities, whose composition and structure are influenced by competition for resources available proportionally to space-whether aboveground (e.g., photosynthetically active radiation) or belowground (e.g., water and nutrients) (Zhang et al., 2011;Hecht et al., 2016). Therefore, w-p temporal trajectories of self-thinning have considerable implications for forest management practices (Ge et al., 2017;Zhang et al., 2019). It is to be noted, however, that p(t) reductions due to ice storms, hurricanes, fires, diseases, or other disturbances are not considered in the w − p relations described by self-thinning.
As shown in Figure 1A, self-thinning does not describe the entire temporal trajectory of w − p and only "kicks in" when the w(t) is large enough for a given initial density to initiate intense resource competition (i.e., space-filling). At the early stage of stand establishment, w(t) may increase rapidly while the density remains at its maximum p(t) = p 0 , where p 0 = p(0) is the initial (or planting) density, until the space-filling requirement is reached.
Quantitative studies on a possible occurrence of universal w − p scaling emerged from data in the early 1930s in forestry (Reineke, 1933). Using measurements collected in many overstocked forests in California (USA), presumed to be experiencing self-thinning, p was empirically linked to the mean diameter at breast height D using where C R is a species specific constant. This equation states that plant density decreases as size increases across forest stands. Equation (1) is commonly referred to as Reineke's rule or Reineke's stand density index in forestry. Contrary to initial expectations by Reineke, the coefficient −1.605 appeared invariant across many species (12 out of 14 studied), age and environmental conditions. Thus, Reineke concluded that determining density of stocking in even-aged stands using Equation (1) has the advantages of freedom from correlation with age and site quality, and thus offers simplicity and general applicability.
Reineke's rule can be recast as a power-law of the form p = e C R D −1.605 . It is evident that when linking w to D using a power-law expression derived from allometry, Reineke's rule can also be formulated as a relation between plant biomass and crowding instead of plant size and crowding. Yoda (1963) and others later popularized similar power-law expressions extending the range downward from mature forest stands to seedlings of herbaceous plants, where both w(t) and p(t) are time-dependent (Figure 2). Equation 2 is hereafter referred to as Yoda's self-thinning rule when setting the exponent α = 3/2. The generality of this expression and the limited variability of the exponent α imply that annual and perennial crops, herbaceous plants, and trees are expected to respond to crowding in a surprisingly similar manner (White and Harper, 1970;Gorham, 1979;Antonovics and Levin, 1980). Moreover, density manipulation experiments seem to yield α = 3/2, yet C varies (Dean and Long, 1985). An α ≈ 3/2 was reported even in mixtures of Sinapis alba and Lepidium sativum, sown together at high densities, after having undergone collective self-thinning as described elsewhere (Bazzaz and Harper, 1976) (see Table S1). Nevertheless, some variability in α has been found (typically between 1 to 3/2), based on both theoretical arguments (reviewed in section 2) and empirical evidence, motivating this review.
FIGURE 1 | Conceptual representations of (A) self-thinning and (B,C) the constant final yield rule. (A) Self-thinning is initiated after crowding occurs, resulting in decreasing plant density p(t) with time t. (B) A high initial planting density p 0 is compensated for by poor growth conditions for each individual resulting in a small individual biomass w at harvest after an initial growth time period (B-3). (C) Conversely, a small p 0 allows for improved growth conditions leading to higher w at harvest after the same growth time period (C-3). The choice of a woody species in column (A) and a herbaceous species in columns (B,C) is to highlight the wide applicability of the self-thinning and constant final yield rules. Symbols: n p refers to the number of plants in a plot, A s is the plot ground area, D is the stem diameter, h is the plant height, p = n p /A s is the plant density, w is the mean individual biomass, and the product wp is the total biomass per unit ground area at time t. Subscripts refer to the different scenarios shown.
considered (Shinozaki and Kira, 1956;Holliday, 1960;Willey and Heath, 1969;Watkinson, 1980) 1 where y is a function of p 0 across stands at a fixed time t after sowing. In Equation (3), C f ,1 and C f ,2 are species-specific parameters that change with time (Kikuzawa, 1999). When intense resource competition has had long enough time to appreciably affect stand structure (t >> 0 where t = 0 is sowing time), the constant final yield rule dictates that C f ,1 p −1 0 becomes small compared to C f ,2 and y(p 0 |t) = y c (t) = C −1 f ,2 (t) (Figure 1). When density-driven mortality or self-thinning is absent (i.e., p(t) = p 0 ), this rule leads to an exponent −1 between w(p 0 |t) and p 0 , as it can be shown by multiplying both sides of Equation (3) by p 0 , and recalling that y(t) = w(t)p(t). This gives a relation between mean individual plant mass at steady-state and initial density, FIGURE 2 | Linkages between dynamic equations for individual plant biomass (w) and density (p; Equation 5) for the self-thinning rule (w = Cp −α ; left panels) and the constant final yield rule (y c = wp; right panels). (A,B) Illustration of the time evolution of w (green curves: mean biomass; green curves with shaded area: probability densities of individual biomass in size-structured models). (C,D) Illustration of the time evolution of p. Combining the dynamic equations in subplots (A-C) and (B-D), respectively, leads to a single equation describing trajectories in the w − p phase space (Equation 6): (E) when starting from a single initial density, different experiments and arguments show that α could vary from 1 to 3/2, or (F) multiple stands starting from different sowing densities achieve the same biomass per unit area y at a fixed point in time after sowing. Brackets right of the subplots: how constraints are imposed on the functions g 1(.) , g 2(.) , and g 3(.) of Equations (5) and (6); constraint numbering refers to the model categories discussed in the text (see also Figure 3). Gray circles and labels correspond to the scenarios depicted in Figure 1.
or w(p 0 |t) = C f ,2 (t) p −1 0 when C f ,1 is negligible compared to C f ,2 p 0 as noted earlier. Equation (4) describes how w varies with p 0 at a fixed time after sowing (comparing different stands) Frontiers in Forests and Global Change | www.frontiersin.org whereas Equation (2) describes temporal changes in w on the same stand. In contrast to Equation (4), Equation (3) also applies with the presence of mortality. In other words, the constant final yield rule and density-driven mortality, or self-thinning, are not mutually exclusive. Constant final yield has received experimental support in crops (Shinozaki and Kira, 1956;Holliday, 1960) and in a number of tree species in forest stands (Pacala and Weiner, 1991;Xue and Hagihara, 1998;Kikuzawa, 1999) (Table S1).

Interpreting Self-Thinning Exponents
A range of α values can be derived from contrasting theoretical arguments suggesting that α = 3/2 (for the self-thinning rule) and α = 1 (for constant, time-independent yield y = f (t)) are not universal values, but rather that α may be contextdependent. On the one hand, if indeed different constraints lead to specific exponents, it might be possible to infer which processes shape forest development from observed α values. For example, different types of plant-plant interactions lead to contrasting exponents in individual-based models (as shown later in section 2.8). In turn, knowledge of the constraints at play would allow predicting how α shifts in response to changing conditions (e.g., climate and land use). On the other hand, as shown in the following, several arguments lead to similar α values, making it difficult to establish which one is dominant (Table 1).
An obvious question to pursue is how the exponent α reflects constraints or mechanisms controlling competition among mono-specific plants. The multiple mechanisms covered in this review are summarized in Table 1 and relations between them are featured in Figure 3. That multiple mechanisms can result in the same α is not new (Pickard, 1983). What is original here is the establishment of a link between the constraint(s) on competition, the transient dynamics leading to power-law relations between w and p, and the numerical value of α. Extensions proposed here are distinguished from published arguments in subsections labeled "Extended Analysis" (some of which are expanded in the Supplementary Material). This effort is motivated both by the lack of a synthesis of the mathematical FIGURE 3 | The various mechanisms leading to the power-law relation w(t) = C[p(t)] −α reviewed here. The mechanisms are grouped based on how growth and mortality are treated and whether time, age class, and resource competition among individuals are explicitly considered. Note: in all temporally dynamic mechanisms, mortality is explicit. However, an explicit account of growth does not preclude an implicit account of mortality.
foundations of the self-thinning and constant final yield rules, and by the need to quantify how the scaling exponents may vary under future conditions, with implications for agricultural and forest management.

THEORY
Additional definitions are now introduced: l is a characteristic dimension of the plant, the one most sensitive to growth, V is the whole plant volume (product of projected crown area and height), ρ is the whole plant density, i.e., the individual plant mass over the entire individual plant volume, and s is the mean ground area covered by an individual plant or tree. From these definitions, w = ρV and crown radius is defined as r = √ s. The symbols and definitions are listed in Table 2.
Mechanistic studies, the subject of this review, typically begin with the carbon balance of the individual plant, where the carbon gains and costs as well as their constraints must be considered. Mortality and its associated effects on stand density must also be parameterized. Hence, these studies lead to a dynamical system coupling the individual scale (e.g., the single plant weight) with the plot or ecosystem scale (e.g., the density). This two-scale system may be represented by the general expression dw dt = g 1 (p, w); dp dt = g 2 (p, w), where g 1 (.) and g 2 (.) are functions that do not explicitly vary in time (i.e., the system of equations is autonomous) and must be determined from structural, hydraulic, energetic, and physiological constraints on competition. In this two-equation system, t can be eliminated to yield dw dp = g 1 (p, w) g 2 (p, w) = g 3 (p, w).

Mechanism 1: Dimensional Analysis and Allometric Constraints
The "Principle of Similitude" is a statement about the dimensional consistency of any mathematical expression relating   physical quantities to each other, such as mass, length, and time as described in Equation (5) (Spencer, 1864). It states, simply, that terms on both sides of an equation describing a physical state need to have the same dimension. Although evident, its consequences, first pointed out by Fourier (in 1822), allow for significant results to be derived (Lemons, 2018). The "Principle of Similitude" is now invoked in the context of w-p powerlaw relations.

Extended Analysis: Applying the Principle of Similitude
The dimensions needed to describe dw/dp are mass ( An expression for g 3 (.) is sought by inspecting a list of variables it might depend on such as w, p, ρ, s, r, and then combining these variables in groups that preserve the dimensions of dw/dp. The analysis is focused only on the period where self-thinning occurs, i.e., p(t) < p 0 , not the entire trajectory linking w to p at all times. Self-thinning only commences when the length scales associated with plant position in a stand (but not necessarily plant height) are related to p. This means that p must be retained to carry [L] into the equation for g 3 (.). If g 3 (.) is assumed to be independent of w, then the only mass unit available in this list of variables is ρ.
Matching the units on the most left-hand side to the most righthand side requires λ 2 = 1, and 2λ 1 + 3λ 2 = −2, or λ 1 = −5/2. Hence, g 3 (ρ, p) ∝ ρ p −5/2 (a power-law as expected from such dimensional analysis). Replacing the proportionality symbol with a dimensionless constant α ′ 1 results in Provided ρ is not impacted by p, although it can vary in time, the ordinary differential Equation (8) can be solved to yield where A d is an integration constant that must be determined from other considerations. This argument apparently recovers Yoda's rule without any explicit considerations to p declining with increasing t as necessary for self-thinning. However, assuming that the length scales are all related to p implicitly means that crowding has occurred. Likewise, if w replaces ρ, then dimensional considerations alone result in where α ′ 2 is a dimensionless constant. Solving this equation leads to w = A d p α ′ 2 , which again is a power-law. In this case, dimensional analysis fails to determine the numerical value of the exponent α ′ 2 , but it still predicts a power-law relation between w and p. Clearly, the choice of variables impacting g 3 (.) or the constraints imposed on it affects the value of the exponent α. For example, if the constraint is a constant total mass in time (i.e., dy/dt = 0), then it directly follows that That is, the constant total mass in time acts as a new constraint, allowing the determination of α ′ 2 = −1 and the achievement of a constant yield y c .
The main constraint on the outcome of competition may be a constant energy (or limiting resource) per unit area supplied by the environment R env . When this constant (in time) supply satisfies the ecosystem metabolic demands per unit area, R env = p R p , where R p is basal metabolic rate per individual plant, then This system yields p ∼ R −1 p . In metabolic theory, R p is uniquely determined by w and the temperature of the environment (Brown et al., 2004). Employing Kleiber's law (Kleiber, 1932) and inserting this result into p ∼ R −1 p , directly recovers the exponent α = 4/3 (i.e., the metabolic argument).

Allometry and Growth Habits as Constraints
Self-thinning is initiated when packing is achieved: the ground area is entirely covered by the plants or trees as discussed elsewhere (Miyanishi et al., 1979). It is emphasized that the probability that some local densities will achieve packing before the majority does is neglected, because p is a property that pertains to the whole ground area. If plant growth is threedimensional (i.e., height and crown diameter are increasing proportionally with increasing individual biomass) and ρ is constant, then V ∼ l 3 , w ∼ ρl 3 , s ∼ l 2 . Self-thinning occurs when s ∼ p −1 . When this point is reached, l ∼ s 1/2 ∼ p −1/2 . Thus, all length scales are now linked to p as foreshadowed earlier. Yoda's rule is directly recovered by noting that w ∼ ρl 3 ∼ ρp −3/2 , which is the key result in Equation (9). This argument assumes that the increment of plant size is isometric and proportional in all three dimensions (Miyanishi et al., 1979). Linking w to l, and all length scales to p, is akin to setting g 3 (.) of Equation (6) to uniquely depend on the density (p) over the course of self-thinning.
Other growth habits may now be analyzed, and two limiting cases are illustrated: prostrate ground cover plants (i.e., 2-D growth) to etiolated seedlings (i.e., 1-D growth). To place these growth habits in a framework that employs allometric scaling, it is assumed that V ∼ l m and s ∼ l n . Hence, w ∼ ρl m and at the incipient point of self-thinning, the condition s ∼ p −1 must be maintained. These assumptions lead to where m = 3 (i.e., 3-D growth) and n = 2 recovers Yoda's rule, and C is a proportionality constant. For prostrate ground cover plants, m = 2 and n = 2 resulting in w ∼ p −1 . For etiolated seedlings, the cross-sectional area is assumed constant and growth only occurs in the vertical (a race to harvest light). Hence, m = 1 and, with mean ground area covered by a plant being constant, requires s to be constant and thus n = 0. Hence, w ∼ p −1/0 yields infinite exponent, or stated differently, no self-thinning is to be expected (Miyanishi et al., 1979).

Extended Analysis: Reineke's vs. Yoda's Rules
Having covered the growth habits, it is now instructive to distinguish between a vertical dimension (canopy height h) and a horizontal dimension (canopy radius r), which allow recasting Equation (14) as where w = ρ(r 2 h). When Equation (15) is combined with an allometric expression linking h to r of the form h ∼ r β a , and when relating r to stem diameter D of the form r ∼ D β b (Enquist et al., 2009), the outcome is Comparing Equation (16) with Reineke's Equation (1) suggests that the exponents of the height-to-canopy radius (β a ) and canopy radius-to-stem diameter (β b ) allometries must be constrained by For m = 3 and n = 2 (i.e., Yoda's rule), β b (2 + β a ) ≈ 2.4. The immediate consequence of this combination is that V ∼ (r 2 h) ∼ D 2.4 . Conversion of D to a characteristic scale l ∼ r, β b and not V ∼ l 3 . Notably, scaling relations discussed elsewhere (Enquist et al., 2009) provide β a = 1.14 and β b = 0.684 so that β b (2 + β a ) ≈ 2.2, lower than the 2.4 value obtained above. Therefore, only for particular choices of the coefficient in Reineke's Equation (1) and of the scaling exponents β a and β b , can the relation V ∼ l 3 be recovered. In fact, this corroborates the point made by the incisive analysis of Weller (1987) that Yoda's exponent α = 3 2 may be the exception rather than the rule.
As the stand becomes crowded, more individuals are suppressed. Acclimation allows suppressed individuals to survive longer by decreasing the carbon investment in diameter relative to height and maintaining smaller crowns closer to the top of the canopy. These adjustments decrease β b of the entire stand. A reduction in β b is expected to reduce the slope of the p − w scaling (e.g., m/n in Equation 17). The reduction in crown size of suppressed individuals reduces the wind-induced drag force, allowing these trees to maintain structural integrity despite the lower taper.

Mechanism 2: Structural and Biomechanical Constraints
Relations between height and diameter can be derived to further constrain allometric scaling based on self-buckling or structural considerations. The key observation regarding self-buckling for trees is that the critical height for buckling (h crit ) plotted against tree base diameter (D) follows a power-law of the form h crit ∼ D 2/3 for nearly every American tree species (McMahon, 1973) (Table S1). This h − D scaling was consistent with the incipient point of buckling of tall columns (due to their own weight) given by the elastic buckling criterion (or Euler-Greenhill formula) where E is the modulus of elasticity, ρ W is the density of wood. If Equation (18) describes an allometric scaling law, then for 3-D growth, l ∼ D, V ∼ D 2 h crit ∼ D 8/3 , and s ∼ D 2 . Employing once more the allometric scaling framework presented previously, w ∼ (ρ)l m=8/3 , and s ∼ l n=2 . Inserting these estimates of m and n in Equation (14) result in (McMahon, 1973) This scaling is intermediate between a steady state biomass per individual (α = 1) and Yoda's rule (α = 3/2). Connections between the aforementioned scaling law in Equation (19) and metabolic arguments (i.e., Kleiber's law) have already been noted (McMahon, 1973). However, the scaling law in Equation (19) can also be derived without resorting to self-buckling, using a variant of the growth-hydraulic constraint (Niklas and Spatz, 2004), as well as metabolic constraints, as described later on. Additional implications of self-buckling are explored in the Supplementary Material.

Mechanism 3: Metabolic Limitations
Metabolic arguments and their connection to the exponent α have been popularized by the work of West, Enquist, and Brown (West et al., 1997;Enquist et al., 1998;Brown et al., 2004). To sustain a total biomass per unit ground area wp requires a rate of energy (or other limiting resource) supply from the environment per unit area of at least R E = pR p , where R p is the metabolic rate per plant (energy or resource use per time per plant). In living organisms, the basal metabolic rate R p varies with size and is given by Kleiber's law (Equation 13) (Kleiber, 1932;Banavar et al., 1999). Hence, R E = pR p ∼ pw 3/4 .

A Steady State Resource Balance
The case of a limiting essential resource is first considered. When the external environmental supply of this resource (= E env ) is balanced, or fully exploited, by the stand (or ecosystem) metabolic demand R E , then resulting in w ∼ p −4/3 when E env is set to a constant (e.g., a given annual shortwave radiation or precipitation rate). For all practical purposes, Equation (20) is an equilibrium argument (constant resource supply) with a constraint shaping g 1 (.) in Equation (5) at a given supply E env . It is also interesting to note that Equation (20) suggests a link between the constant C in Equation (2) and environmental conditions, as C ∼ E 4/3 env . The debate on the difference between the 4/3 and the 3/2 exponents are highlighted as they offer a new perspective on links between the scaling exponent α and the constraints. The α = 4/3 exponent has experimental support when average "mature" plant weight is plotted against p 0 for different species spanning nine orders of magnitude by weight (Enquist et al., 1998) and appears consistent with a number of spatially explicit simulation studies discussed elsewhere (Deng et al., 2012) (Table S1). Such an inter-species comparison, however, fundamentally differs from plotting w(t) against p(t) for a single stand across time (Yoda, 1963).

Extended Analysis: Constraints on the Trans-location Network Distribution
It has been argued that distributed trans-location networks evolved from a need for effective connectivity with increased size (i.e., analogous to the economy of scale in microeconomics). Distributed trans-location networks occur in biological systems (including respiratory networks) and in inanimate systems alike (e.g., river networks). The flow rate in an arbitrary translocation network can be derived as a function of the local connectivity as discussed elsewhere (Banavar et al., 1999). For the problem at hand, this trans-location network may represent the phloem, where metabolic products derived from photosynthesis (mainly carbohydrates) are being translocated from leaves, or the xylem, where water is transported to the leaves. In this network derivation, a moving fluid volume filling the network is assumed to be V f . The V f scales with the product of the number of links in the network and the distance between nodes. In a D i -dimensional space-filling network (i.e., a network that can deliver fluid to all the domain), the number of links is proportional to l D i n , where l n is a linear scale of the network. The distance among links is also proportional to l n . Hence, V f ∼ l D i n × l n = l D i +1 n , or l n ∼ . For a constant fluid density, w ∼ V f and R p ∼ l D i n . It directly follows that the metabolic rate for an individual is given by For R E = pR p = E env = constant, it follows that R p ∼ p −1 as before. As a result, Equation (21) is used to obtain For D i = 3, the 3/4 metabolic scaling exponent is recovered from Equation (21), and α = 4/3 is now recovered from Equation (22) in a manner that is compatible with Kleiber's law without resorting to critical height and self-buckling. Interestingly, the analysis here also suggests that Yoda's 3/2 scaling exponent is recovered for D i = 2 (i.e., R p ∼ w 2/3 ). A 2-D translocation network may be incompatible with Yoda's original assumption of proportional growth in all three dimensions. This incompatibility is one of the salient features of the aforementioned controversy surrounding the 4/3 vs. the 3/2 self-thinning scaling exponent.

Mechanism 4: Hydraulic Constraints on Growth
In addition to structural and energy supply constraints discussed as mechanisms 2 and 3, a hydraulic constraint can be formulated by imposing a steady-state transpiration rate from the roots to the leaves. This constraint may be viewed as a "networkon-network" supply constraint. There are three networks that must be coordinated: a root network that must harvest water and nutrients from the soil, a xylem network that must deliver water to leaves, and distributed end-nodes for water loss through leaves. It is assumed that these three networks are sufficiently coordinated so that no severe "bottleneck" in one network routinely impairs the function of the other two networks (Thompson and Katul, 2012;Huang et al., 2018). Based on this view, a simplified version of a growth-hydraulic constraint (Niklas and Spatz, 2004) is now reviewed. In this mechanism, it is assumed that the annual increment of dry matter per plant scales (i) linearly with standing leaf biomass w L that provides metabolic products and (ii) with w 3/4 as in Kleiber's law. Hence, equating these two assumptions results in where k 0 and k 1 denote allometric constants. With w defined by the sum of leaf, stem, and root mass (i.e., w = w L + w S + w R ) results in The hydraulic component of this argument is framed as follows (Niklas and Spatz, 2004): the amount of water absorbed by roots per unit time must pass through stems, experience a phase transition and then exit through the stomata distributed on leaf surfaces. Because this amount of water loss is conserved throughout the plant (i.e., no storage or capacitive effects on time scales commensurate with stand development), w L must scale with the hydraulically functional cross-sectional area of stems and roots (sapwood area). The key assumption is that the sapwood area is proportional to the square of the stem diameter (i.e., D 2 ). The assumption need not imply that the diameter of the water transporting vessels is proportional to D, but that D reflects the total number of vessels of fixed diameter. Viewed from this perspective, this assumption may also be interpreted as another expression of the so-called da Vinci rule, or the pipe flow model of Shinozaki (Shinozaki et al., 1964;Horn, 2000), and leads to w L = k 2 D 2 . Substituting w L in Equation 24 and rearranging the terms lead to Two additional assumptions are required (Niklas and Spatz, 2004): an allometric relation between root and stem biomass (w R = k 3 w S ) and a relation between stem biomass and stem volume (i.e., w S = k 4 (D 2 h)), where k 4 = ρ W , but the notation of Niklas and Spatz (2004) is maintained in the following. Hence, Equation (25) can be formulated as from which it follows that h ∝ D 2/3 and w S ∝ D 8/3 . Upon comparison with the Euler-Greenhill formula (Equation 18), the same h ∝ D 2/3 scaling has been recovered from metabolic and hydraulic constraints acting in concert (i.e., in coordination), not from mechanical limits on tree height, nor from energy supply by the environment. Combining these outcomes with w ∝ (D 2 h) = D 8/3 , s ∼ D 2 , and s ∼ p −1 (or D ∝ p −1/2 ) at the point where selfthinning commences, recovers the metabolic formulation w ∼ p −4/3 . Here, geometric packing (i.e., s ∼ p −1 ) leading to selfthinning is necessary to arrive at α = 4/3, which was not the case in the metabolic arguments.

Extended Analysis: An Alternative Hydraulic Link to Stem Diameter
The aforementioned arguments may be generalized to include other linkages between sapwood area and stem diameter. One such linkage is the so-called Hess-Murray law that predicts the optimal blood vessel tapering in cardiovascular systems. This linkage leads to w L ∝ D 3 (Murray, 1926;McCulloh et al., 2003) instead of D 2 . Starting again from Equation (24), the aforementioned argument leads to or h ∼ D 2/3 for small D only, not for any D as it was the case for the Da Vinci rule. For intermediate or large D, h ∼ g 11 D 2/3 −g 22 D (g 11 and g 22 are constants linked to k 0 , k 1 , k 2 , k 3 , k 4 ), which does not exhibit a unique exponent provided D < (g 11 /g 22 ) 3 . The connection between the da Vinci rule (along with the pipe flow model) and water transport has been the subject of debate outside the scope of the present work (Bohrer et al., 2005), with some arguing that the da Vinci rule is compatible with structural, rather than water transport theories (Eloy, 2011).

Mechanism 5: Spatial Averaging Arguments
This approach explicitly considers that stands generally comprise individuals of different sizes, even in even-aged mono-cultures, owing to small genetic variability as well as variations in site micro-environmental factors, impacting growth potential and access to resources. It is thus necessary to consider the effect of spatial averaging over individuals within the crop or stand area A s . By definition, p = n p /A s where n p is the number of individual plants within area A s . Also, the arithmetic mean weight of all individuals within A s is defined as where w i is the weight of each individual plant. Equation (28) can be rearranged to yield (Roderick and Barnes, 2004) It was suggested that over an extended life span, the total stand biomass dynamics eventually reaches a steady-state such as in the experiments of Shinozaki and Kira (1956) on soybean, a herbaceous species, where mortality was absent (Table S1). If such steady-state conditions are attained within a single stand, then where K c is a constant carrying capacity determined by the available resources supporting maximum biomass per unit area. Equation (30) implies that α = 1 because w = K c p −1 as long as K c is constant. The derivation of Equation (30) makes no assumption about p 0 , p(t) or w(t), or that y(t) follows logistic growth as in the competition-density effect (Shinozaki and Kira, 1956;Xue and Hagihara, 1998). Equation (30) was also shown to apply for a pine stand (Xue and Hagihara, 1998). For prostrate ground cover plant growth, the emergent scaling law was already shown to be w(t) ∼ p(t) −1 using an entirely different set of assumptions. Evidently, the α = 1 scaling exponent can be recovered from multiple mechanisms. It is demonstrated next that (i) w(t) ∼ p(t) −1 may still reflect a correct minimum exponent under weak self-thinning style competition and (ii) novel links can be established between the newly derived exponent here and other "conservative" ratios describing stand carbon dynamics.

Extended Analysis: Recovering the α = −1 Exponent From a Dynamical System
The previous argument can be extended by relaxing the assumption of steady state, showing that the same result is obtained in a more general case. As a point of departure from prior work (Roderick and Barnes, 2004), the α = 1 exponent is now analyzed using the framework of Equation (5). To facilitate this analysis, the total stand weight W T = n p w is assumed to vary logistically in time using (Verhulst, 1838) where r c is the intrinsic growth rate. This assumption has been used in the original work of Shinozaki and Kira (1956) at the individual level and generalized by others at the stand level (e.g., Xue and Hagihara, 1998). Such assumption is equivalent to prescribing g 1 (.) and g 2 (.) of Equation (5). Instead of analyzing the dynamics at an equilibrium point W T /A s = K c being constant, it is instructive to explore the transient dynamics where n p = pA s begins to decline in time. This type of competition is intended to resemble some but not all aspects of self-thinning (i.e., being a transient and operating when dn p /dt < 0) while maintaining a density-dependent logistic form for total biomass (instead of constant K c ) used by Shinozaki and Kira (1956). In particular, we ask under what conditions such a "stylized competition" remains compatible with scaling laws associated with a steady state yield or the self-thinning rule (or intermediates). A minimal model describing the n p decline is dn p dt = A s dp dt = A s g 2 (p) = −α m n p , where α m is a mortality inverse time constant. Equation (32) specifies the reduction in the number of plants through mortality as proportional to the number of plants n p thus making n p an exponential function of time. Again, viewed from the perspective of Equations (5), these approximations are equivalent to specifying g 2 (p) and g 1 (w, p) via Equation (31) when recalling that W T = n p w. By eliminating time t in Equations (31) and (32) (as before, to obtain Equation 6), an ordinary differential equation describing the variations of w with n p can be explicitly derived, The general solution of Equation (33) is given by where C s is an integration constant. Noting again that p −1 = A s /n p , Equation (34) can be expressed as For Y 0 p (r c /α m ) of the order of unity, no unique scaling exponent exists, although any power-law approximation to this solution must yield exponents exceeding unity in magnitude, which is the sought result. This finding offers an explanation as to why the exponent α varies between 1 and 2 across many data sets a priori conditioned on dn p /dt < 0 (i.e., when mortality begins to play a role).

Extended Analysis: the Effects of Invariant Carbon-Use Efficiency on Self-Thinning
The quantity r c /α m reflects the ratio of two time scales: one associated with net carbon gain of an individual plant (1/r c ) and another associated with its mortality (1/α m ). The time scale for carbon gain may be associated with the net primary productivity (NPP) of an individual plant, so that r c ∼ NPP/w (Thurner et al., 2016). In self-thinning stands where carbon loss in respiration is not compensated by photosynthesis in highly suppressed individuals (under light competition), it may be (simplistically) assumed that carbon starvation is the causal mechanism of mortality. Thus, the mortality time scale is associated with autotrophic respiration R A , so that α m ∼ R A /w and where NPP=CUE×GPP, GPP is the gross primary productivity, R A =GPP-NPP=(1-CUE) GPP, and CUE is the plant carbon use efficiency (0 < CUE < 1). Therefore, this link between r c /α m and CUE offers a new perspective about α and carbon use efficiency; i.e., α = 1 + (r c /α m ) = (1 − CUE) −1 . This estimate of α is expected to be an upper limit, because α m is likely to be underestimated when mortality time scale is estimated from R A . The value of plant CUE typically ranges between 0.4 and 0.8 depending on species, plant age, and growing conditions, with values even lower than 0.4 in mature trees and generally higher values in rapidly growing crop species (Manzoni et al., 2018). For an intermediate CUE=0.47 (typical in forests;Waring et al., 1998), large scaling exponents are obtained, w ∼ p −1.88 as already foreshadowed. For relatively inefficient plants with CUE=1/3, w ∼ p −3/2 , recovering Yoda's rule. The metabolic argument w ∼ p −4/3 can only be recovered for CUE=1/4. This argument is prone to large uncertainties due to both the qualitative link between parameters r c and α m , and CUE, and the uncertainties in CUE estimates. Nevertheless, plants that are more effective in converting resources to biomass are expected to exhibit steeper w −p scaling relations, a conjecture to be explored in future studies.
Up to this point, it was assumed that at the individual plant scale, the entire biomass captured in w is alive and contributes to respiration. However, for a preset total biomass, lower initial density may lead to greater live crown ratio at the incipient point of self-thinning. Hanging onto large branches at the bottom of long crowns contributes little to annual photosynthesis (Oren et al., 1986), but requires investment in maintaining active sapwood, cambium, and phloem. Thus, the initial planting density can play a role in determining the fraction of live to total biomass at the start of self-thinning. At that point, despite similarities in stand density, mean total individual tree biomass (Peet and Christensen, 1987), and leaf area (Dean and Long, 1985), stands characterized by individuals with a higher fraction of live to total biomass may exhibit higher whole-tree respiration rates per unit of leaf area and, therefore, reduced CUE and α.

Mechanism 6: Dynamical Systems Theories for Plant Carbon Balance
A number of approaches have been proposed that recover the self-thinning rule from a mechanistic representation of the plant carbon balance (Hozumi, 1977;Pickard, 1983;Hara, 1984;Perry, 1984;Pahor, 1985;Voit, 1988). Common to all these approaches is the so-called von Bertalanffy equation (von Bertalanffy, 1957;Perry, 1984) or a variant of it that applies to individual plants as discussed in Figure 3. Using the framework of Equation (5), this equation represents g 1 (w, p) as where a ag is the fraction of photosynthesis allocated to biomass, LAP is the leaf area of an individual plant, assumed to vary with w, P m is the photosynthetic rate per unit leaf area, varying with p (e.g.,due to light competition), and k m is the rate of maintenance respiration (whereas mortality is described by Equation 32). The overarching assumption of von Bertalanffy equation is that resource acquisition must traverse a limiting surface area (here LAP; scales as ∼ l 2 ) whereas respiratory and maintenance costs vary with plant size (or mass, w; scales as ∼ l 3 ). Variants to Equation (37) include complex expressions for photosynthetic gains, respiratory losses, connections between P m and p (such connections are the subject of spatially explicit models discussed later), and the partitioning of w into metabolically active and inactive parts. The goal of this section is not to review all of them but to offer links between the von Bertalanffy equation and the general framework set in Equation (5). Equation (37) is coupled to Equation (32) after eliminating time t and substituting n p = A s p to yield Mechanistic models link LAP to w using allometric rules and P m to p assuming that increases in p reduces the main resource driving photosynthesis such as photosynthetically active radiation (Perry, 1984). For example, LAP and P m may be expressed as Here, m a ≈ 0.81 and C f ≈ 0.011 when LAP is treated as all sided (in m 2 ) and w is expressed in grams (determined for a wide range of species), whereas B p = 4.61 and m b was varied as a control parameter (plausible values for most species in Perry, 1984). The representation in Equation (39) preserves the autonomous nature of Equation (38) thereby linking the phase-space of the p − w trajectories directly to model parameters. It also provides a complete description of g 3 (w, p) in Equation (6). However, a unique power-law solution of the form w = Cp −α is not apparent even though model calculations suggest an approximate power-law with exponent α = 1.0 − 1.8 for plausible parameter combinations. We now seek to clarify the connection between the von Bertalanffy equation and the exponent α for certain approximations revising the mathematical form of g 3 (.).

Extended Analysis: Power-Laws From the von Bertalanffy Equation
To extract power-law features from the von Bertalanffy equation and place them in the framework of Equation (5), it is assumed again that individual GPP= a ag P m LAP = R A /(1 − CUE) and R A ≈ k m w, resulting in an estimate of a ag P m LAP = k m w/(1 − CUE). Hence, Equation (38) reduces to The solution to Equation (40) is now a power-law of the form Yoda's rule is recovered when (k m /α m ) = (3/2)(CUE −1 − 1) whereas the metabolic exponent is recovered when k m /α m = (4/3)(CUE −1 − 1). Because CUE ≈ 0.5 (Waring et al., 1998;Manzoni et al., 2018), this analysis leads to a unique relation between respiratory and mortality time scales, and the exponent α given by k m /α m ≈ 3/2, 4/3, or 1. That is, the exponent α may be related to the ratio of the two aforementioned time scales.

Phase Space Trajectories Constraints on α
Dynamical systems theory has been used to explore selfthinning empirically by modeling the w − p time-course in crowded plant populations (Hara, 1984). The dynamical system can be expressed in terms of relative quantities, namely (relative) mortality rate (i.e., p −1 dp dt ) and relative growth rate (i.e., w −1 dw dt ). Among the choices of the functions linking p and w to relative mortality and growth rate, generalized Gompertz functions are empirically well-supported (Hara, 1984;Tsoularis and Wallace, 2002). With such choices, the dynamical systems theory can establish explicit dependencies of the empirical coefficients and the exponent α and plausibility constraints. Such a plausibility constraint is the imposition that equilibrium points are stable fixed points (as expected in selfthinning). A key result is that several combinations of the empirical parameters of the Gompertz function lead to exponents commensurate with the "rule of constant yield" or Yoda's thinning rule, while other exponents are possible with other empirical parameter combinations. The details are illustrated in the Supplementary Material.

Mechanism 7: Size Distribution Arguments
The self-thinning rule can also be obtained by following the temporal evolution of a population of individuals characterized by a certain size, which is interpreted as a stochastic variable. Without loss of generality, stem diameter D can be considered as the relevant size and can be linked to plant height and mass using allometric relations. For size-structured populations, it can be shown that the distribution of individuals of size D conditioned at time t, p D (D|t), is determined by the von Foerster equation (von Foerster, 1959;Hara, 1988;Kohyama, 1992;Strigul et al., 2008) where G is the growth rate (i.e., G = dD/dt) and µ is a mortality rate applied to plant density. In addition, a boundary condition p D (D 0 |t = 0) (i.e., where D 0 is the diameter at birth) must be specified. In principle, the self-thinning and constant yield laws could be obtained from the solution p D (D|t) of Equation (42) for specific choices of the functions G and µ, and the allometric relations between D and w. Here, a simplified approach is followed using the perfect crown plasticity rationale by Strigul et al. (2008) though by no means is this approach unique (Kohyama, 1992). As before, the focus is on a monospecific, even-aged stand with negligible mortality until canopy closure, a constant growth rate (so that D (t) = D 0 + Gt), small but finite initial stem diameter and diameter variance (D 0 ≪ Gt), and constant allometric coefficient linking canopy area per unit ground area to stem diameter (here a c =canopy area per unit ground area over D 2 ). With these conditions and assumptions, integrating the distribution p D (D|t) over all initial sizes, the total canopy area per unit ground area is calculated as, where (D 0 + Gt) 2 ≈ (Gt) 2 is the canopy area per unit area and p 0 is, as before, the initial plant density. When canopy closure occurs, the canopy area per unit ground area reaches 1. Hence, the canopy closure time t * can be calculated as Upon canopy closure (t > t * ), plant growth must adjust to maintain a closed canopy as time progresses, which requires lowering plant density through the death of suppressed, shaded plants according to, These constraints allow finding the scaling relation between plant biomass and density in the two regimes -before (t < t * ) and after (t > t * ) canopy closure. For t < t * , plant biomass w ∼ D 2 h, where h is the plant height as before. However, neither D nor h depend on plant density because they only depend on time before canopy closure. As a consequence, plant biomass w ∼ p 0 when t < t * . In contrast, for t > t * , plant biomass depends on plant density because after canopy closure Equation (45) yields, where isometric scaling of height and diameter (i.e., h ∼ D) was assumed to recover Yoda's rule (last term).

Mechanism 8: Neighborhood Interaction Arguments
As a bridge to the general framework in Equation (5), the equations specifying g 1 (w i ) for an individual i must now include interaction terms with adjacent individuals to explicitly account for competition. Upon specifying mortality and solving w i for each individual, the solution yields the mean biomass w and g 2 (p) by aggregating over all surviving individuals (i.e., the standscale). Hence, w(t) − p(t) trajectories are constructed thereby allowing the determination of α. The previously discussed carbon balance approaches only accounted for competition indirectly by varying the average individual's photosynthetic rate with p. Also, size-structured population approaches accounted for interactions among individuals implicitly. Individual-based models (Aikman and Watkinson, 1980;Westoby, 1982;Hara, 1988;Thomas and Weiner, 1989;Adler, 1996;Li et al., 2000;Stoll et al., 2002;Strigul et al., 2008;Chu et al., 2010;Coomes et al., 2011;Deng et al., 2012;Rivoire and Le Moguedec, 2012;Rüger and Condit, 2012;Lin et al., 2013) are often characterized as either spatially explicit, where plant spatial coordinates are specified, or spatially implicit, where only the zone of influence of each plant is tracked assuming equal spacing among individuals. Such models recover the 3/2 or 4/3 exponents for a wide range of mortality conditions or metabolic thresholds, while others exhibit greater sensitivity to competition between adjacent plants. These models follow a continuum of competition modes bounded by two limiting cases: size asymmetric competition where the largest plants acquire all the resources in overlapping areas to size symmetric competition where resources in overlapping areas are divided equally among interacting individuals regardless of their size (Weiner, 1990;DeMalach et al., 2016). Obviously, the degree of competition among individuals increases in all such models when the plot area A s available for growth is diminished. These models can recover increased variability, skewness, or bi-modality in the histograms of individual plant biomass w i as self-thinning is initiated at the stand level. A parsimonious, spatially implicit model is now considered to explore how different competition modes, initial densities, and experimental durations result in different α values. While some spatially explicit, more complex models are more realistic, the spatially implicit model explored here strikes a balance between simplicity and the ability to grasp all the proposed power-law exponents.

Competition and Mortality in Spatially Implicit Models
In this model, the growth rate of an individual plant i is assumed to be (Aikman and Watkinson, 1980) where a i and b i are constants for a given stand, reflecting growth rate per unit area and the need for more resources as individual plant biomass increases, b i depends on the maximum individual biomass w max , and s i measures the space occupied by plant i, which is linked to its size by a prescribed allometric relation where k g is a constant relating the area or zone of influence s to plant weight w. The 2/3 exponent is derived from dimensional considerations for isometric growth as discussed in section 2.1. The function f (s i ) encodes all of the spatial competition on the growth rate of plant i. To represent the space limitation and the two end-members of symmetric vs. asymmetric size-based competition, f (s i ) was represented as (Aikman and Watkinson, 1980) f where the term j s j A s describes the space availability for resource acquisition (i.e., a measure of crowding) and s s i measures the relative size of plant i compared to the mean size s. The two exponents φ 1 and φ 2 describe the importance of each mode of competition, representing respectively the roles of crowding and size asymmetry. The plot size A s sets the spatial domain for competition. The initial number of uniformly distributed plants within A s defines p 0 . By varying φ 1 and φ 2 , various modes of competition can be explored and their effect on α tracked. Mortality of plant i occurs when its carbon balance first becomes negative (i.e., dw i /dt < 0). Needless to say, mortality need not occur when dw i /dt < 0 (at least not on short time scales) though a negative carbon balance at the individual level implies a progressive competitive disadvantage and an increasing likelihood of mortality. Two related issues are addressed: The effect of f (s i ) on (i) the value of the self-thinning exponent α and (ii) the emergence of constant final yield when varying p 0 across multiple stands, waiting for a fixed duration, and observing w and p at each stand separately as shown in Figure 1 (Weiner and Freckleton, 2010).
Because growth and mortality in Equations (47) and (48) are proportional to powers of biomass (w i ) without distinguishing live and dead parts, this model is more appropriate for herbaceous species rather than forests. The individual tree biomass in high density forests may consist of a considerable proportion of dead biomass, reducing respiration costs. To avert this complexity, large initial densities and growth rates are used as is the case in crops. In fact, the range of parameter values used here (Table S3) are within the range used in Aikman and Watkinson (1980) and which were shown to agree with stand structure observations in even-aged monoculture competition experiments (Ford, 1975).

Effects of Competition Type on α and the Emergence of Constant Final Yield
For the first set of model runs, the power-law relation between individual biomass w(p 0 |T p ) at a fixed time after sowing T p and initial density p 0 is examined, where T p is the integration period of the simulation (Figure 4). A constant integration period of T p = 50 days is maintained for these runs during which no mortality occurs as is the case in the seminal work of Shinozaki and Kira (1956). Here, where the subscript CD stands for competition-density. The model runs here compare different plots at different p 0 and at a fixed period after sowing. Clearly, α CD = 1 corresponds to the constant final yield rule as in Equation (3). In Figure 4A, for small p 0 , normalized biomass per individual at the end of the simulation w(p 0 |T p )/w max appears to be insensitive to p 0 . As p 0 increases, variations in α CD occur depending on choices made about φ 1 and φ 2 . Figure 4B shows that at relatively low crowding exponent (φ 1 = 5) and relatively large size asymmetry exponent (φ 2 = 5), α CD = 1, corresponding to invariant biomass per ground area y(p 0 |T p ) = y c regardless of the initial sowing density. This is therefore a manifestation of the constant final yield rule but not of self-thinning since mortality is absent. As the crowding exponent becomes large (φ 1 → 10), α CD becomes bounded between 4/3 and 3/2, and insensitive to variations in the size asymmetry exponent φ 2 . These cases are compatible with neither the constant final yield rule nor the self-thinning rule. The temporal patterns of y(t) = w(t)p(t) associated with various choices of p 0 and φ 2 are shown in Figure 5 for different FIGURE 4 | (A) Each line represents modeled normalized biomass per individual for multiple simulations of different initial density p 0 after a fixed integration period of 50 days. The lines correspond to the three competition scenarios indicated in the legend. The corresponding scaling exponent α CD is displayed next to each plotted line. (B) Variation in α CD driven by different competition scenarios for crops (as determined by high initial plant densities and growth rates). Model parameters are found in Table S3.
plant properties and a longer integration period of T p = 150 days to assess the robustness of the results (see Table S2). The longer period allows for the presence of mortality whose onset in time is depicted using circles in Figure 5 (p(t) < p 0 ). Here, an intermediate crowding exponent φ 1 = 10 is kept as a constant. Biomass per ground area reaches an equilibrium that is sensitive to p 0 for low φ 2 values of 3 and 5 (Figures 5A,B). This does not conform to the constant final yield rule. For the highly size asymmetric mode of competition set by φ 2 = 7, the steady state biomass becomes independent of p 0 (Figure 5C), consistent with the constant final yield rule as presented by Xue and Hagihara (1998) when density-driven mortality occurs. The fact that a constant biomass is achieved for large φ 2 underscores that the FIGURE 5 | Modeled temporal variations of biomass per unit ground area with time for three competition types with less and more prominent size asymmetry effect (respectively smaller and larger φ 2 ). (A) At low φ 2 , model runs with different initial densities p 0 (m −2 ) result in divergent final yields. Mortality only occurs in the low density plots as indicated by the circles. The final densities [p(150) in m −2 ] as indicated by the color-coded numbers on their corresponding curves are not very different from p 0 . (B) At medium size asymmetry, final yields are closer than in (A), but differences still occur after 150 days. Mortality is present in all four plots. (C) At high size asymmetry, final yields are the same regardless of p 0 90 days into the simulation which corresponds to the constant final yield rule. All plots reach the same low final density of 15 m −2 . Model parameters are found in Table S3. Circles designate the onset of mortality in time and color-coded numbers are the final densities at the end of the 150 days corresponding to each simulation. phenomenon of the constant final yield only applies to certain types of plants competing for certain limited resources.
Self-thinning is shown in Figure 6 and Figure S1. The differences between experiments conducted at a single stand experiencing self-thinning sampled through time (Equation 2), and multiple stands with varying p 0 at a fixed period after sowing (Equations 4 and 50) is seen by comparing Figure 6B and Figure 4B. The discrepancy in the contour plots underscores the fact that the meaning of the scaling exponent α and α CD  Table S3.
is not equivalent. As earlier noted, one of the highly cited critiques on the universality of the self-thinning exponent was an empirical analysis by Weller (1987). Weller noted that when analyzing multiple stands with different p 0 , exponents differing from those tracking self-thinning in a single stand were obtained. The differences between these two experimental setups have been discussed elsewhere (Weiner and Freckleton, 2010), but are quantified here using the same spatially implicit model (as well as a range of φ 1 and φ 2 ).
In Figure 6, it is seen that the stronger the influence of space availability on competition (φ 1 ), the steeper the power-law relation between w(t) and p(t). The effect of size asymmetry on the α is more nuanced (Figure 6B). It is therefore seen that Yoda's definition of self-thinning, where α = 3/2 is achieved only when competition for space is high (φ 1 > 15) and size asymmetry is moderate to high (φ 2 > 3; Figure 6B). The green dashed line of Figure 6A where α = 1.06 is equivalent to the scenario in Figure 5C where the biomass per ground area is invariant with respect to p 0 . Figure 6B shows that all 3 aforementioned exponents (α ≈ 1, 4/3, 3/2) can be recovered from the same zonal model, depending on choices of φ 1 and φ 2 .

IMPACTS ON FOREST MANAGEMENT, FUTURE OUTLOOK, AND CONCLUSIONS
Competition for resources among same-species individuals sharing the same resource niche can be as complex as interactions among individuals of different species (Perry et al., 2008). That such competition among individuals of the same species results in power-law relations between the mean weight of an individual w and plant density p remains scientifically challenging to explain. Yet, such power-law relations are appealing to agricultural and forestry practitioners and have routinely been used in crop and forest management. In this context, mortality is only due to resource competition between individual plants, neglecting mortality due to external factors such as ice storms, hurricanes, forest fires, extended droughts, insect outbreaks, and human thinning of forests for management. As such, they set an "upper bound" on mean individual size for a given stocking (or planting) density (Luyssaert et al., 2011). In the case of crop management, initial planting density emerges as a key determinant of individual biomass and elapsed time when the steady state yield is reached; whereas in the case of forest management the w − p trajectories serve as a guide to when, and how much and often stands must be thinned to either maximize profitability or to mitigate hazards such as forest fires, insect outbreaks, or drought-induced mortality.
As already alluded to by Reineke (1933), forest density management utilizes size-density indices because they are presumably independent of site quality and stand age and the self-thinning line has taken center stage in determining management regimes (Begin et al., 2001). Such a presumably time-invariant power-law relation between w and p enables forest managers to also compare levels of growing stock regardless of differences in site quality or stand age. A particular set of management objectives resulting in an ideal p value can be projected forward or backward in time to a different development stage using the aforementioned w − p trajectories if the powerlaw exponent is known (e.g., α = 3/2). The self-thinning rule is a particularly powerful tool in combination (or as a part of) growth models to inform managers when the stands reach a particular management regime (Landsberg and Waring, 1997).
This review has focused on the many hypothetical mechanisms generating power-law relations between w and p due to the constraints imposed on resource competition in monospecific plots. Depending on the resource constraints (e.g., structural, allometric, hydraulic, supply of energy) and the type of competition imposed, multiple arguments suggest that the exponent of the power-law solution to dw/dp = g 3 (w, p) converges toward one of the three α values: 1, 4/3, 3/2 ( Table 1). The different α values reflect the numerous environmental influences and physiological factors and the degree of asymmetry of the competitive interaction (e.g., light interception is dominated by the tallest trees; Craine and Dybzinski, 2013). Therefore, for foresters aiming to optimize productivity, they should manage tree density based on the specific resource constraints shaping inter-plant competition. Often, the cultivation strategy maximizing plant density also minimizes resource availability such as soil water, i.e., the tragedy of the commons (Hardin, 1968). However, forest managers may be able to avoid this risky strategy by balancing resource use and plant density. This may be increasingly relevant in a changing climate where frequent and extended droughts are becoming a reality in many parts of the world. If storm intensities and inter-arrival times change in relatively short time-scales, then rooting profiles that successfully harnessed soil water in the past might become less effective (Farooq et al., 2009).
The theoretical results presented can be used to generate hypotheses on what controls α, to be tested in specific experiments or simulation studies. For example, species characterized by contrasting growth patterns or hydraulic traits could be grown under the same conditions to test predicted patterns of α. Similarly, trends in α could be assessed along climatic and edaphic gradients to test predicted deviations from the 3/2 or 4/3 values. The focus was purposely restricted to monospecific stands, but self-thinning also occurs in diverse communities, though niche complementarity and facilitation effects can become important drivers of the plant mass-density relations (Loreau and Hector, 2001). It is possible that denser communities containing a greater number of small individuals (belonging to more than one species) emerge when these effects are at play compared to monospecific stands. To tackle this problem, models describing a multitude of species (or functional traits), or capturing differences across individuals, should be used, which we expect will require a broader range of scaling exponents as the communities become more diverse.
The finding that the self-thinning exponent is not invariant has several consequences for forest managers designing their thinning regimes based on an invariant self-thinning rule (Drew and Flewelling, 1977). Time variance has been attributed to stand aging, canopy closure and environmental change, including increasing aridity in many parts of the world. Several amendments to the size-density indices presented in this text have been proposed elsewhere to take these effects into account (Zeide, 2001;Ge et al., 2017;Aguirre et al., 2018;Bravo-Oviedo et al., 2018;Zhang et al., 2019). Forest managers may expect a given self-thinning slope based on data from space-for-time substitution and, thus, set their thinning or harvest operations based on this expectation (Drew and Flewelling, 1979). For example, the self-thinning rule can be used to characterize "reference" conditions and, based on that, define indicators of the degree of land use intensity (Luyssaert et al., 2011). Such indicators quantify how far a given stand is from either a pristine forest or a stand following the self-thinning rule. However, as shown here, the shape of the self-thinning rule may vary depending on growth conditions and therefore indicators based on this curve may be sensitive to the chosen exponent and intercept.
A line of inquiry of increasing relevance to crops and plantation forestry alike is the effect of environmental change (e.g., elevated atmospheric CO 2 , or air temperature) on C or α. Do the w − p trajectories remain the same or are they altered with these changes? Does elevated atmospheric CO 2 simply speed up the trajectories in the w − p phase plane? Does α remain constant but C likely to vary due to changes in leaf area index or other ecosystem properties? These are but a few of the questions that could motivate future work (Brunet-Navarro et al., 2016;Jump et al., 2017;Bravo-Oviedo et al., 2018). The implications of these answers to forest management cannot be overstated. If C or α vary under future conditions, management tactics will need to be adjusted accordingly. Even if the parameters of the self-thinning rule do not change, higher CO 2 and air temperature may promote growth (in absence of other limitations), resulting in faster movements along the w − p trajectory. This alone would require adjusting the thinning schedule. Modeling studies suggest that maintaining under-stocked stands (below the self-thinning curve) by more frequent or intense thinning could compensate negative impacts of future environmental conditions on the tree C balance (Collalti et al., 2018). However, for new thinning approaches to be effective, they will need to be based on a selfthinning rule that accounts for future growth conditions. For example, Equation (41) suggests that lower values of α can be expected if autotrophic respiration increases more than GPP in a warmer world, or if stands become nutrient-limited or age faster, resulting in lower CUE (Collalti et al., 2018). The use of such power-law expressions in forest management map onto the famous quote by the great Russian physicist Lev Landau: Money is in the exponent. And exponent needs to be calculated precisely.
Power law relations between measures of biomass and density have been the subject of over a century of experimentation and theoretical analysis. They not only describe biomass development as a function of density for a single stand but also steady-state biomass as a function of maximum density for species ranging nine orders of magnitude by weight. The ubiquity and relative invariance of these power law relations makes them a research breeding ground to uncover the underlying mechanisms of interplant competition and develop effective management strategies for forests and croplands increasingly suffering from aridity in a changing climate.

DATA AVAILABILITY STATEMENT
All sources of data used in this study are included in the article/ Supplementary Material.

AUTHOR CONTRIBUTIONS
GK, RO, and SM conceived the initial ideas and design of the manuscript. AM and SM designed all the figures and illustrations and were responsible for ideas in the extended analysis. ML and GV contributed to the literature review and checked derivations and dimensional consistencies. AM and GK developed and conducted all model simulations. All authors contributed to the writing of the manuscript and provided their consent to publish this work.