The size-structured PDE

theory

NotePrerequisites

This page builds on The big picture and The maths, gently. Read those first if the notation here is unfamiliar.

The idea first

Tracking a forest plant by plant would quickly become unmanageable. A size-structured model takes a different approach: it describes how many plants occur at each size and follows the resulting distribution through time. Growth moves plants towards larger sizes, mortality removes them, and recruitment adds new plants at the smallest size.

Competition links those processes. In the canopy-light strategies, the plants in a patch create a vertical light profile, and that profile changes their growth, mortality, and reproduction. Disturbance adds another level: a landscape contains patches of different ages, each with its own size distribution. The metacommunity is obtained by combining those patches.

This page develops the equations in that same order: individual life histories, one competing patch, and then a landscape of patches. The formulation used in plant was developed in Falster, Brännström, Dieckmann, & Westoby (2011) and Falster, Brännström, Westoby, & Dieckmann (2015), following general principles from De Roos, Taljapurkar, & Caswell (1997), Kohyama (1993), and Moorcroft, Hurtt, & Pacala (2001).

The companion page on adaptive dynamics builds on this machinery to define invasion fitness and evolutionary endpoints.

The argument proceeds in layers. We first follow one plant through its life, then combine many plants into a size distribution within one patch, and finally combine patches of different ages into a metacommunity. The appendices derive the two transport equations and show how to convert a density from one size coordinate to another.

NoteHow to read size and environment on this page

The notation follows the FF16 and TF24 strategies: \(H\) is plant height and \(E_a\) includes the canopy light profile. The demographic transport equation is more general than that notation. K93 uses stem diameter as its biological size coordinate and cumulative basal area as its competitive environment, while TF24 additionally couples physiological rates to soil water. In each case the strategy supplies the appropriate rates to the same demographic bookkeeping.

Table 1 provides a list of names and definitions used throughout this page.

Table 1: Variable names and definitions in the demographic model of the plant package.
Symbol Unit Description
Patch state variables
\(K\) Number of species in the metacommunity
\(a, a^{\prime}\) yr Patch age (time since disturbance)
\(a_0\) yr Patch age when plant germinates
\(E_a\) Profile of canopy openness within a patch of age \(a\)
\(E_a(z)\) Canopy openness at height \(z\) within a patch of age \(a\)
Plant state variables
\(x\) Vector of traits for a species
\(H\) m Height of a plant
\(H_0(x)\) m Height of a seedling with traits \(x\) after germination
\(z\) m Height in canopy
\(A_{\rm l}(H)\) m\(^2\) Leaf area of a plant with height \(H\)
\(Q(z, H)\) Fraction of leaf area above height \(z\) for a plant with height \(H\)
Abundance measures
\(P(a)\) yr\(^{-1}\) Frequency-density of patches of age \(a\)
\(Y_{x}\) m\(^{-2}\) yr\(^{-1}\) Global birth rate for a species with traits \(x\)
\(N(H, x, a)\) m\(^{-1}\) m\(^{-2}\) Density of plants at height \(H\) per unit ground area for given traits \(x\) and patch age \(a\)
Demographic rates
\(g(x, H, E_a)\) m yr\(^{-1}\) Height growth rate of a plant with traits \(x\) and height \(H\) in the light environment \(E_a\) in a patch of age \(a\)
\(f(x, H, E_a)\) yr\(^{-1}\) Seed production rate of a plant with traits \(x\) and height \(H\) in the light environment \(E_a\) in a patch of age \(a\)
\(d(x, H, E_a)\) yr\(^{-1}\) Instantaneous mortality rate of a plant with traits \(x\) and height \(H\) in the light environment \(E_a\) in a patch of age \(a\)
\(d_{\rm P}(a)\) yr\(^{-1}\) Instantaneous disturbance rate of a patch of age \(a\)
Demographic outcomes
\(H(x, a_0, a)\) m Height of a plant with traits \(x\) in a patch of age \(a\) having germinated at patch age \(a_0\)
\(S_{\rm D}\) Probability a seed survives dispersal
\(S_{\rm G} (x, H_0, E_{{\rm a}0})\) Probability the seed of a plant with traits \(x\) germinates successfully at height \(H_0\) in the light environment \(E_{{\rm a}0}\) in a patch of age \(a_0\)
\(S_{\rm I} (x, a_0, a)\) Probability a plant survives from patch age \(a_0\) to patch age \(a\)
\(S_{\rm P} ( a_0, a)\) Probability a patch remains undisturbed from patch age \(a_0\) to patch age \(a\)
\(\tilde{R}(x, a_0, a)\) Cumulative seed output for plants with traits \(x\) from patch age \(a_0\) to patch age \(a\)
\(R\left(x^\prime, x\right)\) Basic reproduction ratio for a mutant plant with traits \(x^\prime\) growing in a metapopulation of resident plants with traits \(x\)
Miscellaneous constants
\(k_{\rm I}\) Light extinction coefficient
\(\bar{a}\) yr Mean interval between patch disturbances

Metacommunity dynamics

The demographic equations need a physiological sub-model to supply the biology. For a plant with traits \(x\), size \(H\), and environmental state \(E_a\), that sub-model returns its growth, mortality, and fecundity rates. In the canopy-light formulation it also returns the leaf area the plant contributes to shading. Specifically, the physiological model is responsible for calculating the following variables from Table 1: \(A_{\rm l}(H)\), \(Q(z, H)\), \(H_0(x)\), \(g(x, H, E_a)\), \(f(x, H, E_a)\), \(d(x, H, E_a)\), and \(S_{\rm G} (x, H_0, E_{{\rm a}0})\). All other variables in Table 1 are calculated by the demographic model using the equations below. This division is important: changing the physiological strategy changes the supplied rates, while the demographic bookkeeping stays the same.

Individual plants

We first follow one plant. Its traits are \(x\), its size is represented by height \(H\), and \(a\) is the current age of the patch in which it grows. The plant germinated earlier, when the patch was age \(a_0\), at an initial height \(H_0(x)\).

The light environment \(E\) is a vertical profile through the patch rather than a single number. Because that profile changes as the patch develops, we write it as \(E_a\). The functions \(g(x, H, E_a)\), \(f(x, H, E_a)\), and \(d(x, H, E_a)\) are respectively the plant’s instantaneous growth, fecundity, and mortality rates.

Height at patch age \(a\) is initial height plus all growth accumulated since germination: \[ H(x, a_0, a) = H_0(x) + \int_{a_0}^{a} g(x, H(x, a_0, a^\prime), E_{a^\prime}) \, {\rm d}a^\prime \tag{1}\]

Survival has two parts: the probability of surviving germination, followed by survival against the mortality hazard accumulated as the plant grows: \[ S_{\rm I} (x, a_0, a) = S_{\rm G} (x, H_0, E_{{\rm a}0}) \, \exp\left(- \int_{a_0}^{a} d(x, H(x, a_0, a^\prime), E_{a^\prime}) \, {\rm d} a^\prime \right) \tag{2}\]

Here \(S_{\rm G} (x, H_0, E_{{\rm a}0})\) is the probability that a seed germinates successfully. Finally, cumulative seed output adds fecundity over time, but discounts seed production by the probability that the plant is still alive: \[ \tilde{R}(x, a_0, a) = \int_{a_0}^{a} f(x, H(x, a_0, a^\prime), E_{a^\prime}) \, S_{\rm I} (x, a_0, a^\prime) \, {\rm d} a^\prime \tag{3}\]

Thus Eqs. 13 are three running totals along the same life history: accumulated growth, accumulated mortality risk, and accumulated reproduction. Their notation allows every rate to change non-linearly with traits, plant size, light, and patch age.

Patches of competing plants

We now move from one plant to a patch. At patch age \(a\), \(N(H | x, a)\) is the size-density distribution of plants with traits \(x\) and height \(H\). It is a density with respect to height: integrating it over a height interval gives the abundance of plants in that interval. A finite patch can instead be represented by discrete individuals, but for a sufficiently large patch the model treats \(N\) as a continuous function.

Each plant still follows Eqs. 13, but two feedbacks now matter. Existing plants alter the shared light profile \(E_a\), and new seedlings continually enter at the lower end of the size distribution.

For the FF16 and TF24 canopy-light formulation, plants interact by shading one another. TF24’s physiological rates also respond to the soil-water state. Following standard biophysical principles, we let the canopy openness \(E_a(z)\) at height \(z\) in a patch of age \(a\) decline exponentially with the total amount of leaf area above \(z\), \[ E_a(z) = \exp \left(-k_{\rm I} \sum_{i = 1}^{K} \int_{0}^{\infty} A_{\rm l}(H) \, Q(z, H) \, N(H | x_i, a) \, {\rm d}H \right), \tag{4}\] where \(A_{\rm l}(H)\) is the total leaf area of a plant with height \(H\), \(Q(z, H)\) is the fraction of this leaf area situated above height \(z\) for plants of height \(H\), \(k_{\rm I}\) is the light extinction coefficient, and \(K\) is the number of species in the metacommunity.

Read Equation 4 from the inside out. For species \(i\), the integral totals leaf area above \(z\) across all plant heights; the sum then combines the \(K\) species. Multiplying by \(k_{\rm I}\) converts leaf area into optical depth, and the exponential converts optical depth into the fraction of light that remains. The result is a single patch light profile that every plant both contributes to and experiences.

Assuming patches are sufficiently large, the dynamics of \(N\) can be modelled deterministically via the following partial differential equation (PDE) (Kohyama, 1993; De Roos et al., 1997; Moorcroft et al., 2001), \[ \frac{\partial}{\partial a} N(H | x, a) = - d(x, H, E_a) \, N(H | x, a) - \frac{\partial}{\partial H} \left[g(x, H, E_a) \, N(H | x, a)\right]. \tag{5}\] See Section 3.2 for the derivation of this PDE.

The two terms on the right-hand side are a population ledger. The first removes plants through mortality. The second describes the flow of plants along the height axis as they grow: plants leave one height interval and enter the next. It is written as the derivative of the flux \(gN\), so changes in both growth rate and density are handled together.

Equation 5 also needs two boundary conditions. The first supplies recruits at the smallest height \(H_0(x)\). It links the flux of seedlings entering the size-density distribution to the rate \(Y_{x}\) at which seeds arrive in the patch, \[ N(H_0 | x, a_0) = \left\{ \begin{array}{ll} \frac{Y_{x} \, S_{\rm G} (x, H_0, E_{{\rm a}0}) }{ g(x, H_0, E_{{\rm a}0}) } & \textrm{if } g(x, H_0, E_{{\rm a}0}) > 0 \\ 0 & \textrm{otherwise.} \end{array} \right. \tag{6}\] The function \(S_{\rm G} (x, H_0, E_{{\rm a}0})\) denotes survival through germination and must be chosen such that \(S_{\rm G} (x, H_0, E_{{\rm a}0}) / g(x, H_0, E_{{\rm a}0}) \rightarrow 0\) as \(g(x, H_0, E_{{\rm a}0}) \rightarrow 0\), to ensure a smooth decline in initial density as conditions deteriorate (Falster et al., 2011).

The division by \(g\) converts a flux of successful recruits per unit time into a density per unit height. The special case for \(g \leq 0\) prevents seedlings from accumulating at the boundary when conditions do not allow them to grow into the size distribution.

The second boundary condition specifies the state of a patch immediately after disturbance, at \(a = 0\). Throughout, we consider situations starting with an empty patch, \[ N\left(H|x,0\right) = 0, \tag{7}\] although non–zero distributions could alternatively be specified (e.g. Kohyama, 1993; Moorcroft et al., 2001).

Age structure of patches

We now zoom out from one patch to a landscape containing patches of many ages. Here \(a\) is time since the last disturbance, and \(P(a)\) is the frequency-density of patches of age \(a\). Integrating \(P(a)\) over an age range gives the fraction of landscape in that range. The function \(d_{\rm P}(a)\) is the instantaneous disturbance rate, or hazard, for a patch of age \(a\).

We first focus on an equilibrium patch-age distribution. See Section 3.1 for the time-dependent derivation. At equilibrium, patches flow steadily toward older ages and leave each age class when disturbed (McKendrick, 1926; Foerster, 1959): \[ \frac{{\rm d}}{{\rm d} a} P(a) = - d_{\rm P}(a) \, P(a) , \tag{8}\] with the boundary condition \[ P(0) = \int_0^\infty d_{\rm P}(a) \, P(a) \, {\rm d} a. \tag{9}\]

The boundary condition says that the flow of newly disturbed, age-zero patches must equal the total disturbance flow out of all older ages.

The probability that a patch remains undisturbed from patch age \(a_0\) to patch age \(a\) is then given by \[ S_{\rm P} (a_0, a) = \exp\left( - \int_{a_0}^{a} d_{\rm P}(a^\prime) \, {\rm d} a^\prime \right). \tag{10}\]

This has the same survival form as the individual-plant equation: cumulative hazard appears inside a negative exponential. The difference is that the event is patch disturbance rather than plant death.

These equations lead to an equilibrium distribution of patch ages, \[ P(a) = P(0) S_{\rm P} (0, a), \tag{11}\] where \[ P(0) = \frac1{\int_0^\infty S_{\rm P} (0, a) {\rm d}a} \tag{12}\] is the average disturbance frequency of a patch and, at the same time, the frequency-density of patches of age \(0\).

In plain language, the abundance of patches of age \(a\) equals the rate at which new patches are created, \(P(0)\), multiplied by the chance that a patch survives undisturbed to age \(a\). The normalising expression in Equation 12 makes the total area under \(P(a)\) equal to one.

The default approach in plant is to assume that \(d_{\rm P}(a)\) increases linearly with patch age, which leads to a Weibull distribution for \(P\), specified by a single parameter \(\bar{a}\) measuring the mean interval between disturbances (see Section 3.1.2 for details).

Trait-, size-, and patch-structured metacommunities

The one-patch calculation scales up to a metacommunity under three main assumptions:

  • disturbances such as fires, storms, or landslides reset whole patches;
  • plants compete within a patch, while competition between neighbouring patches is negligible; and
  • dispersal connects the patches strongly enough that the same seed pool can recolonise them.

With one species this is a metapopulation; with several it is a metacommunity. Eqs. 512 describe its demographic dynamics.

The disperser pool links the patch and landscape scales. For each strategy \(x\), its landscape-wide seed rain \(Y_x\) is the rate at which viable seeds are produced across all patches, \[ Y_x = \int_0^{\infty} \int_0^{\infty} P(a) \, S_{\rm D} \, f(x, H, E_a) \, N(H | x, a) \, {\rm d} H \, {\rm d} a, \tag{13}\] where \(S_{\rm D}\) is the survival probability of seeds during dispersal.

Again, read the double integral from the inside out. At a given patch age \(a\), the inner integral adds fecundity over every plant height, weighted by plant density. The outer integral averages that patch output across patch ages, weighted by how common each age is. \(S_{\rm D}\) retains only seeds that survive dispersal.

At metacommunity equilibrium, the seed-arrival rate \(Y_x\) is constant. If every disturbance also leaves the same empty initial state, all patches then follow the same developmental trajectory. Their only difference is how long ago they were disturbed. This is why one patch trajectory, combined with the patch-age distribution, can represent the whole metacommunity.

To model the temporal dynamics of an archetypal patch, we need a value for \(Y_x\). But that value is both an input to the patch through Equation 6 and an output of the metacommunity through Equation 13. At demographic equilibrium, the two must match. The numerical challenge is therefore to solve Eqs. 613 simultaneously for all strategies \(x\) in the metacommunity.

Emergent properties of metacommunities

Once \(N(H | x, a)\) and \(P(a)\) are known, stand- and landscape-level outputs are weighted integrals. A hat denotes an average over the entire metacommunity, while a bar below denotes an average across patch ages at a specified plant height.

The average density of individuals with traits \(x\) per unit ground area across the metacommunity is \[ \hat{N}(x) = \int_{0}^{\infty} \int_{0}^{\infty} P(a) \, N(H | x, a) \, {\rm d}a \, {\rm d}H, \tag{14}\] and the average size-density of plants with height \(H\) is \[ \bar{N}(x, H) = \int_{0}^{\infty}P(a) \, N(H | x, a) \, {\rm d}a. \tag{15}\]

Averages for other individual-level quantities use the same weighting. Denote by \(w(x, H, E_a)\) a quantity of interest, either a demographic rate (growth, mortality) or a state (plant height, leaf area, light environment), the average of \(w\) for plants with height \(H\) and traits \(x\) is \[ \bar{w}(x, H) = \frac1{\bar{N}(x, H)}\int_{0}^{\infty}P(a) \, N(H | x, a) \, w(x, H, E_a) \, {\rm d}a. \tag{16}\]

The numerator adds \(w\) across patch ages, weighted by the density of plants to which it applies; dividing by \(\bar{N}(x, H)\) turns that total into an average for plants at height \(H\).

The average of \(w\) across all individuals of the species is \[ \hat{w}(x) = \frac1{\hat{N}(x) }\int_{0}^{\infty} \int_{0}^{\infty}P(a) \, N(H | x, a) \, w(x, H, E_a) \, {\rm d}a \, {\rm d}H. \tag{17}\]

Here the same calculation also integrates over height and divides by the total metacommunity density \(\hat{N}(x)\).

When calculating the average mortality rate, one must decide whether mortality due to patch disturbance is included. Non-disturbance mortality is obtained by setting \(w(x, H, E_a) = d(x, H, E_a)\), while the total mortality due to growth processes and disturbance is obtained by setting \(w(x, H, E_a) = d(x, H, E_a) + d_{\rm P}(a) S_{\rm P}(0, a).\)

Sometimes the desired output is a total rather than a per-plant average. We can integrate over the size-density distribution to obtain an aggregate feature of the vegetation within a patch, \[ W(a) = \sum_{i = 1}^{N} \int_{0}^{\infty} N(H | x_i, a)\, w(x_i, H, E_a) \, {\rm d}H, \tag{18}\] and across the entire metacommunity, \[ \hat{W} = \int_{0}^{\infty} P(a) \, W(a) \, {\rm d}a. \tag{19}\]

For the analysis of invasion fitness and evolutionary endpoints that builds on these dynamics, see the adaptive dynamics page.

Appendices

The main text states the model in its working form. The appendices below show where the transport equations come from and explain one subtle property of any size density: its numerical value depends on the size coordinate in which it is expressed.

Derivation of PDE describing age-structured dynamics

We now distinguish two clocks. Patch age \(a\) is the time since a patch’s last disturbance, while calendar time \(t\) describes change in the whole landscape. Let \(P(a, t)\) be the frequency-density of patches of age \(a\) at time \(t\), and let \(d_{\rm P}(a, t)\) be their instantaneous disturbance rate. Its dynamics are given by \[ \frac{\partial}{\partial t} P(a, t) = -\frac{\partial}{\partial a} P(a, t)-d_{\rm P}(a, t) P(a, t), \tag{20}\]

The first term on the right moves undisturbed patches toward older ages. The second removes patches of age \(a\) when disturbance resets them to age zero. Those resets supply the boundary condition \[ P(0, t) = \int^{\infty}_{0}d_{\rm P}(a, t) P(a, t) \, {\rm d}a. \tag{21}\]

The frequency-density of patches of age \(a < a^{\prime}\) is given by \(\int_{0}^{a^{\prime}}P(a, t) \, {\rm d}a\), with \(\int_{0}^{\infty} P(a, t) \, {\rm d}a = 1\). If \(\frac{\partial}{\partial t}d_{\rm P}(a, t) = 0\), then \(P(a)\) will approach an equilibrium solution given by \[ P(a) = P(0) \, S_{\rm P}(0, a), \tag{22}\] where \[ S_{\rm P}(0, a) = \exp \left( - \int_{0}^{a} d_{\rm P}(a^\prime) \, {\rm d}a^\prime\right) \tag{23}\] is the probability that a patch remains undisturbed for duration \(a\), and \[ P(0) = \frac1{ \int_{0}^{\infty}S_{\rm P}(0, a) \, {\rm d}a} \tag{24}\] is the frequency-density of patches of age \(0\). The rate of disturbance for patches of age \(a\) is given by \(\frac{\partial (1-S_{\rm P}(0, a))}{\partial a} = - \frac{\partial S_{\rm P}(0, a)}{\partial a}\), while the expected lifetime of patches is \(- \int_0^\infty a \frac{\partial}{\partial a} S_{\rm P}(0, a) \, {\rm d} a = \int_0^\infty S_{\rm P}(0, a) \, {\rm d} a = \frac1{P(0)}\) (first step made using integration by parts).

Thus \(1/P(0)\) has two equivalent interpretations at equilibrium: it is the mean time between disturbances and the reciprocal of the rate at which age-zero patches appear.

An equilibrium distribution of patch ages may be achieved under a variety of conditions, for example, if \(d_{\rm P}(a, t)\) depends on patch age \(a\) but not on time \(t\). The rate of disturbance may also depend on features of the vegetation in the patch, rather than on patch age directly, in which case an equilibrium distribution of patch ages can still arise, provided the vegetation is also assumed to be at equilibrium.

Exponential distribution

If the rate \(d_{\rm P}\) of patch disturbance is constant with respect to patch age, then the rates at which patches of age \(a\) are disturbed follow an exponential distribution, \(-\partial S_{\rm P}(0, a)/ \partial a = d_{\rm P} \, \exp(-d_{\rm P} a)\). The distribution of patch ages is then given by \[ S_{\rm P}(0, a) = \exp\left(-d_{\rm P} a\right), \, P(0) = d_{\rm P}. \tag{25}\]

A constant hazard is memoryless: surviving to an old age does not change a patch’s instantaneous chance of disturbance.

Weibull distribution

If the rate of patch disturbance changes with patch age according to \(d_{\rm P}(a) = \lambda \psi a^{\psi-1}\), then the rates at which patches of age \(a\) are disturbed follow a Weibull distribution, \(-\partial S_{\rm P}(0, a)/ \partial a = \lambda \psi a^{\psi -1}e^{-\lambda a^\psi}\). \(\psi>1\) implies that the probability of disturbance increases with patch age, while \(\psi<1\) implies that it decreases with patch age. For \(\psi = 1\), we obtain the exponential distribution, a special case of the Weibull distribution. The Weibull distribution results in the following distribution of patch ages, \[ S_{\rm P}(0, a) = \exp(-\lambda a^\psi), \, P(0) = \frac{\psi \lambda^{\frac1{\psi}}}{\Gamma\left(\frac1{\psi}\right)}, \tag{26}\] where \(\Gamma(x) = \int_{0}^{\infty}e^{-t}t^{x-1} \, dt\) is the gamma function. We can also specify this distribution by the mean disturbance interval \(\bar{a} = \frac1{P(0)}\). From this, we can calculate the relevant value for \(\lambda = \left(\frac{\Gamma\left(\frac1{\psi}\right)}{\psi \bar{a}}\right)^{\psi}\).

The default in plant is to assume \(\psi=2\), such that \(d_{\rm P}\) increases as a linear function of patch age. The distribution of patch ages is then specified by a single parameter, \(\bar{a}\).

The shape parameter \(\psi\) controls how disturbance risk changes with age; \(\lambda\) then sets the time scale. Fixing \(\psi=2\) leaves only the mean disturbance interval to specify.

Derivation of PDE describing size-structured dynamics

The size-structured PDE is a conservation equation: change inside a thin height slice equals plants growing in, minus plants growing out, minus deaths. We now write that ledger explicitly, following De Roos et al. (1997). To keep the notation readable, this subsection suppresses the dependence on traits \(x\).

Assuming that all rates are constant within a thin interval \(\Delta H\), the total number of individuals within the interval spanned by \([H - 0.5\Delta H, H + 0.5\Delta H)\) is \(N(H, a)\Delta H\). The flux of individuals in and out of this size interval can be expressed as \[ \begin{array}{ll} &g(H - 0.5 \Delta H, a) \, N(H - 0.5 \Delta H, a) - g(H + 0.5 \Delta H, a) \, N(H + 0.5 \Delta H, a) \\ & - d (H, a) \, N(H, a)\Delta H\\ \end{array}. \tag{27}\] The first two terms describe the flux in and out of the considered size interval through growth, while the last term describes losses through mortality. The change in the density of individuals within this size interval during a time step \(\Delta a\) is thus \[ \begin{array}{ll} N(H, a + \Delta a)\Delta H - N(H, a)\Delta H = &g(H - 0.5 \Delta H, a) \, N(H - 0.5 \Delta H, a)\Delta a \\ & - g(H + 0.5 \Delta H, a) \, N(H + 0.5 \Delta H, a)\Delta a\\& - d (H, a) \, N(H, a)\Delta H\Delta a. \end{array} \tag{28}\] By rearranging, we obtain \[ \begin{array}{ll} \frac{N(H, a + \Delta a) - N(H, a)}{\Delta a} = & - d (H, a) \, N(H, a) \\ & - \frac{g(H + 0.5 \Delta H, a) \, N(H + 0.5 \Delta H, a) - g(H - 0.5 \Delta H, a) \, N(H - 0.5 \Delta H, a)}{\Delta H}. \end{array} \tag{29}\] The left-hand side above corresponds to the derivative of \(N\) as \(\Delta a\to 0\). For thin slices, \(\Delta H \to 0\), this yields \[ \frac{\partial}{\partial a} N(H, a) = - d (H, a) \, N(H, a) - \frac{\partial}{\partial H} (g(H, a) \, N(H, a)). \tag{30}\]

This is the same PDE used in the main text: mortality removes density locally, while the derivative of \(gN\) moves density through height space.

To complete the model, this PDE must be supplemented with boundary conditions that specify the density at the lower end of heights for all \(a\), as well as the initial distribution \(N(H,0)\). The former is derived by integrating the PDE with respect to \(H\) over the interval \((H_{0}, H_{\infty} )\), yielding \[ \frac{\partial}{\partial a} \int _{H_{0} }^{H_{\infty}}N(H, a) \, {\rm d} H = g(H_{0} , a) \, N(H_{0} , a) - g(H_{\infty}, a) \, N(H_{\infty}, a) - \int _{H_{0} }^{H_{\infty}}d (H, a) \, N(H, a) \, {\rm d} H. \tag{31}\] The left-hand side above is the rate of change of the total density of individuals in the population. On the right, the first two terms are fluxes through the lower and upper size boundaries, and the integral is the total death rate. We assume that \(N(H_{\infty}, a) = 0\). The incoming flux \(g(H_{0} , a) \, N(H_{0} , a)\) must therefore equal the population’s total birth rate \(B\), yielding the boundary condition \[ g(H_{0} , a) \, N(H_{0} , a) = B. \tag{32}\]

Converting density from one size unit to another

The numerical value of a density depends on the unit along its horizontal axis. For example, a density per kilogram is not numerically identical to a density per metre, even when both describe the same plants. What must stay unchanged is the number of individuals represented by a corresponding interval.

Suppose size-density is expressed using size \(M\), but we want to express it using size \(H\). We require a one-to-one conversion \(H = \hat{H}(M)\). Equality of abundance over matching intervals then requires the following: \[ \int_{M_1}^{M_2} N(M | x, a) \, {\rm d}M = \int_{\hat{H}(M_1)}^{\hat{H}(M_2)} N^\prime(H | x, a) \, {\rm d}H . \tag{33}\] For very small size intervals, this equation is equivalent to \[ \left(M_2 - M_1 \right) \, N(M_1 | x, a) = \left( \hat{H}(M_2) - \hat{H}(M_1)\right) \, N^\prime(\hat{H}(M_1) | x, a). \tag{34}\] Rearranging gives \[ N^\prime(\hat{H}(M_1) | x, a) = N(M_1 | x, a) \, \frac{M_2 - M_1}{\hat{H}(M_2) - \hat{H}(M_1)}. \tag{35}\] Noting that the second term on the right-hand side is simply the definition of \(\frac{{\rm d} M}{{\rm d} H}\) evaluated at \(M_1\), we have \[ N^\prime(H | x, a) = N(M | x, a) \, \frac{{\rm d} M}{{\rm d} H}. \tag{36}\]

The derivative \(\mathrm{d}M/\mathrm{d}H\) is the Jacobian of the change of variables. It stretches or compresses the density exactly enough that integrating it still returns the same number of plants.

References

De Roos, A. M., Taljapurkar, S., & Caswell, H. (1997). A gentle introduction to physiologically structured population models. In Structured population models in marine, terrestrial and fresh-water systems (pp. 119–204). New York: Chapman & Hall.
Falster, D. S., Brännström, Å., Dieckmann, U., & Westoby, M. (2011). Influence of four major plant traits on average height, leaf-area cover, net primary productivity, and biomass density in single-species forests: A theoretical investigation. Journal of Ecology, 99(1), 148–164. doi:10.1111/j.1365-2745.2010.01735.x
Falster, D. S., Brännström, Åke, Westoby, M., & Dieckmann, U. (2015). Multi-trait eco-evolutionary dynamics explain niche diversity and evolved neutrality in forests. bioRxiv, 014605. doi:10.1101/014605
Foerster, H. V. (1959). Some remarks on changing populations. In F. Stohlman (Ed.), The kinetics of cellular proliferation (pp. 382–407). Grune et Stratton, New York.
Kohyama, T. (1993). Size-structured tree populations in gap-dynamic forest: The forest architecture hypothesis for the stable coexistence of species. Journal of Ecology, 81(1), 131–143.
McKendrick, A. G. (1926). Applications of mathematics to medical problems. Proceedings of the Edinburgh Mathematical Society, 44, 1–34.
Moorcroft, P. R., Hurtt, G. C., & Pacala, S. W. (2001). A method for scaling vegetation dynamics: The Ecosystem Demography model (ED). Ecological Monographs, 71(4), 557–586. doi:10.1890/0012-9615(2001)071[0557:AMFSVD]2.0.CO;2