APPENDIX 1.
DISCUSSION OF ANALYSIS OF DYNAMICAL SYSTEM MODELS

The resilience of ecological systems has been assessed through the bifuraction analysis of dynamical systems (Rinaldi and Scheffer 2000). In this appendix, I demonstrate this approach via the analysis of a simple model of the fire-mediated relationships between longleaf pine and hardwood. An Excel version of this model and these analyses is available in Appendix 2.

A simple model

This model represents a forest site as being composed of a mix of hardwood species and longleaf pine, so that the proportion of hardwood and the proportion of longleaf in combination comprise the entire forest at that site. That is HW + LL = 1.0, where HW is the proportion of the site that is hardwood, and LL is the proportion of the site that is longleaf pine. In addition, HW and LL must both be greater than equal to zero and less than equal to one.

Because hardwood invades longleaf pine in the absence of fire, I model the dynamic of hardwood. Let X = HW, then LL = 1 - X, where the dynamics of hardwood are determined by the balance between the growth of hardwood and its destruction by fire so that dX/dt = growth - destruction.

The increased prevalence of hardwood species in longleaf pine stands is modeled as a logistic growth process:

growth = r*X*(1-X),
(A1.1)

where r is the growth rate of hardwood.

The growth of hardwood is inhibited by fire. Mortality produced by fire depends on the amount of hardwood (X), but the ability of fire to spread depends on longleaf pine. Longleaf pine (1 - X) facilitates the spread of fire. Sites dominated by longleaf pine are easier to burn, and the increased presence of longleaf pine eases the ability of fire to spread across the landscape. Consequently, the destruction of hardwood by fire was represented as proportional to the density of hardwood multiplied by the density of longleaf pine cubed, or

destruction = f*X*(1 - X)3,
(A1.2)

where f is the rate at which fire destroys hardwood.

Therefore,

dX/dt = r*X*(1-X) - f*X*(1 - X)3.
(A1.3)

A version of this model is shown in Fig. A1.1. This figure illustrates that, although the mortality due to fire is greater than the growth rate of hardwood at low densities of hardwood, as the density of hardwood increases, the hardwood growth rate exceeds its fire-produced mortality, leading to a net growth of hardwood. In this example, there are three densities of hardwood at which change in hardwood density is 0, which means that these points are stable.

Stability analysis

Stable states of this system are mixes of hardwood and longleaf pine that do not change over time. In the model, this corresponds to values of X at which dX/dt = 0. Equation A1.3 can therefore be rewritten as

dX/dt = X*(1 - X)*[r - f* (1 - X)2].
(A1.4)

Letting dX/dt = 0 reveals that X = 0 (all longleaf pine and no hardwood) and X = 1 (all hardwood and no longleaf pine) are both roots of dX/dt. The remainder of the equation could also equal zero, with r - f*(1 - X)2 = 0. This leaves two more roots: X = 1 ± sqrt( r/f). Because r and f both have positive values, X = 1 + sqrt(r/f) will be greater than zero and therefore invalid. If r > f, then 1 - sqrt(r/f) < 0 will also be invalid.

To assess the stability of these roots, I check the derivative of dX/dt with respect to X at each of the roots. This indicates the slope of dX/dt at a value of X. If the slope is negative, then an increase in X leads to a negative value of dX/dt, and a decrease in X leads to a positive value of dX/dt, both of which move X back toward the root. Therefore, if d2X/dtdX < 0, then a root is stable; if d2X/dtdX > 0, then it is unstable. In that case,

d2X/dtdX = r*(1 - 2*X) + f*(1 - X)2*(4*X - 1).
(A1.5)

At X = 0, d2X/dtdX = r - f. Therefore, if r < f, then X = 0 is stable. At X = 1, d2X/dtdX = -r. If r > 0, then X = 1 is stable. Because r will always be greater than zero, X = 1 is always stable.

At X = 1 - sqrt(r/f), the root is stable if 2r*[1 - sqrt(r/f)] < 0. Presuming that this root exists, i.e., that r < f, then this root will always be unstable, because, if the root exists, then r must be positive and 1-sqrt( r/f) must be between 0 and 1. Therefore, either two roots exist, X = 0, which is unstable, and X = 1, which is stable, or there are three roots: X = 0 and X = 1, which are both stable and separated by the unstable root X = 1 - sqrt(r/f).

To determine system behavior further requires considering specific values of r and f. In the absence of fire, hardwoods invade longleaf pine sites relatively quickly, so that even sites in which longleaf is relatively dominant (X = 0.1) will be half dominated by hardwoods within a decade. Hardwoods generally take more than a century to completely dominate a site. On the basis of this information, r was set to 0.075. The effect of fire on hardwood depends on two factors: the frequency of fire and the ability of fire to destroy hardwood. Consequently, f = proportion of hardwood killed by a fire/fire frequency. The proportion of hardwood killed by a fire was set to 1/2, and fire frequency was allowed to vary.

Bifurcation analysis of the system is performed by varying the parameters r and fire frequency and observing how the roots of the system change. This process is illustrated in Fig. A1.2, which shows how variation in fire frequency changes the roots on the system, and in Fig. A1.3, which shows how variation in the growth rate of hardwood, r, changes the behavior of the system.

The stable roots at X = 0 and X = 1 represent the system states that the system moves toward over time. The state that a particular state in the system moves toward depends on its position. States that are between the unstable state and a stable state are attracted to that stable state. For example, if the system is at X = 0.2 in Fig. A1.4 and fire frequency = 3, then the system will be attracted toward the longleaf state.

The unstable state divides values of X attracted from one stable state from those attracted to another. Another way of conceptualizing this division is that the distance to the unstable state from a given point, e.g., X = 0.2, represents how much disturbance a state can absorb before its resilience is exceeded, and it switches from being attracted from one stable state to another. By calculating how much change is required to cause a system state to be attracted to a different stable state, the resilience of that system state can be calculated. Resilience as measured by this distance is shown in Fig. A1.4 for fire frequency = 1. Visually, it is clear that the resilience of the longleaf pine stable state decreases as fire frequency (fires/year) increases, until the state loses all resilience and becomes unstable at slightly below 1 fire/7 yr. Similarly, in Fig. A1.5, as the hardwood growth rate increases, the resilience of the longleaf state decreases before vanishing completely.

Resilience can be visualized as a landscape by subtracting the integral of Eq. A1.3 from a constant. This method works, because Eq. A1.3 represents the forces (dX/dt) acting on the system state X. Given that these forces correspond to the slopes of the stability landscape, integrating them produces a stability landscape itself (Fig. A1.5). However, to ensure that stable points are minima rather than maxima, the integral must be multiplied by -1. Subtracting the values of this integral from a constant is done merely for graphical convenience. This adjustment does not matter, because the potential energy axis only relates different portions of the curve to one another; the values themselves have no intrinsic meaning.