Second derivatives of population growth rate: calculation and applications
Hal Caswell
INTRODUCTION
The rate of increase of a population is a function of the age- or stage-specific rates of survival, reproduction, growth, development, etc. (collectively, the “vital rates”). Sensitivity analysis explores this dependence by calculating the change in the rate of increase that would result from a change in any of the vital rates. Mathematically, this is done by computing the derivatives of the rate of increase with respect to each of the vital rates (Hamilton 1966, Demetrius 1969, Keyfitz 1971, Caswell 1978). Sensitivity analysis has become standard practice in demographic analysis. A nonexhaustive sample of recent applications, to a variety of taxa, includes McFadden (1991), Kalisz and McPeek (1992), Sivertown et al. (1992), McDonald (1993), Nault and Gagnon (1993), O’Connor (1993), Ang and DeWreede (1993), and Svensson et al. (1993). Because the rate of increase measures fitness in age-classified populations (Charlesworth 1980, 1993, Lande 1982), sensitivity analysis is also important in life history theory, where it measures the selection gradient on the vital rates (Lande 1982; for extensive discussions see Roff 1993, Stearns 1992).
The widespread use of sensitivity calculations raises the question of how the sensitivities themselves respond to changes in the vital rates. Stearns (1992:34) calls this the “situational sensitivity” problem: “The strength of selection acting on a given life history trait depends on the kind of life history in which the trait occurs. Sensitivity is situational.” Its solution requires calculation of the sensitivities of the sensitivities, i.e., of the second derivatives of the rate of increase with respect to the vital rates. This paper derives an easily computed formula for these second derivatives, and shows how to apply it to several important biological problems.
DERIVATION
Notation
I will use [x.sup.T] to denote the transpose of the vector x, and [x.sup.*] to denote the complex conjugate transpose. The complex conjugate of any number y will be denoted [Mathematical Expression Omitted]. The scalar product of two column vectors x and y is defined by = [y.sup.*]x.
Population growth is described by a population projection matrix A. I will use superscripts to distinguish the eigenvalues and eigenvectors, and subscripts to denote the elements of matrices and vectors. Thus the eigenvalues [[Lambda].sup.(i)] and the right and left eigenvectors [w.sup.(i)] and [v.sup.(i)], respectively, of A satisfy
A[w.sup.(i)] = [[Lambda].sup.(i)][w.sup.(i)] (1)
[v.sup.(i)*]A = [[Lambda].sup.(i)][v.sup.(i)*]. (2)
The eigenvalues are numbered so that [[Lambda].sup.(i)] is the dominant eigenvalue; it gives the population rate of increase. The eigenvectors can be scaled as desired; I will assume that they have been scaled so that
= [[Delta].sub.ij], (3)
where [[Delta].sub.ij] is the Kronecker delta
[[Delta].sub.ij] {0 if i [not equal to] j {1 if i = j. (4)
For technical reasons, the matrix A must be diagonalizable. A sufficient condition for this is that the eigenvalues of A are distinct, as they almost always are in demographic applications. If you encounter one of the unusual exceptions, I suggest checking to make sure that the matrix is diagonalizable.
Second derivatives
As is now well known (Caswell 1978), the sensitivity of any eigenvalue [[Lambda].sub.(m)] to changes in [a.sub.ij] is given by
[Mathematical Expression Omitted].
These sensitivities of each eigenvalue can be conveniently calculated and displayed in a sensitivity matrix
[Mathematical Expression Omitted]
whose (i, j) entry is
[Mathematical Expression Omitted].
Now let us focus attention on [[Lambda].sup.(1)], which is the rate of increase. Differentiating Eq. 5 again with respect to [a.sub.kl] gives
[Mathematical Expression Omitted].
(It is not necessary to include the scalar product that is often written in the denominator of Eq. 5; the derivative of this term with respect to [a.sub.kl] can be shown to be zero.)
From Eq. 8 it is clear that calculating the second derivative of [Lambda] requires the first derivatives of the eigenvectors [w.sup.(l)] and [Mathematical Expression Omitted]. These are given in Caswell (1980, 1989a) (Equation 6.93 in Caswell (1989a) contains a typographical error; the [Lambda] terms in the denominator of that formula should be replaced with their complex conjugates):
[Mathematical Expression Omitted]
[Mathematical Expression Omitted].
Substituting Eqs. 9 and 10 into Eq. 8 yields
[Mathematical Expression Omitted].
This can be rewritten in terms of the entries of the sensitivity matrices [S.sup.(m)]:
[Mathematical Expression Omitted].
Since the sensitivity matrices [S.sup.(m)] can be calculated directly from the eigenvectors, Eq. 12 is easy to evaluate, given only software capable of calculating all the eigenvalues and eigenvectors of a matrix. The Appendix gives a simple MATLAB routine for the calculation.
Note that the terms [[Lambda].sup.(l)] – [[Lambda].sup.(m)] appearing in the denominator of Eq. 12 mean that the formula would not work if [[Lambda].sup.(l)] were a repeated eigenvalue, since at least one of the denominators would be zero. This problem will not arise for calculations on the dominant eigenvalue, which is guaranteed to be distinct by the Perron-Frobenius theorem, but might occur if the formula (Eq. 12) is applied to one of the other eigenvalues, if it is repeated or nearly repeated.
As far as I know, Eqs. 11 and 12 for the second derivatives are new, although they are certainly implied by results in Desoer (1967). They have certainly never appeared in the ecological or demographic literature.
Completely different approaches to the second derivatives of eigenvalues have been taken by Cohen (1978) and Deutsch and Neumann (1984). Cohen derived the second derivatives of [Lambda] with respect to the diagonal elements of A, i.e., derivatives of the form [Mathematical Expression Omitted]. He proved that these second derivatives are nonnegative, so that [Lambda] is a convex function (a convex function has a positive second derivative; a concave function has a negative second derivative) of the diagonal elements of A. However, his approach does not generalize to the other second derivatives. Deutsch and Neumann (1984) derived formulae for all the second derivatives of [[Lambda].sup.(l)]. Kirkland and Neumann (1994) used these results to study the second derivatives of [[Lambda].sup.(l)] with respect to the entries in the first row and on the sub-diagonal of an age-classified Leslie matrix, i.e., derivatives of the form [Mathematical Expression Omitted]. They proved several theorems about the convexity or concavity of these derivatives (discussed below in comparison with a numerical example). The methods of Deutsch and Neumann (1984) and Kirkland and Neumann (1994) are completely different from those presented here. They rely on the “group inverse,” calculation of which is notoriously difficult (Campbell and Meyer 1979, but see Anstreicher and Rothblum 1987 for a promising simple approach). I have found that the results of Eq. 12 agree with those of the Deutsch and Neumann method.
An example
As an example, I calculated the second derivatives of the rate of increase for an age-classified population, using the matrix for the human population of the United States in 1965 (Keyfitz and Flieger 1968: 177, projection first row). This matrix, based on 5-yr age classes, has survival probabilities [P.sub.i] on the subdiagonal and fertilities [F.sub.i] in the first row, and zeros elsewhere. The second derivatives are plotted in three figures [ILLUSTRATION FOR FIGURE 1-3 OMITTED] showing the second derivatives with respect to pairs of fertilities, pairs of survival probabilities, and survival and fertility.
The second derivatives [[Delta].sup.2][Lambda]/[Delta][F.sub.i][Delta][F.sub.j] are shown in Fig. 1. (I am dropping the superscript on [[Lambda].sup.(l)] where no confusion seems likely to ensue.) Note that [Mathematical Expression Omitted], as Cohen’s (1978) results say it must be. The rate of increase [Lambda] is a convex function of [F.sub.i] for the first three age classes, and a concave function of [F.sub.i] for later ages.
These results have implications for life history theory. The classical explanation for the evolution of reproductive senescence is that the selective gradient [Delta][Lambda]/[Delta][F.sub.j] declines with age. The steepness of this decline can be measured by the slope
[b.sub.j] = [Delta][Lambda]/[Delta][F.sub.j+1] – [Delta][Lambda]/[Delta][F.sub.j]. (13)
How is this slope affected by a change in fertility at some other age i? The derivative of the slope with respect to [F.sub.i] is given by
[Delta][b.sub.j]/[Delta][F.sub.i] = [[Delta].sup.2][Lambda]/[Delta][F.sub.i][Delta][F.sub.j+1] – [[Delta].sup.2][Lambda]/[Delta][F.sub.i][Delta][F.sub.j]. (14)
In Fig. 1, the derivative [[Delta].sup.2][Lambda]/[Delta][F.sub.i][Delta][F.sub.j] decreases with both i and j; this means that [Delta]b/[Delta][F.sub.i] [less than] 0. Thus, increasing fertility at any age makes the selective pressure on fertility fall off more rapidly with age.
The second derivatives [[Delta].sup.2][Lambda]/[Delta][P.sub.i][Delta][P.sub.j] show a different pattern [ILLUSTRATION FOR FIGURE 2 OMITTED]. Population growth rate is a concave function of [P.sub.i] (i.e., [Mathematical Expression Omitted]), especially at young ages. At old ages, [Mathematical Expression Omitted] is very small, i.e., changes in survival at later ages have little or no effect on the sensitivity of [Lambda] to further changes in survival. At early ages, the survival second derivatives show an interesting alternation in sign. An increase in [P.sub.i] reduces [Delta][Lambda]/[Delta][P.sub.i], but increases [Delta][Lambda]/[Delta][P.sub.j] for j [not equal to] i. That is, as survival at one age increases, further increases in survival at that age become less important, and increases in survival at other ages become more important.
The results for the cross partial derivatives [[Delta].sup.2][Lambda]/[Delta][F.sub.i][Delta][P.sub.j] are shown in Fig. 3. Consider the curve labeled j = 1. It shows that an increase in [P.sub.l] increases the sensitivity to fertility through age 7, and reduces it for ages after 7. The curve labeled j = 5 shows that an increase in [P.sub.5] reduces the sensitivity to fertility changes for ages 1-5, increases it dramatically for age 6, and to a lesser degree for all subsequent ages. On the other hand, looking at the curves at, e.g., the point i = 4 shows that an increase in fertility at age 4 reduces the sensitivity to survival at age 4 and all subsequent ages, and increases the sensitivity to survival at ages 1-3.
This pattern reflects the way that survival and fertility combine to influence [Lambda]. A change in fertility at one age will have its biggest effect when individuals survive to realize that fertility, so survival at earlier ages becomes more important. Similarly, an increase in survival at any age will have its biggest effect when the age class to which the individual survives has increased fertility. The relative importance of fertility at earlier and much later ages is reduced.
I have found the same patterns in age-classified data from a number of other species, including long-tailed mannakins (Chiroxiphia linearis; McDonald 1993) and the polychaete Streblospio benedicti (Levin et al. 1987). Kirkland and Neumann (1994) prove that some of these patterns are general to all age-classified Leslie matrices. They show that the second derivatives [Mathematical Expression Omitted] decrease with age, becoming negative no later than halfway through the life cycle. The second derivatives with respect to survival, [Mathematical Expression Omitted] may be positive and decrease with age for early ages, but this seems to be difficult to produce and was not seen in any of the examples examined here. More commonly, the derivatives [Mathematical Expression Omitted], in which case they must increase (or at least not decrease) with age.
These results suggest that it may be possible to produce useful generalizations about the way that life histories mold the age-specific patterns of sensitivity.
APPLICATIONS
The second derivatives of [Lambda] have many potential ecological applications. This section describes four of them: the perturbation analysis of elasticity patterns (which may be useful in conservation and management), the sensitivity analysis of stochastic growth rates, the inclusion of higher order terms in the decomposition analysis of life table response experiments, and the analysis of nonlinear and correlational selection. More detailed explorations of these applications will be presented elsewhere.
Perturbation analysis of elasticities
The elasticities, or proportional sensitivities, of [Lambda] to changes in the entries of A are given by
[e.sub.ij] = [a.sub.ij]/[Lambda] [Delta][Lambda]/[Delta][a.sub.ij]. (15)
The elasticities sum to one, and can be interpreted as giving the contribution of the corresponding matrix element to the rate of increase (DeKroon et al. 1986). Elasticity analysis has been applied in evaluating management strategies in conservation biology: all other things being equal, management efforts should focus on vital rates with large elasticities (e.g., Crouse et al. 1987, Olesiuk et al. 1990, Brault and Caswell 1993, Crowder et al. 1994, Doak et al. 1994, Heppell et al. 1994). Thus it is of interest to know how the elasticities, and the corresponding management decisions, would respond to changes in the [a.sub.ij].
The sensitivity of the elasticity to changes in the matrix elements can be calculated by differentiating Eq. 15:
[Delta][e.sub.ij]/[Delta][a.sub.kl] = [a.sub.ij]/[Lambda] [[Delta].sup.2][Lambda]/[Delta][a.sub.ij][Delta][a.sub.kl] – [a.sub.ij]/[[Lambda].sup.2] [Delta][Lambda]/[Delta][a.sub.ij] [Delta][Lambda]/[Delta][a.sub.kl] + [[Delta].sub.ik][[Delta].sub.jl]/[Lambda] [Delta][Lambda]/[Delta][a.sub.ij], (16)
where [[Delta].sub.ik] and [[Delta].sub.jl] given by Eq. 4.
As an example, consider one of the population projection matrices derived by Doak et al. (1994) for the threatened desert tortoise (Gopherus agassizii). Their model includes eight size classes, with a projection matrix
[Mathematical Expression Omitted]
The corresponding elasticity matrix is
[Mathematical Expression Omitted]
By far the most important coefficients to [Lambda] are those on the main diagonal, which together account for 69.9% of [Lambda]; the growth probabilities on the subdiagonal account for 25.8% and the fertilities. in the first row for only 4.3%. The largest elasticity in the life cycle ([e.sub.77] = 0.279) is that for the probability of surviving and remaining in stage 7. Doak et al. (1994) conclude from this that “any management policy having a substantial impact on the survival of [stage 7] females is likely to produce the most dramatic change in population trends.”
A sensitivity analysis of [e.sub.77] can shed light on why it is so large, i.e., on the aspects of the life history that make this transition so important for the desert tortoise. Let D be defined as a matrix whose (i, j) entry is [d.sub.ij] = [Delta][e.sub.77]/[Delta][a.sub.ij]. Applying Eq. 16 gives
[Mathematical Expression Omitted]
Changes in fertility (the first row) have little impact on [e.sub.77], although what effect they have is negative. Increasing the probability of remaining in stage 7 and the probability of growing into stage 8 would have the biggest effects (positive and negative, respectively) on [e.sub.77]. Reducing all of the diagonal and subdiagonal terms for the other stages would also increase [e.sub.77]. Thus [e.sub.77] will be large when growth into stage 7 is slow, as is growth out of stage 7, while the probability of surviving and remaining in stage 7 is high.
This analysis can be extended. The diagonal and sub-diagonal elements in size-classified models like Eq. 17 conflate the probabilities of growth and survival. These effects can be separated by parameterizing the model in terms of the stage-specific survival probability [[Sigma].sub.i] and the stage-specific growth probability [[Gamma].sub.i]. Denote the diagonal and subdiagonal elements of Eq. 17 by [P.sub.i] and [G.sub.i], respectively; these terms can be written
[P.sub.i] = [[Sigma].sub.i](1 – [[Gamma].sub.i]) (20)
[G.sub.i] = [[Sigma].sub.i][[Gamma].sub.i]. (21)
The elasticities of [[Sigma].sub.i] give the proportional effects on [Lambda] of a change in survival probability alone, without the confounding effect of growth.
Let e([[Sigma].sub.i]) denote the elasticity of [Lambda] with respect to [[Sigma].sub.i]. Doak et al. (1994) found that e([[Sigma].sub.7]) was the largest of all the survival elasticities for the desert tortoise. Thus a sensitivity analysis of e([[Sigma].sub.7]) would be useful in seeing how the life history molds this elasticity. This calculation is easy to do, because for models of this form,
[Delta][Lambda]/[Delta][[Sigma].sub.i] = (1 – [[Gamma].sub.i]) [Delta][Lambda]/[Delta][P.sub.i] + [[Gamma].sub.i] [Delta][Lambda]/[Delta][G.sub.i]. (22)
Multiplying both sides by [[Sigma].sub.i]/[Lambda] and using Eqs. 20 and 21, it is not hard to show that
e ([[Sigma].sub.i]) = [e.sub.ii] + [e.sub.i+1,i], (23)
which means that it is easy to calculate the sensitivity of e([[Sigma].sub.i]) to changes in the entries of A:
[Delta]e([[Sigma].sub.i])/[Delta][a.sub.kl] = [Delta][e.sub.ii]/[Delta][a.sub.kl] + [Delta][e.sub.i+1,i]/[Delta][a.sub.kl] (24)
where the right-hand side is calculated using Eq. 16.
The results for [Delta]e([[Sigma].sub.7])/[Delta][a.sub.ij] in the case of the desert tortoise are very similar to the results for [Delta][e.sub.77]/[Delta][a.sub.ij] given in D and are not presented here. Again, slow growth and low survival in preceding stages combined with high survival and a long duration in stage 7 yield a high elasticity for survival in that stage.
These results for the desert tortoise are remarkably similar to those obtained from the loggerhead sea turtle (Caretta caretta) model of Crouse et al. (1987) and Crowder et al. (1994). The diagonal, subdiagonal, and first-row elasticities of the Crowder et al. (1994) model sum to 71.9, 22.3, and 5.8% of [Lambda], respectively. The sensitivities of these elasticities to changes in the matrix entries are similar to the pattern in the desert tortoise (H. Caswell, unpublished). The elasticity analysis of the loggerhead sea turtle model (Crouse et al. 1987) was influential in the 1989 decision to require the use of turtle excluder devices (TEDs) on shrimp trawls. The fact that the elasticity patterns in these two species respond in the same way to the structure of the life cycle raises the hope that there may be general explanations of how life histories determine elasticity patterns. Such explanations would be useful in life history theory (cf. Silvertown et al. 1992) as well as conservation biology.
Sensitivity of stochastic growth rate
Population growth in a stochastic environment can be described by a product of random matrices
n(t) = [A.sub.t-1][A.sub.t-2] … [A.sub.0]n(0), (25)
where the [A.sub.t] are generated by a stochastic process. Under a wide range of conditions on this process, such a population is characterized by a stochastic growth rate [[Lambda].sub.s] satisfying
[Mathematical Expression Omitted]
where N(t) is population size (Cohen 1976, 1979, Tuljapurkar and Orzack 1980; see Tuljapurkar 1989, 1990a for a thorough discussion). This quantity gives the almost certain long-run rate of increase of any realization of the population process, and is the appropriate measure of fitness in a stochastic environment.
Tuljapurkar (1982) derived a useful approximation for log [[Lambda].sub.s] for cases where the environment produces relatively small amounts of variability around a mean set of vital rates. Suppose that the matrix [A.sub.t] = A + [B.sub.t], where A is the mean matrix and [B.sub.t] is a random Xmatrix with mean zero. Suppose that the mean matrix A has a dominant eigenvalue [Lambda], and that the covariances of the entries of [B.sub.t] are given by C(ij, kl). For the special case of temporally independent variation, Tuljapurkar’s approximation is
[Mathematical Expression Omitted],
where C(ij, kl) is the covariance of [a.sub.ij] and [a.sub.kl]. The approximation is quite accurate, even for moderate levels of variation in the vital rates (e.g., Orzack and Tuljapurkar 1989).
Tuljapurkar’s approximation to the stochastic growth rate has been used to study the evolution of life history characteristics (e.g., Orzack 1985, Orzack and Tuljapurkar 1989, Tuljapurkar 1990b). In such studies, it would be useful to know how log [[Lambda].sub.s] changes in response to changes in the mean population matrix A; this is a stochastic analog of the familiar deterministic sensitivity analysis. (Tuljapurkar [1990a] gives a numerical procedure to calculate these sensitivities in cases where log [[Lambda].sub.s] is calculated without using Approximation 27; the results of the two formulae will be compared elsewhere.) Because log [[Lambda].sub.s] depends on the first derivatives of [Lambda], calculating its sensitivity requires the second derivatives of [Lambda] with respect to the entries of A.
Denoting the sensitivity matrix of [Lambda] as S, differentiating Approximation 27 yields
[Mathematical Expression Omitted].
The sensitivity of log [[Lambda].sub.s] depends on [Lambda], on the first and second derivatives of [Lambda], and on the covariation of the vital rates. When there is no environmental variation, so that C = 0, then log [[Lambda].sub.s] = log [Lambda], and from Eq. 28 it follows that
[Delta] log [[Lambda].sub.s]/[Delta][a.sub.mn] = [s.sub.mn]/[Lambda]
(29)
as it should.
The calculation of Eq. 28 can be simplified as follows. Define vec A as the vector obtained from the matrix A by stacking its columns one above the other. Write the derivative of the sensitivity matrix S with respect to [a.sub.mn] as a matrix T. Then
[Delta] log [[Lambda].sub.s]/[Delta][a.sub.mn]
= [S.sub.mn]/[Lambda][1 + [(vec S).sup.T] C vec S/[[Lambda].sup.2]]
– [[(vec T).sup.T] vec S + [(vec S).sup.T] C vec T/2[[Lambda].sup.2]]. (30)
These results permit a sensitivity analysis of the stochastic growth rate, using the small noise approximation, for the temporally independent case. When the environment is autocorrelated, there is an additional term in the approximation, involving the autocorrelations of the vital rates and the sensitivities of all the eigenvalues of A, rather than just the dominant eigen-value. The derivative of this term with respect to [a.sub.mn] can also be calculated, involving the second derivatives of all the eigenvalues; this problem will be explored elsewhere.
Higher order terms in LTRE analysis
Life table response experiments (LTREs) are manipulative experiments or comparative observations in which the dependent variable is a complete set of vital rates (loosely speaking, a life table; Caswell 1989b, 1996). Such experiments are often summarized by using the rate of increase [Lambda] to integrate treatment effects on survival and reproduction throughout the life cycle.
Treatment effects on [Lambda] can be decomposed into contributions from the effects on each of the vital rates (Caswell 1989b); this decomposition makes it possible to pinpoint where in the life history the treatment has its greatest impact. The decomposition uses a first-order linear approximation to the effect on [Lambda]. Consider the simplest case of two treatments, producing rates of increase that differ by [Delta][Lambda]. This effect can be written
[Delta][Lambda] [approximately equal to] [summation over ij] [Delta][a.sub.ij] [Delta][Lambda]/[Delta][a.sub.ij], (31)
where [Delta][a.sub.ij] is the treatment effect on the entry [a.sub.ij] of the population matrix, and the partial derivatives are evaluated at an appropriate point midway between the two matrices. Each term in the summation gives the contribution of the treatment effects on one vital rate to the overall effect on [Lambda]. (More complicated factorial fixed effects designs are analyzed in Caswell 1989a, b).
In a number of applications of the method (Levin et al. 1987, 1996, Caswell 1989b, Levin and Huggett 1990, Walls et al. 1991, Silva et al. 1991, Brault and Caswell 1993) Approximation 31 has been found to be quite accurate. However, in some cases it may fail. This means that second and higher order terms are too important to be neglected. The second-order terms can be included by using the second derivatives:
[Mathematical Expression Omitted].
Comparing the accuracy of the first- and second-order approximations can give an indication of the importance of the second derivatives. Specific second-order contributions will identify important nonlinearities in the response of [Lambda] to the treatments.
Nonlinear selection
Lande and Arnold (1983) characterized multivariate natural selection in terms of the directional selection gradients [[Beta].sub.i]. These gradients are measured by the partial regressions of fitness on trait i, holding all other traits constant. The selection gradient is given locally by the first partial derivative of fitness with respect to trait i. Lande (1982) applied this theory to selection on life history traits in age-structured populations, measuring fitness by [Lambda] (or r = log [Lambda]) and the selection gradient by the first derivative of [Lambda] with respect to age-specific vital rates.
Phillips and Arnold (1989) later extended this theory to the second order, treating fitness as a quadratic function of the traits (or as a function that can be approximated by a quadratic function). In the univariate case, this “nonlinear selection” includes the traditional categories of “stabilizing” and “disruptive” selection. Nonlinear selection is characterized by the quadratic selection gradient [[Gamma].sub.ij], which is locally the second derivative of fitness with respect to traits i and j. Thus the second derivatives [[Delta].sup.2][Lambda]/[Delta][a.sub.ij][Delta][a.sub.kl] give the quadratic selection gradients on [a.sub.ij] and [a.sub.kl].
Nonlinear selection changes the variance of the trait distribution, beyond the changes that which may be caused by linear selection. Phillips and Arnold (1989) refer to the case where [Mathematical Expression Omitted] as concave selection; the result is a decrease in variance, analogous to stabilizing selection in the univariate case. Convex selection, where [Mathematical Expression Omitted], increases the variance. The second derivatives with respect to two different traits produces correlational selection; if [[Delta].sup.2][Lambda]/[Delta][a.sub.ij][Delta][a.sub.kl] [less than] 0 ([greater than] 0) selection leads to a decrease (increase) in the correlation between [a.sub.ij] and [a.sub.kl].
The results of the age-classified example in Figs. 13 can be interpreted in terms of nonlinear selection on survival probability and fertility:
There is convex selection on fertility in the first three age classes, and concave selection thereafter. This implies that selection will increase the variance in fertility at young ages and reduce it at later ages. The results of Kirkland and Neumann (1994) show that this pattern is general to all age-classified models.
Selection on survival probability is concave, especially at young ages. Thus selection will reduce the variance in survival probability at young ages.
There is negative correlational selection on survival at any age and fertility at the same or earlier ages. Such selection should lead to negative genetic correlations between fertility and subsequent survival, an expected pattern of “costs of reproduction.”
There is positive correlational selection on survival at any age and fertility for the following few age classes.
There is negative correlational selection on survival at early ages and fertility at much later ages.
These are the first numerical results on nonlinear selection gradients on life history traits. My limited exploration of age-classified models for different species, and the results of Kirkland and Neumann (1994), suggest that some of these patterns may be general properties of age-classified life cycles. They deserve more investigation.
DISCUSSION
The first derivatives of [Lambda] with respect to changes in the vital rates have found a place as a standard method in demographic analysis. Their properties are well understood and they are familiar in both theoretical analyses and practical applications. The second derivative formulae presented here have a number of potential applications, but little or no accumulated experience. The formulae themselves are not revealing, and, except for the results of Cohen (1978) and of Kirkland and Neumann (1994) on convexity, little is known analytically about the second derivatives in general.
One misconception bears mention here. Real demographic data are never perfect. I have heard the suggestion that the sophistication of this analysis outstrips the possible accuracy of ecological data. This is not so. If demographic data are too noisy to be trusted, the place to stop the analysis is before calculating [Lambda]. The sensitivities of [Lambda] are implicit in the calculation of itself. Indeed, one of the potential uses of the sensitivities is to evaluate how reliable the estimate of [Lambda] is, given the uncertainty in the data. Similarly, if sensitivities or elasticities are reported, the ingredients of the second-derivative calculations are implicitly present, and the second derivatives can be used to evaluate the reliability of the sensitivities and elasticities.
The examples of applications of the second derivatives in this paper are typical of the uses of [Lambda] – as a measure of population success, as a measure of fitness, as a measure of response to environmental factors, and as a component of stochastic models – but are not the only possible ones.
ACKNOWLEDGMENTS
I thank Joel Cohen, Solange Brault, Selina Heppel, Steve Kirkland, and several anonymous reviewers for comments. This research was supported by NSF Grant DEB 9211945, EPA Grant R818408-01-0, and ONR Grant URIP N00014-92-J-1527. Woods Hole Oceanographic Institution Contribution 8608.
LITERATURE CITED
Ang, P. O., and R. E. DeWreede. 1993. Simulation and analysis of a Fucus distichus (Phaeophyceae, Fucales) population. Marine Ecology Progress Series 93:253-265.
Anstreicher, K. M., and U. G. Rothblum. 1987. Using Gauss-Jordan elimination to compute the index, generalized null-spaces, and Drazin inverse. Linear Algebra and its Applications 85:221-239.
Brault, S., and H. Caswell. 1993. Pod-specific demography of killer whales (Orcinus orca). Ecology 74:1444-1454.
Campbell, S. L., and C. D. Meyer, Jr. 1979. Generalized inverses of linear transformations. Pitman, London, England.
Caswell, H. 1978. A general formula for the sensitivity of population growth rate to changes in life history parameters. Theoretical Population Biology 14:215-230.
—–. 1980. On the equivalence of maximizing reproductive value and maximizing fitness. Ecology 61:19-24.
—–. 1989a. Matrix population models: construction, analysis, and interpretation. Sinauer, Sunderland, Massachusetts, USA.
—–. 1989b. The analysis of life table response experiments. I. Decomposition of treatment effects on population growth rate. Ecological Modelling 46:221-237.
—–. 1996. Demography meets ecotoxicology: untangling the population level effects of toxic substances. Pages 255-292 in M. Newman and C. Jagoe, editors. Ecotoxicology: a hierarchical treatment. Lewis, Boca Raton, Florida, USA.
Charlesworth, B. 1980. Evolution in age-structured populations. Cambridge University Press, Cambridge, England.
—–. 1993. Natural selection on multivariate traits in age-structured populations. Proceedings of the Royal Society, London, B 251:47-52.
Cohen, J. E. 1976. Ergodicity of age structure in populations with Markovian vital rates. I. Countable states. Journal of the American Statistical Association 71:335-339.
—–. 1978. Derivatives of the spectral radius as a function of non-negative matrix elements. Mathematical Proceedings of the Cambridge Philosophical Society 83:183-190.
—–. 1979. Ergodic theorems in demography. Bulletin of the American Mathematical Society 1:275-295.
Crouse, D. T., L. B. Crowder, and H. Caswell. 1987. A stage-based population model for loggerhead sea turtles and implications for conservation. Ecology 68:1412-1423.
Crowder, L. B., D. T Crouse, S. S. Heppell, and T. H. Martin. 1994. Predicting the impact of turtle excluder devices on loggerhead sea turtle populations. Ecological Applications 4:437-445.
DeKroon, H., A. Plaisier, J. van Groenendael, and H. Caswell. 1986. Elasticity: the relative contribution of demographic parameters to population growth rate. Ecology 67:1427-1431.
Demetrius, L. 1969. The sensitivity of population growth rate to perturbations in the life cycle components. Mathematical Biosciences 4:129-136.
Desoer, C. A. 1967. Perturbations of eigenvalues and eigen-vectors of a network. Pages 8-11 in the Fifth Annual Allerton Conference on Circuit and System Theory. University of Illinois, Urbana, Illinois, USA.
Deutsch, E., and M. Neumann. 1984. Derivatives of the Perron root of an essentially nonnegative matrix and the group inverse of an M-matrix. Journal of Mathematical Analysis and Applications 102:1-29.
Doak, D., P. Kareiva, and B. Klepetka. 1994. Modeling population viability for the desert tortoise in the western Mojave desert. Ecological Applications 4:446-460.
Hamilton, W. D. 1966. The moulding of senescence by natural selection. Journal of Theoretical Biology 12:12-45.
Heppell, S. S., J. R. Walters, and L. B. Crowden 1994. Evaluating management alternatives for red-cockaded woodpeckers: a modeling approach. Journal of Wildlife Management 58:479-487.
Horn, R. A., and C. A. Johnson. 1985. Matrix analysis. Cambridge University Press, Cambridge, England.
Kalisz, S., and M. A. McPeek. 1992. Demography of an age-structured annual: resampled projection matrices, elasticity analyses, and seed bank effects. Ecology 73:1082-1093.
Keyfitz, N. 1971. Linkages of intrinsic to age-specific rates. Journal of the American Statistical Association 66:275-281.
Keyfitz, N., and W. Flieger. 1968. World population. University of Chicago Press, Chicago, Illinois, USA.
Kirkland, S. J., and M. Neumann. 1994. Convexity and concavity of the Perron root and vector of Leslie matrices with applications to a population model. SIAM Journal of Matrix Analysis and Applications 15:1092-1107.
Lande, R. 1982. A quantitative genetic theory of life history evolution. Ecology 63:607-615.
Lande, R., and S. J. Arnold. The measurement of selection on correlated characters. Evolution 37:1210-1226.
Levin, L. A., H. Caswell, T. Bridges, D. Cabrera, G. Plaia, and C. DiBacco. 1996. Demographic responses of estuarine polychaetes to sewage, algal, and hydrocarbon contaminants. Ecological Applications, in press.
Levin, L. A., H. Caswell, K. D. DePatra, and E. L. Creed. 1987. The life table consequences of larval development mode: an intraspecific comparison of planktotrophy and lecithotrophy. Ecology 68:1877-1886.
Levin, L. A., and D. V. Huggett. 1990. Implications of alternative reproductive modes for seasonality and demography in an estuarine polychaete. Ecology 71:2191-2208.
McDonald, D. B. 1993. Demographic consequences of sexual selection in the long-tailed manakin. Behavioral Ecology 4:297-309.
McFadden, C. S. 1991. A comparative demographic analysis of clonal reproduction in a temperate soft coral. Ecology 72:1849-1866.
Nault, A., and D. Gagnon. 1993. Ramet demography of Allium tricoccum, a spring ephemeral, perennial forest herb. Journal of Ecology 81:101-119.
O’Connor, T. G. 1993. The influence of rainfall and grazing on the demography of some African savanna grasses: a matrix modelling approach. Journal of Applied Ecology 30:119-132.
Olesiuk, P. F., M. A. Bigg, and G. M. Ellis. 1990. Life history and population dynamics of resident killer whales (Orcinus orca) in the coastal waters of British Columbia and Washington State. Report of the International Whaling Commission, Special Issue 12:209-243.
Orzack, S. H. 1985. Population dynamics in variable environments V. The genetics of homeostasis revisited. American Naturalist 125:550-572.
Orzack, S. H., and S. Tuljapurkar. 1989. Population dynamics in variable environments. VII. The demography and genetics of iteroparity. American Naturalist 133:901-923.
Phillips, P. C., and S. J. Arnold. 1989. Visualizing multivariate selection. Evolution 43:1209-1222.
Roff, D. A. 1992. The evolution of life histories. Chapman and Hall, New York, New York, USA.
Silva, J. G., J. Raventos, H. Caswell, and M. C. Trevisan. 1991. Population responses to fire in a tropical savanna grass Andropogon semiberbis: a matrix model approach. Journal of Ecology 79:345-356.
Silvertown, J., M. Franco, and K. McConway. 1992. A demographic interpretation of Grime’s triangle. Functional Ecology 6:130-136.
Stearns, S. C. 1992. The evolution of life histories. Oxford University Press, Oxford, England.
Svensson, B. M., B. A. Carlsson, P. S. Karlsson, and K. O. Nordell. 1993. Comparative long-term demography of three species of Pinguicula. Journal of Ecology 81:635-645.
Tuljapurkar, S. D. 1982. Population dynamics in variable environments. III. Evolutionary dynamics of r-selection. Theoretical Population Biology 21:141-165.
—–. 1989. An uncertain life: demography in random environments. Theoretical Population Biology 35:227-294.
—–. 1990a. Population dynamics in variable environments. Springer-Verlag, New York, New York, USA.
—–. 1990b. Delayed reproduction and fitness in variable environments. Proceedings of the National Academy of Sciences (USA) 87:1139-1143.
Tuljapurkar, S. D., and S. H. Orzack. 1980. Population dynamics in variable environments. I. Long-run growth rates and extinction. Theoretical Population Biology 18:314-342.
Walls, M., H. Caswell, and M. Ketola. 1991. Demographic consequences of Chaoborus-induced defenses in Daphnia pulex: a sensitivity analysis. Oecologia 87:43-50.
Wilkinson J. H. 1965. The algebraic eigenvalue problem. Oxford University Press, Oxford, England.
APPENDIX
A MATLAB PROGRAM FOR SECOND DERIVATIVES
The following simple program takes as arguments a matrix A and a pair of indices (k, l). It returns a matrix whose (i, j) entry is [[Delta].sup.2][[Lambda].sup.(1)]/[Delta][a.sub.ij][Delta][a.sub. kl].
For information on MATLAB, contact The MathWorks, 24 Prime Park Way, Natick MA 01760.
function d2= secder(a,K,L)
% function d2=secder(a,K,L)
%
% returns a matrix d2 whose (i,j) entry is the second
% derivative of the dominant eigenvalue of the matrix a % with respect to a(i,j) and a(K,L) % n=length(a(1,:));
% find eigenvalues and right eigenvectors [w,d] =eig(a); lambda=diag(d); imax= find(lambda= = max(lambda));
% arrange eigenvalues so that dominant eigenvalue is first if imax[greater than]1 & imax [less than]n
lambda = lambda([imax 1:imax- 1 imax + 1:n]); w=[w(:,imax) w(:,1:imax-1) w(:,imax + 1:n)];
elseif imax= =n
lambda=lambda([n 1:n-1]); w=[w(:,n) w(:,1:n-1)]; end % if
% find left eigenvectors
v=inv(conj(w));
% generate sensitivity matrices; store as columns of svec() for i=1:n
vi= v(i,:).’; %note use of .’ for nonconjugate transpose
wi=w(:,i);
senmat=conj([vi).sup.*]wi.’; %note use of .’ for nonconjugate %transpose
svec(:, i) = senmat(:);
end %for i
% next calculate the second derivatives
sl = svec(:,1);
% scale the sensitivities by lambda(1)-lambda(m)
for m=2:n
scalesens(:,m- 1) = svec (:,m)/(lambda(1)-lambda(m));
end % for m
vecsum = sum(scalesens.’);
% now find the array of second derivatives
for i=1:n
for j=1:n
% the following two indices find the necessary entries in the
% vectorized sensitivity matrices x1=(L-1)*n+1 + (i-1);
x2=(j-1)*n+1 + (K-1);
d2(i,j)=s1(x1)*vecsum(x2) + s1(x2)*vecsum(x1);
end % for j
end % for i
% remove imaginary parts, which are zero for a real dominant
% eigenvalue
d2=real(d2),
end
COPYRIGHT 1996 Ecological Society of America
COPYRIGHT 2004 Gale Group
