UFDC Home  Search all Groups  UF Institutional Repository  UF Institutional Repository  UF Theses & Dissertations   Help 
Material Information
Thesis/Dissertation Information
Subjects
Notes
Record Information

Full Text 
PAGE 1 1 FOUR PERMUTATIONS OF GROWTH, REPRODUCTION, AND SURVIVAL By MOLLIE ELIZABETH BROOKS A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2012 PAGE 2 2 2012 Mollie Elizabeth Brooks PAGE 3 3 To my parents and Grandma Brooks PAGE 4 4 ACKNOWLEDGMENTS I thank my advisor s and committee Ben Bolker, Mary Christman, Bob Holt, Colette St. Mary, and Emilio Bruna for great advice. M y collaborators Yoh Iwasa, Emilio Bruna, and Mike McCoy have contributed to the work in this dissertation Jose Ponciano Jake Ferguson Brett Melbourne, and Sebastian Schreiber helped with stimulating discussions. I thank teachers Karen Donathan, Dave Fournier, John Sibert, Mark Maunder, and Steve Martel. I thank the National Science Foundation (NSF) Quantitative Ecology, Evolution, and Environment Integrative Graduate Education and Research Traineeship NSF East Asia and Pacific Summer Institutes, and the Automatic Differentiation Model Builder Foundation for funding. PAGE 5 5 TABLE OF CONTENTS page ACKNOWLEDGMENTS ................................ ................................ ................................ ... 4 LIST OF TABLES ................................ ................................ ................................ ............. 8 LIST OF FIGURES ................................ ................................ ................................ ........... 9 ABSTRACT ................................ ................................ ................................ .................... 10 CHAPTER 1 INTRODUCTORY REMARKS ................................ ................................ ................. 12 2 SIZE DEPENDENT SEX CHANGE CAN BE THE ESS WITHOUT ANY SIZE ADV ANTAGE OF REPRODUCTION WHEN MORTALITY IS SIZE DEPENDENT ................................ ................................ ................................ ........... 14 Introduction ................................ ................................ ................................ .............. 14 Motivational Case of Pollen Limitation ................................ ................................ ..... 16 Model ................................ ................................ ................................ ....................... 17 Dynamic Programming ................................ ................................ ...................... 20 Direct Search for the Optimal Solution ................................ .............................. 24 Parameter Dependence ................................ ................................ .................... 25 Results and Discussion ................................ ................................ ............................ 26 Sex Change in the Absence of Size Dependent Reproduction ......................... 26 Patterns of Sex Ratio ................................ ................................ ......................... 29 Unimportant Parameters ................................ ................................ ................... 30 Model Generality ................................ ................................ ................................ ...... 31 3 DETECTING VARIATION IN GROWTH RATES WITHOUT MARKING INDIVIDUALS ................................ ................................ ................................ .......... 36 Introduction ................................ ................................ ................................ .............. 36 Methods ................................ ................................ ................................ ................... 40 Derivation ................................ ................................ ................................ .......... 41 Simulating Growth Data ................................ ................................ ..................... 42 Model Fitting ................................ ................................ ................................ ...... 42 Procedures ................................ ................................ ................................ ........ 43 Repeated Measures ................................ ................................ .......................... 44 Size Dependent Mortality Simulations ................................ ............................... 44 Results and Discussion ................................ ................................ ............................ 45 Compari son of Variance Pattern and Repeated Measures Approaches .......... 46 Size Dependent Mortality Causes Bias ................................ ............................. 46 Conclusion ................................ ................................ ................................ ............... 47 PAGE 6 6 4 DOES FRAGMENTATION ALTER VITAL RATES AND INCREASE HETEROGENeiTY? ................................ ................................ ................................ 52 Introduction ................................ ................................ ................................ .............. 52 Effects of Fragmentation ................................ ................................ ................... 54 Hierarchical Models ................................ ................................ ........................... 55 Aims ................................ ................................ ................................ ................... 55 Methods ................................ ................................ ................................ ................... 56 Data Description ................................ ................................ ................................ 56 Model Formulations ................................ ................................ ........................... 56 Growth ................................ ................................ ................................ ........ 57 Reproduction ................................ ................................ ............................... 58 Survival ................................ ................................ ................................ ....... 59 Parameter Estimation and Uncertainty ................................ .............................. 59 Allowing variance to differ by habitat ................................ .......................... 59 Confidence interval estimation ................................ ................................ .... 60 Maximum likelihood estimation and MCMC Sampling using AD Model Builder ................................ ................................ ................................ ...... 61 Model selection ................................ ................................ ........................... 63 Assessing Spatial Patterns ................................ ................................ ................ 63 Ignoring Demographic Heterogeneity ................................ ................................ 64 Correlation Among Vital Rates ................................ ................................ .......... 64 Results and Discussion ................................ ................................ ............................ 65 Fixed Effects ................................ ................................ ................................ ...... 65 Random Effects ................................ ................................ ................................ 66 Ignoring Demographic Heterogeneity ................................ ................................ 67 Correlation of Vital Rate Random Deviates ................................ ....................... 67 Spatial Pattern ................................ ................................ ................................ ... 68 Conclusion ................................ ................................ ................................ ............... 69 5 SCALING UP VARIATION ................................ ................................ ....................... 76 Introduction ................................ ................................ ................................ .............. 76 Study System ................................ ................................ ................................ .... 78 Predicted Effects of Heterogeneity ................................ ................................ .... 79 Changes to the mean population growth rate ................................ ............. 79 Changes to the variance of the population growth rate .............................. 80 Methods ................................ ................................ ................................ ................... 81 Data Description ................................ ................................ ................................ 81 Parameter Estimation ................................ ................................ ........................ 81 Covariance Among Vital Rates ................................ ................................ .......... 82 Population Growth Simulations ................................ ................................ ......... 83 Decomposing Effects of Heterogeneity on Population Growth .......................... 85 Scaling Up Variation Among Individuals in Inflorescence Rates ....................... 88 S imulations of Fragmentation ................................ ................................ ............ 90 Results and Discussion ................................ ................................ ............................ 90 Effects Of Heterogeneity On Average Population Growth ................................ 90 PAGE 7 7 Effects Of Heterogeneity On Stochastic Population Growth ............................. 91 Effects Of Heterogeneity On Variance In Population Growth ............................ 92 Inflorescence Production ................................ ................................ ................... 93 Transition From Continuous Forest To Fragment ................................ ............. 94 Conclusions ................................ ................................ ................................ ............. 94 6 CONCLUDING REMARKS ................................ ................................ .................... 100 APPENDIX A DYNAMIC PROGRAMMING ................................ ................................ ................. 101 B MAXIMIZTION ................................ ................................ ................................ ....... 103 C SURVIVORSHIP ................................ ................................ ................................ .... 105 D DIRECT OPTIMIZATION ................................ ................................ ....................... 107 E DERRIVATION OF STATISTICAL GROWTH MODEL ................................ ......... 111 F PARAMETER ESTIMATES ................................ ................................ ................... 112 G ZERO INFLATED POISSON DEMOGRAPHIC VARIANCE ................................ 113 H DISTRIBUTION OF DEVIATES ................................ ................................ ............. 115 LIST OF REFERENCES ................................ ................................ .............................. 117 BIOGRAPHICAL SKETCH ................................ ................................ ........................... 126 PAGE 8 8 LIST OF TABLES Table page 2 1 Parameters with their biological meanings, units, and default values. ................ 34 2 2 Elasticities of ESS life history traits with respect to growth, mortality, and reproduction rate para meters ................................ ................................ .............. 35 5 1 Simulated heterogeneity situations. ................................ ................................ ..... 99 5 2 Baseline log geometric growth rate and effects of heterogeneity in continuous forest and fragments.. ................................ ................................ .......................... 99 F 1 Parameters used for both habitats. ................................ ................................ ... 112 F 2 Habitat specific parameter values. ................................ ................................ .... 112 PAGE 9 9 LIST OF FIGURES Figure page 2 1 ESS growth schedule. ................................ ................................ ......................... 33 2 2 Survivorship during the ESS growth schedule.. ................................ .................. 33 3 1 S imulated growth of individuals. ................................ ................................ .......... 48 3 2 Estimates of growth autocorrelation ................................ ................................ .... 49 3 3 Number of individuals needed to detect positive gro wth autocorrelation. ........... 50 3 4 Power comparison with repeated measures approaches.. ................................ 50 3 5 Size dependent mortality causes bias. ................................ ................................ 51 4 1 P osterior probability density of having different vital rates at a given plant size .. 71 4 2 Changes in plant size (shoot number) relative to stasis (black line) (a), inflorescence production (b), and survival (c). ................................ ..................... 72 4 3 Each panel represents variation quantified by fitting a GLMM to a different vital rate: ................................ ................................ ................................ .............. 73 4 4 Inflorescence rate estimated with mixed effects models with a random effect of individual ................................ ................................ ................................ .......... 74 4 5 Variation in vital rates through time was fairly synchronous in forest (red circles) and fragments (blue triangles). ................................ ............................... 75 5 1 Geometric growth rates. ................................ ................................ ...................... 96 5 2 Inflorescence per individual was averaged across repeated simulations of the same population ................................ ................................ ................................ .. 97 5 3 Predicted and observed variance among replicates of inflorescence production per individual ................................ ................................ ...................... 97 5 4 Transition from continuous forest to fragment,. ................................ ................... 98 5 5 Population structure 30 years after fragmentation. ................................ .............. 98 H 1 Squared value s of log growth rate minus mean log growth rate ...................... 116 PAGE 10 10 Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirem ents for the Degree of Doctor of Philosophy FOUR PERMUTATIONS OF GROWTH, REPRODUCTION, AND SURVIVAL By Mollie Elizabeth Brooks August 2012 Chair: Benjamin Michael Bolker Cochair: Mary Christman Major: Zoology This dissertation contains four studies in quantitative ecology. The first study deals with the optimal life history strategy of flowering perennial plants. Previous work has shown that it is optimal for plants to be non reproductive when small, produce m ale flowers at a medium size, and produce female flowers when large, with the assumption that reproductive rate increases with size. However, reproductive rate does not always increase with size. The study presented here shows that size dependent mortality can select for sex change in the same way as size dependent reproduction. The second study develops a statistical method for detecting differences in growth rate among individuals. Detecting individual variation in growth is important because (1) it could indicate a life history tradeoff and (2) fast growers contribute disproportionately to population growth rates. Detecting among individual differences in growth rates is straightforward when individuals are marked and measured repeatedly. However, marking individuals is not always ethical or feasible. We used maximum likelihood estimation to fit a statistical model to simulated data to test our method. For a given sample size, more information can be gained by marking individuals, but our method PAGE 11 11 increases the range of experimental designs for detecting among individual variation in growth rates. The third and fourth studies deal with demographic data collected from Heliconia accuminata a perennial plant that grows in the Amazon understory. The Am azon is being fragmented by development by humans cutting roads and farms into the forest. The plants rates of growth, survival, and reproductive effort were modeled using generalized linear mixed models; uncertainty was estimated using Markov chain Monte Carlo sampling. The study indicates that, certain environmental conditions increase mortality and reproduction Fragmentation may also be reducing temporal heterogeneity in vital rates. Individual based models were used to predict how heterogeneity affect s the population growth rate. A ll forms of heterogeneity increased the growth rates of our simulated populations PAGE 12 12 CHAPTER 1 INTRODUCT ORY REMARKS My research focuses on life history evolution and population dynamics. They might seem disparate, but they are quantified by a common currency and inherently connected. Life history involves the timing of events and rates of a species' ontogeny. Natural selection optimizes those events and rates to maximize fitness within certain constraints. During a n individual's ontogeny, energy is allocated in different amounts toward growth, reproduction, and survival. Fitness the currency of natural selection, can be measured as either the lifetime production of offspring or the population growth rate depending on specifics o f the dynamics (Mylius and Diekmann, 1995) A population's growth rate is fundamentally a sum of the life histories of the individuals in the population. My research is both theoretical and applied In general, I use statistical models to quantify patterns and use life history theory to find potential explanations for patterns. However, I keep in mind that not every pattern is adaptive. Statistical tools I develop can be used by researchers to get the most information from their data and design informative experiments. My work on population modeling can help conservationists understand how heterogeneity can affect population viability. My work on life history evolution enables us to anticipate how species might adapt to environmental changes. Fitness is frequently discussed as a function of survival and reproduction because these rates directly translate to the currency of fitness life, but my work additionally highl ights the importance of growth. Growth is important because survival probabilities and reproductive rates are usually higher for larger individuals. PAGE 13 13 In many cases, individuals in a population follow life histories that are similar enough to be treated equa lly in population models. However, substantial variation among individuals is often observed in at least one vital rate. Part of my research deals with detecting this variation and predicting how it affects population dynamics. PAGE 14 14 CHAPTER 2 SIZE DEPENDENT SEX CHANGE CAN BE THE ESS WITHOUT ANY SIZE ADVANTAGE OF REPRODUCTION WHEN MORTALITY IS SIZE DEPENDENT Introduction Life history theory predicts reproductive characteristics based on optimal resource allocation with respect to tradeoffs among growth, maint enance, and reproduction (Williams, 1966; Gadgil and Bossert, 1970; Charnov, 1982; Iwasa, 1991; Iwasa et al., 2000) Growth indirectly increases fitness because larger individuals often have better survival rates and may have more resources for reproduction. Reproductive individuals have reduced growth rates because acquired resources are allocated toward reproduction instead of growth or maintenance (Wilson et al., 1994; Obeso, 2002) The sexes often differ in their energetic cost of reproduction. Gametes produced by females are larger and typically require more energetic investment. In sessile organisms, such as plants and corals, which do not require territories or elaborate mati ng behaviors, reproduction as a female typically costs more than reproduction as a male in terms of energy expenditure and reduced growth rate (Zimmerman, 1991; Burd, 1995; Loya and Sakai, 2008) Sex allocation theory deals with allocating resources towards male and femal e reproduction in a way that maximizes lifetime reproductive success (or fitness) in relation to these growth, maintenance, reproduction tradeoffs (Charnov, 1982; Knight et al., 2005) In the past, the evolution of size dependent sex change has been studied theoretically using size advantage, mortality advantage, and growth rate advantage models (Warner et al., 1975; Charnov, 1982; Iwasa, 1991; Nishizawa et al., 2005) These models predict the optimal size of sex change based on differences between the PAGE 15 15 sexes in various vital rates that influence fecundity while accounting for the benefit of being the rarer sex. Each model assumes that fecundity increases with body size for at least one of the sexes. In the size advantage model (Warner et al., 1975; Policansky, 1981; Charnov, 1982; Iwasa 1991; Wright and Barrett, 1999) fecundity increases with size more rapidly for one sex; thus, it is optimal to start out as the sex with shallower gain in fecundity rate and later switch to the sex with a steeper gain rate. In the mortality advantage model, fecundity increases with size in the same way for both sexes, but mortality is different between sexes. It is optimal to be the sex with lower mortality risk first, and then switch sex (Charnov, 1982; Iwasa, 1 991) In the growth rate advantage model, fecundity increases evenly for the sexes, but they differ in their growth rates, so it is optimal to be the fast growing sex first (Charnov, 1982; Iwasa, 1991) For math ematical simplicity, all of these previous models mentioned above assumed that mortality risk is independent of body size, but size dependent mortality should not change their general results. However, it may change specific results and the method we prese nt here to incorporate size dependent mortality may improve previous models' predictive ability. In this paper, we extend the current models of sex change evolution to cases when reproductive rate is perfectly independent of size for both sexes and morali ty decreases with size. This is similar to the growth rate advantage model except that the fitness advantage of growth in size is given not by a higher reproductive rate but by a lower mortality risk for larger individuals. Fitness is measured as lifetime reproduction and depends on both reproduction rate and longevity. Under the assumption of size independent reproductive rate, longevity is the most influential part of this combination. PAGE 16 16 Thus, the directionality of sex change evolves to optimize longevity. With size dependent mortality, staying small longer is more risky; hence, longevity is optimized when individuals quickly outgrow small sizes. Outgrowing small sizes occurs faster when an individual is the faster growing sex (male) first. The dynamic pro gramming calculation, similar to that in Iwasa (1991), is applied. We show that the ESS life history schedule is to be non reproductive at small size, male at intermediate sizes and to be female at large sizes. Motivational Case of Pollen Limitation Thus far, most theories on sex change evolution have assumed that the reproductive success of at least one sex increases with size, but this is not always observed. Typically, female fecundity is limited by size dependent resources and not by access to mates; t his can cause reproductive success to be size dependent. In plants, maximum potential female fecundity depends on the number of ova and the amount of resources available for seed maturation, and these generally increase with plant size. However, realized f ecundity is also constrained by the amount of pollen received (Wilson et al., 1994) If larger plants are not more attractive to pollinators, then reproductive success might not increase with size. Seed set in low pollen systems could be independent of plant size. Stochastic pollen supplies select for females to have more ova per flower than will be pollinated in most years, because an overabundance of ova allows females to ta ke advantage of good pollen years (Burd, 1995) However, if a plant never experiences a good pollen year or if the pollinator population is permanently reduced, then the plant's over allocation toward ova will not improve its fecundity. Under perma nently reduced pollen levels, selection should act to evolve floral traits that increase pollinator attraction (Knight et al., 2005) However, increasing pollinator PAGE 17 17 attraction may be hindered due to genetic or environmental constraints (e.g. if plant and/or pollinator populations are at low density). Number of ova per plant may increase with size, but if pollination is below this level, then seed set may be independent of size and other life history traits, such as size at sex change, may be selected for in the absence of size dependent fecundity. There is some empiric al evidence that pollen limitation causes size independent fecundity. Seed set was observed to be independent of size in two year studies of four small, pollen limited Jack in the Pulpit ( Arisaema triphyllum ) populations (8 to 34 females per population, wi th a mean of 25; Bierzychudek, 1982; C. Heckel, pers. comm.). Similar results were found in a single year study of a population of 38 female Arisaema serratum (Nishizawa et al., 2005) Seed set was also size independent in a single year study of a pollen limited population of 31 Trillium grandiflorum (Wright and Barrett, 1999) In these studies, the number of ova per plant increased significantly with plant size, but due to pollen limitation, seed set did not. These species all allocate resources to male and female reproduction with the ratio changing with size, but only the Arisaema species shows complete switching of sex. These observations moti vated us to ask the theoretical question of whether sex change is still the ESS when reproductive success is independent of size. Model We consider a n individual plant whose reproductive activity is a function of age. Let t be the plant's age. f t ( ) and m t ( ) are reproductive activity functions as a female and male, respectively. If a plant is non flowering at age t then f t ( ) = m t ( ) = 0 ; if it is fully active as a female, then f t ( ) = 1 and m t ( ) = 0 ; if it is fully active as a male, then PAGE 18 18 f t ( ) = 0 and m t ( ) = 1 In general a plant can be in a mixture of these three states and its activity functions satisfy: f t ( ) # 0 m t ( ) # 0 and f t ( ) + m t ( ) # 1 ( 2 1) Let f be the rate of seed production by a fully active female. We assume that f is independent of size because the amount of pollen that females receive is below a level that would cause female seed production to b e resource limited. Let m be the rate of pollen production by a fully active male. The number of seeds fertilized by this pollen should depend on the sex ratio in the population. The expected number of offspring to be sired by a male is mR where R is the ratio of total seeds produced to total pollen production in the population. By this definition, R is also equal to the number of females times seed production rate divided by the number of males and their pollen production rate. If R =1, then each pollen grain fertilizes an ova and produces a seed; this would be an extreme. Smaller values of R indicate that not all of the pollen is being used (typically because of pollinator limitation rather than overproduction of pollen). In this way, m and R are closely coupled parameters. The expected l ifetime reproductive success (RS) of an individual plant at age 0 is : = f f ( t ) + m R m ( t ) ( ) S t ( ) 0 d t ( 2 2) where S t ( ) is survivorship to age t : the probability for a young plant of age 0 to survive until age t Eq. ( 2 2) is equivalent to fitness. We assume that mortality risk is smaller for larger individuals. To be specific, we assume that instantaneous mortality is the sum of two terms: size independent mortality PAGE 19 19 and the size dependent mortality z 0 L t ( ) the latter being inversely proportional to size L t ( ) Then S t ( ) = e + z 0 L x ( ) ( ) dx 0 t # ( 2 3) The size L t ( ) is correlated with the plant's stored energy (e.g. mass, leaf area) and it increases with age as described by the following differential equation: dL dt = a + bL c # m d # f ( 2 4) a is the basic rate of increase in size (Table 2 1). We assume that growth rate increases with size ( b >0), because larger plants have more leaf area and roots with whic h they acquire resources. Growth rate is fastest when there is no reproductive activity ( f = m = 0 ). When reproduction takes place, growth is slower. The cost of being a male and being a female are c and d, respectively, expressed in terms of the reduction of growth rate ( c and d are positive). We assume that being a female is more costly than being a male ( c < d ). Let L 0 be the size of a newly germinated plant. In this paper, we focus on t he case in which both female and male RS are independent of the plant size. Hence f and m are independent of L We search for the optimal reproductive schedule f t ( ) m t ( ) for t > 0 { } that maximizes the lifetime RS Eq. ( 2 2) under the constraint of Eq. ( 2 1 ) and ( 2 4). This is a typical problem of dynamic optimization (Bellman, 1957) In addition, the sex ratio in the population R is determined by the distribution of plant size. If all the plants in the population experience the same environment, we can calculate the evolutionarily stable strategy of reproductive schedule. R is chosen to sa tisfy the following equality: PAGE 20 20 f f ( x ) e # + z 0 L ( y ) $ % & ( ) dy 0 x 0 + dx = mR m ( x ) e # + z 0 L ( y ) $ % & ( ) dy 0 x 0 + dx ( 2 5) which indicates that total female RS summed over all the individuals in the population is exactly the same as total male RS. This relation holds because each offspring has one mother and one father (pollen donor), and they give exactly the same contribution of nuclear genes to the offspring (Fisher 1930) Dynamic P rogramming To solve this dynamic optimization problem, we first consider the following "reproductive value", which indicates the e xpected total RS of an individual from time t until the end of its life: V ( L ( t )) = max f f ( x ) + mR m ( x ) ( ) e # + z 0 L ( y ) $ % & ( ) dy t x t + dx / 0 1 2 ( 2 6) where the function to be integrated, with respect to x from t to infinity, is the reproductive rate at time x multiplied by the survivorship until time x The symbol of "max" implies that the plant should choose the reproductive activity to maximize the lifetime RS. This quantity is a function of current plant size L t ( ) We consider a short tim e interval of length t and separate the integral in Eq. ( 2 1) into the term within this interval and the rest of time (from t + t to infinity). Eq. 2 1 can be rearranged to give the following basic equation (Appendix A ): 0 = max f m f # d dV dL $ % & ( ) f + mR # c dV dL $ % & ( ) m + a + bL ( ) dV dL # + z 0 L $ % & ( ) V ( L ) + / ( 2 7) The max symbol implies that f and m are chosen to maximize the terms in the braces under the constraint of Eq. ( 2 1). PAGE 21 21 Three kinds of intervals: The optimal strategy for the plant should satisfy Eqs. ( 2 7) and ( 2 1). We search for a growth schedule that satisfies these constraints. Since the sum of these terms is a linear function of f and m the optimal choice is always at the extremes, su ch as f m ( ) =(0, 0), or (1, 0) or (0, 1). Here, we consider three different intervals in which [I] the plant does no reproduction (vegetative growth period), [II] the plant is fully active as a male (pollen donor), and [III] the plant is fully active as a female (or fruit producer). Details are ex plained in Appendix B. [I] Vegetative growth phase In this interval, the optimal choice satisfies: f t ( ) = 0 and m t ( ) = 0 Reproductive activity as a female or a male would reduce the growth rate too much, and its benefit to th e plant is smaller than the cost of reducing growth rate. In the vegetative phase, reproductive value changes with size L as (from B 1c): dV dL = + z 0 L a + bL V L ( ) ( 2 8) We define the boundary between vegetative and male phases to be at size L 1 without making any assumptions about the order of the ph ases. At this boundary (from (B 1b) and (B 2b)), dV dL L 1 ( ) = mR c ( 2 9) holds. [II] Male phase PAGE 22 22 In this interval, the optimal choice satisfies: f t ( ) = 0 and m t ( ) = 1 Reproduction as a male is more valuable than vegetative growth or reproduction as a female. Reproductive value changes with size L according to the following equation (from B 2c): dV dL = + z 0 L # $ % & V ( mR a + bL ( c ( 2 10) We define the boundary between male and female phases to be at size L 2 without assuming anything about the o rder. At this boundary (from (B 2b ) and (B 3b)), dV dL L 2 ( ) = f mR d c ( 2 11) holds. [III] Female phase In this interval, the optimal choice satisfies: f t ( ) = 1 and m t ( ) = 0 In this interval, reproduction as a female is more valuable than vegetative growth or reproduction as a male. Reproductive value changes with size L accordi ng to the following equation: dV dL = + z 0 L # $ % & V ( f a + bL ( d ( 2 12) The boundary sizes L 1 and L 2 as well as Eqs. ( 2 7, 2 8, 2 9, 2 10, 2 11, and 2 12) arose from maximizing Eq. ( 2 7) under the constraint of Eq. ( 2 1) without making any assumptions about the order of the phases. These equations make up the dynamic programming conditions and determine the ESS. However, we need candidate strategies to check against the dynamic programming conditions. We examined can didates of the ESS life history schedule that are vegetative (nonreproductive) at small sizes, male at intermediate sizes, and female at large sizes because this pattern PAGE 23 23 is observed in sex changing plants. Let L 0 be the initial size of an individual, L 1 be the size of switching from vegetative phase to male phase, and L 2 be the size of switching from male phase to female phase. Given starting conditions L 0 and V 0 it would be possible to integrate Eq. ( 2 8) and find L 1 according to Eq. ( 2 9), then use V ( L 1 ) and L 1 as starting conditions to integrate Eq. ( 2 10) and find L 2 according to Eq. ( 2 11). However, beca use V L 0 ( ) = V 0 is defined as the expected total reproductive value for an individual of age 0, it must satisfy the following equality: V 0 = mR m t ( ) e # + z 0 L ( x ) $ % & ( ) dx 0 t 0 + dt + f f t ( ) e # + z 0 L ( x ) $ % & ( ) dx 0 t 0 + dt ( 2 13) In Appendix C, we calculate the survivorship to size L and using this expression, we can rewrite Eq. ( 2 13) as follows: V 0 = a + bL 0 a + bL 1 # $ % & b ( z 0 a L 0 L 1 # $ % & z 0 a ) mR a + bL 1 ( c a + bL ( c # $ % & b ( z 0 a ( c L 1 L # $ % & z 0 a ( c dL a + bL ( c L 1 L 2 + + a + bL 1 ( c a + bL 2 ( c # $ % & b ( z 0 a ( c L 1 L 2 # $ % & z 0 a ( c f a + bL 2 ( d a + bL ( d # $ % & b ( z 0 a ( d L 2 L # $ % & z 0 a ( d dL a + bL ( d L 2 / 0 1 1 1 ( 2 14) The second integral in the right hand side converges to a finite value because > 0 If we can find a set of values: V 0 R L 1 and L 2 that satisfy all the equations above, the reproductive schedule is the optimal solution and the ESS. PAGE 24 24 Direct S earch fo r the Optimal S olution We tried different starting values of V 0 and R then calculated new starting values based on the results of L 1 and L 2 using ( 2 14) and ( 2 5). We attempted to find an ESS solution using Eq. ( 2 9) and ( 2 11). However, it was difficult to numerically converge on the ESS parameter sets by this method b ecause these conditions simultaneously depend on and determine V 0 and R Instead, we found a candidate for the optimal solution within a limited life history pattern that has two critical sizes ( L 1 and L 2 ). The expected fitness of an individual can be expressed as a function of these two sizes V 0 = L 1 L 2 ( ) We derive the conditions for L 1 L 2 ( ) to be maximized by taking the derivative with respect to the life history parameters of i nterest. In addition, the sex ratio R should satisfy Eq. ( 2 5). In simplifying this system of equations, it is useful to introduce a quantity s R = mR f which is equal to the total number of females divided by total number of males in the pop ulation. Using these conditions we can determine the ESS value of critical sizes, L 1 and L 2 and the sex ratio, R It was much easier to converge on a solution to these conditions compared to Eq. ( 2 9) and ( 2 11) be cause these conditions are independent of V 0 Solutions to the direct optimization criteria were found using the function uniroot in the statistical programming language R version 2.9.0. The range of sex ratios for uniroot to search had to be chosen carefully because if the lower limit of search values was too small, then there was no solution to Eq. (D 11) and the algorithm would stop before finding an ESS solution to the system. This means that (in a population not follow ing the ESS) when sex ratios are severely male skewed, there is no solution for L 1 PAGE 25 25 and thus, an optimal individual life history contains no male stage. The details of the direct optimization procedure are explained in Appendix D and cod e is available online. After the solution of two critical sizes ( L 1 and L 2 ) was determined, we then numerically confirmed that the solution from this direct optimization method satisfies the conditions obtained from dynamic optimization, Eq. ( 2 9) and ( 2 11). Calculations were done using function lsoda from package odesolve (Setzer, 2011) With approval of the dynamic optimization conditions, we can claim that the solution is the ESS among all life history schedules of the form of Eq. ( 2 1 ) t he life history pattern with two switching sizes ( L 1 and L 2 ) has a higher fitness of any schedule under the constraint of Eq. ( 2 1). P arameter D ependence To determine how the model's behavior changes with parameters, we performed sensitivity analys e s. We first chose a set of default parameters as shown in Table 2 1 The relative values of certain groups of parameters are more relevant than absolute values. The results of the model are general as long as some dimensionless ratios fit the model assumptions. Our default values were not based on data, but on the general model assumptions. Growth parameters were chosen so that a plant that w as reproductive at initial size would have a negative growth rate ( c/(a+bL 0 ) >1 and d/(a+bL 0 ) >1), with female being more costly than male ( d/c >1). In addition to our assumptions, the ratio of size dependent growth, b and size independent mortality, became another constraint of the model. This dimensionless ratio is contained in the survivorship equations Eq. (B 11a, b, c) as explained in Appendix B We were unable to find the ESS solutions for /b greater than 0.9009 (see discussion below). A growt h PAGE 26 26 schedule with the ESS for the default parameter values is illustrated in Fig. 2 1 and survivorship is in Fig 2 2. We assessed the dependence of the ESS reproductive schedule on parameters included in the model by adjusting each parameter separately. We found the ESS for values ten and twenty percent above and below the default values. Then we calculated elasticity as the slope of the log log plot with five points per parameter centered at the default value. Each log log plot was visually inspected for l inearity. Regressions were done using the function lm in R. Results and Discussion Sex C hange in the Absence of Size Dependent Reproduction We set out to determine whether sex change can be the ESS when individuals' reproductive rates are size independent The method we created to answer this question involved two steps: first we chose a candidate life history and then we confirmed that the candidate is the optimal out of any reproductive schedule. We found a candidate for the optimal schedule by directly maximizing lifetime reproductive success (Appendix D). Although the candidate was selected from a limited class of schedules, we confirmed that it was the optimum of any reproductive schedule satisfying Eq s ( 2 1) ( 2 5) and ( 2 6) This was determined by numerically confirming that the candidate met the conditions given by the dynamic programming calculation (Appendix B). By this calculation, we know that the candidate life history is the ESS. Thus, we conclude that sex change can be the ESS even when there is no fecundity advantage to being large. Th e analys e s support our conjecture that size dependent sex change can be the ESS even when the reproductive success rate is perfectly independent of size in both PAGE 27 27 sexes, which is observe d in terrestrial plants engaging in sequential hermaphroditism ( see citations in pollen limitation section ). The growth rate advantage of males combines with size dependent mortality to create a fitness advantage for switching from male to female even when there is no size advantage to reproduction. This acts to maintain the observed sex changing pattern during periods when the population is pollen limited; or it could even explain the origin of sex change. This finding highlights that patterns of reprodu ction, growth, and maintenance all have significant effects on fitness and selection on life histories. Thus, it is important to consider their effects. In light of this finding, incorporating size dependent mortality into other models of size dependent s ex allocation may improve their accuracy. This will require field biologists to monitor mortality in addition to other vital rates. These results support the finding of (Day and Aarssen, 1997) that mortality is important in size dependent sex allocation even when mortality rate is independent of sex. Day and Aarssen assu me that mortality risk decreases with size as we do here. However, other assumptions of their model differ from the one presented here: they ignore an individual's accumulated growth from year to year, assume that potential pollen and ova production increa se with size, assume it takes longer to mature fertilized ova than to release pollen, and assume there is a non linear tradeoff between pollen and ova production. By assuming that the pollen ova tradeoff is concave down, they focus their work on simultaneo us hermaphrodites (Charnov et al., 1976) Despite these differing assumptions, the general resu lt holds size dependent mortality is important. This could explain why some empirical study results do not fit predictions of simpler models; perhaps size dependent mortality and growth need to be accounted for in PAGE 28 28 addition to reproductive rates. The effe cts of size dependent mortality on sex change elucidated by our model should hold generally, even when reproductive rates are size dependent. Switch sizes ( L 1 and L 2 ) were strongly affected by the respective costs of male and female reproduction ( c and d )(Table 2 2) Higher growth costs mean that individuals have a higher risk of mortality per size increment. Thus, it is best to have more leaf area and root mass with th e potential for more resource acquisition before paying the growth cost of pollen and ova production. The two mortality terms had opposite effects on the optimal switch sizes. Size independent mortality promotes earlier switching. This result has been obs erved in other studies because with an uncertain future, there may not be many future opportunities to reproduce; so current reproduction is a better bet than the future (McNa mara et al., 2004) Size dependent mortality promotes delayed reproduction to more quickly outgrow the risky sizes. The ratio /b had strong effects on the optimal switch sizes because, after mortality and growth rates were translated into survivorship, this ratio was in the exponent of most terms of the survivorship equations (Appendix C). Thus, it affects the second derivative (i.e. concavity) of survivorship with respect to size (Werner and Gilliam, 1984; Filin and Ovadia, 2007) Larger values of /b reduce survivorship in all stages. When the popula tion is constrained to be at it s ESS sex ratio (total fertilization by males equals total seed production (Eq. ( 2 5)) then the range of sizes in which male reproduction is optimal ( L 2 L 1 ) log linearly decreases as /b increases. This is because survivorship is a key component of Eq. ( 2 5) and is strongly affected by /b PAGE 29 29 In contrast to the male female switch, the timing of maturation ( L 1 ) is more strongly influenced by size dependent mortality risk, which increases the hazard only for small plants. Small plants can outgrow this hazardous size range more quickly by delaying reproduction and the rate at which they outgrow it depends partly on the size dependent growth rate b Because L 2 indirectly depends on L 1 (see above paragraph), it has a non linear relationship with /b on a log scale Thus, and b had non linear effects on the optimal size to be come female, L 2 and survivorship to become female on a log scale. Patterns of Sex Ratio Sex ratio ( s R ) was strongly affected by the growth costs of being male and female ( c and d ). When the cost of being male increases, fewer males ma ke it through that stage and the sex ratio becomes less male skewed (higher s R ). Similarly, when the cost of being female increases, females do not survive as long and the sex ratio becomes more male skewed (lower s R ). The sex ratio (number of females per male) becomes more male skewed as /b increases because survivorship is reduced throughout the female phase proportionally more relative to the male phase. This holds as long as the population follows the sequence: non reproductive to male to female, even when the switch sizes are not the ESS. We confirmed thi s numerically by solving Eq. (D 13) for s R and inputting different values of and b without resetting L 1 and L 2 These resu lts offer an explanation for why plants may have more females per male in higher quality sites (empirical studies reviewed in (Bierzychudek and Eckhart, 1988) ; recent example in (Heckel et al., 2010) ). As Day and Aarssen (1997) suppose, PAGE 30 30 higher quality sites may have lower mortality: or z 0 Our model shows that a reduction in or z 0 would result in more females per male (Table 2 2). Also, a higher quality site may have higher nutrient availability and thus higher size dependent growth rates, b. In our model, this scenario also leads to more females per male. Note, that the sex ratio is still always biased toward the first sex (ma le) as Charnov and Bull (1989) predict (see also: (Frank and Swingland, 1988; Charnov, 1989; Allsop and West, 2004) ). Although the results we present in Table 2 are under the constraint of t he population being at the ESS, they also hold for genotypes adapted to follow the ESS of a site with differing quality than the one t hey are in (tested using Eq. (D 13) as stated above). For example, if genotypes that had adapted to follow the ESS of one site found themselves in a lower quality site (by either gene flow or environmental change) higher mortality rates would reduce survivorship to the female stage or, alternatively, slower growth rates would cause fewer individuals to reach the genetically d etermined size of sex change. Thus either growth or mortality differences could cause there to be fewer females per male in a lower quality site. These patterns should hold even when reproductive rates are size dependent. Thus, to determine the cause of di ffering sex ratios among populations, field biologists should measure factors that affect either growth or mortality. Unimportant Parameters A few parameters had no effect on the ESS. Initial plant size had no effect on optimal switch sizes because L 0 is not part of the criteria for optimal switching: Eqs. ( 2 8), ( 2 9), ( 2 10), ( 2 11), (D 11), (D 12), and (D 13). Optimal switch sizes depend only on current and future mortality and growth, not on previous conditions. Increasing the pollen production rate o f males only increases the amount of excess pollen (lower R ), PAGE 31 31 because there is already more pollen produced than reaches and fertilizes ova. Without an increase in seed production (i.e. pollen delivery and fertilization), pollen produ ction has no value. Because the population's sex ratio, s R is a function of gamete production rates and the population's seed per pollen ratio, R the changes cancel each other out and the ESS sex ratio is unaffected. Model Generality This type of model predicts sex change dependent only on size; age is only a means of accounting for growth and future expectations. Because sex change is optimized at a particular size threshold, females that shrink below this threshold due to damage should switch back to being male; the same holds for the vegative to male switch size. We do not account for any possible sensing of other individuals so this would not apply to organisms such as certain fish species that change their sex due to d ominance hierarchies and mate availability (Kuwamura and Nakashima, 1998; Munday, 2002) This same set of assumptions may apply to fre e spawning coral populations that are sperm limited. Coral populations at low density in highly turbid water can experience severe sperm limitation (Yund, 2000) and environmental variables may erase the size dependence of gamete production when it comes to reproductive success. This possibility is something to consider in future studies of sex changing coral such as mushroom stony corals (Loya and Sakai, 2008) T he model does have some limitations. First, the model i s formalized in a continuous time version under a constant environment in which there is no seasonality. In a seasonal environment, reproduction, mor t ality, and growth may vary seasonally, and a discrete time model may be more realistic. In the case of Jac k in the Pulpits and PAGE 32 32 many other herbaceous perennials, mortality risk is highest in the winter, growth occurs in summer, and sex is determined once per year (Bierzychudek, 1984) Second we did not consid er potential density dependence explicitly in the analysis because pollen and sperm limited populations are likely to be at low density (but see (Mylius and Diekmann, 1 995) for density dependence in life history optimization). Third, we assumed that there is n either fitness cost of sex change nor selfing benefit to hermaphrodites, and the results of the model might be different when these assumptions are modified. The exten s ion of the analysis considering these modifications would be a good theme of future theoretical study of sex change evolution of terrestrial plants and free spawning corals PAGE 33 33 Fig 2 1 ESS growth schedule Parameters are in Table (2 1 ) L 1 is the ESS size to begin reproducing as a male and L 2 is the ESS size to switch from male to female. Fig. 2 2 Survivorship during the ESS growth schedule Parameters are the same as in Fig. 2 1 (see Table 2 1). L 1 is the ESS size to begin reproducing as a male and L 2 is the ESS size to switch from male to female. PAGE 34 34 Table 2 1 Parameters with their biological meanings, units, and default values. parameter biological meaning units default model input a base growth rate size/time 0.1 b size dependent growth rate 1/time 0.1 c growth cost of male reproduction size/time 0.4 d growth cost of female reproduction size/time 0.5 size independent mortality rate 1/time 0.05 z 0 size dependent mortality rate size/time 0.2 L 0 initial size of organism size 1 m gamete production rate of males RS/time 1 f offspring production rate of females RS/time 1 model output (ESS) L 1 size to begin reproduction as male size 4.17 L 2 size to switch from male to female size 6.25 s R females per male in population unitless 0.90 R offspring per male gamete in population unitless 0.90 PAGE 35 35 Table 2 2 E lasticities of ESS life history traits with respect to growth, mortality, and reproduction rate parameters calculated as the slope at the default value of the log log plot (see text for details). L 1 is the size to begin reproducing as a male and L 2 is the size to switch from male to female. s R is the number of females per male in the population. Survivorship to maturity is the probability that a plant of size L 0 will survive to be male, and survivorship to be female is the probability that a plant of size L 0 will survive to become female. M is the expe cted RS while male given that a plant survives to become male: mR *Eq. (C. 2 ). F is the expected RS while female given that a plant survives to become female: f *Eq. (C. 3 ). R is the population's offspring per pollen ratio. Entries marked represent nonlinear relationships. L 1 L 2 s R s urvivorship to maturity s urvivorship to be female M F R a 0.17 0.05 0.01 0.63 0.6 0.02 0.01 0.01 b 0.64 0.24 0.21 1.4 7 1.2 6 0.19 0.39 0.2 1 c 0.76 0.14 0.49 0.60 0.1 9 0.47 0.06 0.49 d 0.28 0.8 0.4 7 0.22 0.4 8 0.1 0.15 0.4 7 0.37 1.33 0.21 0.18 0.02 1. 2 1.4 0.21 z 0 0.13 0.12 0.03 1.0 5 1.1 7 0.3 5 0.2 3 0.03 L 0 0 0 0 1.26 1.26 0 0 0 m 0 0 0 0 0 0 0 1 f 0 0 0 0 0 1 1 1 PAGE 36 36 CHAPTER 3 DETECTING VARIATION IN GROWTH RATES WITH OUT MARKING INDIVIDU ALS Introduction Ecologists and evolutionary biologists have long been interested in growth in body size. Studies of growth typically focus on differences among means of popu lations or treatment groups, striving for low variability around the mean to increase statistical power and v ariation within groups is often treated as noise obscuring the phenomena of interest. However, ecological studies are increasingly incorporating a mong individual variation into analyses as either a treatment or a response variable (Pfister and Stevens, 2002; Benedetti Cecchi, 2003; Inouye, 2005; Peacor and Pfister, 2006; Peacor et al., 2007) T hese studies have shown that variation among individuals is itself the result of important biological processes and that population dynamics are sensitive to among individual variation ( Morris and Doak, 2002; Pfister and Wang, 2005; de Valpine, 2009; Zuidema et al., 2009; Bolnick et al., 2011) Th ese studies have also highlighted challenges for both quantifying and explaining the mechanisms behind observed variation among individuals i n a population or cohort (Pfister and Stevens, 2002; Peacor and Pfister, 2006) Here, we present a new method that will permit researchers to separate among and within individual variation in growth rates based on data from individuals within cohorts measured several times over the course of their ontogeny, but where (in contrast to existing methods) individual organisms need not be marked. By expanding the range of organisms and experimental designs where among individual variation can be estimated this method will enable researchers to better understand the sources of variation in body size data. PAGE 37 37 Many ecological and evolutionary processes depend on body size (Weiner, 1985; Brown et al., 2004; Weitz and Levin, 2006) For example, anurans and other organisms with complex life histories must attain a minimum size before they metamorphose and leave the water (Werner and Gilliam, 1984) Because the body size at which individuals undergo life history transitions is thought to be correlated wit h fitness (Alford and Harris, 1988; Altwegg and Reyer, 2003; Marshall and Keough, 2005) much life history theory has focused on predicting the size and timing of these transitions (Wilbur and Collins, 1973; Werner and Gilliam, 1984; Rowe and Ludwig, 1991) While t hese models largely ignore among individual variation, a variety of neutral and adaptive processes can maintain among i ndividual variation and modulate the expected results For example, prey may be forced to choose between a proactive strategy, foraging and growing normally in the presence of predators, and a reactive strategy, reducing foraging to avoid predation but g rowing more slowly as a result. In unpredictable habitats such as vernal pools, reactive individuals trade off predation resistance against a greater risk of mortality from other factors such as desiccation from pond drying. Such growth mortality risk trad eoffs can lead to flat or bimodal fitness curves that maintain variation in populations (Mangel and Stamps, 2001) Variation in size among individuals also m odifies population dynamics (Persson et al., 1998; Zuidema and Franco, 2001; Zuidema et al., 2009) For example, among individual variation in development al rates changes the amplitude and periodicity of population cycles in host parasitoid models (Wearing et al., 2004) Neglecting among individual variation in surv ival probability leads to overestimation of extinction risk in population viability analyses, due to the variance reduction effect, (Fox and Kendall, PAGE 38 38 2002) and underestimation of the population's asymptotic growth rate (Kendall et al., 2011) ; similarly, population viability analyses that neglect variance in fecundity among individuals ma y either over or underestimate the population growth rate depending on the shape of the mean variance relationship (Fox and Kendall, 2002) and the direction of covariance between parent and offspring fecundity (Kendall et al., 2011) Because survival probabilities and fecundity rates are closely linked to individual body size, variation in individual growth rate may drive changes in demography; thus incorporating growth autocorrelation in models may allow more accurate predictions of population size structure (Pfister and Wang, 2005) Three main growth processes lead to growth depensation (increased size v ariation within a cohort through time): within individual variation in growth rate, among individual variation in growth rate (i.e., positive growth autocorrelation), and size dependent growth (Fig. 3 1a ) (Pfister and Stevens, 2002) Within individual variation in growth ra te occurs when environmental heterogeneity causes growth rates to vary within individuals in ways that are uncorrelated through time. For the purposes of this paper, we take the pattern of growth depensation caused by within individual variation to be the null expectation, and see if we can detect particular deviations from it. Among individual variation in growth rate, or positive growth autocorrelation, is defined as positive temporal correlation in the growth rate of individuals Many ecological processe s can generate positive growth autocorrelation. The proactive and reactive behavior types discussed above generate permanent autocorrelation (autocorrelation that applies throughout the entire life stage), as individuals consistently PAGE 39 39 express the same behav ior pattern (Coleman and Wilson, 1998; Sih et al., 2004) In tree populations, variation in liana load generates permanent autocorrelat ion (Zu idema and Franco, 2001) while growth autocorrelation driven by extra light availability near treefall gaps is temporary, acting only until an individual near the gap grows up to fill it. In this paper, we focus on growth autocorrelation that persists throughout the time period of interest, although our methods could in principle be adapted to detect temporary autocorrelati on. Size dependent growth, where larger individu als have higher expected growth rates, can result from size dependent gape limitation or size dependent range s ize in animals. In plants, size dependent growth often results from size dependent resource uptake and asymmetric competition (Weiner, 1985) While size dependent growth is impor tant, and has frequently been suggested as a mechanism of positive growth autocorrelation (Peacor et al., 2007) we do not focus on it here although our methods could be extended to detect it, or (with sufficient available data) to partition growth depensation into components arising from size dependent growth and from positive growth autoco rrelation These three classes of mechanisms lead to different patterns of variation among individuals in a cohort through time (Fig 3 1 b ). Because size dependent and autocorrelated growth persist thro ugh time, they lead to larger growth depensation than within individual variation in growth rates. Because methods to separate the contributions of all three mechanisms acting simultaneously would be both complicated and data hungry, we focus on the rela tive contribution to growth depensation of within and among individual variation in growth. PAGE 40 40 We thus assume that individual growth rates are independent of size, or equivalently that a cohort's mean body size grows linearly through time. Although this assu mption may seem restrictive, many organisms grow approximately linearly in size over some window in their ontogeny (Ricklefs 1967, Laird et al. 1965) More generally, our method will apply whenever body size data can be transformed to be a linear function of time: in particular, if organisms grow exponentially with time (a common pattern early in ontogeny : (LAIRD et al., 1965; RICKLEFS, 1967; RICKLEFS and Cullen, 1973) ), then log transforming the data will make our method applicable. In the past, teasing apart the relative importance of within and among individual variation in growth for growth depensation has required scientists to mark individuals and follow each individual's growth pattern (Brienen et al., 2006; Peacor and Pfister, 2006) or to create distinct size classes in a starting cohort and monitor the intermixing of size classes (Peacor et al., 2007) In many systems, however, neither of these approaches is feasible. Marking i ndividuals can bias results (McCarthy and Parris, 2004; Saraux et al., 2011) and requires extra time that limits the scope of studies (Lavine et al., 2002) Lavine et al. describe a method to estimate seedling mortality without marking individuals, using only observations of the numbers of old and new seedlings through time. Similarly, our method for quantifying among individual variation is observational, fitting a model of the expected changes in variance through time to repeated measures of a cohort's variance. Methods We first derive the equat ions for the changes in cohort variance over time as a function of average growth rate, total variance in growth rate, and level of growth autocorrelation (non technical readers can safely skip this section). We then discuss our PAGE 41 41 protocol for simulating coh ort growth dynamics to test the statistical power of our approach and summarize the practical aspects of the estimation procedure for researchers interested in applying the method to their own data. We compare our method to standard repeated measures metho ds that are available only when individuals are marked. Finally, we add size dependent mortality to the data simulations and describe its effects on parameter estimates. Derivation Suppose that individuals in a cohort grow linearly with mean growth rate g per time step t E ach individual, i consistently deviates from this average growth rate by i a normal deviate with mean 0 and variance 2 # g 2 (Assuming normality is convenient for statistical inference on the parameters, but the derivation depends only on the mean and variance of this and other values.) At each time step, each individual 's growth rate also has an uncorrelated deviation i t with mean 0 and variance 1 # 2 ( ) $ g 2 % t Then an individual's size changing through time can be modeled as S i ( t + t ) = S i ( t ) + t g + t # i + # i t ( 3 1 ) Modeling an individual's growth in this way is equivalent to using two normal distributions with un related variances. The advantage of using these two variance components is that g 2 and 2 are meaningful parameters. T he cohort's size variance increase s quadratically through time when 2 >0 : s 2 t + # t ( ) = s 2 t ( ) + # t 2 $ 2 g 2 + # t ( 1 % $ 2 ) g 2 + 2 t # t $ 2 g 2 ( 3 2) See Appendix E for a more detailed derivation. PAGE 42 42 Without loss of generality, we scale the units of t so that t ranges from 0 to 1 during the period of observation, so that g 2 i s the total increase in cohort variance during the period of o bservation ; set average growth rate g (which does not appear in the variance equation, Eq. ( 3 2 )) to 1 ; and set the initial variance 0 2 to 1, equivalent to setting the units of size this also redefines g 2 as the relative increase in variance over the observation period. Simulating Growth D ata To evaluate the performance of our statistical model, we simulated growth data in the R statistical programming environment (version 2.12.2) and fit the variance equation to the simulated data using maximum likelihood estimation. The increase in variance, g 2 was set to 4, 16, 32, 64, 128, and 256. This range of variances includes those observed in empirical studies on tadpoles (MWM cCoy unpublished, Peacor and Pfister 2006). The growth autocorrelation parameter ranged from 0 to 0.9 by increments of 0.1. Observations were sampled at n t =8, 16, or 32 evenly spaced time steps. Growth autocorrelation was simulated by assigning each individual a normally distrib uted mean growth rate with mean g and variance 2 # g 2 t where t =1/( n t 1). The within individual variation was simulated for each individual at each time step by choosing random deviates with mean 0 and variance ( 1 # 2 ) $ g 2 % t We simulated observations of body size for 8, 16, 32, 64, 128, 256, and 512 individuals in a cohort and ran 1000 replicates for each combination of parameters. Model Fitting Because our model assumes only process and not measurement error we fit the parameters by step ahead prediction, equivalent for a normal response to fitting the PAGE 43 43 between step changes in variance as independent and normally distributed values with mean (Bolker, 2008) : s 2 t + t ( ) s 2 t ( ) = t 2 # 2 g 2 + t ( 1 # 2 ) g 2 + 2 t t # 2 g 2 ( 3 3 ) We fit this equation to the data using maximum likelihood estimation in A utomatic D ifferentiation Model Builder (ADMB) via the R package R2admb version 0.7.4/r65 to get estimates of 2 and g 2 (Fournier et al., 2011; Bolker and Skaug, 2012) The same maximum likelihood estimation could be done using the package bbmle in R, but ADMB uses a more robust algorithm, automatic differentiation (Bolker, 2012) We bounded to be between 0.001 and 0.999 with initial value 0.5. We started g 2 at the value of the observed increase in variance and bounded g 2 between 0.01 and 1.5 times the observed increase in variance. We computed 95% likelihood profile confidence intervals on 2 We calculated the coverage for each combinati on of parameter values as the proportion of simulations in which confidence intervals contained the true value of 2 Procedures A step by step protocol for quantifying growth autocorrelation with our method is as follows: (1) C onfirm that the mean growth rates of the cohort are roughly independent of the mean body size (or equivalently that growth is approximately linear), transforming the data (e.g. by taking logarithms) if necessary. (2) C alculate the cohort's variance at each time step and take the diff erences to find the change in variance at each time step. (3) Use maximum likelihood equation to estimate 2 from the data and PAGE 44 44 the equatio n for change in size variation Eq. ( 3 3 ) Use likelihood profiling to f ind 95% confidence intervals for 2 Repeated Measures When it is possible to mark individuals, more traditional repeated measures analyses can be used to estimate growth autocorrelation To compare our method to repeated measures methods we fit a linear mixed model (LMM) to individual growth rates w ith a random effect of individual, using the same simulated data described above We fit the model using lmer from the R package lme4 version 0.999375 40 W e estimated 2 from the variance of the random effect of individual and the residual variance of the fitted model: e s t i m a t e 2 = n t # 1 ( ) $ i n d i v i d u a l 2 n t # 1 ( ) $ i n d i v i d u a l 2 + $ r e s i d u a l 2 ( 3 4 ) To test if individual growth rates varied significantly (equivalent to testing the null hypothesis 2 =0 ), we did an ex act restricted likelihood ratio test on the random effect of individual using the function exactRLRT from the R package RLRsim version 2.0 10 (Scheipl et al., 2008) Size Dependent Mortality Simulations Our model assumes that all individuals survive throughout the experim ent. However, this assumption may be violated in experimental and especially in observational studies The worst case scenario is when individual mortality rates depend on size; we tested our method's performance in this scenario, specifically assuming tha t smaller individuals have a higher mortality rate. Each individual survived according to a Bernoulli trial at each time step, with a probability equal to a logistic function of its size at time t scaled by the duration of the time step: PAGE 45 45 t 1 + e x p # r S i t ( ) # x 0 ( ) ( ) ( 3 5 ) We used r = 0.4 and x 0 = 10, 9, 8, and 7. For each value of x 0 and combination of parameters specified above, we estimated 2 for 1000 replicate simulation s. Results and Discussion Sampling more individuals both improved point esti mates and narrowed confidence intervals ( Fig 3 2 ). Sampling more time points (in the range 8 to 32) had negligible effect on point estimates of 2 Fewer time points gave slightly narrower (undercovering) confidence intervals (most on the order of 0.01 but up to 0.04 narrower for 8 time points compared to 32 time points). In simulations with 8 time points, confidence intervals contained the true value of 2 slightly less than 95% of the time; their coverage leveled off at 90% for more than 50 individuals sampled. Whenever fewer than 50 individuals were sampled, 95% confidence intervals contained the true value of 2 less than 95% of the time (i.e. undercoverage); the worst case was 64% coverage for 8 individuals sampled at 8 time points with 2 = 0.81. Th us, we recommend sampling a minimum of 50 individuals on more than 8 occasions. Simulations with greater increase in variation, g 2 had slightly less biased point estimates of 2 (most on the order of 0.01 but up to 0.04 less bias for simulations with g 2 =256 compared to 4). Increase s in variation, g 2 had no effect on confidence interval width. Our simulation results can be used to guide exp erimental designs for detecting growth autocorrelation in cohorts of unmarked individuals. Preliminary growth data from pilot lab or field studies, or data from the literature, can be used to guess an PAGE 46 46 approximate 2 Given this informatio n, researchers can use Figs. 3 2 and 3 3 to make decisions about feasible precision and necessary sample sizes. At a minimum, researchers will want to confirm whether observed growth depensation could be the result of growth autocorrelation (i.e., to test the null hypothes is that 2 =0 versus the alterative that 2 >0). The number of measured individuals needed for 80% power to detect 2 greater than zero depends strongly on the true value of 2 (Fig. 3 3 ). For example, if the true value of 2 is 0.64 then only 30 individuals are needed for 80% power (although with fewer than 50 individuals, estimates of 2 may be biased : Fig. 3 4 ). If the true value of 2 is 0.36 and g 2 is 16, then 290 individuals are needed for 80% power For the simulated experiments, only values of 2 greater than 0.36 were ever distinguishable from zero with 80% power, regardless of sample size Comparison of Variance Pattern and Repeated Measures A pproaches When individuals of the study species can be marked, traditional repeated measures analyses can be used to estimate growth autocorrelation As more individuals are sampled, both variance pattern and repeated measures methods approach 100% power, although repeated measures power is always higher and i ncreases more rapidly (Fig. 3 4 ). Nevertheless, variance pattern power never has more than 30% less power in the range of scenarios we examined. Both methods are more powerful when detecting larger true values of 2 Size Dependent Mortality Causes Bias In simulations incorporating size dependent mortality estimates of 2 were biased downwards The strength of the bias increased with the proportion of individuals that died (Fig 3 5 ) Estimates of g 2 were also biased downward with inc reasing magnitude PAGE 47 47 as a larger proportion of individuals died. Because smaller individuals were selectively removed the cohort's variance increased by less than the nominal amount g 2 With less final variance, the cohort's change in var iance through time followed a more linear, less quadratic pattern than the model predicts without mortality. Mortality rates were higher in simulations with larger values of g 2 and hence bias increased, because the larger cohort varian ce resulted in more individuals falling within the high mortality size range determined by Eq. ( 3 5 ) For a given value of x 0 realized mortality varied greatly, and realized mortality is a better predictor of bias, as well as being more directly related to ecological information that would be available to empirical researchers. Conclusion Detecting autocorrelated growth patterns previously required the marking of individuals, but marking is logistically infeasible in many ecological systems. Here, we have shown that, with a large enough sample, one can detect growth autocorrelation observationally by analyzing the patterns of increasing variance over time. This allows researchers to make a choice of how to best allocate their effort: sampling more individu als without marking them (a cheaper and faster design) or marking fewer individuals (a design for better information). We recommend applying this method to experimental designs where mortality can be kept below 10% over the course of sampling. It can be ap plied to open populations as long as individuals of the same age can be clearly identified. Better quantification of the patterns and genesis of size variation will help ecologists improve predictions of population dynamics as well as their basic understan ding of ecological systems. PAGE 48 48 Fig. 3 1 We simulated growth of individuals with the same growth parameters except for the assumptions: (panels from left to right) uncorrelated variation among individuals, autocorrelated variation i.e. variati on among individuals, and positively size dependent growth. (a) Patterns of growth rate Growth rates of five individuals are represented in each graph by five separate lines. (b) Patterns of size variation We replicated 2000 cohorts of 50 individuals and plotted the average cohort's variance (solid line) and the 2.5% and 97.5% quantiles (grey ribbons). a b PAGE 49 49 Fig. 3 2 Estimates of growth autocorrelation 2 Estimates of 2 (solid lines) true values of 2 (dashed horizontal lines) and 95% confidence intervals (gray ribbons), averaged over 1000 replicates of each parameter combination. Number of time points n t =16. The true value of 2 increases moving across the panel PAGE 50 50 Fig. 3 3 Number of individuals ne eded to detect positive growth autocorrelation The line represents the minimum number of individuals in a cohort needed to statistically detect that 2 is greater than 0 at least 80 percent of the time. These are based on the 95% confidence intervals of 1000 simulations with g 2 = 16 and n t = 16 (see Simulating growth data section for details). When 2 < 0.36, more than 512 individuals are needed, beyond the range we simulated. Fig. 3 4 Power comparison with repeated measures approaches. Our method for detecting growth autocorrelation (solid line) is less powerful than doing an exact restricted likelihood ratio test on a linear mixed model fit to data on marked individuals (dotted line). The true value of 2 increases moving ac ross the panels. PAGE 51 51 Fig. 3 5 Size dependent mortality causes bias. Estimates were biased downward with increasing magnitude as a larger proportion of individuals died in simulations with size dependent survival following Eq. 5. True values of 2 increase across the panels (horizontal dotted lines). Solid and dashed lines represent different parameterizations of increases in total variance. PAGE 52 52 CHAPTER 4 DOES FRAGMENTATION A LTE R VITAL RATES AND IN CREASE HETEROGEN EI TY? Introduction Environmental heterogeneity has important consequences for populations. When an organism's vital rates (growth, survival, and reproduction) are influenced by environmental factors, then environmental heterogeneity can affect population dynamics and extinct ion risk. The effects of environmental heterogeneity are scale dependent. Small scale spatial heterogeneity can increase variation among individuals and can either increase or decrease the population's mean outcome due to Jensen's inequality (Ruel and Ayres, 1999) Spatia l heterogeneity on the landscape scale reduces extinction risk if offspring settle in environments similar to their parents (Fox, 2005) and allows populations to attain higher densities due to a process known as habitat association in which more individuals become concentrated in the higher quality patches (Bolker, 2003) Temporal heterogeneity in demographic rates tends to reduce a population's long term growth rate (Boyce, 1977; Tuljapurkar and Orzack, 1980) In this study, we quantify heterogeneit y on multiple scales and multiple habitats: yearly, individual, and landscape scales, in fragmented (F) and continuous forest (CF) habitats. We quantify heterogeneity in the vital rates growth, reproduction, and survival of plants, which we expect to b e indicative of their environmental quality. We are specifically interested in stochastic heterogeneity rather than periodic heterogeneity or gradients. As previous work has shown, this variability will have important effects on populations' probability of persistence Although our study species is not in foreseeable PAGE 53 53 risk of extinction, variation in its vital rates are a good indicator of environmental heterogeneity and this will have implications for other species. Variation among individuals in their expe cted demographic rates is also called demographic heterogeneity; which differs from demographic stochasticity that can occur even when all individuals have the same trait and same expected rate (Kendall et al., 2011) Demographic heterogeneity is a form of phenotypic trait variation and can be the result of both genetics an d environment. Here we are interested in quantifying demographic heterogeneity because of its link to small scale environmental heterogeneity. Recent studies have tried to decompose the relative strength of the consequences of variation on multiple scales. Melbourne and Hastings (2008) examined the relative contribution of demographic heterogeneity, demographic stochasticity, and environmental stochasticity to extinction risk. They found that demographic heterogeneity increases the risk of extinction in den sity dependent populations more than environmental stochasticity, but if it's left out of the statistical model, variation due to demographic heterogeneity is erroneously lumped in with environmental variation. Jongejans and Kroon (Jongejans and De Kroon, 2005) and Ramula et al. (Ramula et al., 2009) compare the relative variation of population growth rate in space versu s time; results differed. Jongejans et al. (Jongejans et al., 2010) quantified the net variance contribution of vital rates to variance in population growth at the regional, site, and yearly temporal scale; however with only 2 years of data temporal variability cann ot be well described. PAGE 54 54 Effects of F ragmentation Fragmentation is a pervasive form of landscape change. Conditions in fragments can be highly variable (e.g., from the edge to the interior), so different individuals in the same fragment can be exposed to diff erent environmental conditions. Fragments can also differ from each other because of disturbance events such as fires and blowdowns. Fragmentation alters forest hydrology and wind patterns, thereby creating more environmental heterogeneity both temporally and spatially Clearings promote desiccation up to 2.7 km into nearby fragments and in turn, increase spatial heterogeneity (Malcolm, 1998; Didham and Lawton, 1999; Cochrane and Laurance, 2008; Briant et al., 2010) Because clearings have less evapotransporation and rainfall interception by vegetation streams respond more quickly to precipitation events which promotes localized flooding in the wet season and reduced stream flow in the dry season ( unpublished Trancoso 2008). These more extreme hydrological events will spill over into forest fragments increasing temporal stochasticity and spatial heterogeneity. Additional edge effects include increased windshear and wind turbulence that increase rates of tree mortality and damage (Saunders et al., 1991; Ferreira and Laurance, 1997; Laura nce, 1997) Tree mortality and damage can create gaps in the forest canopy that are temporally stochastic and increase small scale spatial heterogeneity. In addition to the abiotic effects of fragmentation discussed above, there is evidence that popula tion, community, and landscape dynamics have higher frequency or amplitude in fragments (i.e. hyperdynamic) (Laurance, 2002; Laurance et al., 2011) For example, butterfly species richness is more spatially and temporally variable in f ragments (Leidner et al., 2010) This could be an indicator of environmental PAGE 55 55 heterogeneity, either abiotic factors or biotic interactions. However, another study found that fragmentation reduced temporal turnov er among beetle species (Davies et al., 2001) Hierarchical M odels Hierarchical models are useful tools for quantifying variability without observing the underlying pro cesses (Clark, 2003) They allow us to quantify the variance on multiple scales without overfitting because they enable us to assume that the variation at each scale has some underlying distribu tion. Mixed models with multiple random effects are the tool we use here. Random effects are a way to account for how individuals, plots and years may differ from the average pattern. They allow us to quantify how much heterogeneity there is at each level of replication in our long term demographic data set. To quantify the effect of fragmentation on the vital rates of our study species we used generalized linear mixed models (GLMMs). GLMMs allow us to model count data and quantify variation while account ing for correlations and lack of independence (Bolker et al., 2009) We analyze repeated measures of individuals within plots and years to quantify variation at these three scales. Aims Here we used hierarchical models to quantify variation in growth rates of Heliconia acuminata an under story herb in the Amazon. We addressed the following questions: 1) Do es fragmentation alter demographic rates on average? 2) Does fragmentation increase spatial and temporal variation in demographic rates? PAGE 56 56 Methods Data Description Individual plants within plots within habitats were marked and censuses yearly. There are four plots one half hectare in size each located within a different one hectare experimental fr agment There are six plots one half hectare in size located within continuous forest. All plan ts within plots were censused and new seedlings were added to the census each year for twelve years. Seedlings became adults in their second year. Here, we only analyze observations of adult plants. A total of 1364 adults in fragments and 4873 in continuou s forest were marked. Plant size was recorded as the number of shoots for each plant. Reproductive effort was monitored by recording the number of inflorescences for each plant. The census also included survival. We excluded plants with zero shoots because this class was only recorded in 2008 and 2009. More details can be found in Gagnon et al. (2011). Details of the Biological Dynamics of Forest Fragments Project, which this data set was part of can be found in Laurance et al. (2011). Model Formulations We assumed that plant size and habitat (fragment versus continuous forest) were good predictors of plant vital rates because these were previously shown to be good predictors of demographic rates in our study species (Gagnon et al., 2011) We initially formulated all models by examining graphs of the data with smoothed means (Wickham, 2009) For each habitat, we plotted shoot number on the x axis. To quantify variation among individuals, plots, and years, our models included random effects of individual plant, plot, and year. PAGE 57 57 Growth The average number of shoots in the next year S t+1 followed a linear relationship with number of shoots in the current year, S t A linear model of the form (1) is an autoregressive model in which the absolute value of must be less than one to keep S from going to infinity. The slope of our smoothed mean was less than one so it fits this criterion. This agrees with previous work on H. accuminata showing that small plants grow on average, but plants with more than four shoots tend to shrink on average (Gagnon et al., 2011) If growth was deterministic and size was continuous, this would make size at time t be S t = S 0 t 1 + 1 t 1 and size would asymptote at 1 In our non deterministic model, plants below the asymptotic size will probably grow and those above will probably shrink. We included Gaussian random effects on the term An individual with a higher would have an above average asymptote or maximum plant si ze as could be the result of a nutrient rich patch of soil. All plants within a plot with a higher would have higher asymptotes. A year with a higher would increase the expected number of shoots in the following year for all plants in all plots with the absolute amount of increase being the same for all sizes. This could be caused by better weather, e.g. longer rainy season. Because shoot number is not continuous but a count, the potential distributions are Poisson and negative binomial. Although it i s typical to fit these distributions with log link functions (meaning that the response increases exponentially as a function of the predictors), we fit our distributions with identity links (meaning that the response S t + 1 = S t + PAGE 58 58 increases linearly with the predictors as in ) The mean of the Poisson distribution has to be a positive number, so fitting this model with an identity link can be delicate (Marschner, 2010) For this particular data set, we did not have to constrain the parameters in any way or use any special algorithms to find the maximum likelihood estimates We allowed both par ameters slope and intercept is to vary by habitat. Reproduction Plants with five or fewer shoots typically produced no inflorescences, but occasionally produced up to four inflorescences. The smoothed mean number of inflorescences accelerated as plan t size increased, but many large plants still produced no inflorescences. A typical form for an accelerating mean as a function of the predictors is e a S + b i.e. a log link. As with shoot number, inflorescences are discrete counts so their distribution is either Poisson or negative binomial. However, inflorescence number may be zero inflated meaning that there are more zeros than would be predicted by the statistical distribution. Therefore, we fit zero inflated versions of the Poisson and n egative binomial in addition to the standard versions. We included plant size as predictors for both parts of the model zero inflation and conditional mean and allowed all parameters to vary by habitat. We used information criteria AIC to determine if zero inflated models were better fits to the data. We included random effects in the count part of the mode l only, not on zero inflation. Gaussian r andom effects were applied to the intercept term, b. Plants with larger values of b produce more inflores cence s than average plants and the difference increases with size. Plants within plots with larger values of b produce more S t + 1 = S t + PAGE 59 59 inflorescences than average plots and the difference is more extreme for larger plants in that plot. Because a random effect of indi vidual has an influence with strength dependent on the individual's size and not all individuals were measured at all sizes, there was potential for plants measured at small sizes to have larger random intercepts than plants measured at large sizes. T o mak e sure this was not an issue, w e visually accessed the model fits by plotting each individual's size versus random deviate overlaid with a smoothed mean. Survival Survival was more probable for larger plants. Plants with ten or more shoots always survived. Because probability is constrained between zero and one, it is standard to use a logit link. We modeled survival probability as a function of shoot number: logit 1 ( !" We allowed both parameters to vary by habitat. We included Gaussian random effect s on the intercept, d. Units of replication (individuals, plots, years) with larger values of d have higher survival probabilities at smaller sizes and reach the maximum value, one, at smaller sizes. The effect of a larger d diminishes at larger sizes as u nits with average values of d also reach the maximum. Parameter Estimation and Uncertainty Allowing variance to differ by habitat Our goal was to compare the variability in fragments versus continuous forest, so we needed separate random effect variances f or the two habitats. Therefore, we separated our data by habitat type and fit separate models to the data subsets. It would have been possible in ADMB to fit the full model with habitat as a predictor of the fixed effects and the random effect variance; ho wever, it would have been more difficult without giving us more information. PAGE 60 60 The subsets contained four and six experimental plots each, too few levels for a random effect; so we modeled plot as a fixed effect (Bolker et al., 2009) We fit separate intercepts for each plot and calculated the variation among plots as the variation among their separate intercepts. This is not equivalent to a random effect because it lacks shrinkage; however it serves the same purpose of allowing us to quantify variation. F or simplicity, we will use "intercept" to refer to the mean of the plot intercepts ; random effects of plot would be approximate to the deviations from this mean but with shrinkage that pushes them closer to the mean and d epends on sample size Before subsetting the data, we c entered plant size to the mean across habitats. Confidence interval estimation The least controversial method for testing if parameters in GLMMs differ from zero or differ from one another is to estimate the probability distribution of the parameters ; thi s avoids the controversy of having to degrees of freedom necessary for hypothesis testing with F tests or t tests (Baayen et al., 2008) Confidence intervals and inference statistics can be computed using the probability distribution of the parameters. P robability distribution s can be estimated using either Markov chain Monte Carlo (MCMC) sampling or parametric bootstrapping. Here, we use MCMC sampling to estimate the poste rior probability density of the parameters. We used flat priors for all parameters. We fit random effect variances on the log scale to ensure the variances were positive; this meant that their priors were flat on the log scale. We visually assessed the tr ace plots of the samples to check if the chains had run long enough and increased the samples as needed. We ran all models for 40 000 steps and saved 1000 sample s ; additionally, w e ran the survival model for 80,000 steps These samples were insufficient fo r achieving an effective sample size of 1000, so all PAGE 61 61 models should be run longer or using a different MCMC algorithm before publishing the results in peer reviewed literature. In addition to quantifying heterogeneity in vital rates, we wanted to test if th ere were significant differences between habitats in the vital rate after controlling for heterogeneity. Therefore, using our MCMC samples from the models fit to the habitats separately, we calculated the posterior probability of having a greater vital rat e in either of the two habitats. For each habitat separately, for each MCMC sample, we use the estimated fixed effect to predict the expected vital rate for a plant of a certain size. We compared predictions for plants with one to fifteen shoots because fi fteen was the maximum number of shoots observed in the fragments. Then for each number of shoots, we calculated the posterior probability as the proportion of MCMC samples that predict a greater vital rate in one habitat or the other. For a n error rate of .05 in a two tailed test, we used .975 as a significance cutoff. Maximum likelihood e stimation and MCMC Sampling using AD Model Builder Only a subset of our GLMMs could have been fit using standard R packages like glmer MCMCglmm and glmmADMB For our growth model, glmer would fit Poisson with identity link, but not a negative binomial version. Currently, the identity link is not available in MCMCglmm and is not well tested in glmmADMB For our growth model, zero inflation decreases with plant size, so we needed to have predictors on that part of the model; this is possible only without random effects using the function zeroinfl in the R package pscl version 1.04.1 (Zeileis et al., 2008) A zero inflated mixed model can be fit in glmmADMB but not with predictors on the zero inflation (Skaug et al., 2012) Our PAGE 62 62 survival model is a standard form and could be fit with any of the three packages mentioned above. We wrote our own code to fit these models using Automatic Differentiation Model Builder (ADMB) (Fournier et al., 2011) ADMB is an open source program that facilitates maximum likelihood estimation in a very flexible way. It can es timate random effects using Lap lace approximation or Gauss Herm ite quadrature, but we used Lap lace (Skaug and Fournier, 2006) ADMB can do MCMC sampling using either random walk Metropolis Hastings or hybrid MCMC; we used Metropolis Hastings. By default, ADMB starts MCMC sam pling from the maximum likelihood estimate so no burn in period is needed We fit all of the final models in ADMB for consistency and to facilitate MCMC sampling. Consistency is important when using AIC for model selection because some programs, such as g lmer do not include constants when calculating the joint negative log likelihood making them incomparable to other programs that do include the constant We organized data in R and read it into ADMB using the R package R2admb version 0.7.5.2/r73 (Bolker and Skaug, 2012) We also used R2admb to extract maximum likelihood estimates and MCMC samples from the model fits. We used the package coda to analyze MCMC results (Plummer et al., 2006) Initial values for all parameters were first estimated using functions in R for each habitat separately. For our growth models, we used glmer in the R package lme4 version 0.999375 42 to fit a Poisson GLMM wit h an identity link (Bates et al., 2011) B ecause glmer estimated zero variance among individuals in growth, we initialized this to 10 on the log scale in our ADMB model We initialized our zero inflated inflorescence models with estimates from PAGE 63 63 a zero inflated negative binomial GLM using the funct ion zeroinfl in the R package pscl (Zeileis et al., 2008 ) We fit our survival models using glmer from R package lme4 (Bates et al., 2011) Model s election In our growth model we had to choose between the Poisson and negative binomial distributions. We also had that distinction for our reproduction model in addition to determining if the data had m ore zeros than the distribution (i.e. zero inflation). We fit all model versions to the entire data set. We allowed slopes, intercepts, and variance parameters to vary by habitat. The formulation of the negative binomial we used had mean variance and no covariates on k. We bounded k to be between 0.001 and 1000. In both the growth model and reproduction model, k went to 999.9 indicating that there was no overdispersion relative to the Poisson. Also, graphs of the mean variance relationship in our data did not look overdispersed. For the reproduction model, we used information criteria (AIC) to assess if the zero inflated Poisson was a better fit than the Poisson. Assessing Spatial P atterns We calculated conditional residuals for all models, bo th mixed effects models and fixed effect only versions. To check for spatial patterns and edge effect s in the fragments, we plotted the spatial locations of the residuals as well as the random effect estimates. We also calculated viariograms using the func tion Variogram from the R package nlme version 3.1 103 PAGE 64 64 Ignoring Demographic H eterogeneity We fit each of our models without random effects of individual, to compare with the full models and test if this variation would be lumped in with environmental variation as observed by Melbourne and Hastings (2008). Correlation Among Vital R ates We estimated the correlation among vital rates of our random effect s with in each unit of replication (year, plot and individual). The advantage of looking at random devia tes rather than the average vital rates for subsets of the data is that deviates account for differences in the size distribution of the population. A population's size distribution could be changing through time and space and an individual's size affects its vital rates. Looking at the random deviates allows us to independently examine the effects of temporal and spatial environmental heterogeneity or demographic heterogeneity among individuals while controlling for size structure. Our main objective in th is study was not to quantify the correlation among vital rates. If we wanted to get the best estimate of the relationships among vital rates, then we should have modeled a random effect's influence on the three viral rates as a multivariate normal distribu tion instead of three separate normal distributions. This method is known to give better estimates of the correlations, but similar accuracy for point estimates as the methods used here (Evans and Holsinger 2012). The combined model would be quite complex given that we have three random effects influencing the outcomes of three vital rates. The random effect structure would be three multivariate normal distributions, one for each unit of replication, and each multivariate normal distribution would have mean zero and a three by three variance covariance matrix to characterize the relationships between the vital rate responses. This model could be PAGE 65 6 5 implemented in ADMB, but compared to the three separate models (one for each vital rate), it would take longer to converge and to do MCMC sampling for. Also, the three separate models can be run in parallel. Results and Discussion Fixed Effects Growth is slower in the fragments than in the continuous forest for plants with more than five shoots (Fig s. 4 1a 4 2a ) Su rvival rates are lower in fragments than in the continuous forest for plants with more than three shoots (Fig s. 4 1c 4 2c ). The best model of inflorescence production was zero inflated Poisson compared to Poisson ( AIC 499.6). Inflorescence production is significantly higher in fragments than in the continuous forest for plants with more than ten shoots (Fig s. 4 1b 4 2b ). Higher reproduction at a given size in fragments could be the result of evolution because higher mortality, as observed in fragments, can select for earlier maturation (Cohen, 1971; Kozlowski, 1992; Jacquemyn et al., 2012). This is also a potential explanation for why growth rates are reduced for plants in fragments: they are allocating a greater proportion of their resources to reprodu ction compared to similar sized plants in continuous forest. However, for such a long lived perennial plant it seems unlikely for evolution to have occurred during thirty years since fragmentation; we discuss other potential explanations below. Our growth model is equivalent to a continuous form of the growth model previously done by Gagnon et al 2011 Our growth estimates agreed with those of Gagnon et al with the single exception that our estimate of the growth rate of plants with one shoot was 2% abov e their upper confidence limit PAGE 66 66 Random E ffects Reproduction was the most variable vital rate (Fig 4 3b ). Reproduction varied greatly among individuals, approximately equivalently in continuous and fragmented forest. Among year variation in reproduction w as marginally greater in continuous forest than in fragments (. 924 posterior probability). W e know that variation can change the mean outcome of a process. We consider the expected value, across individuals, of inflorescences per individual [ ( 1 p ) ] equal to ( 1 p ) c o r ( p ) v a r ( p ) v a r ( ) The correlation between p and is always negative because p is a decreasing function of shoot number and is an increasing function of shoot number. Therefore, heterogeneity in either parameter will increase the average rate of inflorescence in the population and therefore increase the population growth rate. Growth was the least variable vital rate (Fig 4 3 ). There was marginally mo re variation in growth rate among fragment plots, than among continuous forest plots (.9 25 posterior probability). There was marginally more variation among years in continuous forest than in fragments (.9 35 posterior probability). There was little to no d etectable variation among individuals in growth rate; this agrees with the results of simulations we conducted and with findings of Gagnon et al. (2011) It is difficult to detect a random effect of individual growth when each individual in measured at mos t eleven times. Also, our method was not designed to detect variation among individuals that changes through time, the kind that would arise from temporary canopy gaps that quickly become overgrown. PAGE 67 67 Survival varied approximately equally among plots and am ong years. Variation among individuals was undetectable in the fragmented forest, and only marginally detectable in the continuous forest (Fig 4 3 ). For individual variation in survival in continuous forest, the mean of our MCMC samples differed from the MLE even after running 80,000 samples. We suspect that there is little to no information in the data to inform this parameter. Ignoring Demographic H eterogeneity In most models, the random effect of individual was close to zero so ignoring individual variation had little effect on the model. Reproduction was the exception where estimated variation among individuals was significantly greater than zero. As expected d ropping the random effect of individual caused the random effect of plot to increase in both models (! 0.01 CF 0.03 F ) and caused the random effect of year to increase in the fragments (! 0.09) with a negligible change to random effect of year in the continuous forest (! 0.005). These results agree with the findings of Melbourne and Hastings (2008): ignoring variation among individuals causes variation to be lumped in with the other sources of variation. However, only a small portion of the variation among individuals was transferred to the other random effects so ov erall variation was under estimated in models that ignore individual variation (Fig 4 4). Correlation of Vital Rate Random Deviates Variation in vital rates was fairly synchronous between fragments and continuous forest (Fig 4 5). The random deviates showed strong within year correlation between habitat types: growth 0.7, reproduction 0.5, and s urvival 0.6. There was some synchronicity among vital rates. Growth and reproduction were positively correlated within years with 0.6 co rrelation in continuous forest, and much less in fragments with PAGE 68 68 0.0 8 correlation. Growth and survival were positively correlated within years in continuous forest 0.5, but slightly negatively correlated in fragments 0.00 2 Reproduction and survival were n egatively correlated within years 0.3 in continuous forest and 0.5 in fragments. This indicates that years with high reproductive rates have reduced survival rates, especially in fragments. Within plots, reproduction and growth were positively correlate d (0.7 in forest, 0.9 in fragments). As with years, reproduction and survival were negatively correlated within plots ( 0.6 in forest, 0.3 in fragments). Growth and survival were positively correlated in fragments (0.2), but slightly negatively correlated in continuous forest ( 0.01). Within individuals, growth and reproduction wer e positively correlated (0.1 CF and F ); growth and survival were negat ively correlated ( 0.2 CF, 0.03 F ); survival and reproduction were slightly negatively correlated ( 0.06 CF 0.05 F ). At all levels of replication reproduction and survival were negatively correlated. At the individual level, this could be due to a life history trade off. However, the correlations were stronger at the plot and year level indicating an environ mentally driven process. An environment that is bad for survival could be good for reproduction or a plant could up regulate reproduction when it senses an environment that is bad for survival. This warrants more investigation using field experiments. Spat ial P attern Our spatial graphs and variograms showed no patterns in residuals nor individual random effects. This indicates that if spatial heterogeneity within plots was affecting the vital rates, then it was on too fine of a scale for us to detect. PAGE 69 69 Con clusion Fragmentation caused plants to have lower survival, lower growth, and higher reproductive rates. If fragmentation increases mortality due to tree falls, air temperatures, or soil moisture then natural selection would cause plants to begin allocatin g more energy to reproduction than growth. So the reduction in growth rates could be in indirect effect of fragmentation. This pattern was also seen in the random effects: survival and reproduction random deviates were negatively correlated at all levels. This indicates that it could be a plastic response of the plants to some detectable environmental cue such as soil moisture. Cohen (1971) postulated for annual plants that if growth is terminated at a low soil moisture level, then soil moisture would be a good indicator of the optimal time to reproduce within a season. In future work, this optimization problem could be extended to perennial plants and coupled with experiments. Alternatively, the same environmental variable could be directly affecting both vital rates (i.e. increasing mortality and increasing reproduction). There was significant temporal heterogeneity and landscape level heterogeneity in all vital rates. There was strong demographic heterogeneity in reproduction. Individuals in good microenv ironments may be allocating all additional resources to reproduction. Although there were no significant differences in the heterogeneity between continuous forest and fragments, there were some interesting patterns. In general, it is difficult to detect a difference in variance estimates because statistical power is low. There seemed to be more individual variation in reproduction in fragments. There seemed to be more temporal variation in all three vital rates in the continuous forest. This could lead to greater variability in the population growth rate in continuous forest and therefore a reduced long term growth rate. PAGE 70 70 Overall, vital rates were higher in continuous forest, but because they were temporally more variable, the long term growth rates of the populations could in theory be equal across habitats. PAGE 71 71 Fig. 4 1 The solid line is the posterior probability density of having different vital rates at a given plant size. The dashed lines are the cutoffs for a two tailed significance test with .05 error rate: probabilities above .975 indicate that the vital rate is significantly greater in fragments and values below.025 indicate that the vital rate is significantly greater in continuous forest. Each panel is a separate vital rate modeled with separa te GLMMs for each habitat. PAGE 72 72 a b c Fig. 4 2 Changes in plant size (shoot number) re lative to stasis (black line) (a), inflorescence production (b), and survival (c ) were modeled separately for habitats of continuous forest (red) and fragmented forest ( blue) using GLMMs with fixed effect of shoot number (x axis) and random effects of individual, plot, and year. The fixed effect predictions plus and minus two standard deviations of year are plotted to show the effect of greater variance in continuous fore st. Circles represent observations with size representing the number (n) of observations at that point. PAGE 73 73 Fig. 4 3 Each panel represents variation quantified by fitting a GLMM to a different vital rate: (from top to bottom) changes in shoot number, inflorescence production, and survival. The y axis is the estimated standard deviation in random intercepts among units of replication divided by the overall intercept estimate for the same model. The color represents continuous forest (red circles) or fra gmented forest ( blue triangles). The points represent the mean of MCMC samples with bars for the 95% CI. Growth and reproduction models were run for 40,000 MCMC steps; the survival model was run for 80,000 steps. PAGE 74 74 Fig. 4 4 Inflorescence rate estimated with mixed effects models with a random effect of individual (circular points) and without (triangular points) a random effect of individual. The random effect is on the x axis. The color represents fragmented forest ( blue ) or continuous forest (red). The maximum likelihood estimates are plotted. PAGE 75 75 Fig. 4 5 Variation in vital rates through time was fairly synchronous in forest (red circles) and fragments (blue triangles), but the different vital rates (different panels) showed both positive (growth and reproduction) and negative (reproduction and survival) correlation. PAGE 76 76 CHAPTER 5 SCALIN G UP VARIATION Introduction Demographers have long known that variation has consequences for population growth. Predictions from stoch astic models can greatly differ from those of deterministic models. Temporal variation increases the chance of extinction (Kendall, 1949) Landscape level spatial variation helps populations persist (Bolker, 2003) Demographic stochasticity (i.e. random variation in the fates of individuals) increases extinction risk of small populations (Morris and Doak, 2002) Demographic heterogeneity is a component of demographic variance defined as variation among individuals in traits that affect demographic rates (Fox et al., 2006; Melbourne and Hastings, 2008; Kendall et al., 2011) These traits could be observable (e.g. size, sex, habitat) or they could be unobservable (e.g. small scale spa tial environmental heterogeneity in plants, territory quality in animals). Demographic heterogeneity differs from demographic stochasticity because demographic stochasticity is variation in the fates of individuals due to random chance (Roughgarden, 1975) Like demographic stochasticity but unlike other forms of environmental stochasticity, the effect of demographic heterogeneity is strongest for small populations (Melbourne and Hastings, 2008) However, unlike demographic stochasticity, demographic heterogeneity can reduce extinction risk. Some studies look at the effects of demographic heterogeneity while others use matrix models or integral projection models to examine the effects of temporal and spatial heterogeneity. It is not simple to include all three in struc tured population models. Demographic heterogeneity in growth was inco rporated into matrix models PAGE 77 77 (Pfister and Wang, 2005; Zuidema et al., 2009) and an integral project ion model (de Valpine 2009) ; its inclusion improved predictions of the population's size structure, growth rate, and individuals' fates. It would be challenging to expand the capabilities of either matricies or integral projection models such that they can incorporate continuously distributed demographic heterogeneity combined with spat ial and temporal heterogeneity. Variation among individual reproductive rates can be modeled as a multi type branching process (Lloyd Smith et al., 2005) but the problem becomes intractable with environmental and demographic heterogeneity. Individual based models are thus far the only framework that can incorporate variation on multiple scales into a model of population growth, so that is the tool we use here (Conner and White, 1999; Buckley et al., 2003) Melbourne and Hastings (2008) were able to tease apart the effects of demographic stochasticity, demographic heterogeneity, and environmental stochasticity in a population of red four beetles Their lab raised populations had non overlapping generations. In their simple system, they were able to derive a tractable form for the expected variation in the number of individuals in the next generation as a function of the different sources of heterogeneity and population size They found that extinction risk could only be correctly estimated if demographic heterogeneity was included in the a nalyses in addition to environmental and demographic stochasticity. The study presented here differs from that of Melbourne and Hastings in that we model a hermaphroditic organism with overlapping generations and we ignore density dependence. Also, the siz e structure of our population is known which accounts for most of the demographic heterogeneity. In this study, we are primarily interested in the PAGE 78 78 part of demographic heterogeneity caused by unobservable factors because, just like other components of demog raphic stochasticity, it has consequences for population viability and yet it is rarely incorporated in population models. Repeated measures of individuals' vital rates are typically needed to quantify individual variation in vital rates due to unobservab le factors. This type of data is commonly collected and analyzed in a repeated measures statistical framework, but among individual variation is not typically the result of int erest. In our previous work (Chapter 4 ) we quantified variation among individual s in vital rates using generalized linear mixed models. We included random effects of year, experimental plot, and individual to quantify variation on those scales. We found that there was substantial variation among individuals in reproductive effort, mor e variation than among years or among plots. Here we quantify the effect it will have on population growth relative to other sources of stochasticity. Study System Heliconia acuminata is a long lived perennial plant common to the understory of the Amazon. Experimental fragments were isolated in the mid 1980s at the Biological Dynamics of Forest Fragments Project (described in Laurence et al. 2002). Heliconia acuminata have almost no a sexual reproduction. Average densities at our study site were 5870 +/ 858 (s d ) per continuous forest hectare and 1696 +/ 210 (s d ) per fragment hectare. Dispersal is low (Uriarte et al., 2011) PAGE 79 79 Predicted Effects of H eterogeneity Changes to the mean population growth rate Demographic heterogeneity should incr ease the population growth rate. Demographic heterogeneity in survival increases the population growth rate because individuals with below average survival will die and the remaining individuals have above average survival, a process called cohort selection (Kendall et al., 2011) When inflor escence is distributed as zero inflated Poisson as it is in the case of Heliconia demographic heterogeneity in inflorescence rate should increase the average infloresce nce rate of the population (Chapter 4 ). Although we did not observe individual variatio n in growth rates of Heliconia we would expect it to increase the population growth rate because fast growers can contribute disproportionately to population growth (i.e. the fast growth effect) (Pfister and Wang, 2005; Zuidema et al., 2009) Because dispersal is low for Heliconia acuminata we expect landscape level spatial variation to increase the population growth rate through time. Locations with higher survival and reproductive rates will contain an increasing proportion of the population and thus the population's mean survival and reproduction will increase. With low distance dispersal, new individuals are born into locations that are similar to their parent. Th us, landscape level spatial environmental variation and its effect on an individual's vital rates can also be interpreted as heritability in those vital rates (Fox, 2005) Kendall et al. (2011) found that if there is a positive correlation between the reproductive rates of parents and offspring (such as that caused by a shared environment) then heterogeneity in reproductive rates will increase the population's growth rate. Here, we model density independent population growth, but this princi ple PAGE 80 80 also applies to density dependent population growth and is sometimes known as habitat selection (Bolker, 2003) Changes to the variance of the population growth rate A population's stochastic growth rate (geometric mean) is more informative than its arithmetic mean growth rate be cause the stochastic growth rate reflects how variation reduces the geometric mean. So we also quantify how heterogeneity affects the variation around the mean growth rates. Demographic heterogeneity in reproduction can either increase or decrease the vari ance in the population's reproduction depending on the distribution of the reproductive rate (Kendall and Fox, 2003 ) Inflorescence production in Heliconia acuminata was previously shown to be zero inflated Poisson so we estimate the expected variance of that distribution in the methods section. We expect the effects of spatial heterogeneity at the landscape scale to have effects similar to demographic heterogeneity b ecause both forms of heterogeneity are consistent through time. We expect temporal heterogeneity in all vital rates to increase variance in the population growth rate. Here, we use data from our real system to quantify heterogeneity at different scales and to parameterize an individual based model simulating the population dynamics in situations with different combinations of heterogeneity. W e analyze results of our simulations to answer the following questions: (1) How does heterogeneity at different sca les affect the mean growth rate of a population? (2) How much variation in the population growth rate is attributable to each scale of heterogeneity? (3) Does individual variation in reproduction scale up to increase variability in the population's growth and thus lower its stochastic growth rate? (4) Do observed differences in PAGE 81 81 demographic rates explain the change in population structure since fragmentation occurred 30 years ago? Methods Data D escription We had sufficient data to quantify demographic, spat ial, and temporal heterogene i t y in adult vital rates. We quantified spatial and temporal heterogen eit y in seedling survival. Seed germination was habitat specific. All other parameters were measured in continuous forest only. We use parameter estimates of adult vital rates modeled in chapter 3. The number of flowers per inflorescence was observed on 57 inflorescences from multiple plants in 1997 in continuous forest only. The proportion of flowers becoming fruits was observed on 22 plants in 1997 in continu ous forest only. The number of seeds per fruit was observed in 873 fruits taken from continuous forest habitat in 1998. Seed germination probabilities, in both continuous forest and fragments, were published in (Bruna, 2002) Seedling survival and growth were monitored in the same data set as adult demography for 12 years in ten plots in both habitats. Parameter E stimation Seedling survival was modeled separately for each habitat type using logistic regression using glmer in R. Plot was the only fixed effect and year was the only random effect. Plot was used as a fixed effect because there were too few to fit it as a random effect (Bolker et al., 2009) .We calculated the landscape level spatial heterogeneity of this process as the standard deviation of the plot specific intercepts. Number of flowers per inflorescence was modeled as negative binomial because AIC scores indicated it was a better distribution than the Poisson distribution ( PAGE 82 82 AIC=106). We used the function glm.nb in the MASS package in R. The proportion of flowers that became fruits on a plant was modeled as binomial (logit link) with a random effect of plant using glmer; AIC scores indicated that this was better than a model without a random effect of plant ( AIC=83). However, because individuals whose fruit production was measured do not correspond to those in the long term study, we cannot say if this random effect of individual would be persistent through time or if it could be explained by covariates. Each fruit produced up to three seeds; the probabilities of each outcome from zero to three seeds were calculated directly because they were underdispersed compared to the Poisson distribution. The distribution of the number of shoots for new adults was qu antified by subsetting the long term dataset to look at only seedlings; then we directly calculated the probability that a seedling would have zero, one, two, three, or four shoots in the following year. Covariance Among Vital R ates We calculated variance covariance matrices for the adult vital rates at each level of replication: individual, plot, and year. Random deviates from the mean vital rates were estimated using a normal distri bution with mean zero (chapter 4 ). Thus, we calculated separate variance covariance matrices for individuals, years, and plots. The diagonal elements of the matricies contain the variance in each vital rate. The off diagonal elements contain the covariance among vital rates at that unit of replication. A more accurate way of es timating the covariance of vital rates would be to fit one combined model with the three vital rates as a multivariate response dependent on the underlying variance covariance structure we secondarily estimated here (Evans and Holsinger, 2012) However, this would have added substantial complexity to our model fitting procedures and estimation of covariance was not a primary goal of this study. PAGE 83 83 Population Growth Simulations We simulated population growth using parameters from the vital rate models described above. When applicable, we used values specific to either habitat type: continuous or fragmented forest. We initialized each population with 10 individuals and assigned each individual three random dev iates for the intercepts of each vital rate: growth, survival, and conditional mean inflorescence number. The zero inflation part of our inflorescence model did not include heterogeneity other than covariates (chapter 3). We simulated the deviates from a m ultivariate normal distribution with mean zero and variance covariance matrix calculated as described in the section Covariance among vital rates and using mvrnorm in the package MASS in R. They should be approximately normally distributed because random e ffects were fit using a normal distribution with mean zero. We assigned the starting number of shoots for each individual by resampling with replacement from the observed shoot number of adults in the habitat. We located each individual randomly in one of five plots. Using the same method as for individuals, f or each simulated plot and year, we created random deviates from the overall intercepts using multivariate normal distributions with mean zero and variance covariance matrix calculated as descr ibed in above. To project the population forward in time, we simulated survival, growth, and reproduction repeatedly. We iterated the population (storing its size each year) for as many years as we expected it to take for the population to grow to 100,000 individuals based on preliminary estimates of lambda. We chose 100,000 individuals as a target population size because it was few enough to store in 20GB of RAM and many enough to see a drastic reduction in demographic stochasticity. We ran 100 replicates To PAGE 84 84 estimate demographic stochasticity at smaller population sizes, we ran another 10 replicates for 10 years each, starting with a single individual. All simulations included fixed effects. We included random effects (RE) of individual, year, and plot in 6 combinations (Table 5 1 ). Therefore, in situations 1 through 4, all plots were identical, but in situations 5 and 6, each replicate had a different combination of five randomly generated plots. So this random effect increased variation among replicates but not among years within a replicate. In each year, we calculated each adult's probability of survival, expected number of inflorescence s and expected number of shoots using the models fit to our real data. All parameter values depended on the habitat type being simulated except for models of flowers, fruits, and seeds. We included appropriate random effects for each situation Individual survival was simulated as a binary process. We removed deceased individuals from the population. We simulat ed the number of shoots for an individual in the next year as a Poisson deviate with mean equal to the expected number of shoots. Simulating reproduction required multiple steps (from inflorescence to new adult in 2 years): plants produce a number of inflo rescences; inflorescences contain an average of 23 flowers; each flower can become one fruit; each fruit can produce up to 3 seeds; each seed can germinate to become a seedling; each seedling has a probability of surviving its first year and becoming an ad ult. We simulated an individual's total flowers as the sum of samples from a negative binomial distribution, one sample for each inflorescence. We gave each individual a random deviate to include as the random effect of individual on the probability of fl owers becoming fruits. Then that individual's PAGE 85 85 fruits were simulated as binomial with n equal to that individual's simulated number of flowers. Then we summed the fruits produced by the population and treated them equally. For each fruit, we simulated a uni form random number and based on that number, gave the fruit between zero and three seeds with cutoffs determined by our model fit to seed data. We summed the population's seeds in a year and used that as n in a binomial random number with habitat specific probability of germination. Seedling survival through the first year was a binomial random number with n equal to the number of seeds that germinated. We stored the number of new adults to be added to the population two years from the year when the inflore scences were produced. New adults were added to the population based on the simulation from the previous two years because it takes two years for a seed to become an adult; for this reason, we ignored the first three years of all simulations when analyzing results. New adults were assigned random deviates for their three vital rates. We simulated shoot number by drawing a uniform random number between zero and one and using cutoffs based on the observed number of shoots that seedlings had the following year Decomposing Effects of H ete rogeneity on Population G rowth For our analyses, we pooled the geometric growth rates across all simulated replicates and all years except the first three (before recruitment), and kept habitats and situations separate. We cal culated geometric population growth rate as the population size in the next year divided by population size in the current year, l t = n t+1 /n t We kept track of their corresponding current population sizes, n t because variation in population growth rates due to demographic stochasticity should decay with 1/ n t The log population growth rates are normally distributed (at time goes to infinity) so we analyzed these values (Tuljapurkar and Orzack, 1980; Morris and Doak, 2002) We'll refer to the PAGE 86 86 log of the population growth rate as r. For each situation in each habitat separately, we calcula ted the mean of the log growth rates, r Heterogeneity at each scale has effects on the mean and variance of the growth rate as described in the introduction section. To enable testing for these effects, we created predictor variables of the random effects included in each situation with elements equal to zero or one ac cording to the columns of T able 5 1 Variance in the population growth rate can be decomposed as the sum of the effects of environmental stochasticity (temporal Y 2 and spatial P 2 ), individual variation (demographic heterogeneity) and demographic stochasticity both proportional to the inverse of the number of individuals in the population, I 2 ( 1 / n ) and D 2 ( 1 / n ) (Morris and Doak, 2002; Melbourne and Hastings, 2008) D emographic heterogeneity cause s population reproduction per y ear contributing to n t+1 to vary on the scale n t {Kendall:2003vu} thus var( n t+1 ) depends on n t but in calculating the geometric growth rate at a given population size var ( n t+1 / n t ) depends on 1 / n t because var ( X / a ) = var ( X) /a 2 r 2 = Y 2 + P 2 + I 2 ( 1 / n ) + D 2 ( 1 / n ) For each observed log growth rate, we calculated the squared deviation from the mean for that particular situation r r ( ) 2 which allowed us to estimate the effect sizes of the three sources of heterogeneity on variance r 2 We could not directly calculate the variance for each treatment because demographic stochasticity depends on population size and for many population sizes, there was only one observation. The distribution of r r ( ) 2 is related to, but not equal to the chi squared d istribution with one degree of PAGE 87 87 freedom (see appendix H for details ). Because of the unusual distribution, we fit the model in ADMB. We combined the values for all simulations but kept habitats separate. The model included separate intercepts for variance due to temporal heterogeneity and landscape heterogeneity, a slope on the inverse number of individuals (1/ n ) for demographic stochasticity and an interaction term between individual heterogeneity and (1/ n ). We calculated likelihood profiles of our estimates to get confidence intervals. To estimate the stochastic population growth rate, we found the geometric mean of the geometric growth rates. The geometric mean of the population's geometric growth rate is defined as 1/ 123 (...) t Gt lllll = and on the log scale l o g ( l G ) = 1 t l o g ( l 1 ) + l o g ( l 2 ) + + l o g ( l t ) ( ) = 1 t r 1 + r 2 + + r t ( ) This allows us to use linear regression to estimate log( l G ) (Dennis et al., 1991; Morris and Doak, 2002) The standard method developed by Dennis et al. is to fit a linear model with only an intercept, but in our case, we allow that in tercept to change depending on the sources of heterogeneity included in the simulation. Our simulations with different types of heterogeneity had different variances around r as estimated in the previous section. To account for the unequal variances, we u sed fits of our variance model to predict the variance associated with each observation of r We used the inverse of each estimated standard deviation, r 1 as weights to do weighted least squares regression on r using the function lm in R For the mean, we included predictors for each type of heterogeneity included in the model. To account for the fact that plot heterogeneity could be interpreted as heritability and lead to evolution through time, we fit versions with and without an inte raction between plot heterogeneity and PAGE 88 88 year. We used information criteria to see which was a better fit to the data. We used t tests to test if the log of the geometric growth rate was greater than zero and to test if the contrasts differed from zero. Scaling Up Variation Among Individuals in I nfloresce nce R ates We expected that variation among individuals in their reproductive rate would scale up and increase variation in the population's reproductive rate, and therefore in its growth rate (Kendall an d Fox 2003, L l oyd Smith et al. 2005). Kendall and Fox (2003) showed that variance among individuals in some trait affecting reproduction will affect year to year variation in the population's total reproduction in a manner dependent on the second derivativ e of the demographic variance in reproduction evaluated at the mean trait value. Specifically, v a r ( F ) N V D m [ ] ( ) + 1 2 V D # # m [ ] ( ) + 1 $ % & ( ) v a r ( m ) + 1 3 V D # # # m [ ] ( ) s k e w ( m ) + + / where var(F) is the variance among years in the population's total reproductive rate, E[] denotes the expected value, V D is the demograp hic variance (i.e. the variance in reproduction caused by demographic stochasticity determined by the statistical distribution), and m is the trait that determines reproductive rate of an individual. The reproductive trait we examine here is the number of inflorescences produced in a year. In our case, as in many cases, inflorescence production is zero inflated Poisson. Therefore, the relevant trait m is actually a function of two variables : (1 p ) and demographic variance is V D ( p ) = ( 1 p ) + p ( 1 p ) 2 We expanded the calculations of Kendall and Fox (2003) so that the equation could apply to cases with two parameters (see appendix B for details). PAGE 89 89 v a r ( F ) N V D p [ ] [ ] ( ) + 1 2 v a r ( p ) # 2 V D # p 2 + 2 c o v ( p ) # 2 V D # p # + v a r ( ) # 2 V D # 2 $ % & ( ) + v a r m ( ) + h o t + / Here, we've used a bar to denote the expectation. The var( m ) term in the equation above is the variance among individuals in their expected reproductive rate. For simplicity, we will ignore higher order terms ( h.o.t. ). More specifically, for the zero inflated Poisson distribution the variance in the population's inflorescen ce production is v a r ( F ) N ( 1 p ) + p ( 1 p ) 2 + v a r ( ( 1 p ) ) + 1 2 2 2 v a r ( p ) + 2 1 + 2 ( 1 2 p ) ( ) c o v ( p ) + 2 p ( 1 p ) v a r ( ) # $ % & ( ) + + + (5 1) We simulated inflorescence production using parameters fit to our data set. This enabled us to numerically verify the above calculations and to estimate the scale of variance caused by this process in a real system. For both continuous forest and fragmented forest, we simulated inflorescence production for populations of size 1,000, 10,000, 100,000 and 1,000,000. We initialized each population by resampling (with replacement) from the observed distribution of shoots for that specific habitat type. We simulated the dynamics with and without random effects of individual; in the situation without a random effect, individuals still varied in their shoot number. We assigned each individual a random intercept for the conditional mean of its inflorescence rate, Based on the individuals' shoot number and random intercept (when applicable), we calculated their values of and p. For each replicate, for each individual, we samp led from a zero inflated Poisson distribution with that individual's parameters, and p; then we summed across individuals to get the total number of inflorescences. We did 1000 replicates for each population and resampled to get new p opulations 10 times. For each PAGE 90 90 population, we calculated the average and p as well as the variance of and p and their covariance; we also calculated the variance among individuals in (1 p) so that these quantities could be plugged into Eq. (5 1) to check if it was a good predictor of the variation in the population's total inflorescent number. We calculated the variation in total inflorescence number across the 1000 replicates for a given po pulation. Simulations of F ragmentation To check if it was possible that the demographic structure of populations in continuous forest could be transformed into that of forest fragments within the timespan of the fragmentation experiment, we simulated this transition. We started with populations representative of the continuous forest and simulated their dynamics for 30 years with rates dependent on models fit to fragmented forest. We initialized each population with 1000 individuals. Each individual was ass igned a number of shoots by resampling with replacement from the distribution of number of shoots observed in adults in continuous forest between 1998 and 2009. Simulations were done in the same way as described in the previous section, but we kept track o f the distribution of number of shoots through time. Results and Discussion Parameter estimates can be found in appendix F and chapter 4. Effects Of H eterogeneity On Average Population Growth All forms of heterogeneity increased the average value of the geometric growth rate, l. In continuous forest, the average was 1.23 (0.034 SE); demographic heterogeneity increased the average by 0.05 (0.036 SE); temporal heterogeneity increased the average by 0.10 (0.043 SE); landscape level spatial heterogeneity increased the average by 0.09(0.046 SE). In fragments, the average was 1.12 (0.005 PAGE 91 91 SE); demographic heterogeneity increased the average by 0.003(0.005 SE); temporal heterogeneity increased the avera ge by 0.02 (0.006 SE); landscape level spatial heterogeneity increased the average by 0.05 (0.007 SE). These results indicate that Jensen's inequality was acting at all levels to increase the population's mean growth rate. This may have been caused by the log link in our equation of inflorescence production. A log link in Poisson regression implies that the mean of the Poisson is an accelerating function of the predictors. Accelerating functions with variation in the predictors have higher averages than th e function evaluated at the average predictor value (Ruel and Ayres, 1999). More variation in the predictor increases the mean function value. In our model, random effects of individual, year, and plot were predictors and therefore they increased the mean. Effects Of Heterogeneity On Stochastic Population Growth Demographic heterogeneity and temporal heterogeneity increased the stochastic geometric growth rate in both habitats. Demographic heterogeneity increased the stochastic geometric growth rate by 1.1 % in continuous forest and .8% in fragments. Temporal heterogeneity increased the stochastic geometric growth rate by 1% in continuous forest and .4% in fragments. The strong effect of demographic heterogeneity could only have been caused by heterogeneity in inflorescence production because that was the only vital rate with significant variation among individuals. The best models in both habitats indicated that landscape level spatial heterogeneity reduced the average growth rate in early years ( 4.2% CF; 1.1% F), but that the growth rate increased through time (.4% per year CF; .1% per year F) ( AIC= 48 CF and 95 F). This model had a slope on year only in situations that included spatial heterogeneity; the slope accounts for the potential for a greater pro portion of individuals PAGE 92 92 to be located in plots with higher random deviates through time. In contrast, the model without a slope on year indicated that landscape level spatial heterogeneity increased the average growth rate by 3% in both habitats. Effects O f Heterogeneity On V ariance In Population Growth All types of heterogeneity increased variation in the log population growth rate. Demographic stochasticity had the largest influence in small populations more so in continuous forest (5.7; 95%CI: 5.49, 5.90 ) than in fragments (2.5; 95%CI: 2.43, 2.57), but its influence decays as population size increases because these parameter estimates are slopes on the inverse population size. A population living at average density in one hectare would experience varianc e due to demographic stochasticity on the order of the other sources of heterogeneity: 0.0008 in continuous forest or 0.001 in fragments. Demographic heterogeneity, which also depends on population size, was also stronger in continuous forest (1.44; 95%CI: 1.12, 1.77) than in fragments (.82; 95%CI: .70, .93), but in a typical population in one hectare the direction of strengths would be reversed (.00025 CF, .00048 F) because the population sizes differ. If we consider large populations, in fragments, landsc ape level spatial heterogeneity was the most influential source of variation (0.010; 95%CI: 0.0094, 0.011) followed by temporal heterogeneity (0.0063; 95%CI: 0.0060, 0.0065). In continuous forest, temporal heterogeneity was the most influential source of v ariation (0.046; 95%CI: 0.043, 0.049) followed by landscape level spatial heterogeneity (0.035; 95%CI: 0.029, 0.040). These values might only have meaning in their relation to one another. At a given population size, demographic stochasticity had twice as much influence in continuous forest than in fragmented forest. In our models of growth and reproductive effort, variance is an in creasing function of the mean; l arger individuals with higher PAGE 93 93 predicted means also have more variance in their fates (both in the statistical model and in the raw data). I ndividual sizes in continuous forest span nearly twice the range as in fragments (maximum 24 versus 15 shoots) so there is more demographic stochasticity The influence of temporal heterogeneity was an order of magnitude larger in continuous forest than in fragments; this can be attributed to greater temporal heterogeneity in all three vital rates in continuou s forest (see chapter 4 for parameter estimates). The temporal heterogeneity observed in each demographi c rate was on the same order in continuous forest and fragments; so there must have been some synergy that increased the order of the differences at the population level. Inflorescence Production Demographic heterogeneity increased inflorescence production in our simulations, confirming the effect we expected to see from Jensen's inequality. For populations of size 10 3 and 10 4 variation among populations caused overlap in the mean inflorescence rate with and without demographic heterogeneity, but for larger populations there was clear separation (Fig 5 2 ). Observed variance in the population's inflorescence production was near the approximat ion calculated using Eq. ( 5 1 ) (Fig 5 3). On average, without demographic heterogeneity, the variance amon g replicates in inflorescence rate per individual was 0.023 in continuous forest and 0.0009 less in fragments. Demographic heterogeneity increased the variance by 0.0014. To scale up from the variance per individual to the entire population, multiply by th e population size. Most of the variation among individuals was due to differences in their size. PAGE 94 94 Transition From Continuous Forest To Fragment In our simulation of the transition from the stage structure of the continuous forest to the stage structure of the fragments, the average population converged to the steady stage structure quickly: within the first 10 years when ignoring heterogeneity and after 30 years wh en including heterogeneity (Fig 5 4 ). The stage distribution quickly changed as a result of t he growth dynamics: plants from continuous forest that were far above the asymptotic size defined by the fragment parameters were highly likely to shrink (see chapter 4 for details). Compared to observations of the real fragments in 2008 and 2009, our simu lations after 30 years under predicted the number of plants with 0, 2, or 3 shoots, and over predicted the number of p lants with 6 to 8 shoots (Fig 5 5 ). The version without heterogeneity over predicted the number of seedlings, and the version with hetero geneity varied greatly in the 95% CI of seedlings. Conclusions The transition to the stage structure of a fragmented forest occurred quickly in our simulations. This supports the idea that the differences in stage structure observed between habitats, 30 y ears after fragmentation, are caused by the changes in vital rates due by fragmentation. Heterogeneity increased mean growth rates through Jensen's inequality. Stochastic growth rates were lower than the means due to increased variability, but higher than stochastic growth rates in simulations without heterogeneity. Thus, heterogeneity had a net positive effect. In all situations, the stochastic growth rate was higher in continuous forest than in fragments. Thus, as man y other studies have shown, fragmentat ion has a negative effect on population growth. Our results are unique PAGE 95 95 because we've shown this pattern holds even after accounting for three scales of heterogeneity PAGE 96 96 Fig. 5 1 Geometric growth rates were calculated from simulations of population growth in two habitats: continuous forest and fragments (upper and lower rows of panels respectively). Each column of panels represents a different situation. All situations include demographic stochasticity; panels marked I' include demographic heterogen eity; panels marked Y' include temporal heterogeneity; panels marked P' include landscape level spatial heterogeneity. The blue line is a smoothed condition mean. PAGE 97 97 Fig. 5 2 Inflorescence per individual was averaged across repeated simulations of the same population. Each point is a different resampling of the observed population structure in each habitat. Simulations were done with (blue triangles) and without (red circles) demographic heterogeneity. Even without demographic heterogeneity, individual s varied in their sizes (more so in continuous forest) which determined their inflorescence production. Fig. 5 3 Predicted and observed variance among replicates of inflorescence production per individual are near the one to one line. Variance is similar between habitats or slightly higher in continuous forest (left panel) and slightly higher with demographic stochasticity (blue triangles) than without (red circles). PAGE 98 98 Fig. 5 4 Transition from continuous fo rest to fragment, years (represented by color): 1 to 35, was simulated by resampling from adult states observed in the continuous forest and projecting forward using demographic rates observed in fragments. Seedlings are the left most point on the x axis. The pink line represents the observed structure of the fragments in 2008 and 2009. The first panel was simulated with only demographic stochasticity; the second panel additionally incorporated heterogeneity among years, plots, and individuals. Fig. 5 5 Population structure 30 years after fragmentation occurred was simulated by resampling from adult states observed in the continuous forest and projecting forward using demographic rates observed in fragments. Grey ribbons represent the 95% quantiles among replicates. Seedlings are the left most point on the x axis. The pink line represents the observed structure of the fragments in 2008 and 2009. The left panel was simulated with only demographic stochasticity; the right panel additionally incorporated het erogeneity among years, plots, and individuals. PAGE 99 99 Table 5 1 Simulated heterogeneity situations. S ituation RE individual RE year RE plot 1 2 3 4 5 6 Table 5 2 Baseline log geometric growth rate and effects of heterogeneity in continuous forest and fragments. The effects were fit as offsets using weighted least squares regression weighted by the inverse of the estimated standard deviation. All estimates were significantly different from zero (P values of t tests were all less than .05). C F e st F e st E ffect W itho ut y ear:RE plot W ith y ear:RE plot W ithout y ear:RE plot W ith y ear:RE plot B aseline 0.132(0.002) 0.132(0.003) 0.089(0.001) 0.089(0.001) RE year 0.010(0.004) 0.010(0.004) 0 .00 4(0.001) 0 .00 4(0.001) RE individual 0.011(0.003) 0.011(0.003) 0 .00 8(0.001) 0.008(0.001) RE plot 0.029(0.006) 0 043(0.012) 0.031(0.002) 0 .011 (0.005) Y ear:RE plot 0 .00 4(0.000) 0 .001 (0.000) PAGE 100 100 CHAPTER 6 CONCLU DING REMARKS This work has made contributions to the fields of sex change evolution, statistics, stochastic demography, and fragmentation ecology. The life history strategy of changing sex depending on size can be evolutionarily advantageous even when reproduction does not increase with size. Previous theories were based on the assumption that reproductive rate increases with size even though this is not always observed in nature. We have created a method that allows researchers to detect among individual variation in growth even when they are unable to mark individuals. Marking individuals in ethically question able and demanding of resources. Our work increases the variety of systems where demographic heterogeneity can be detected. We have shown that heterogeneity can increase the stochastic growth rate by increasing the average growth rate. PAGE 101 101 APPENDIX A DYNAMIC PROGRAMMING We consider a short time interval of length t and separate the integral in Eq. ( 2 6) into the term within this interval and the rest of time (from t + t to infinity): V ( L ( t )) = max f f ( t ) + mR m ( t ) ( ) # t + max f f ( x ) + mR m ( x ) ( ) e $ + z 0 L ( t ) % & ( ) # t e $ + z 0 L ( y ) % & ( ) dy t + # t x + t + # t + dx / 0 / 1 2 / 3 / / / 0 / / 1 2 / / 3 / / (A 1) where the first symbol of maximization is related to the choice of the reproductive activity during the short time interval t and the second symbol for maximization is the optimal choice of reprod uctive activity after t + t Eq. (A 1) can be rewritten as = max f ( t ), m ( t ) f f ( t ) + mR m ( t ) ( ) # t + e $ + z 0 L ( t ) % & ( ) # t max f f ( x ) + mR m ( x ) ( ) e $ + z 0 L ( y ) % & ( ) dy t + # t x + t + # t + dx / 0 / 1 2 / 3 / / / 0 / / 1 2 / / 3 / / Now, we note that the integral in the left hand side is expressed in terms of Eq. ( 2 6). It is rewritten as = max f ( t ), m ( t ) f f ( t ) + mR m ( t ) ( ) # t + e $ + z 0 L ( t ) % & ( ) # t V L t + # t ( ) ( ) + / 0 1 = max f ( t ), m ( t ) f f ( t ) + mR m ( t ) ( ) # t + e $ + z 0 L ( t ) % & ( ) # t V L t ( ) ( ) + dV dL dL dt L t ( ) f t ( ) m t ( ) ( ) # t % & ( ) + / 0 1 where the Taylor expansion of V L t + t ( ) ( ) is calculated and expressed in terms of the reproductive activities that affect growth during the short time interval t We adopted: L t + t ( ) = L t ( ) + dL dt L t ( ) # f t ( ) # m t ( ) ( ) t PAGE 102 102 Now we note e + z 0 L ( t ) # $ % & ( ) t 1 + z 0 L ( t ) # $ % & ( ) t and expand the right hand side according to the order of t After removing V L t ( ) ( ) from both sides, we divide both sides by t and ignore terms of orders smaller than t : 0 = max f ( t ), m ( t ) f f ( t ) + mR m ( t ) ( ) # + z 0 L t ( ) $ % & ( ) V L t ( ) ( ) + dV dL g L t ( ) f t ( ) m t ( ) ( ) + / 0 (A 2) Now we specify how growth rate is reduced by reproductive activity as Eq. ( 2 4), where growth rate during vegetative growth ( m = f = 0 ) is a linear function of the current size: dL dt L t ( ) 0 0 ( ) = a + bL ( t ) When the pla nt is reproductive, growth rate is reduced. Using Eq. ( 2 4), Eq. (A 2) is rewritten as Eq. ( 2 7). PAGE 103 103 APPENDIX B MAXIMIZTION Maximization in Eq. ( 2 7) under the linear constraint Eq. ( 2 1) is achieved by comparing three quantities: f d dV dL ( ) mR c dV dL ( ) and 0. If f d dV dL ( ) is the largest among the three, the optimal strategy is to be a female ( m f ( ) = 0 1 ( ) ). If mR c dV dL ( ) is the largest, the optimal strategy is to be a male ( m f ( ) = 1 0 ( ) ). If both quantities are negative, the optimal strategy is to stay immature ( m f ( ) = 0 0 ( ) ). These implications arise directly from the maximization of Eq. ( 2 7) and require no assumption about the order of the phases. We will not assume anythin g about the order of the phases until Appendix C and D. We have three intervals as follows: [I] Vegetative growth phase : In this interval, the optimal choice satisfies: f t ( ) = 0 and m t ( ) = 0 (B 1a) which implies (fr om Eq. 2 7) f d dV dL < 0 and mR c dV dL < 0 (B 1b) and 0 = a + bL ( ) dV dL + z 0 L # $ % & ( V L ( ) (B 1c) holds In this interval reproductive activity as a female or a male would reduce the growth rate too much, and it is less beneficial to the plant. [II] Male reproductive phase: In this interval, the optimal choice satisfies: f t ( ) = 0 and m t ( ) = 1 (B 2a) which implies (from Eq. 7) f d dV dL < mR c dV dL and mR c dV dL > 0 (B 2b) PAGE 104 104 and 0 = mR + a + bL ( ) c ( ) dV dL + z 0 L # $ % & ( V L ( ) (B 2c) holds. In this interval, reproductively active as a male is more valuable than being vegetative growth or being reproductively active as a female. [III] Female reproductive phase: In this interval, the optimal choice satisfies: f t ( ) = 1 and m t ( ) = 0 (B 3a) which implies (from Eq. ( 2 7)) f d dV dL > mR c dV dL and f d dV dL > 0 (B 3b) and 0 = f + a + bL ( ) d ( ) dV dL + z 0 L # $ % & ( V L ( ) (B 3c) holds. In this interval, reproductively active as a female is more valuable than vegetative growth or being reproductively active as a male. PAGE 105 105 APPEND IX C SURVIVORSHIP Now we consider the survivorship of an individual to a given size L The mortality per unit time is z 0 L irrespective of the reproductive activity. However growth rate depends on reproductive activity so mortality per unit increase in size differs depending on reproductive state. Specifically, mortality per unit increase in size would be lowest for a vegetatively growing plant, and higher for a male plant an d the highest for a female plant. Considering these, we can derive the following results for survivorship. Denote the survivorship of an individual from size L 0 to L by S L 0 ; L ( ) For L 0 < L < L 1 we have S L 0 ; L ( ) = exp + z 0 L ( ) dL a + bL L 0 L # $ % & & ( ) ) = a + bL 0 a + bL + / b z 0 a L 0 L + / z 0 a (C 1) For L 1 < L < L 2 we have S L 0 ; L ( ) = S L 0 ; L 1 ( ) exp + z 0 L ( ) dL a + bL c L 1 L # $ % & & ( ) ) = a + bL 0 a + bL 1 + / b z 0 a L 0 L 1 + / z 0 a a + bL 1 c a + bL c + / b z 0 a c L 1 L + / z 0 a c (C 2) For L > L 2 we have S L 0 ; L ( ) = S L 0 ; L 2 ( ) exp + z 0 L ( ) dL a + bL d L 2 L # $ % & & ( ) ) = a + bL 0 a + bL 1 + / b z 0 a L 0 L 1 + / z 0 a a + bL 1 c a + bL 2 c + / b z 0 a c L 1 L 2 + / z 0 a c a + bL 2 d a + bL d + / b z 0 a d L 2 L + / z 0 a d (C 3) Using these expressions, we can rewrite Eq. ( 2 13) as Eq. ( 2 14). PAGE 106 106 Note that when we change the variable of integration from age t to size L we must include the inverse of the growth rate because dt = dL a + bL d ( ) in the female phase and dt = dL a + bL c ( ) in the male phase. The right hand side of Eq. ( 2 14) gi ves the expected lifetime RS of the schedule. It must be equal to V 0 as shown in Eq. ( 2 14). However for an arbitrary choice of V 0 Eq. ( 2 14) does not hold and we need to find a suitable choice of V 0 that satisfies Eq. ( 2 14). The procedure is described in the text. PAGE 107 107 APPENDIX D DIRECT OPTIMIZATION Here we develop the direct calculation of the optimal condition. We express the lifetime RS as a function of two critical sizes ( L 1 and L 2 ), then calculate the derivatives with respect to them and set them equal to 0. Finally, we combine these conditions with the ESS sex ratio condition: total male RS is equal to total female RS for the population. Let S I L 0 L ( ) be the survivorship from size L 0 to size L, during which the plant grows vegetatively Let S II L 1 L ( ) be the survivorship from size L 1 to size L, during which the plant grows vegetat ively Let S III L 2 L ( ) be the survivorship from size L 2 to size L, during which the plant grows vegetatively We have the following expressions: S I L 0 L ( ) = exp + z 0 L # $ % & ( dL a + bL L 0 L ) + / / for L 0 < L < L 1 (D 1a) S II L 1 L ( ) = exp + z 0 L # $ % & ( dL a + bL c L 1 L ) + / / for L 1 < L < L 2 (D 1b) S III L 2 L ( ) = exp + z 0 L # $ % & ( dL a + bL d L 2 L ) + / / for L 2 < L (D 1c) Hence we have L 0 S I L 0 L ( ) = + z 0 L 0 # $ % & ( 1 a + bL 0 S I L 0 L ( ) (D 2a) L S I L 0 L ( ) = # + z 0 L $ % & ( ) 1 a + bL S I L 0 L ( ) (D 2b) PAGE 108 108 L 1 S II L 1 L ( ) = + z 0 L 1 # $ % & ( 1 a + bL 1 ) c S II L 1 L ( ) (D 2c) L S II L 1 L ( ) = # + z 0 L $ % & ( ) 1 a + bL # c S II L 1 L ( ) (D 2d) L 2 S III L 2 L ( ) = + z 0 L 2 # $ % & ( 1 a + bL 2 ) d S III L 2 L ( ) (D 2e) L S III L 2 L ( ) = # + z 0 L $ % & ( ) 1 a + bL # d S III L 2 L ( ) (D 2f) Now we define the lifetime RS as L 1 L 2 ( ) = S I L 0 L 1 ( ) mR S II L 1 L ( ) dL a + bL # c L 1 L 2 $ % & + fS II L 1 L 2 ( ) S III L 2 L ( ) dL a + bL # d L 2 ( $ ) + + (D 3) The derivatives with respect to two critical sizes are: "# L 1 = $ + z 0 L 1 % & ( ) 1 a + bL 1 # + + z 0 L 1 % & ( ) 1 a + bL 1 $ c # + mRS I L 0 L 1 ( ) $ 1 ( ) 1 a + bL 1 $ c = + z 0 L 1 % & ( ) 1 a + bL 1 $ c $ 1 a + bL 1 % & ( ) # $ mR S I L 0 L 1 ( ) a + bL 1 $ c "# L 2 = S I L 0 L 1 ( ) mRS II L 1 L 2 ( ) 1 a + bL 2 $ c $ fS II L 1 L 2 ( ) 1 a + bL 2 $ d % & $ f + z 0 L 2 ( ) + 1 a + bL 2 $ c S II L 1 L 2 ( ) S III L 2 L ( ) dL a + bL $ d L 2 / + f + z 0 L 2 ( ) + 1 a + bL 2 $ d S II L 1 L 2 ( ) S III L 2 L ( ) dL a + bL $ d L 2 / 0 1 2 2 Setting these equal to 0, we have PAGE 109 109 0 = + z 0 L 1 # $ % & c a + bL 1 ( ) mRS I L 0 L 1 ( ) (D 4a) 0 = f a + bL 2 c ( ) mR a + bL 2 ^ d ( ) ( ) + f d c ( ) + z 0 L 2 # $ % & ( S III L 2 L ( ) dL a + bL d L 2 ) (D 4b) In addition to these optimization conditions, we have the condition for the equality of male/female contributions. mR S II L 1 L ( ) dL a + bL c L 1 L 2 # = fS II L 1 L 2 ( ) S III L 2 L ( ) dL a + bL d L 2 $ # (D 5) From these four relatio ns (D 1), (D 4a), (D 4b) and (D 5), we can determine three unknown variables: L 1 L 2 R, and Define the following quantities: A L 1 L 2 ( ) = S II L 1 L ( ) dL a + bL c L 1 L 2 # (D 6a) B L 2 ( ) = S III L 2 L ( ) dL a + bL d L 2 # $ (D 6b) Then = S I L 0 L 1 ( ) mRA [ + fS II L 1 L 2 ( ) B ] (D 7) + z 0 L 1 # $ % & c a + bL 1 ( = mRS I L 0 L 1 ( ) (D 8) f a + bL 2 c ( ) mR a + bL 2 d ( ) ( ) = f d c ( ) + z 0 L 2 # $ % & ( B (D 9) mRA = fS II L 1 L 2 ( ) B (D 10) We remove and use s R = mR f We can derive the following three equations: PAGE 110 110 + z 0 L 1 # $ % & c a + bL 1 s R A [ + S II L 1 L 2 ( ) B ] ( s R = 0 (D 11) a + bL 2 c ( ) s R a + bL 2 d ( ) ( ) d c ( ) + z 0 L 2 # $ % & ( B = 0 (D 12) s R A S II L 1 L 2 ( ) B = 0 (D 13) The three unknowns are two critical sizes ( L 1 L 2 ) and sex ratio s R = mR f First we set sex ratio s R = mR f as a trial val ue, and find the root of Eq. (D 12) to obtain L 2 Then we find the root Eq. (D 11) to obtain L 1 T hen examine the value of Eq. (D 13), and change s R Searching for the solution to all three unknowns was somewhat automated using the uniroot function in the statistical programming language R. PAGE 111 111 APPENDIX E DERRIVATION OF STATI STIC AL GROWTH MODEL Assume individual i 's size through time follows: S i ( t + t ) = S i ( t ) + t g + t # i + # i t where i has mean 0 and among individual variance 2 # g 2 and i t has mean 0 and variance through time 1 # 2 ( ) $ g 2 % t Let denote the expected value across individuals. 2 t + t ( ) = S i t + t ( ) ( ) 2 S i t + t ( ) 2 = S i t ( ) + g t + # i t + # i t ( ) 2 S i t ( ) + g t ( ) 2 = S i t ( ) 2 + 2 S i t ( ) g t + 2 S i t ( ) # i t + g 2 t 2 + # i t ( ) 2 + # i t ( ) 2 S i t ( ) 2 2 S i t ( ) g t g 2 t 2 = S i t ( ) 2 S i t ( ) 2 + 2 t S i t ( ) # i + # i t ( ) 2 + # i t ( ) 2 = S i t ( ) 2 S i t ( ) 2 + 2 t S i 0 ( ) + g t + # i t + # i $ $ = 0 t # $ % & ( # i + # i t ( ) 2 + # i t ( ) 2 = S i t ( ) 2 S i t ( ) 2 + 2 t # i 2 t + # i t ( ) 2 + # i t ( ) 2 = 2 t ( ) + 2 t % 2 g 2 t + t 2 % 2 g 2 + ( 1 % 2 ) g 2 t PAGE 112 112 APPENDIX F PARAMETER ESTIMATES Table F 1 Parameters used for both habitats. Parameter Source Value Mean flowers per inflorescence Bruna unpublished data 22.625 Flowers overdisp par (NB size) Bruna unpublished data 7.56 Logit prob flower producing a fruit Bruna unpublished data 0.9677 Individual stdev logit prob flower to fruit Bruna unpublished data 0.86292 Prob 0 seeds per fruit Bruna unpublished data 0.00458 2 Prob 1 seed per fruit Bruna unpublished data 0.365407 Prob 2 seeds per fruit Bruna unpublished data 0.355097 Prob 3 seeds per fruit Bruna unpublished data 0.274914 Prob seedling has 1 shoot next year 12 year Bruna dataset 0.79969 3 Prob seedling has 2 shoots next year 12 year Bruna dataset 0.16034 9 Prob seedling has 3 shoots next year 12 year Bruna dataset 0.02920 1 Prob seedling has 4 shoots next year 12 year Bruna dataset 0.010758 Table F 2. Habitat specific parameter values. Parameter Source CF value F value Prob seed germination Bruna 2002 Oecologia 0.055 0.02 Logit Prob seedling survival 12 year Bruna dataset 0.9664731 1.618192 Plot stdev logit seedling survival 12 year Bruna dataset 0.5469587 0.2034674 Year stdev logit seedling survival 12 year Bruna dataset 0.28094 0.30733 PAGE 113 113 APPENDIX G ZERO INFLATED POISSON DEM OGRAPHIC VARIANCE Here we show the calculations to determine the expected variation in the population's inflorescence rate. v a r ( F ) = N [ V D ( m ) ] + v a r ( m ) ( ) [ V D ( m ) ] = V D ( p ) f ( p ) f ( ) # p # V D ( p ) $ V D ( p ) + ( p % p ) # V D # p + ( % ) # V D # & ( ) + + 1 2 ( p % p ) 2 # 2 V D # p 2 + 2 ( p % p ) ( % ) # 2 V D # p # + ( % ) 2 # 2 V D # 2 & ( ) + + 1 3 ( p % p ) 3 # 3 V D # p 3 + 3 ( p % p ) 2 ( % ) # 3 V D # p 2 # + 3 ( p % p ) ( % ) 2 # 3 V D # p # 2 + ( % ) 3 # 3 V D # 3 & ( ) + Then by substituting the above Taylor expansion of V D (p, ) into the formula for E[V D (m)] and integrating, we get [ V D ( m ) ] V D ( p ) + 1 2 v a r ( p ) # 2 V D # p 2 + 2 c o v ( p ) # 2 V D # p # + v a r ( ) # 2 V D # 2 $ % & ( ) + 1 3 s k e w ( p ) # 3 V D # p 3 + s k e w ( ) # 3 V D # 3 + 3 ( p p ) 2 ( ) # 3 V D # p 2 # + ( p p ) ( ) 2 # 3 V D # p # 2 + / 0 f ( p  ) f ( ) # p # 1 1 $ % & & & & ( ) ) ) ) Then more specifically for the zero inflated Poisson: PAGE 114 114 V D ( p ) = ( 1 p ) + p ( 1 p ) 2 V D p = + ( 1 2 p ) 2 2 V D p 2 = 2 2 V D = ( 1 p ) + 2 p ( 1 p ) 2 V D 2 = 2 p ( 1 p ) 2 V D p = 1 + 2 ( 1 2 p ) 3 V D 2 p = 2 ( 1 2 p ) 3 V D p 2 = 4 3 V D 3 = 0 3 V D p 3 = 0 PAGE 115 115 APPENDIX H DISTRIBUTION OF DEVIATES Let x = r r ( ) 2 r 2 ~ d f = 1 2 so the probability density of x is 1 2 1 / 2 ( 1 / 2 ) x 1 / 2 1 e x / 2 = 1 2 1 / 2 ( 1 / 2 ) x 1 / 2 e x / 2 And because it is a probability density, then by definition 1 2 1 / 2 ( 1 / 2 ) x 1 / 2 e x / 2 d x = 1 Let y = r r ( ) 2 = x r 2 Then d y = d x r 2 1 2 1 / 2 ( 1 / 2 ) y r 2 # $ % & 1 / 2 e y 2 r 2 # $ $ % & d y 0 ( ) = 1 2 1 / 2 ( 1 / 2 ) x ( ) 1 / 2 e x / 2 0 ( ) r 2 d x = r 2 Therefore, to get the probability density of y to integrate to 1, we have to divide by r 2 The probability density of y is 1 r 2 1 2 1 / 2 ( 1 / 2 ) x 1 / 2 e x / 2 = 1 2 1 / 2 ( 1 / 2 ) r 2 y r 2 # $ % & 1 / 2 e y 2 r 2 # $ $ % & In our case, r 2 is a linear function of covariates related to sources of heterogeneity and the number of individuals in the population. Define model matrix M coefficient vector and error vector such that r 2 = + # Then we can calculate the negative log likelihood of each observed value of y For observation i the negative log likelihood is l o g 1 2 1 / 2 ( 1 / 2 ) # [ ] i y i # [ ] i $ % & & ( ) ) 1 / 2 e y 2 # [ ] i $ % & & ( ) ) $ % & & & & & & ( ) ) ) ) ) ) = 5 l o g ( 2 ) + l o g ( ( 5 ) ) + 5 l o g ( # [ ] i ) + 5 l o g ( y i ) + y i 2 # [ ] i where [ ] i is the i th element of the vector formed by multiplying the model matrix by the coefficient vector. We use this formula to fit the model r 2 = + # to our observed values of y = r r ( ) 2 by minimizing the negative log likelihood in ADMB We assume that the errors are i.i.d. We had to restrict elements of to be greater than 1e 6 to avoid PAGE 116 116 taking the log of a negative number. We penalized the negative log likelihood if the optimizer went into regions of parameter space that would cause negative values. Fig. H 1 S quared values of log growth rate minus mean log growth rate. The mean was taken for each situation separately. Each column of panels represents a different situation. All situations include demographic stochasticity; panels marked I' include demographic heterogeneity; panels marked Y' include temporal heterogeneity; panels marked P' include landscape level spatial heterogeneity. PAGE 117 117 LIST OF REFERENCES Alford, R.A., Harris, R.N., 1988. Effects of larval growth history on anuran metamorphosis. Am Nat 91 106. Allsop, D., West, S., 2004. Sex ratio evolution in sex changing animals. Evolution 1019 1027. Altwegg, R., Reyer, H. U., 2003. Patterns of natural selection on size at metamorphosis in water frogs. Evolution 57, 872 882. Baayen, R.H., Davidson, D.J., Bates, D.M., 2008. Mixed effects modeling with crossed random effects for subjects and items. Journal of Memory and Language 59, 390 412. Bates, D., Maechler, M., Bolker, B.M., 2011. lme4: Linear mixe d effects models using S4 classes. Bellman, R., 1957. The theory of games. Benedetti Cecchi, L., 2003. The importance of the variance around the mean effect size of ecological processes. Ecology 84, 2335 2346. Bierzychudek, P., 1984. Determinants of gen der in jack in the pulpit: the influence of plant size and reproductive history. Oecologia 65, 14 18. Bierzychudek, P., Eckhart, V., 1988. Spatial segregation of the sexes of dioecious plants. Am Nat 132, 34. Bolker, B., 2003. Combining endogenous and ex ogenous spatial variability in analytical population models. Theoretical Population Biology 64, 255 270. Bolker, B.M., 2008. Ecological models and data in R. Princeton Univ Pr. Bolker, B.M., 2012. bbmle: Tools for general maximum likelihood estimation. Bolker, B.M., Brooks, M.E., Clark, C.J., Geange, S.W., Poulsen, J.R., Stevens, M.H.H., White, J. S.S., 2009. Generalized linear mixed models: a practical guide for ecology and evolution. Trends in Ecology & Evolution 24, 127 135. Bolker, B.M., Skaug, H., 2012. R2admb: ADMB to R interface functions. Bolnick, D.I., Amarasekare, P., Arajo, M.S., Brger, R., Levine, J.M., Novak, M., Rudolf, V.H.W., Schreiber, S.J., Urban, M.C., Vasseur, D.A., 2011. Why intraspecific trait variation matters in community ecology. Trends in Ecology & Evolution 26, 183 192. PAGE 118 118 Boyce, M.S., 1977. Population growth with stochastic fluctuations in the life table. Theoretical Population Biology 12, 366 373. Briant, G., Gond, V., Laurance, S.G.W., 2010. Habitat fragmentation and the desiccation of forest canopies: a case study from eastern Amazonia. Biological Conservation 143, 2763 2769. Brienen, R., Zuidema, P.A., During, H., 2006. Autocorrelated growth of tropical forest t rees: unraveling patterns and quantifying consequences. Forest Ecology and Management 237, 179 190. Brown, J.H., Gilliam, J.F., Allen, A.P., Savage, V.M., West, G.B., 2004. Toward a Metabolic Theory of Ecology. Ecology 85, 1771 1789. Bruna, E., 2002. Eff ects of forest fragmentation on Heliconia acuminata seedling recruitment in central Amazonia. Oecologia 132, 235 243. Buckley, Y.M., Briese, D.T., Rees, M., 2003. Demography and management of the invasive plant species Hypericum perforatum II. Constructi on and use of an individual based model to predict population dynamics and the effects of management strategies. J Appl Ecol 40, 494 507. Burd, M., 1995. Ovule packaging in stochastic pollination and fertilization environments. Evolution 100 109. Charnov E., Bull, J., Smith, J., 1976. Why be an hermaphrodite? 263, 125 126. Charnov, E.L., 1982. The theory of sex allocation. Monogr Popul Biol 18, 1 355. Charnov, E.L., 1989. Evolution of the breeding sex ratio under partial sex change. Evolution 43, 1559 1561. Clark, J., 2003. Uncertainty and variability in demography and population growth: a hierarchical approach. Ecology 84, 1370 1381. Cochrane, M.A., Laurance, W.F., 2008. Synergisms among fire, land use, and climate change in the Amazon. Ambio 37, 522 527. Coleman, K., Wilson, D., 1998. Shyness and boldness in pumpkinseed sunfish: individual differences are context specific. Animal Behaviour 56, 927 936. Conner, M.M., White, G.C., 1999. Effects of individual heterogeneity in estimating the persistence of small populations. Natural Resource Modeling 12, 109 127. PAGE 119 119 Davies, K., Melbourne, B., Margules, C., 2001. Effects of within and between patch processes on community dynamics in a fragmentation experiment. Ecology 82, 1830 1846. Day, T., Aarssen, L., 1997. A time commitment hypothesis for size dependent gender allocation. Evolution 51, 988 993. de Valpine, P., 2009. Stochastic development in biologically structured population models. Ecology 90, 2889 2901. Dennis, B., Munholland, P., Sco tt, J., 1991. Estimation of growth and extinction parameters for endangered species. Ecol Monogr 61, 115 143. Didham, R.K., Lawton, J.H., 1999. Edge Structure Determines the Magnitude of Changes in Microclimate and Vegetation Structu re in Tropical Forest Fragments Biotropica 31, 17 30. Evans, M.E.K., Holsinger, K.E., 2012. Estimating covariation between vital rates: A simulation study of connected vs. separate generalized linear mixed models (GLMMs). Theoretical Population Biology 1 8. Ferreira, L.V., L aurance, W.F., 1997. Effects of Forest Fragmentation on Mortality and Damage of Selected Trees in Central Amazonia. Conserv Biol 11, 797 801. Filin, I., Ovadia, O., 2007. Individual Size Variation and Population Stability in a Seasonal Environment: A Disc rete Time Model and Its Calibration Using Grasshoppers. The American Naturalist 170, 719 733. Fournier, D., Skaug, H., Ancheta, J., Ianelli, J.N., Magnusson, A., Maunder, M., Nielsen, A., Sibert, J., 2011. AD Model Builder: using automatic differentiation for statistical inference of highly parameterized complex nonlinear models. Optimization Methods and Software 1 17. Fox, G., 2005. Extinction risk of heterogeneous populations. Ecology 86, 1191 1198. Fox, G.A., Kendall, B.E., 2002. Demographic stochasticity and the variance reduction effect. Ecology 83, 1928 1934. Fox, G.A., Kendall, B.E., Fitzpatrick, J.W., Woolfenden, G.E., 2006. Consequences of heterogeneity in survival probability in a population of Florida scrub jays. Journal of Animal Eco logy, Heterogeneity in survival of Florida scrub jays 75, 921 927. Frank, S.A., Swingland, I.R., 1988. Sex ratio under conditional sex expression. Journal of Theoretical Biology 135, 415 418. PAGE 120 120 Gadgil, M., Bossert, W., 1970. Life historical consequences o f natural selection. The American Naturalist 104, 1. Gagnon, P.R., Bruna, E.M., Rubim, P., Darrigo, M.R., Littell, R.C., Uriarte, M., Kress, W.J., 2011. Growth of an understory herb is chronically reduced in Amazonian forest fragments. Biological Conserva tion 144, 830 835. Heckel, C., Bourg, N., McShea, W., Kalisz, S., 2010. Nonconsumptive effects of a generalist ungulate herbivore drive decline of unpalatable forest herbs. Ecology 91, 319 326. Inouye, B.D., 2005. The importance of the variance around th e mean effect size of ecological processes: comment. Ecology 86, 262 265. Iwasa, Y., 1991. Sex change evolution and cost of reproduction. Behavioral Ecology 2, 56 68. Iwasa, Y., Levin, S., Iwasa, Y., 2000. Dynamic optimization of plant growth. Evolutiona ry Ecology Research 2, 437 455. Jongejans, E., De Kroon, H., 2005. Space versus time variation in the population dynamics of three co occurring perennial herbs. J Ecology 93, 681 692. Jongejans, E., Jorritsma Wienk, L.D., Becker, U., Dostl, P., Mildn, M., De Kroon, H., 2010. Region versus site variation in the population dynamics of three short lived perennials. J Ecology 98, 279 289. Kendall, B., Fox, G., 2003. Unstructured individual variation and demographic stochasticity. Conserv Biol 17, 1170 1172. Kendall, B.E., Fox, G.A., Fujiwara, M., Nogeire, T.M., 2011. Demographic heterogeneity, cohort selection, and population growth. Ecology 92, 1985 1993. Kendall, D.G., 1949. Stochastic processes and population growth. Journal of the Royal Statistical Society. Series B (Methodological) 11, 230 282. Knight, T.M., Steets, J.A., Vamosi, J.C., Mazer, S.J., Burd, M., Campbell, D.R., Dudash, M.R., Johnston, M.O., Mitchell, R.J., Ashman, T. L., 2005. Pollen Limitation of Plant Reprodu ction: Pattern and Process. Annu. Rev. Ecol. Evol. Syst 36, 467 497. Kuwamura, T., Nakashima, Y., 1998. New aspects of sex change among reef fishes: recent studies in Japan. Environmental Biology of Fishes 52, 125 135. Laird, A., Tyler, S., Barton, A., 1 965. Dynamics of Normal Growth. Growth 29, 233. PAGE 121 121 Laurance, W., 2002. Hyperdynamism in fragmented habitats. Journal of Vegetation Science 13, 595 602. Laurance, W., Camargo, J., Luizo, R., Laurance, S., Pimm, S., Bruna, E., Stouffer, P., Bruce Williamson G., Bentez Malvido, J., Vasconcelos, H., 2011. The fate of Amazonian forest fragments: A 32 year investigation. Biological conservation 144, 56 67. Laurance, W.F., 1997. Biomass Collapse in Amazonian Forest Fragments. Science 278, 1117 1118. Lavine, M ., Beckage, B., Clark, J., 2002. Statistical modeling of seedling mortality. Journal of agricultural, biological, and environmental statistics 7, 21 41. Leidner, A.K., Haddad, N.M., Lovejoy, T.E., 2010. Does Tropical Forest Fragmentation Increase Long Ter m Variability of Butterfly Communities? PLoS ONE 5, e9534. Lloyd Smith, J.O., Schreiber, S.J., Kopp, P.E., Getz, W.M., 2005. Superspreading and the effect of individual variation on disease emergence. Nature 438, 355 359. Loya, Y., Sakai, K., 2008. Bidir ectional sex change in mushroom stony corals. Proceedings of the Royal Society B: Biological Sciences 275, 2335 2343. Malcolm, J.R., 1998. A Model of Conductive Heat Flow in Forest Edges and Fragmented Landscapes. Climatic Change 39, 487 502. Mangel, M., Stamps, J., 2001. Trade offs between growth and mortality and the maintenance of individual variation in growth. Evolutionary Ecology Research 3, 583 593. Marschner, I.C., 2010. Stable Computation of Maximum Likelihood Estimates in Identity Li nk Poisson Regression. Journal of Computational and Graphical Statistics 19, 666 683. Marshall, D.J., Keough, M.J., 2005. Offspring size effects in the marine environment: A field test for a colonial invertebrate. Austral Ecology 30, 275 280. McCarthy, M ., Parris, K., 2004. Clarifying the effect of toe clipping on frogs with Bayesian statistics. J Appl Ecol 41, 780 786. McNamara, J., Welham, R., Houston, A., Daan, S., Tinbergen, J., 2004. The effects of background mortality on optimal reproduction in a s easonal environment. Theoretical Population Biology 65, 361 372. Melbourne, B.A., Hastings, A., 2008. Extinction risk depends strongly on factors contributing to stochasticity. Nature 454, 100 103. PAGE 122 122 Morris, W.F., Doak, D.F., 2002. Quantitative conservatio n biology, theory and practice of population viability analysis. Sinauer Associates Inc. Munday, P.L., 2002. Bi directional sex change: testing the growth rate advantage model. Behavioral Ecology and Sociobiology 52, 247 254. Mylius, S., Diekmann, O., 19 95. On evolutionarily stable life histories, optimization and the need to be specific about density dependence. Oikos 74, 218 224. Nishizawa, T., Watano, Y., Kinoshita, E., Kawahara, T., Ueda, K., 2005. Pollen movement in a natural population of Arisaema serratum (Araceae ), a plant with a pitfall trap flower pollination system 1. Am J Bot 92, 1114 1123. Obeso, J.R., 2002. The costs of reproduction in plants. New Phytologist 155, 321 348. Peacor, S., Schiesari, L., Werner, E., 2007. Mechanisms of nonlethal predator effect on cohort size variation: ecological and evolutionary implications. Ecology 88, 1536 1547. Peacor, S.D., Pfister, C.A., 2006. Experimental and model analyses of the effects of competition on individual size variati on in wood frog ( Rana sylvatica ) tadpoles. J Anim Ecol 75, 990 999. Persson, L., Leonardsson, K., de Roos, A.M., Gyllenberg, M., Christensen, B., 1998. Ontogenetic scaling of foraging rates and the dynamics of a size structured consumer resource model. Th eoretical Population Biology 54, 270 293. Pfister, C.A., Stevens, F.R., 2002. The genesis of size variability in plants and animals. Ecology 83, 59 72. Pfister, C.A., Wang, M., 2005. Beyond size: matrix projection models for populations where size is an incomplete descriptor. Ecology 86, 2673 2683. Plummer, M., Best, N., Cowles, K., Vines, K., 2006. CODA: Convergence Diagnosis and Output Analysis for MCMC. R News 6, 7 11. Policansky, D., 1981. Sex choice and the size advantage model in jack in the pulpi t ( Arisaema triphyllum ). Proceedings of the National Academy of Sciences 78, 1306 1308. Ramula, S., Dinntz, P., Lehtil, K., 2009. Spatial data replacing temporal data in population viability analyses: An empirical investigation for plants. Basic and App lied Ecology 10, 401 410. PAGE 123 123 Ricklefs R., 1967. A Graphical Method of Fitting Equations to Growth Curves. Ecology 48, 978 983. Ricklefs R., Cullen, J., 1973. Embryonic Growth of Green Iguana Iguana Iguana. Copeia 296 305. Roughgarden, J., 1975. A simple model for population dynamics in stochastic environments. Am Nat 713 736. Rowe, L., Ludwig, D., 1991. Size and Timing of Metamorphosis in Complex Life Cycles: Time Constraints and Variation. Ecology 72, 413. Ruel, J., Ayres, M., 1999. Jensen's inequality predicts effects of environmental variation. Trends in Ecology & Evolution 14, 361 366. Saraux, C., Bohec, C.L., Durant, J.M., Viblanc, V.A., Gauthier Clerc, M., Beaune, D., Park, Y. H., Yoccoz, N.G., Stenseth, N.C., Maho, Y.L., 2011. Reliability of flipper banded penguins as indicators of climate change. Nature 469, 203 206. Saunders, D., Hobbs, R., Margules, C., 1991. Biological Consequences Of Ecosystem Fragmentation A Review. Co nserv Biol 5, 18 32. Scheipl, F., Greven, S., Kuechenhoff, H., 2008. Size and power of tests for a zero random effect variance or polynomial regression in additive and linear mixed models. Computational statistics & data analysis 52, 3283 3299. Setzer, R .W., 2011. odesolve: Solvers for Ordinary Differential Equations. Sih, A., Bell, A., Johnson, J., 2004. Behavioral syndromes: an ecological and evolutionary overview. Trends in Ecology & Evolution 19, 372 378. Skaug, H., Fournier, D., 2006. Automatic app roximation of the marginal likelihood in non Gaussian hierarchical models. Computational statistics & data analysis 51, 699 709. Skaug, H., Fournier, D., Nielsen, A., Magnusson, A., Bolker, B.M., 2012. glmmADMB: Generalized Linear Mixed Models using AD Mo del Builder. Tuljapurkar, S., Orzack, S., 1980. Population dynamics in variable environments I. Long run growth rates and extinction. Theoretical Population Biology 18, 314 342. Uriarte, M., Anciaes, M., da Silva, M.T.B., Rubim, P., Johnson, E., Bruna, E .M., 2011. Disentangling the drivers of reduced long distance seed dispersal by birds in an experimentally fragmented landscape. Ecology 92, 924 937. PAGE 124 124 Warner, R.R., Robertson, D.R., Leigh, E.G., 1975. Sex change and sexual selection. Science 190, 633 638. Wearing, H., Rohani, P., Cameron, T., Sait, S., 2004. The dynamical consequences of developmental variability and demographic stochasticity for host parasitoid interactions. Am Nat 164, 543 558. Weiner, J., 1985. Size hierarchies in experimental populat ions of annual plants. Ecology 743 752. Weitz, J.S., Levin, S.A., 2006. Size and scaling of predator prey dynamics. Ecology Letters 9, 548 557. Werner, E., Gilliam, J.F., 1984. The ontogenetic niche and species interactions in size structured populations Annual review of ecology and systematics 15, 393 425. Wickham, H., 2009. Ggplot2. Wilbur, H.M., Collins, J.P., 1973. Ecological Aspects of Amphibian Metamorphosis: Nonnormal distributions of competitive ability reflect selection for facultative metamor phosis. Science 182, 1305 1314. Williams, G.C., 1966. Natural selection, the costs of reproduction, and a refinement of Lack's principle. Am Nat 100, 687 690. Wilson, P., Thomson, J.D., Stanton, M.L., Rigney, L.P., 1994. Beyond floral Batemania: gender b iases in selection for pollination success. Am Nat 283 296. Wright, S.I., Barrett, S.C.H., 1999. Size dependent gender modification in a hermaphroditic perennial herb. Proceedings of the Royal Society B: Biological Sciences 266, 225 232. Yund, P., 2000. How severe is sperm limitation in natural populations of marine free spawners? Trends in Ecology & Evolution 15, 10 13. Zeileis, A., Kleiber, C., Jackman, S., 2008. Regression models for count data in R. Journal of Statistical Software 27, 1 25. Zimmerman, J.K., 1991. Ecological correlates of labile sex expression in the orchid Catasetum viridiflavum Ecology 597 608. Zuidema, P.A., Brienen, R.J.W., During, H.J., Gneralp, B., 2009. Do Persistently Fast Growing Juveniles Contribute Disproportiona tely to Population Growth? A New Analysis Tool for Matrix Models and Its Application to Rainforest Trees. The American Naturalist 174, 709 719. PAGE 125 125 Zuidema, P.A., Franco, M., 2001. Integrating vital rate variability into perturbation analysis: an evaluation f or matrix population models of six plant species. J Ecology 89, 995 1005. PAGE 126 126 BIOGRAPHICAL SKETCH Mollie was born in Charleston, West Virginia in a good school district with excellent public school teachers. They taught her how to think critically, write design experiments, do calculus, and program in C++. It would be great if everyone could have the same opportunity After high school, Mollie went to the University of Pittsburgh. She chose Pitt because she was undecided about her major and it offered bachelor of science (B.S.) degrees in both computer science and mathematics. Mollie graduated with a B.S. in math under the guidance of William J. Layton. Dr. Layton studied numerical methods in addition to fluid dynamics. He taught my first course in math ematical modeling. He advised Mollie to stay in math and continue her studies in graduate school. Mollie started graduate school in the same math department at Pitt witho ut much motivation. Eventually she realized that she wa nted to apply math to ecology She began working as a modeler in Susan Kalisz's lab in the Biological Sciences D epartment. Dr. Kalisz's graduate student Christopher Heckel was studying sex ratios of local Jack in the Pulpit populations. They soon realized that existing models could no t explain the skewed sex ratios he observed; this was the motivation for Mollie's first dissertation chapter. Mollie started graduate school in the Zoology D epartment at the Univers ity of Florida in August 2006. She wanted to learn quantitative ecology f rom Benjamin Bolker. She did not know the specific area of research she wanted to work in, so it appealed to her that the department had many researchers working on interesting topics 