ABSTRACT Title of Dissertation: THE INFLUENCE OF ISOLATION ON THE DYNAMICS OF POPULATIONS AND COMMUNITIES Justin Marco Calabrese, Doctor of Philosophy, 2005 Dissertation directed by: Dr. William F. Fagan Department of Biology Isolation is defined as the separation in time or space of individuals, populations, or of species within a community. Though isolation can be the result of many ecological processes, its role in affecting the structure and dynamics of populations and communities is not often acknowledged directly. For example, spatial heterogeneity is a frequently recognized as a significant ecological factor, but the effects of spatial heterogeneity are manifested through the isolation that heterogeneity imposes on the focal populations or communities. Isolation is an important, but hidden, component of many other ecological theories and frameworks as well. In this dissertation, I explore the role of isolation per se as an organizing theme in ecology by studying the effects of isolation in time and in space on both populations and communities. Chapter 1 explores how isolation in time among individuals in a population may affect the population?s dynamics and risk of extinction. Through a combination of modeling and meta-analysis, Chapter 1 demonstrates that reproductive asynchrony, a form of temporal isolation, can have profound negative effects at the population level in species that feature annual lifecycles. Chapter 2 reviews and synthesizes the literature on habitat connectivity, the inverse of spatial isolation, and lays out a novel framework for organizing and understanding the different metrics used to measure the connectivity. Chapter 3 examines the role of spatial isolation among species in an assemblage of Costa Rican bark beetles in mediating species interactions. The chapter uses a combination of modeling and field-collected observational data to test the hypothesis that isolation among species in this bark beetle assemblage results in a community that behaves neutrally. The studies presented in this dissertation represent a broad sweep of the ways in which the concept of isolation may be applied to better understand the dynamics of populations and communities. Individually, each chapter is an original contribution to the ecology literature. Taken together, these papers demonstrate the power of isolation as an organizing theme in ecology and will hopefully stimulate increased research effort and theoretical development around the concept of isolation. THE INFLUENCE OF ISOLATION ON THE DYNAMICS OF POPULATIONS AND COMMUNITIES by Justin Marco Calabrese Dissertation submitted to the Faculty of the Graduate School of the University of Maryland, College Park in partial fulfillment of the requirements for the degree of Doctor of Philosophy 2005 Advisory Committee: Professor William F. Fagan, Chair Professor Howard V. Cornell Professor Robert F. Denno Professor Maile C. Neel Professor Margaret A. Palmer ? Copyright by Justin Marco Calabrese 2005 ii ACKNOWLEDGEMENTS A great many people have contributed to the work here presented and have facilitated, in one way or another, my completion of this degree. First, I thank my mother Bette, father Chuck, and brother Jesse for all their love and support, and for their unwavering belief in me. I truly could not have done this without them. Next, my girlfriend, Tatjana Schmidt, has been a continual source of love, support and inspiration. She always seems to know when to tell me everything will be alright, and when to tell me to stop complaining and get back to work. I would like to thank my advisor, Bill Fagan for his patience and guidance while I, somewhat randomly, explored many different corners of ecology before settling into my program of study. Bill has definitely led by example, and in so doing, has pushed me to learn and grow even when I didn?t think I could. As I have developed, Bill has become, and hopefully will continue to be, a great colleague and collaborator of mine. I thank my committee from Arizona State University ? John Alcock, Steve Carroll, Stan Faeth and Dave Pearson ? for their guidance in overseeing the early part of my graduate studies, and their belief that I could accomplish what I set out to do. I also wish to thank my committee from the University of Maryland ? Mary Christman, Howard Cornell, Bob Denno, Maile Neel, and Margaret Palmer ? for their support, guidance and flexibility in the latter stages of my degree. Though difficult, switching universities mid-way through my Ph.D. allowed me to benefit greatly from having so many generous and talented scientists as mentors. iii Many other colleagues and friends have contributed to my intellectual and personal development during my graduate studies. I thank my collaborator Larry Kirkendall from the University of Bergen, Norway, for his expert advice and assistance with the Cecropia beetle work, and for providing a consistent check on my theoretical exuberance. Alan and Karen Masters and Deedra McClearn have all been fantastic mentors and have each helped shape the way I think about science and life. Jesse Brunner has been a great friend and has patiently and logically served as a sounding board for my ideas, both good and bad, over the years. Jessamy Rango, Lisa Eby, Craig Aumann and Leslie Ries have all been friends and mentors along the way and have helped me navigate difficult times and difficult decisions. I thank my volunteer assistants for their tireless help, good ideas, and ability to deal with a boss who many times had no idea what he was doing. Julie Asmar was a fantastic field and lab assistant and admirably bore the brunt of my inexperience in the field. Her many practical ideas and reluctance to let me do silly and inefficient things greatly improved the quality of the research. I also thank my Mother, Bette Calabrese, for suspending her disbelief and taking a break from her normal life to come down to Costa Rica and help her son carry sticks through the jungle. Julio Contreras-Martinez provided invaluable help as an onsite lab assistant both while I was doing fieldwork and while I was back in the U.S. sitting in front of a computer. His thoroughness and attention to detail allowed me to expand the scope of my fieldwork and still attend to my other responsibilities back in the States. I also thank my girlfriend, Tatjana Schmidt, iv for her invaluable help, advice and energy in the field and lab. Our relationship not only survived, but was strengthened by two months of working together. Finally, I thank the following organizations for generously funding my research: the Achievement Rewards for College Scientists Foundation, the National Science Foundation, the Organization for Tropical Studies, Sigma Xi, the Center for Insect Science at the University of Arizona, the Cosmos Club Foundation, Arizona State University, and the BEES program at the University of Maryland. v TABLE OF CONTENTS List of Tables .....................................................................................vi List of Figures....................................................................................vii Introduction........................................................................................1 Chapter I: Lost in time, lonely and single: Reproductive asynchrony and the Allee effect...................................................6 Abstract..................................................................................7 Introduction............................................................................8 Methods..................................................................................12 Results....................................................................................19 Discussion..............................................................................24 Chapter II: A comparison-shopper?s guide to connectivity metrics.....................................................................40 Abstract..................................................................................41 Introduction............................................................................42 The Data Dependent Framework...........................................45 Modifications .........................................................................54 Scale Dependence ..................................................................55 Conclusions............................................................................57 Chapter III: Neutral theory as a null hypothesis for species diversity in aggregated arthropod assemblages: a test with Cecropia petiole beetles.......................................................64 Abstract..................................................................................65 Introduction............................................................................67 Methods..................................................................................71 Results....................................................................................80 Discussion..............................................................................82 Literature Cited ..................................................................................95 vi LIST OF TABLES 1.1 Individual versus population-level reproductive periods for selected species ......................................................................31 1.2 Parameters of stretched Beta distribution used to model the effects of reproductive asynchrony and protandry on population persistence..................................................................33 2.1 A summary of the data-dependent classification framework for connectivity metrics.............................................59 3.1 Abundances and spatial distribution statistics for the 19 species present in the sample .......................................................89 vii LIST OF FIGURES 1.1 Schematic diagram of reproductive asynchrony..........................34 1.2 Reproductive overlap as a function of asynchrony for even and peaked maturation distributions....................................35 1.3 Consequences of asynchrony for reproductive success...............36 1.4 Joint effects of population density and the duration of population-level reproductive period on realized population growth rate ...................................................................................37 1.5 Extinction risk profiles for various levels of asynchrony ............38 1.6 Effects of protandry on reproductive overlap and extinction risk in populations with reproductive asynchrony.......................39 2.1 Photos of different types of habitat edges....................................60 2.2 Schematic representation of the three types of connectivity .......61 2.3 A simplification of the spatial scales discussed in this paper, based on a recent landcover classification for Jamaica ...............62 2.4 Schematic representation of the tradeoff between information content and data requirements among connectivity metrics........63 3.1 The maximum likelihood fit of the analytical solution to the species abundance data ................................................................90 3.2 Likelihood profiles and 95% confidence intervals for each parameter......................................................................................91 3.3 Dominance-diversity curves for the simulated neutral model predictions and the beetle data.....................................................92 3.4 Frequency distributions of N obs,i under the joint action of process and sampling error for four species ................................93 3.5 Bootstrapped frequency distribution of the modified V statistic .........................................................................................94 1 INTRODUCTION Many population and community-level processes in ecological systems require close spatial and temporal proximity of individual organisms. At the population level, social interactions such as mating, group foraging and group defense obviously require both spatial and temporal proximity of individuals (Allee 1949, Courchamp et al. 1999, Dennis 2002). In fragmented landscapes, the population dynamics of a species will often depend on how close or far apart its constituent subpopulations are located (Hanski 1991, Hanski and Ovaskainen 2003). In communities, species interactions such as certain forms of competition, mutualism and predation require physical proximity of individuals. Classical ecological models have typically incorporated proximity by assuming that populations and communities occur in well-mixed homogeneous spatial and temporal environments (Nicholson and Bailey 1935, MacAurthur and Levins 1967, Hassell and Comins 1976). For example, simple population growth models such as the logistic assume panmictic populations, where any individual could mate with any other. Classical community-level models such as the Lotka- Volterra competition equations assume that entire communities are well mixed and therefore every individual would be able to interact with any other. Mathematical models necessarily omit much of the detail and complexity inherent in natural systems, and assumptions of spatial and temporal proximity often greatly facilitate the modeling process (Case 2000). When compared to empirical data, such modeling exercises inevitably reveal which details may safely be 2 neglected and which are crucial to the functioning of the system (Hilborn and Mangel 1997). Though much has been learned from homogenous classical models, this modeling approach has often misestimated the importance of many ecological interactions (both intra and interspecific) and has led to qualitatively erroneous predictions about the dynamics of populations and communities (Hanski 1991, Hubbell 2001). For example, classical competition models such as Lotka-Volterra suggest that species that do not meet fairly stringent coexistence criteria will be quickly eliminated from the community, yet in natural systems many species that do not appear to meet the required coexistence criteria persist (Hutchinson 1959, Atkinson and Shorrocks 1981, Hubbell 2001). Many advances in ecology have therefore come by contradicting these original assumptions and introducing heterogeneity into ecological models (Levins 1969, Levins and Culver 1971, Horn and MacArthur 1972, Atkinson and Shorrocks 1981, Chesson 1994). Heterogeneity facilitates isolation of individuals, populations and species, with drastic consequences for the dynamics of populations and communities. For example, population dynamics play out very differently when groups of individuals are isolated spatially as they are in a metapopulation (Levins 1969, Hanski 1991, Hanski and Ovaskainen 2003). The degree to which subpopulations are isolated is a critical parameter in determining the resulting metapopulation dynamics (Tischendorf 2001, Moilanen and Nieminen 2002). However, even when individuals occur together in the same place they may not all overlap in time, thus effective population growth rates maybe lower and a population?s extinction risk may be higher than census data 3 pooled throughout the course of a year or breeding season may suggest (Augsperger 1981). Similarly, spatiotemporal heterogeneity and/or life-history tradeoffs among species in communities may isolate species enough either in time or in space that interspecific interactions such as competition may not be nearly as important as homogeneous models may suggest (Levins and Culver 1971, Horn and MacArthur 1972, Atkinson and Shorrocks 1981, 1984, Ives 1988, 1991, Tilman 1994). Despite the advances ecology has made by considering the causes and consequences of isolation in populations and communities, isolation per se is not generally recognized as an integrating concept in ecology. This lack of recognition is due partly to a focus on spatial and temporal heterogeneity, and partly to the difficulty inherent in defining and measuring isolation. Heterogeneity sets the stage for isolation to occur, but it is often isolation that alters ecological interactions. More emphasis on isolation and especially on how it is defined and measured in different contexts is therefore warranted. In this dissertation, I explore the role of isolation as an organizing theme in ecology by studying three different manifestations of isolation in ecological systems. Chapter 1 uses a combination of modeling and meta-analysis to examine how temporal isolation of individuals within a population affects the population?s dynamics and risk of extinction. The effects of isolation in time have received relatively little attention compared to those of isolation in space. Chapter 1 demonstrates how reproductive asynchrony, a usually advantageous bet hedging strategy in temporally unpredictable environments, can create enough temporal 4 isolation among individuals to have a major impact on population dynamics and extinction risk at low population densities. Furthermore, Chapter 1 proposes a simple measure, the ratio of the average individual breeding period to that of the entire population, as a means to quantify temporal isolation in naturally asynchronous populations. Chapter 2 reviews the literature on how connectivity, the inverse of spatial isolation, is measured for populations that have patchy or fragmented spatial structure. Though this is a large and active literature, little consensus as to how connectivity should be defined or measured has emerged. Instead of attempting to derive a universally applicable definition of connectivity (as previous authors have), Chapter 2 deals with the complexity of quantifying connectivity by proposing a new organizational scheme that focuses on the types of ecological data that are required to compute different connectivity metrics. Chapter 3 explores how isolation among species might affect patterns of species diversity at the community level. This chapter focuses on a field study of an assemblage of Costa Rican bark beetles that breed in the fallen petioles of Cecropia insignis trees, and quantifies both spatial distribution and species diversity patterns in this assemblage. Species in this system are strongly aggregated intraspecifically and are distributed with very little covariance (i.e., different species tend not to cue in on the same resource units), and thus are generally isolated from one another inside the resource units (petioles) where competition occurs. Because strong, intraspecific aggregation has been shown to facilitate the coexistence of species of unequal competitive ability, Chapter 3 tests 5 the hypothesis that interspecific competition may not be an important community structuring force by comparing species diversity patterns in the Cecropia assemblage to the predictions of a neutral model that assumes species are functionally equivalent. All three chapters were written as stand-alone manuscripts, and thus the relevant literatures are reviewed within each chapter. Additionally, the specific concepts and techniques used are presented and explained in detail within each chapter. As in any other science, ecology has sought general themes and integrating concepts that can be used to structure our thinking about, and our study of, natural systems. Though many successful themes have emerged (e.g., spatial heterogeneity, allometric scaling, ecological stoichiometry), the search continues for new ways in which to synthesize theory and data. Isolation has been an implicit feature of many other theoretical frameworks in ecology, especially those dealing with spatial or temporal heterogeneity, but its role as an integrating theme has been relatively unexplored. The studies presented in this dissertation represent a broad sweep of the ways in which the concept of isolation may be applied to better understand the dynamics of populations and communities. Individually, each chapter is an original contribution to the ecology literature. Taken together, these papers demonstrate the power of isolation as an organizing theme in ecology. 6 Chapter 1 Lost in time, lonely and single: Reproductive asynchrony and the Allee effect Published in The American Naturalist 164(1): 25-37, 2004. 7 ABSTRACT Identifying linkages between life history traits and small population processes is essential to effective multispecies conservation. Reproductive asynchrony, which occurs when individuals are reproductively active for only a portion of the population-level breeding period, may provide one such link. Traditionally, reproductive asynchrony has been considered from evolutionary perspectives as an advantageous bet-hedging strategy in temporally unpredictable environments. Here, we explore the dynamic consequences of reproductive asynchrony as a density-dependent life history trait. To examine how asynchrony affects population growth rate and extinction risk, we used a general model of reproductive timing to quantify the temporal overlap of opposite-sex individuals and to simulate population dynamics over a range of initial densities and empirical estimates of reproductive asynchrony. We also considered how protandry, a sexually selected life history strategy that often accompanies asynchrony, modulates the population-level effects of reproductive asynchrony. We found that asynchrony 1) decreases the number of males a female overlaps with, 2) decreases the average probability of mating per male/female pair that does overlap, and 3) leaves some females completely isolated in time. This loss of reproductive potential, which is exacerbated by protandry, reduces population growth rate at low density and can lead to extinction via an Allee effect. Thus reproductive asynchrony and protandry, both of which can be evolutionarily advantageous at higher population densities, may prove detrimental when population density declines. 8 INTRODUCTION Maintaining mate-finding efficiency at low population density is of paramount importance to both individual fitness and population persistence. Reduced mate-finding efficiency at low density can cause an Allee effect, where population growth rate is an increasing function of population density (Allee et al. 1949, McCarthy 1997, Wells et al. 1998). Such inverse density dependence may select for increased mate-finding efficiency by favoring individuals that aggregate spatially or employ more efficient mate-location strategies. However, if traits affecting mate-finding efficiency cannot evolve quickly enough in response to this selection pressure, an Allee effect can translate into a lower critical density (termed the ?Allee threshold?) below which population growth rate becomes negative, dooming the population to extinction. If ? is the finite annual rate of increase under conditions of perfect mate finding, the Allee effect can be demonstrated phenomenologically in the context of a geometric growth equation ))(1( 1 NqNN tt ?= + ? Eq. 1.1 where N is female population density, and q is the proportion of females that go mateless (assumed constant across time for a given density). The population will decline to extinction when ? ? )1( )( ? >Nq . Eq. 1.2 Mate finding is generally considered from a spatial perspective, where concerns about the relative locations of male and female individuals or gametes are the focus (McCarthy 1997, Wells et al. 1998, Groom 1998). In this view, high population density results in higher encounter rates among potential mates. In 9 contrast, temporal variation in effective population size has been largely neglected in considerations of mate-finding efficiency. Nevertheless, the framework provided by Equations 1.1 and 1.2 makes clear that isolation in time could lead to an Allee effect in the same way as isolation in space. Reproductive asynchrony, which occurs when individuals are reproductively active at different times within a larger population-level reproductive period, could cause Allee dynamics by reducing the temporal overlap of potential mates. To see this, assume that the probability that a given female and male mate is proportional to their temporal overlap (where d is their maximum possible temporal overlap), and that for each female, each encounter with a male is an independent event and does not influence her probability of mating with any other male. Assume further that females need only mate once to reproduce fully. The probability that a female does not mate, given encounters with n males is then ? = ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?= n i i d overlap malesnmatingNotP 1 1)|( Eq. 1.3 which can be large when n is small (low density) and approaches 0 when n is large (high density). Though several authors have suggested this possibility (Waldbauer 1978, Augspurger 1981, Bullock and Bawa 1981), the interactive effects of asynchrony and population density on population dynamics and extinction risk have not yet been studied in detail. Reproductive phenology, in general, is frequently under strong natural and/or sexual selection and could influence population dynamics and extinction risk because it is often a key determinant of individual reproductive success (del 10 Castillo and Nunez-Farfan 1999, Satake et al. 2001). To date, both theoretical and empirical studies of within-season reproductive phenology have focused on the selective pressures that favor synchronous or asynchronous reproductive strategies in populations where density is not an issue. Some of these have focused on natural selection acting on asynchrony among individuals (Augspurger 1981, Iwasa and Levin 1995, Ollerton and Diaz 1999, Post et. al 2001, Satake et al. 2001), whereas others have focused on sexual selection for asynchrony between the sexes, usually in the form of protandry (Wiklund and Fagerstr?m 1977). Studies on asynchrony among individuals have found that in a temporally unpredictable, coarse-grained environment, reproductive asynchrony ensures that some individuals of an asynchronous genotype attempt to reproduce at a favorable time during the breeding season each year. Much of the population-level variance in reproductive timing in this type of bet-hedging strategy is the product of alternative phenotypes of a given genotype, and not of a polymorphism for maturation time (Simmons and Johnston 1997, Tammaru et al. 1999). An evolutionarily stable distribution of maturation times can result from a single genotype expressing a variety of phenotypes that mature on different dates (Satake et al. 2001). Such ?coin-flipping? plasticity in reproductive timing maximizes a genotype?s geometric mean fitness over multiple generations (Cooper and Kaplan 1982, Seger and Brockmann 1987, Philippi and Seger 1989, Satake et al. 2001). Furthermore, theory predicts (Iwasa and Levin 1995, Satake et al. 2001) and empirical results confirm (Post et al. 2001) that asynchrony 11 increases with the magnitude of large-scale, temporally unpredictable, environmental disturbance. Protandry, where modal reproductive maturity of males precedes that of females, can be advantageous to males when females mate only once and males must compete for receptive females (e.g., Wiklund and Fagerstr?m 1977, Iwasa et al. 1983, Stephenson and Bertin 1983). Under these conditions, protandry increases a male?s chance of successfully mating, and can therefore be strongly favored via sexual selection. However, by separating male and female modal maturation times, protandry could aggravate temporal separation of potential mates at low density, and thus may intensify any population-level effects of reproductive asynchrony among individuals. For asynchrony to be advantageous, a population?s effective density must remain high enough throughout the breeding season that opposite-sex individuals overlap with one another in time. Reproductive asynchrony thus creates a tension between spreading risk in an unpredictable environment and maintaining enough temporal overlap of potential mates throughout the breeding season to ensure reproductive success (Waldbauer 1978). Here we explore the population-dynamic consequences of reproductive asynchrony as a density-dependent life history trait. We assume that the degree of asynchrony in a population remains constant as population density declines because we found no data in the literature quantifying how heritable variance in reproductive timing might be nor how reproductive asynchrony might evolve in response to rapid changes in population density. Exploring the evolutionary 12 dynamics of reproductive asynchrony as population density changes is a good next step, but it is beyond the scope of the present paper. Instead, we focus on demographics to explore the potential of asynchrony to affect small populations. We use a general model of reproductive timing to quantify the temporal overlap of opposite-sex individuals in a population as a function of asynchrony. Empirical data on reproductive timing from a range of asynchronous species, some of which are also protandrous, allow us to restrict our analyses to biologically relevant levels of asynchrony. We find that reproductive asynchrony among individuals can decrease a population?s growth rate at low densities and induce an Allee effect; even small amounts of protandry can exacerbate these effects. In real systems, the population-level consequences of asynchrony will depend on how responsive traits affecting reproductive phenology are to selection at low population density, with both an increased risk of extinction or increased reproductive synchrony as possible outcomes. In either case, asynchrony among individuals and asynchrony between the sexes, both of which can be strongly favored in high-density conditions, appear to be critical but little studied factors at low density. METHODS Compilation of Empirical Data on Reproductive Timing We conducted a literature search to identify representative species for which the timing of reproductive events has been studied in detail (Table 1.1 in the online appendix). We recorded the duration of the reproductive period at both 13 the individual and population levels. When data were available, we recorded the individual reproductive period for each of the sexes separately. For insects, information on the timing of individual reproductive activity was generally not available. Instead, we assumed that the individual reproductive period was equal to adult lifespan or residence time. Insofar as some individuals may not be capable of reproducing throughout their entire adult life or residence in a population, these data overestimate the length of the individual reproductive period, making our estimates of asynchrony somewhat conservative. We quantified the degree of asynchrony in these species as the ratio of the individual- level reproductive period to the population-level reproductive period. When applicable, we also recorded the extent of protandry in the population. The types of empirical data underlying published reports on species? phenologies vary widely among authors. For example, some studies report the mean or median duration of reproductive activity while others report ranges. Because such differences may affect the accuracy of our estimates of the degrees of asynchrony and protandry in these species, we explicitly report in Table 1.1 the types of data used to characterize species? phenologies. Despite these methodological uncertainties, it is clear that wide discrepancies exist between individual-level and population-level phenologies in many natural populations. Overall, our goal was to use this phenological dataset to constrain our mathematical analyses to a range of realistic levels of asynchrony and protandry, not as a basis for precise, quantitative studies of particular species. 14 Development of a Reproductive Timing Model We focused directly on the effects of an asynchronous life history, and thus intentionally omitted other factors that may affect small populations, such as inbreeding depression and skewed sex ratios. We first addressed asynchrony among individuals and later added protandry. For simplicity, we separated the problems of quantifying the effects of asynchrony on a population?s reproductive potential and quantifying its effects on population dynamics. First, we developed a static model that builds asynchronous populations for a given set of parameters and then records several statistics that quantify the loss of reproductive potential due to asynchrony. We then developed a dynamic model that incorporates geometric population growth, and recorded the probability of extinction due to reproductive asynchrony across replicate populations for each parameter set. Both the static and dynamic models have a stochastic element in that we dealt with random draws of individual phenologies from a larger population of possibilities. We used the stretched Beta distribution (Hastings and Peacock 1975, Morris and Doak 2002) to represent the distribution of times at which individuals within a population become reproductively mature. The probability density function of the stretched Beta distribution is () 11 ),( 1 ,,| ?? ? ? ? ? ? ? ? ? ? ? ? ? ? = ?? ?? ?? M xM M x MB Mvx Eq. 1.4 where ?andv are shape parameters controlling the distribution, M (days) is the total duration of the population?s reproductive maturation period, from day 0 when the first individual becomes reproductively active to day M when the last 15 individual in the population initiates its reproductive activity, and B(? , ? ) is the beta function with parameters ?andv ? ?? ?= 1 0 11 .)1(),( duuuB v ? ?? Eq. 1.5 An advantage of the stretched Beta distribution is that the maturation times can be completely constrained to finite intervals while retaining extensive flexibility in shape (Figure 1.1). With this modeling approach, we assume that species-level traits determine both the length of the population-level reproductive maturation period (M) and how concentrated maturation events are within that period. In both models, we first considered non-protandrous populations where male and female Beta distributions were identical and overlapped completely (asynchrony among individuals). We drew reproductive maturation times of individuals at random for populations with N f females and N m males (where N is population density and N f = N m , thus fixing the sex ratio at 1:1) and male and female Beta distribution parameters v = ? = 1, and then v = ? = 4 (Figure 1.1). When v = ? = 1, the Beta distribution is formally equivalent to the uniform distribution (Hastings and Peacock 1975) and individuals are evenly distributed throughout the population-level maturation period, M. For v = ? > 1, a mid- season peak in maturation exists, and this peak becomes more strongly pronounced with further increases in the governing parameters. Once the reproductive maturity time for an individual was drawn, d, the duration of the individual reproductive period, was added to it to obtain each individual?s reproductive activity period. This process was repeated until the reproductive activity periods of all individuals in the population had been determined. We 16 defined the population-level breeding period, D (days), as the length of time during which a non-zero probability of individual reproductive activity exists, which is D = M + d. Though the total densities of males and females were kept equal each season (N f = N m ), the sex ratio at particular times within a breeding season could fluctuate because individual activity periods were defined by randomly chosen, Beta distributed initiation times. For each set of Beta parameters, we examined reproductive asynchrony by varying the population-level maturation period M across 19 levels between five and 50 days. We fixed d at five days for both males and females, assuming that finite resource stores or accumulated damage (e.g., wing wear in butterflies or flower injury in plants) would constrain individuals? reproductive activity. The assumption of equal male and female d is justified based on the empirical dataset (Table 1.1), but in the Discussion we describe the consequences of relaxing this assumption. We quantified the baseline potential for asynchrony in a population as the d/D ratio, and thus could have obtained similar effects by fixing M and varying d. The d/D values we considered in our analyses were within the range of d/D values observed in the empirical dataset (Table 1.1). For each of the 19 levels of reproductive asynchrony, we built 500 replicate populations for each of 29 population densities ranging from 10 to 150 individuals per unit area. In the static model, we tracked several measures of asynchrony and its effects. First we quantified how reproductive asynchrony reduced temporal overlap of females with males at the population level. Summing individual overlaps across all female-male pairs, we calculated realized ?reproductive 17 overlap? as a proportion of maximum possible overlap (which is calculated as N f * N m *d). Reproductive overlap was then averaged over 500 replicate populations of each parameter combination. Second, we tracked the mean number of males that each female overlapped with in each population and then computed a grand mean across the 500 replicate populations. Third, for each population, we recorded the mean overlap for male/female pairs that had overlap > 0, and then calculated a grand mean for this measure over the 500 replicates. Fourth, we recorded the coefficient of variation (CV) of the total temporal overlap of individual females with males within each population to characterize the degree of variability in overlap among females. These CV?s were then averaged across replicate populations to obtain the mean CV of individual overlap. Fifth, to quantify the most extreme effects of reproductive asynchrony on a population?s reproductive potential, we recorded the mean proportion of females that went mateless, either due to complete temporal isolation from males, or probabilistic failure to mate (q(N) from Eq. 1.1). Once the fraction of mateless females was known for a particular replicate, we could calculate what reduction (if any) in population growth rate would be realized over a one-year interval. Building off of the static model, the dynamic model considers populations that reproduce annually and have discrete, non-overlapping generations, such as annual plants and many insects. For each parameter combination we conducted 500 replicate simulations of a stochastic variant of the simple discrete-time geometric population growth model given in Eq. 1.1. In our dynamic model, q(N) is a stochastic term that varies among years based on randomly chosen 18 phenologies of individual males and females. For each replicate, each year, we probabilistically determined if each female in the population would mate based on her temporal overlap with each of the males. Specifically, for each male-female pair, we defined the probability of mating as the pair?s temporal overlap, in days, divided by d. A female needed only to mate with one male to enter the mated pool; multiple successful matings had no effect on fitness. (This constitutes a conservative assumption because fitness of female insects can increase with multiple matings [e.g., Oh 1979]). We set ? = 1.03 so that in the absence of stochastic effects attributable to reproductive asynchrony, the population would grow at a reasonably fast rate. Notice that no density dependence or predetermined Allee threshold is built into this population growth model: Eq. 1.1 has no functional dependence on density and the parameters d, D, v, and ? are assumed independent of density. The shape of the stochastic function q(N) is a consequence of an interaction between asynchrony and population density. In these simulations, we focused on the dynamic consequences of reproductive asynchrony, recording the proportion of the 500 replicate populations that went extinct within 100 years (beyond which time extinction was unlikely to occur because of our assumption of geometric growth). Finally, we considered the combined effects of asynchrony among individuals and protandry on a population?s reproductive potential and dynamics. To simulate protandry, we manipulated the shape parameters of the male and female Beta distributions such that modal reproductive maturity occurred earlier 19 for males than for females. The mode of the Beta distribution, with M scaled to one day, is (Hastings and Peacock, 1975) .1,, )2( )1( > ++ ? = ? ? vfor v v Mode Eq. 1.6 Rescaling yields the actual extent of protandry in the population for each parameter combination. Because rescaling affects the total number of days of protandry for a given set of shape parameters, we used a numerical direct search routine to identify parameter combinations that yielded exactly two days of protandry for all values of M (Table 1.2). Otherwise, we used the same parameters and the same analyses in both the static and dynamic models as above. This approach allowed us to compare asynchronous populations with and without protandry. RESULTS Empirical Data We found data on the reproductive phenologies of 21 species including 16 butterflies, a bee, a stonefly and three flowering plants (Table 1.1). A wide range of reproductive asynchrony (defined as the d/D ratio) was apparent with the lowest degree being 0.52-0.82 for the self-incompatible perennial Discaria toumatou and the maximum being 0.02-0.05 for the monoecious (though rarely selfing) annual Arum maculatum. Insects spanned a slightly narrower range of d/D from the largely synchronous 0.45 (females) and 0.60 (males) stonefly Megarcys signata to the highly asynchronous 0.04-0.24 (females) and 0.05-0.13 (males) butterfly Mellicta athalia. Based on these data, we restricted our analyses 20 to d/D of 0.33 to 0.09, which falls inside the natural range. Male and female phenologies differed for 11 of the 21 species, but quantitative data were available for only five species. Among this subset of species the degree of protandry ranged from ?1 day (technically protgyny) to 21 days. To be conservative, we used only a two-day separation between the modes of male and female reproductive activity in our modeling. Measures of Reproductive Asynchrony Reproductive asynchrony can decrease the number of mating opportunities in a population by reducing the mean temporal overlap of potential mates. For a given combination of the Beta parameters v and ? , holding the individual d constant while increasing the population D decreases mean reproductive overlap (Figure 1.2). However, when expressed as a proportion of the maximum possible overlap, the reduction in reproductive overlap remains constant across the population densities we considered and appears to be an intrinsic feature of that population?s level of asynchrony (determined by d/D and the parameters of the Beta distribution) (Figure 1.3a). This independence of density arises because both the maximum possible overlap and the realized overlap scale as functions of N f * N m . Although proportional overlap itself is density-independent, the reduction in reproductive potential it causes behaves in a density-dependent manner. Reproductive asynchrony acts in three ways to reduce reproductive potential through effects on temporal overlap. First, asynchrony reduces the mean number of males with which each female in the population overlaps in time (Figure 1.3b). 21 Second, a decrease in total reproductive overlap in the population also decreases the average overlap (relative to the maximum possible pair-wise overlap (i.e., d)) of those male/female pairs that do overlap in time. Because we have defined the probability of mating per male/female pair in the population as the realized proportion of their maximum possible overlap, asynchrony increases the number of females in the population that are mateless due to probabilistic failure to mate. Finally, at low population densities (generally < 20 individuals/unit area, but dependent on d/D) some females are mateless by virtue of complete temporal isolation. These three effects conspire to increase, q(N), the mean proportion of females that fail to reproduce as density declines in an asynchronous population (Figure 1.3c). The coefficient of variation of female overlap with males increases with the level of asynchrony in the population and with decreasing population density (Figure 1.3d). This variability in overlap among females can be considered a form of demographic stochasticity. At population densities where the CV begins to climb rapidly in the static model (Figure 1.3d), populations have already gone extinct in the dynamic model, suggesting that inter-individual variability in overlap modulates, but does not drive, the observed extinction dynamics (see below). The consequences of altering v and ? to manipulate the shape of the reproductive maturation distributions for a given d/D ratio were weak compared to manipulating the d/D ratio for a given distribution shape. When d/D was held constant, spreading individuals more evenly across the breeding season (v, ? = 1) 22 slightly decreased mean reproductive overlap, exacerbating the negative effects of asynchrony relative to the case where individuals were more concentrated in time (v, ? = 4) (Figure 1.3a). Similarly, for a given d/D ratio, spreading individuals more evenly across the breeding season slightly increased the CV of overlap among females (Figure 1.3d). Thus, the shape of the reproductive maturation distribution acts only to modulate the effects of asynchrony determined by the d/D ratio. Effects of Asynchrony on Population Growth Rate and Extinction Risk The proportion of mateless females, q(N), in the population directly affects realized population growth rate, which in turn determines the probability of extinction in the dynamic model. Because increases in D for a fixed d increase mean q(N), increasing D strongly reduces mean realized growth rate over one- year intervals (Figure 1.4a & b) and increases the fraction of replicate populations in decline during a given time step (Figure 1.4c & d). The effects of asynchrony on population dynamics scale nonlinearly with density in that a given increase in D has larger consequences for small populations than for large (to see this note that the contour lines in Figure 1.4 are not parallel). The predominance of d/D over (v, ? ) is also apparent in Figure 1.4. For a given d/D, shifting from a broadly asynchronous reproductive distribution of reproductive activity (v = ? = 1) to a distribution that is quite concentrated in time (v = ? = 4) makes only small changes to the slopes of the contour lines describing realized growth rate (contrast Figure 1.4a & c with 1.4b & d). 23 For the most extreme levels of asynchrony we considered (d/D = 0.09), populations regularly went extinct at total densities (males + females) of 70-80 individuals / unit area (Figure 1.5). Even for minimal levels of asynchrony (d/D = 0.33), populations still regularly went extinct at total densities of 10-20 individuals / unit area. Thus, even acting alone, loss of reproductive potential due to asynchrony among individuals can drive an otherwise-growing population extinct. The shape of the probability of extinction profiles is consistent with the expectation that reproductive asynchrony causes an Allee effect (Figure 1.5). This result is also in agreement with other studies of reduced mating efficiency at low population density, but because our model includes stochasticity in q(N), there was no specific Allee threshold, per se. Instead, in all cases, a population?s probability of extinction transitioned from 0 to 1 over a small range of density (Figure 1.5). The d/D ratio had the strongest effect on a population?s probability of extinction, whereas manipulating the shape parameters of the Beta distributions for a given d/D had small effects on the probability of extinction (Figure 1.5). Effects of Asynchrony on Protandrous Populations Protandry, as expected, exacerbates the negative effects of reproductive asynchrony among individuals by further reducing mean reproductive overlap between potential mates. As an example, consider a population with d/D = 0.2, male v = 3.86, ? = 3.14 and female v = 3.14, ? = 3.86. These parameters result in male and female maturation distributions that are symmetrical with respect to one another and feature two days of protandry (based on the difference in modes of the distributions). Even this minimal level of protandry negatively affects 24 reproductive overlap relative to non-protandrous populations with similar shape parameters (Figure 1.6a). Protandry had comparable effects on populations across a range of d/D values, relative to similar non-protandrous populations (results not shown). For a given population density, the extra reduction in reproductive potential due to protandry increases a population?s risk of extinction compared to a population that is asynchronous but not protandrous (Figure 1.6b). Thus, protandry can act synergistically with asynchrony among individuals to increase the risk of extinction at low population density. DISCUSSION Anthropogenically driven declines of many species have forced ecologists and evolutionary biologists to consider density dependence in the population-level effects of life history traits. Such analyses can both identify life history traits that may influence population persistence (Pimm et al. 1988, Saether 1997, Fagan et al. 2001, Johnson 2002, Green 2003) and highlight selection pressures that can affect life history evolution. Though reproductive asynchrony?either among individuals or between the sexes?can be advantageous at high density, we have demonstrated here that it can have hidden consequences at low population density. We found that biologically realistic levels of reproductive asynchrony (Table 1.1) reduce the reproductive potential of the population by decreasing the temporal overlap of potential mates. Reduced mating efficiency at low population density, regardless of the specific mechanism that causes it, leads to an Allee effect (McCarthy 1997). 25 The Allee effects observed in this study emerge from the interaction of population density with reduced mating efficiency caused by variable reproductive timing among individuals; they do not derive from a predetermined ?Allee threshold? in our population growth model. Instead, reproductive asynchrony itself acts as a mechanism generating the Allee effect. Specifically, a female?s total probability of mating within a breeding season depends on the density of males during her reproductive activity period. Male density at any point during the breeding season, in turn, is affected by both the total male population density and the temporal distribution of male reproductive activity across the breeding season. Reproductive asynchrony therefore satisfies the criterion of inverse density dependence at low population density necessary for the operation of an Allee effect (Courchamp et al. 1999). Accordingly, both realized population growth rate (Figure 1.4) and extinction risk are affected (Figure 1.5). Variability among females in total reproductive overlap due to sampling effects at small population densities can be considered a form of demographic stochasticity. It causes the realized population growth rate contours to be ?messy? (Figure 1.4) and the (0, 1) step function for extinction probability in deterministic Allee effect models to be ?blurred? into a sigmoidal curve that decreases as a function of population density (Figure 1.5) (see also Boukal and Berec 2002). Several authors have noted this ?stochastic blurring? effect in models that explicitly include both Allee effects and demographic stochasticity (Dennis 1989, 2002, Berec et al. 2001). Stochastic Allee effects are characterized by 1) a probability of extinction versus initial population density curve that 26 exhibits an inflection point and 2) a sharp transition in the probability of extinction near that inflection point (Dennis 1989, 2002, Boukal and Berec 2002). These patterns differ markedly from the dynamic characteristics of demographic stochasticity alone, in which the probability of extinction increases smoothly and gradually with decreasing population size (Dennis 2002). Populations suffering from an Allee effect induced by reproductive asynchrony are therefore more likely to exhibit sudden crashes than those suffering from demographic stochasticity per se. Protandry, which separates the modal maturation times of males and females within a population, clearly exacerbates the effects of reproductive asynchrony among individuals, placing populations at greater risk of extinction for a given density (Figure 1.6). Even the minimal degree of protandry we considered (two days) had significant effects on a population?s probability of extinction. Empirical data suggest that protandry can be far more extreme (Table 1.1). For example, the meadow brown butterfly Maniola jurtina, had approximately 21 days of protandry and a d/D ratio between 0.11 and 0.2! It must be recognized however that M. jurtina is a common, and occasionally abundant, species, and it is not clear that such extreme protandry would persist in small populations. Indeed, the density dependence of protandry appears to be quite open as an area of inquiry. The timing of the initiation of male reproductive activity is frequently under strong sexual selection in populations of butterflies (Wiklund and Fagerstr?m 1977, Wiklund and Solbreck 1982, Iwasa et al. 1983), dioecious 27 plants (Purrington 1993, Purrington and Schmitt 1998), and other species (e.g., del Castillo and Nunez-Farfan 1999, Holzapfel and Bradshaw 2002). Consequently, the same kinds of species that feature major discrepancies between individual and population-level reproductive periods (i.e., small d/D ratios) frequently exhibit significant protandry (Table 1.1). Several recent papers have noted the potential influence of certain sexually selected traits on extinction risk, but none, to our knowledge, have dealt with the added risks associated with sexual selection acting on phenology (Doherty et al. 2003, Kokko and Brooks 2003, M?ller 2003). The potential to explore issues like phenology that may differ between males and females is one advantage of working with two-sex models when examining extinction risk (see also Engen et al. 2003) Clearly, a variety of changes to the model, such as making ? larger, making the sex ratio consistently male-biased or lengthening male d relative to female d, will lessen the severity of the loss of reproductive potential caused by asynchrony. For example, consider that in many species, a few individuals will have long individual reproductive periods, while most individuals hover close to the population mean. Adding this kind of inter-individual variability in d would likely decrease the negative effects of asynchrony, but would require additional model complexity relating to individual senescence or limits on the number of matings per male per unit time or per lifetime. Still, our results show that even when all individuals in the population are long-lived (high d/D ratio), asynchrony reduces a population?s growth rate and elevates its extinction risk relative to a synchronously breeding population. In contrast to the above suite of factors that 28 could lessen the effects of asynchrony, any spatial processes that reduce effective population density, such as limited search area or imperfect mate locating ability, will exacerbate the effects of asynchrony. Additional development of this modeling framework is warranted to explore its sensitivity to assumptions we made concerning density independence of the key parameters d and D and the emphasis on annual life cycles. For example, if individual reproductive timing is highly heritable, then D could narrow with decreasing density, as those females that were closely synchronized with the bulk of the male population would be more likely to reproduce. Although numerous experimental studies, especially in plants, have assessed heritability of the date of first reproduction (Matziris 1994, Kelly and Levin 1997, Nikkanen 2001, Tikkanen and Lyytikainen-Saarenmaa 2002), the degree to which variance in reproductive timing is heritable appears little explored. Likewise, the degree to which individual d can evolve in response to selection for more synchronized reproduction at low densities appears worthy of study. Another obvious extension would be to explore the dynamic consequences of reproductive asynchrony in perennial populations. Quantifying the effects of reproductive asynchrony in perennial species would require modifying the population growth model we employed (e.g., shifting the focus to geometric average growth rates per generation). However, the phenomenon seems likely to remain important because, even in perennial species, reproductive asynchrony could reduce an individual?s lifetime reproductive success and alter population-level recruitment patterns. 29 Other authors have noted potential effects of asynchrony on populations of insects (Waldbauer 1978), plants (Primack 1980, Augspurger 1981, Ollerton and Diaz 1999), the maintenance of both plant-pollinator mutualisms (Anstett et al. 1995) and host-parasitoid interactions (Godfray et al. 1994), but none have studied quantitatively the density-dependent effects of asynchrony per se on population growth and extinction risk. Our results demonstrate that reproductive asynchrony may strongly affect low-density populations, particularly when the ratio of the individual-level reproductive period to the population-level reproductive period is less than one-third. Several species in our analysis exhibit population parameters that, in our model, cause considerable decreases in population growth rate and make a population quite vulnerable to extinction at low density (Table 1.1). Acting alone or synergistically with life history traits such as protandry, reproductive asynchrony among individuals can reduce population growth rate and increase extinction risk. The severity of these effects may hinge upon how quickly traits affecting individual reproductive timing can respond to selection for reproductive synchrony at low density. A quick response at low density could serve as a buffer against the negative effects of asynchrony, permitting variable populations to be more asynchronous at high density. In contrast, a slow response might allow the Allee effect to limit the degree of asynchrony that is advantageous in natural populations. It is ironic that reproductive asynchrony and protandry, both of which may be under strong positive selection at high density, may be quite disadvantageous to population persistence at low density. Taken to the extreme, 30 reproductive asynchrony could provide another example of evolutionary suicide. Reproductive asynchrony should therefore be recognized as a mechanism of the Allee effect and be included among the suite of life history characters analyzed when determining a species? extinction risk at low population density. More generally, the consequences of phenological variation among individuals have not received adequate attention in relation to population dynamics and extinction risk. It is clear from the literature that, as commonly used, ?phenology? usually refers to population-level events such as the flight period in butterflies or blooming time in flowering plants. Our results highlight the importance of the distinction between the phenology of individuals and the phenology of populations and outline some of the consequences of this relationship for ecological systems. 31 Table 1.1 Individual versus population-level reproductive periods for selected species. Organism Type Species Individual Reproductive Period, d (days) Population Reproductive Period, D (days) Ratio of Individual to Population Reproductive Period (d / D) Protandry (days) Reference Butterflies Papilio polyxenes ? M: 7.3-12.4 ? 35-44 ' 0.17-0.35 No Data Lederhouse 1983 Leptidea sinapis % F: 8.0-10.6 ? M: 8.2-9.8 ? 46-73' * F: 0.11-0.23 M: 0.11-0.21 0-20 ? * Warren et al. 1986 Brassolis sophorae $ 6.1-11.9 ? 36 ? 0.17-0.33 - 1-13 ? Carvalho and Queiroz 1998 Maniola jurinata % F: 7.6-12.7 ? M: 6.7-8.7 ? 60 ? * F: 0.13-0.21 M: 0.11-0.15 21? * Pollard 1981 Euphydryas editha bayensis % 4? 21-35' 0.11-0.19 No Data Singer and Ehrlich 1979, Cushman et al. 1994 Euphydryas aurinia % F: 8.9 + M: 10.7 + 31 ? F: 0.29 M: 0.35 No Data Wahlberg et al. 2002 Euphydryas maturna % F: 3.3 + M: 13.3 + 35 ? F: 0.09 M: 0.38 No Data Wahlberg et al. 2002 Melitaea cinxia % F: 3.0 + M: 8.2 + 29 ? F: 0.10 M: 0.28 No Data Wahlberg et al. 2002 Melitaea diamina % F: 6.7 + M: 6.0 + 29 ? F: 0.23 M: 0.21 No Data Wahlberg et al. 2002 Melitaea athalia % F: 7.2 + M: 12.5 + 35 ? F: 0.21 M: 0.36 No Data Wahlberg et al. 2002 Mellicta athalia % F: 2.3-10.8 ? M: 2.7-6.0 ? 45-60' * F: 0.04-0.24 M: 0.05-0.13 No Data Warren 1987a, b Proclossiana eunomia % F: 2.2-13.0 ' M: 3.8-11.1 ' 21-35 ' F: 0.06-0.62 M: 0.11-0.53 Yes Schtickzelle et al. 2002 Euphiolotes enoptes % ~ 2-9 ? 29-49 ' 0.04-0.31 Yes Peterson 1995 Icaricia icariodes fenderi % 15 + 28-42' 0.36-0.54 No Data C. Schultz, pers. comm. Lysandra bellargus $ F:10.6-12 ? M: 4.2-9.5 ? 60 ? * F: 0.18-0.20 M: 0.07-0.16 No Data Davis et al. 1958, Pollard and Yates 1993 Lysandra coridon % F: 4.7 + M: 6.6 + 70 ? * F: 0.07 M: 0.09 7? * Davis et al. 1958, Pollard and Yates 1993 Polyommatus icarus $ 5.4 + 49? * 0.11 Yes Dowdeswell et al. 1940, Pollard and Yates 1993 Solitary Bees Amegilla dawsoni % 6.2-9.1 ? 35 ? * 0.18-0.26 Yes Alcock 1996, 1997, 1999 Stoneflies Megarcys F: 18 + 40? * F: 0.45 5 ? * Taylor et al. 32 Footnotes for Table 1.1: Individual reproductive periods for plants are underestimates; they explicitly exclude persistence times of dispersed pollen. F=Female, M=Male. A ?yes? entry in the protandry column indicates that the species is known to be protandrous, but no quantitative estimate was available. A ?No data? entry in the protandry column indicates that we could not find information on whether or not the species was protandrous. Key to footnote symbols: % univoltine $ bivoltine ? multvoltine ~ species also known to bet-hedge across years through variation in the duration of the pupal stage @ monoecious perennial, self-pollination rare in nature # monoecious perennial, self-incompatible ^ dioecious perennial + mean ? mean ? 1 standard deviation ? range of mean from different samples ' range ? single value ? median & mode ? difference between modes of male and female emergence distributions ? range of difference between median male and female emergence dates from different samples * estimated from graph signata % M: 24 + M: 0.60 1998 Flowering Plants Arum maculatum @ 1 & 19.5-40.5 ? 0.02-0.05 Yes Ollerton and Diaz 1999, Sowter 1949 Couratari multiflora ^ 15-60 ' * 195 ? * 0.08-0.31 No Data Lepsch- Cunha and Mori 1999 Discaria toumatou # 17-23 ? 28-33 ' 0.52-0.82 Yes Primack 1980 33 Table 1.2 Parameters of stretched Beta distribution used to model the effects of reproductive asynchrony and protandry on population persistence. Parameters were chosen to obtain a protandrous reproductive activity pattern with a constant two-day separation between modes of symmetrical male and female distributions. Symmetry arises when Male v= Female ? and Male ? = Female v. M (days) Female v Female ? 5 3.900 1.100 8 4.000 2.000 10 3.965 2.310 12 3.970 2.550 15 3.930 2.770 18 3.800 2.840 20 3.785 2.915 22 3.800 3.000 25 3.860 3.140 28 3.950 3.290 30 3.880 3.270 32 3.760 3.200 35 3.810 3.290 38 3.800 3.320 40 3.725 3.275 42 3.730 3.300 45 3.935 3.515 48 3.875 3.485 50 3.745 3.380 34 Figure 1.1 Schematic diagram of reproductive asynchrony. A) The relationship between individual and population-level reproductive periods (d and D, respectively). Horizontal bars represent male and female individual reproductive activity periods, whereas vertical bars demonstrate how one would quantify overlap between individual males and females. B) The stretched Beta distribution is flexible enough to treat situations in which reproductive activity is broadly asynchronous (v = ? = 1.5) or concentrated and highly skewed (v = 1.5, ? = 2.5). C) The stretched Beta distribution can be used to study protandry by generating different distributions for male and female reproductive activity. 35 Figure 1.2 Reproductive overlap as a function of asynchrony for even (v, ? = 1) and peaked (v, ? = 4) maturation distributions. Reproductive asynchrony increases with decreasing d/D ratio. 36 Figure 1.3 Consequences of asynchrony for reproductive success. A) Proportional overlap differs for each level of asynchrony (Beta parameters and d/D ratio), but, for a given combination of parameters, remains constant over the range of population densities considered. B) Mean number of males a female overlaps with in populations characterized by different levels of asynchrony. An increase in reproductive asynchrony decreases the slope of the relationship between mean number of males per female and population density. C) Average proportion of females in the population that are mateless due to reproductive asynchrony. Failure to mate increases sharply with decreasing population density and with increasing degrees of asynchrony in the population. D) Among-female variability in reproductive overlap with males. Variability increases with decreasing population density and with increasing levels of asynchrony. Changes in the Beta parameters, and thus the variance, of maturation distributions for a given d/D had only minor effects on variability. 37 Figure 1.4 Joint effects of population density and the duration of population-level reproductive period on realized population growth rate. Panels A and B provide contours of realized population growth rate (arithmetic mean across 500 replicates). Values < 1.0 correspond to populations that would on average decline due to reproductive asynchrony. Panels C and D provide contours of the proportion of 500 replicate populations with realized population growth rate < 1.0. Panels A and C are for populations with uniform distributions of reproductive activity (v = ? =1) whereas Panels B and D represent populations with a mid- season peak in reproductive activity (v = ? = 4). 38 Figure 1.5 Extinction risk profiles for various levels of asynchrony. A population?s d/D ratio was the main determinant of extinction risk for that population. Changes in the variance of the maturation distributions for a given d/D ratio had minor effects. 39 Figure 1.6 Effects of protandry on reproductive overlap and extinction risk in populations with reproductive asynchrony. A) Reproductive overlap for a population with two days of protandry (male v = 3.86, ? = 3.14, female v = 3.14, ? =3.86) and a population with no protandry (male and female v, ? = 3.14). For both populations, d/D = 0.2. B) Extinction risk profiles for a population with two days of protandry and a population with no protandry (parameters as in A). 40 Chapter 2 A comparison-shopper's guide to connectivity metrics Published in Frontiers in Ecology and the Environment 2(10): 529-536, 2004 41 ABSTRACT Connectivity is an important but inconsistently defined concept in spatial ecology and conservation biology. Theoreticians from various sub disciplines of ecology argue over its definition and measurement, but no consensus has yet emerged. Despite this disagreement, measuring connectivity is an integral part of many resource management plans. A more practical approach to understanding the many connectivity metrics is needed. Instead of focusing on theoretical issues surrounding the concept of connectivity, we describe a data-dependent framework for classifying these metrics. This framework illustrates the data requirements, spatial scales, and information yields of a range of different connectivity measures. By highlighting the costs and benefits associated with using alternative metrics, this framework allows practitioners to make more informed decisions concerning connectivity measurement. 42 INTRODUCTION Dispersal, the movement of individuals among populations, is a critical ecological process (Ims andYoccoz 1997). It can maintain genetic diversity, rescue declining populations, and re-establish extirpated populations. Sufficient movement of individuals between isolated, extinction-prone populations can allow an entire network of populations to persist via metapopulation dynamics (Hanski 1991). As areas of natural habitat are reduced in size and continuity by human activities, the degree to which the remaining fragments are functionally linked by dispersal becomes increasingly important. The strength of those linkages is determined largely by a property known as "connectivity", which, despite its intuitive appeal, is inconsistently defined. At one extreme, metapopulation ecologists argue for a habitat patch-level definition, while at the other, landscape ecologists insist that connectivity is a landscape-scale property (Merriam 1984, Taylor et al. 1993, Tischendorf and Fahrig 2000, Moilanen and Hanski 2001, Tischendorf 2001a, Moilanen and Nieminen 2002). Differences in perspective notwithstanding, theoreticians do agree that connectivity has undeniable effects on many population processes (Wiens 1997, Moilanen and Hanski 2001). It is therefore desirable to quantify connectivity and use these measurements as a basis for decision-making. Currently, many reserve design algorithms factor in some measure of connectivity when weighing alternative plans (Siitonen et al. 2002, 2003, Singleton et al. 2002, Cabeza 2003). Consideration of connectivity during the reserve design process could highlight 43 situations where it really matters. For example, alternative reserve designs that are similar in other factors such as area, habitat quality, and cost may differ greatly in connectivity (Siitonen et al. 2002). This matters because the low-connectivity scenarios may not be able to support viable populations of certain species over long periods of time. Analyses of this sort could also redirect some project resources towards improving the connectivity of a reserve network by building movement corridors or acquiring small, otherwise undesirable habitat patches that act as links between larger patches (Keitt et al. 1997). Reserve designs could therefore include the demographic and genetic benefits of increased connectivity without substantially increasing the cost of the project (e.g., Siitonen et al. 2002). If connectivity is to serve as a guide, at least in part, for conservation decision-making, it clearly matters how it is measured. Unfortunately, the ecological literature is awash with different connectivity metrics. How are land managers and decision makers to efficiently choose between these alternatives, when ecologists cannot even agree on a basic definition of connectivity, let alone how it is best measured? Aside from the theoretical perspectives to which they are tied, these metrics differ in two important regards: the type of data they require and the level of detail they provide. Here, we attempt to cut through some of the confusion surrounding connectivity by developing a classification scheme based on these key differences between metrics. Connectivity depends on the interaction between particular species and the landscapes in which they occur (Schumaker 1996, Wiens 1997, Tischendorf and Fahrig 2000, Moilanen and Hanski 2001). Put another way, a single landscape or 44 habitat patch will possess different degrees of connectivity, depending on the behaviors, habitat preferences, and dispersal abilities of the species being considered (Johnson and Gaines 1985; Figure 2.1). Strategies exist for developing multi-species connectivity metrics (Fagan and Calabrese in press), but here we stick to the standard, single species view. We distinguish three classes of connectivity metrics, based on interactions between focal species and the landscape. Listed in increasing order of detail, they are: structural, potential, and actual connectivity (Figure 2.2). Structural connectivity is derived from physical attributes of the landscape, such as size, shape, and location of habitat patches, but does not factor in dispersal ability (Figure 2.2a). Potential connectivity combines these physical attributes of the landscape with limited information about dispersal ability to predict how connected a given landscape or patch will be for a species (Figure 2.2b). Examples of limited dispersal information include estimates of mobility derived from body size or energy budgets (Cresswell et al. 2000, Porter et al. 2000), or measurements with little spatial detail, such as mean or maximum recapture distances from mark-recapture studies (Clark et al. 2001). Actual connectivity relates to the observation of individuals moving into or out of focal patches, or through a landscape, and thus provides a concrete estimate of the linkages between landscape elements or habitat patches (Figure 2.2c). To facilitate classification of connectivity metrics according to their data- dependence, the various types of data used to estimate connectivity are simplified into six frequently encountered categories (see below). Within each data category, the spatial scales at which the metrics are usually calculated are simplified to four 45 levels: point occurrences, individual habitat patches, landscape classes, and entire landscapes (Figure 2.3). Our approach here is to sketch the relationships between the three types of connectivity described above and the basic data requirements of the various connectivity metrics (Table 2.1). We also discuss the common modifications to many connectivity metrics and the scale-dependence of connectivity. The DATA-DEPENDENT FRAMEWORK Nearest Neighbor Distance: Patch Occupancy Data and Interpatch Distance Field surveys of a species' occupancy pattern in a habitat patch and measurements of the distance to the nearest occupied patch provide a simple, patch-level structural connectivity metric. Interpatch distance is, technically, a patch isolation measure, and connectivity is its inverse. Though simple to obtain, distance to the nearest occupied neighbor is a crude connectivity metric. Moilanen and Nieminen (2002) demonstrated the poor performance of this metric through a meta-analysis of published studies that quantified connectivity, and by using various connectivity metrics to predict colonization events in two detailed empirical butterfly metapopulation datasets. Overall, they found that nearest neighbor measures were less likely to detect a significant effect of connectivity and were more sensitive to sample size than were other, more complex connectivity metrics. Bender et al. (2003) obtained similar results using a computer-simulated dispersal process on both real - derived from a geographic information system (GIS) - and artificially generated landscapes. They found that 46 nearest neighbor distance was consistently the worst or second worst performer of the four proximity indices they studied, and that it performed especially poorly when patch size and shape were varied (Bender et al. 2003). The weak performance of nearest neighbor distance can be attributed to several factors. First, this metric counts only the contribution of the patch nearest to the focal patch, thus ignoring how all other patches affect the connectivity of the focal patch (Bender et al. 2003). Furthermore, in its most basic form, the nearest neighbor measure includes no information about the population size of the focal species in the nearest patch. Finally, no knowledge of the species' dispersal ability is incorporated into the metric. Despite these limitations, the nearest neighbor distance is one of the most commonly used connectivity metrics (Moilanen and Nieminen 2002, Bender et al. 2003). This is most likely due to its simplicity and modest data requirements. Unfortunately, these advantages do not adequately compensate for its limitations. Spatial Pattern Indices: Spatially Explicit Habitat Data Spatially explicit habitat data are often remotely sensed, cover a large area, and are represented in either raster or vector form in a GIS. Spatial pattern indices quantify the number, size, extent, shape, or aspects of the spatial arrangement of landscape elements. The use of these indices as connectivity metrics relies on the assumption that the spatial patterns these indices quantify actually affect species' ability to move through the landscape. Examples of spatial pattern metrics include number of patches, patch area, core area, patch perimeter, contagion, perimeter-area ratio, shape index, fractal dimension, and patch 47 cohesion (Haines-Young and Chopping 1996, Schumaker1996). The increasing availability of this type of data and software packages such as Fragstats (McGarigal et al. 2002) make the metrics in this category relatively easy to calculate. Although spatial pattern indices are sometimes assumed to represent actual connectivity, we consider them estimators of structural connectivity because they do not incorporate dispersal data. The lack of dispersal data does not, however, preclude the possibility that these indices could show predictable relationships with actual connectivity. There has been little empirical research regarding this possibility, but several simulation-modeling studies have explored the relationships between spatial pattern indices and dispersal success. For example, Schumaker (1996) demonstrated that shape index and patch cohesion were the best predictors of dispersal success, while fractal dimension, number of patches, patch area, core area, patch perimeter, contagion, and perimeter-area ratio were, at best, weakly correlated with dispersal success. Similarly, Tischendorf (2001b) found that, while some spatial pattern indices were strongly correlated with simulated dispersal success, 68% of the statistical relationships between the 26 metrics and three measures of dispersal success considered were inconsistent when landscape structure and dispersal behavior were varied. The simulation results therefore suggest that relationships between spatial pattern indices and dispersal success might not generalize well across landscapes or species. A potential advantage of spatial pattern indices is that they could be used to quickly characterize connectivity for large areas. However, the weak or 48 inconsistent relationships between spatial pattern indices and dispersal success suggest that further research is required before these indices can be relied upon to estimate actual connectivity. The lack of empirical work in this area only underscores this point. As several authors have noted (Schumaker 1996, Tischendorf 2001b, Fortin et al. 2003), focusing on the relationships between the spatial pattern that these metrics quantify and the underlying ecological processes that influence connectivity, such as demographics, dispersal, and behavior, may be the most effective way to develop these metrics further. Scale-Area Slope: Point- or Grid-Based Occurrence Data Another approach to quantifying structural connectivity can be used when records of species' spatial occurrences are available, but the locations of actual habitat patches are unknown. Datasets fitting into this category include those assembled from museum records or long-term surveys of species presence or absence, where patch boundaries are not known or may have changed since the data were collected. This approach builds from individual occurrences of a species to a landscape-level connectivity metric known as the "scale-area slope". Both point data, where considerable spatial detail is available, and grid data, where spatial descriptions are less precise, can be used to estimate structural connectivity based on the slope of a scale-area curve (Kunin 1998, Fagan et al. 2002). Scale-area slopes are derived by dividing a landscape into a series of equal-sized grid cells at several map resolutions, with a fixed number of fine- resolution cells inside each coarser-resolution cell. Presence or absence of the focal species in each cell at each resolution is determined and the map area 49 occupied by the species (assuming a cell with at least one incidence record is occupied) is plotted against grid cell size at each map resolution. Scale-area slope is then estimated via power-law regression. Steep scale-area slopes characterize species that have fragmented distributions, whereas shallow slopes identify species with less fragmented (i.e., more contiguous) spatial occurrences. A shallow (i.e., numerically small) scale-area slope would therefore be associated with higher structural connectivity. The use of the scale-area slope as a connectivity metric assumes that proximity is the major determinant of the connectivity among occurrences. Such an assumption is clearly justified in certain circumstances. For example, Fagan et al. (2002) demonstrated that for Sonoran Desert fishes, species that were historically distributed more compactly (i.e., species with shallow scale-area slopes) were at a distinct advantage when it came to weathering the ensuing decades of anthropogenic alterations to their habitats and landscape. In contrast, species with steep scale-area slopes, whose distributions were more fragmented historically, were at greater risk of local extinction. Despite this promising result, the relationships between scale-area slope and various measures of actual connectivity have not yet been established. Although scale-area approaches do not provide a direct linkage between connectivity and dispersal, the techniques can help to identify the spatial scales over which processes affecting connectivity are most important. 50 Graph-Theoretic Measures: Spatially Explicit Habitat Data with Dispersal Data Graph-theoretic measures combine spatially explicit habitat data derived from a GIS with data acquired from independent studies on the dispersal biology of species. Inclusion of species-specific dispersal data represents a substantial increase in data requirements, but allows these metrics to go beyond structural connectivity and address potential connectivity. In their most basic form, graph- theoretic approaches entail making a mathematical "graph" of a network of habitat patches for a species that incorporates information on the spatial arrangement of patches as well as patch attributes (Cantwell and Forman 1993, Keitt et al. 1997, Bunn et al. 2000, Urban and Keitt 2001). The graph is simply a means of summarizing the spatial relationships between landscape elements in a concise way. Next, potential connections between all pairwise combinations of habitat patches are established by considering the dispersal ability of the focal species. If the distance between a given pair of patches is less than or equal to the measure of dispersal ability used, the patches are considered connected. Measures of dispersal ability typically include a fixed critical dispersal distance (Keitt et al.1997, D'Eon et al. 2002) or a random draw from a dispersal kernel. A fixed critical distance represents the distance after which a species' probability of dispersal is assumed to decline rapidly (van Langvelde 2000), while a dispersal kernel is a function describing the relationship between dispersal distance and a species' probability of dispersal (e.g., Kot et al. 1996, Havel et al. 2002). These potential connections are depicted on the graph as lines ("edges" in the terminology of graph theory) drawn between each pair of connected patches. 51 After establishing pairwise connections, graph-theoretic approaches scale up to consider the connectivity of the entire patch network (landscape level), using metrics including correlation length, distance to cluster edge, number of graph components, and diameter of the largest graph component (Keitt et al. 1997, Urban and Keitt 2001, D'Eon et al. 2002). These metrics are different ways of quantifying how connected the graph is overall. For example, a graph that had one large cluster of interconnected patches would be considered to have higher connectivity than a graph that had several small, isolated clusters of interconnected patches. An advantage of these methods is that graph operations that simulate the destruction of habitat patches or dispersal corridors can be used to rank habitat patches by their contributions to landscape-level connectivity (Keitt et al. 1997). The graph-theoretic approach could therefore allow land managers to make decisions based on which patches are most critical to landscape connectivity. Buffer Radius and Incidence Function Metrics: Spatially Explicit Patch Occupancy, Patch Area, and Dispersal Data Spatially explicit patch occupancy data are usually obtained by directly sampling habitat patches for a species of interest and spatially referencing patch locations. With such data, one can calculate buffer radius or incidence function measures (see below) of patch-level, potential connectivity, depending on assumptions about the dispersal biology of the species in question. These metrics incorporate patch occupancy information, usually for a large number of patches. Such data allow the potential contribution of each patch to be assessed by its occupancy status as well as by proxies for population size, such as area, if it is 52 occupied. The result of this extra information is that these indices can give a more detailed estimate of patch-level potential connectivity than other metrics. For buffer radius measures, patch-occupancy data for all patches that lie within a fixed distance, or "buffer radius", of the focal patch are required. The connectivity of a patch is a function of the number and areas of all occupied patches that lie within the buffer radius. Though buffer radii are often arbitrarily selected, Moilanen and Nieminen (2002) have shown that the performance of these measures is sensitive to the buffer radius chosen, suggesting that incorporation of even the most basic dispersal information could substantially improve the performance of these metrics. A similar set of connectivity metrics derive from the incidence function metapopulation model (IFM) (Hanski 1994, Hanski et al. 1996). These measures require spatially explicit patch-occupancy data for a large number of patches in a metapopulation, and also a dispersal kernel describing how the focal species' probability of dispersal decays with distance. The dispersal kernel can be parameterized either with independent data on the dispersal ability of the focal species, or by model fit to patch-occupancy data. If such data are used to estimate dispersal ability, it is desirable to have more than one year of data to obtain robust parameter estimates (Hanski 2001). The basic IFM connectivity measure essentially sums the potential contribution of all occupied patches in a metapopulation, weighted by area and distance, to the connectivity of a focal patch. 53 Buffer and IFM metrics can still be calculated in the absence of patch- occupancy data, but give a less detailed estimate of potential connectivity. This method, called the "connectivity of landscape elements" (Moilanen and Hanski 2001), is similar to the graph-theoretic approaches described above. When patch-occupancy data are available, buffer radius and IFM measures provide detailed descriptions of patch-level potential connectivity, but do not necessarily scale up to landscape levels. However, if sufficient data are available to parameterize a stochastic patch-occupancy metapopulation model (e.g., the IFM ), one could then calculate the "metapopulation capacity" of the study system (Hanski and Ovaskainen 2000, Ovaskainen and Hanski 2001). Although connectivity is not directly quantified by metapopulation capacity, it may be a more useful quantity than landscape-level connectivity per se because it quantifies a landscape's potential to maintain a viable metapopulation over time. Observed Movement Rates: Individual Movement Data Data on the individual movements of organisms provide the most direct estimate of actual connectivity. Many methods exist for obtaining such data (Ims and Yoccoz 1997), but often these types of studies are too labor intensive to be conducted at even moderately large, let alone landscape, scales. Depending on the taxa in question, detailed tracking of the movement pathways of individual animals via radio telemetry or other methods (Gillis and Krebs 1999, 2000, Turchin 1998), mark-release-recapture studies (Southwood 1978, Sutherland 1996), or mass mark-recapture methods (where individuals do not have a unique marking) may be used. In addition, measurements of patch-level immigration or 54 colonization rates for unmarked animals can, by themselves, serve as a connectivity metric (van Langevelde 2000). This approach is difficult in practice, however, because immigration or colonization rates must be sufficiently high that useful data can be collected over a reasonable period of time. Despite the difficulty, many techniques for the direct measurement of movement can be applied to a variety of taxa, and these methods provide direct information about short-term dispersal. Alternatively, to quantify the extent of past dispersal over long time scales, metrics based on genetic data (e.g., Andreassen and Ims 2001) could be used. Although landscape-level estimates of actual connectivity are possible for wide-ranging species that can be radio tracked (e.g., Florida panthers (Meegan and Maehr 2002) the data-intensive nature of direct measurement methods will generally limit the spatial scales to which they can be applied. Still, in situations where movement data are already available or only a few habitat patches are of interest, quantifying emigration, immigration, or dispersal rates provides a detailed estimate of how well particular patches are connected in a fragmented landscape. MODIFICATIONS For many of the metrics discussed here, additional data, not included in the basic definition of the metric, can be incorporated to improve performance. The most common modification is weighting patch contributions to connectivity by area or some other proxy for population size. Such "area-informed" metrics 55 generally perform better than those that lack area considerations (Moilanen and Nieminen 2002, Bender et al. 2003, Tischendorf et al. 2003). Additionally, parameters that scale patch emigration or immigration according to patch area or population size can be used to capture some aspects of the dispersal behavior of species (Moilanen and Nieminen 2002). For example, for a given habitat area or population size, individuals of different species may not be equally willing to leave or enter habitat patches (Haddad 1999). Another commonly modified component of many connectivity metrics is the definition of interpatch distance. While it is the simple Euclidean distance most often used, other distance measures, such as least-cost movement pathways, can be used when appropriate (Bunn et al. 2000). Alternate movement pathways may be especially important to assess connectivity when landscape features such as rivers or mountains force organisms to disperse along pathways not well described by Euclidean distances (Dunham and Rieman 1999, Fagan 2002). In addition to the quantity of data, the effects of data quality on metric performance should also be considered. Such a discussion is beyond the scope of this paper, but Ruckelshaus et al. (1997, 1999) and Moilanen (2002) provide effective starting points. SCALE DEPENDENCE Two issues of scale dependence arise when considering connectivity. First, on which scale should connectivity be defined? Though several papers have debated this point (Tischendorf and Fahrig 2000, Tischendorf 2001a, 56 Moilanen and Hanski 2001), there is no evidence that connectivity should be limited to a particular spatial scale. This leads to the second issue: connectivity will change with spatial scale. How does one decide which scale is most appropriate for a particular problem? Clearly, the dispersal ability of the species imposes a relevant scale on the landscape (Wiens 1997), but dispersal ability is often unknown or poorly known. In such cases, explicitly calculating connectivity at a series of nested spatial scales and examining how connectivity changes as a function of scale is likely to provide a more robust picture of connectivity for the study area. Many of the approaches to connectivity detailed in this review have, at least to some degree, utilized this method. Tischendorf (2001b) showed that spatial pattern indices were generally better predictors of dispersal success when calculated at the landscape element (class) level than at the landscape level. The scale-area approach of Kunin (1998) and Fagan et al. (2002) is defined by a nested spatial scale methodology, scaling up from individual occurrences to the entire landscape. Similarly, graph theory naturally lends itself to such multi-scale analyses and allows the integration of patch-level and class- or landscape-level connectivity (Keitt et al. 1997, Urban and Keitt 2001). Though not purely a connectivity measure, metapopulation capacity accomplishes a similar scaling-up from patch to class or landscape levels by focusing on how landscape structure, which affects patch-level connectivity, influences population persistence for particular species (Hanski and Ovaskainen 2000, Ovaskainen and Hanski 2001). 57 Such multiscale methodologies could be used to look for connectivity thresholds (Keitt et al. 1997) or to assess the sensitivity of connectivity estimates to assumptions about the dispersal ability of the focal species. CONCLUSIONS Across the different connectivity metrics, a tradeoff exists between information content and data requirements (Figure 2.4). For example, the nearest neighbor measures and spatial pattern indices do not require extensive data to calculate, but provide only a crude estimate of structural connectivity. In contrast, buffer radius and IFM approaches provide very detailed estimates of potential connectivity at the individual patch level, but are extremely data-intensive. Likewise, the direct observation methods provide the only estimates of actual connectivity, but are, again, applicable mainly to small scales and are extremely data-intensive. Given the tradeoff between information content and data requirements, the graph-theoretic approaches may possess the greatest benefit to effort ratio for conservation problems that require characterization of connectivity at relatively large scales. These measures provide a reasonably detailed picture of potential connectivity, but have relatively modest data requirements. When habitat patches cannot be reliably delimited, the scale-area approach might be the only option. However, the relationship between scale-area slopes and actual connectivity needs to be better developed. Unfortunately, no all-purpose method exists for choosing which of the many connectivity metrics to use in addressing real-world problems. Future 58 research will undoubtedly illustrate which of these metrics perform best and which need to be left by the wayside. However, for many urgent conservation decisions, we do not have the luxury of waiting until a consensus is reached. Our goal in developing this classification system was to give non-theoreticians a starting point from which to choose appropriate connectivity metrics. Hopefully, knowledge of data requirements and informational detail, as well as the strengths and weaknesses of different approaches to connectivity, will allow practitioners to invest limited funds and efforts wisely when connectivity is used to evaluate alternative conservation strategies. 59 Table 2.1 A summary of the data-dependent classification framework for connectivity metrics. Connectivity Metrics Type of Connectivity / Level of Detail Habitat- Level Data Species- Level Data Methodology Nearest neighbor distance Structural Nearest neighbor distance Patch occupancy Patch-specific field surveys Spatial pattern indices Structural Spatially explicit None GIS / Remote sensing Scale-area slope Structural None Point- or grid- based occurrences Occurrence databases, presence/absence sampling Graph-theoretic Potential Spatially explicit Dispersal ability GIS / Remote sensing + dispersal studies Buffer radius, IFM Potential Spatially explicit including patch area Patch occupancy and dispersal ability Multi-year patch- specific field surveys or single year patch occupancy study with dispersal study Observed emigration, immigration or dispersal rates Actual Variable, depends on methodology Movement pathways or location- specific dispersal ability Track movement pathways (specific methods depend on study organism), mark-release- recapture studies 60 Figure 2.1 Photos of different types of habitat edges. (a) A pronounced edge in semi-arid grassland habitat of the Chiricahua Mountains, Arizona, induced by different grazing practices. Habitat edges like this represent semi-permeable barriers, disrupting the dispersal behaviors of some species but not others. Interspecific differences in edge responses are one reason why ecologists need to be alert to the species-specific nature of connectivity metrics. (b) A more complex landscape near Wuerzburg, Germany. Different species may have different perceptions about which landscape elements are usable. For example, some may be restricted to the forest fragments while others will move freely through forest as well as vineyards. 61 Figure 2.2 Schematic representations of the three types of connectivity. (a) Structural connectivity depends mainly on physical attributes of landscape elements, such as spatial proximity. Therefore the elements in the left column have higher structural connectivity than those in the right column. (b) Potential connectivity depends on physical attributes, but also on the dispersal ability of focal species. The red and blue bars represent measures of dispersal ability for two hypothetical species. If the distance between patches is greater than this measure of dispersal ability, the patches are not connected. Thus, the landscape on the left is connected for both species while the landscape on the right is connected for the blue species but not for the red species. (c) Actual connectivity is based on observed movement pathways. While factors considered in the other two classes of connectivity metrics certainly influence actual connectivity, movement must be observed or quantified. The left and right columns represent different observed pathways that would not necessarily be predicted by the structural or potential connectivity approaches. Thicker arrows indicate higher movement rates, and thus, higher actual connectivity. 62 Figure 2.3 A simplification of the spatial scales discussed in this paper, based on a recent landcover classification for Jamaica (Evelyn and Camirand 2003). The entire inset represents a 19 054.74-ha landscape scale. Eight landcover classes are represented within the landscape, as described by the legend. For example, closed broadleaf forest, in dark green, represents a single landcover class. An individual patch within the broadleaf forest class is outlined in red and highlighted with a red arrow. The blue dot highlighted by the blue arrow represents a hypothetical point occurrence of a focal species. Land classification provided by the Forestry Department of Jamaica. 63 Figure 2.4 Schematic representation of the tradeoff between information content and data requirements among connectivity metrics. Both information content and data requirements increase going from nearest neighbor measures to actual movement rates. The embellishments to the metrics mentioned in the "modifications" section may alter the position of various metrics in the hierarchy, but in general, the tradeoff between information content and data requirements holds. 64 Chapter 3 Neutral theory as a null hypothesis for species diversity in aggregated arthropod assemblages: a test with Cecropia petiole beetles Coauthored with Lawrence Kirkendall 65 ABSTRACT The study of arthropod assemblages that utilize discrete, ephemeral, patchily distributed resources has focused on how intraspecific aggregation mediates coexistence in the face of asymmetric competition, while paying little attention to species abundance. In contrast, neutral theory has sought to explain species abundance patterns in a range of assemblages by assuming equivalence among species or individuals. Recent work on neutral theory suggests that many non-neutral mechanisms that facilitate coexistence may also lead to ?functional equivalence? among species in an assemblage by minimizing interspecific differences among species. Here, we treat neutral theory as a null hypothesis for species abundance patterns and explore the possibility that intraspecific aggregation decouples interspecific competition from species abundance patterns in an assemblage of wood-boring beetles that breed in the fallen petioles of Cecropia insignis trees in Costa Rica. We use analytical methods to fit the neutral model to species abundance data from the Cecropia assemblage and simulation/randomization methods to account explicitly for both process and sampling errors when considering the model?s goodness of fit. We demonstrate that while species? spatial distributions are consistent with the aggregation model of coexistence, the species abundance distribution of the Cecropia assemblage deviates from the neutral model prediction in the direction of excess dominance. These results suggest that aggregation, though strong, does not allow species to act as functional equivalents, and that a model that incorporates non-neutral 66 mechanisms would be necessary to capture the complete pattern of species abundance in this community. 67 INTRODUCTION Neutral community theory has attracted much recent attention for its ability to predict species diversity patterns in some communities at some spatial scales (Hubbell 1997, 2001, Volkov et. al 2003, Alonso and McKane 2004, Chave 2004). The theory assumes that individuals, regardless of species, within a fixed size community are equivalent in their probabilities of reproducing, dying and migrating (Hubbell 1997, 2001). Local community structure arises as a consequence of random birth, death and migration events; a process termed ecological drift (Hubbell 2001). The neutral theory has been the focus of considerable debate over its assumptions and how accurate a description of the dynamics of natural communities it provides (Abrams 2001, Bell 2001, Fargione 2003, Harte 2003, McGill 2003, Ricklefs 2003, Volkov et al. 2003, Wootton 2005). Despite the controversy, the theory has been successful at providing a dynamic null hypothesis of community structure (Caswell 1976, Bell 2000, Bell 2001, Hubbell 2001, Etienne 2005). In this role, failure to reject the null hypothesis of neutrality does not necessarily imply that the neutral theory is ?true? for the test community. Instead, it suggests that asymmetric species interactions are not required to explain the observed community structure. Deviations from the theory can then be used to help guide future research efforts. Neutral theory in ecology is still at an early stage, and both its domain of applicability and the methods used to test it are still being defined (McGill 2003, Volkov et al. 2003, Alonso and McKane 2004, Chrisholm and Burgman 2004, Etienne and Olff 2004, Hubbell and Borda-de-Agua 2004, Etienne 2005, Wootton 68 2005). To date, most tests of the theory have focused on space-limited assemblages such as closed canopy tree, intertidal fouling and coral reef communities (Hubbell 2001, McGill 2003, Volkov et al. 2003, Wootton 2005). Despite this focus, the analysis of other types of communities could also benefit from using neutral theory as a null hypothesis for community structure. Assemblages of arthropods that breed in discrete, ephemeral, patchily distributed resources (hereafter DEP systems) are one such class of communities. Examples of DEP systems include flies that breed in carrion, fruit or fungus and beetles that breed in dung or wood (Atkinson 1985, Ives 1988, Hanski and Camberfort 1991, Shorrocks and Bingley 1994, Jordal and Kirkendall 1998). In DEP systems, the focus has been on the mechanisms that mediate coexistence of species in the face of potentially strong and asymmetrical resource competition. Relatively little attention has been paid to species abundance patterns in these assemblages (but see Krijer and Sevenster 2001, Warren et al. 2003), which in our view, is due at least partly to the lack of a formal theory of abundance. The persistent focus on the mechanisms that promote coexistence also provides a reasonable conceptual basis for using neutral theory as a null hypothesis in DEP assemblages. The aggregation model of coexistence has been developed to explain how multiple species can coexist in DEP systems with little or no apparent resource partitioning (Atkinson and Shorrocks 1981, Hanski 1981, Atkinson and Shorrocks 1984, Ives and May 1985, Ives 1991, Sevenster 1996, Hartley and Shorrocks 2002). The aggregation model builds from the empirical observation that species within DEP systems tend to have independent, 69 intraspecifically aggregated distributions of individuals (larvae) over resource units (Atkinson and Shorrocks 1981). In other words, individuals tend to occur and compete with conspecifics more frequently than they do with heterospecifics. Because most competition in these assemblages occurs at the larval stage, within resource units, intraspecific aggregation increases the strength of intraspecific competition relative to that of interspecific competition and therefore facilitates the coexistence of unequal competitors under a wide range of conditions (Atkinson and Shorrocks 1981, Atkinson and Shorrocks 1984, Ives 1988, Sevenster 1996, Hartley and Shorrocks 2002). The reason that the aggregation model provides a reasonable conceptual basis for testing neutral theory in DEP systems is that recent work on neutral theory suggests that other coexistence mechanisms such as niche partitioning (Hubbell 2001, Ch. 9) and life-history tradeoffs (Chave et al. 2002, Chave 2004) may facilitate neutral-like dynamics by minimizing the functional differences among species. Thus intraspecific aggregation, a nearly ubiquitous feature of DEP systems, may also facilitate neutral-like dynamics by reducing the importance of asymmetrical species interactions (e.g., Shorrocks et al. 1984). Veech et al. (2003) have made a similar suggestion based on their meta-analysis of species richness patterns in aggregated arthropod assemblages. By assuming that intraspecific aggregation renders species within DEP systems ?ecologically equivalent?, we can then use neutral theory as a null hypothesis and test for the signature of asymmetric competition (i.e., systematic deviations from the predictions under neutrality) in empirical species abundance distributions from DEP assemblages. 70 Here, we test species abundance predictions of Hubbell?s (2001) neutral theory in an assemblage of wood boring beetles that breed in the fallen petioles of Neotropical Cecropia trees (Jordal and Kirkendall 1998), using both analytical and simulation/randomization methods. First we perform a standard aggregation analysis to demonstrate that species? spatial distributions in the Cecropia assemblage are consistent with the aggregation model. Next, we fit Hubbell?s (2001) two-scale, spatially implicit neutral model to the species abundance data. To test the quantitative predictions of any theoretical model against empirical data, the ways in which uncertainty enter in to the data must be taken into account (Hillborn and Mangel, 1997). To do this, we place the problem in a process/sampling error framework (Hillborn and Mangel 1997). First, we use the analytical solution to Hubbell?s neutral model developed by Volkov et al. (2003) to estimate the relevant parameters. Unfortunately, the mean-field approach of Volkov et al. (2003) does not allow estimation of the process error of the underlying stochastic model. To estimate this source of error, we simulate the neutral model with the parameters estimated from fitting the analytical solution. To incorporate sampling error, we note that individual resources (i.e., petioles, fruits, dung pats, etc?) are the usual sampling units in DEP systems, and that individuals of each species tend to be non-randomly distributed over those resource units (i.e., there is intraspecific aggregation). As a result, randomly sampling resource units is different than randomly sampling individuals from a population, and we consider explicitly the error that results from this type of sampling. We then use these two sources of error to estimate the joint (process + 71 sampling) error distribution expected for each species under the neutral model and that species? spatial distribution. We use these joint error distributions to assess the fit of the model to the empirical species abundance data both species-by- species and for the entire abundance distribution. METHODS Study System Fieldwork was conducted in the lowland tropical wet forest at La Selva Biological Station, Heredia Province, Costa Rica. Cecropia trees grow throughout the Neotropics and are widespread pioneers of light gaps, forest edges and disturbed areas. These trees shed leaves throughout the year, and the large (20-120cm), woody petioles provide a continually available but spatially clustered resource base for a guild of bark beetles (Beaver 1979, Jordal and Kirkendall 1998). Beetles breeding in Costa Rican Cecropia petioles comprise 36 known species from two families (Coleoptera: Curculionidae and Cerambycidae) (Jordal and Kirkendall 1998). Cecropia petiole beetles complete the bulk of their lifecycles inside the petioles, emerging after pupation to disperse, and search for fresh petioles to colonize. Individual petioles generally only support one generation of beetles before decomposition makes them inhospitable to beetle larvae. Jordal and Kirkendall (1998) and Jordal (1998) provide detailed descriptions of the basic biology and taxonomy of beetles in this assemblage. 72 Sampling Sixteen Cecropia insignis trees were chosen in the 1168 ha area of primary forest at La Selva in October 2002. Sites were selected to maximize coverage of the entire area of primary forest, within the constraints of accessibility. Petioles taken in the sample represent a subset of the decaying petioles at each tree that were old enough to have been colonized by beetles but not so old that beetles would have already completed development and emigrated. Petioles from which beetles have already emerged can be distinguished because they tend to be partially decomposed and have extensive larval tunneling and exit holes left by emerging beetles. Of the suitable petioles around each tree, a sample of 10 was randomly selected from the ground or the vegetation within 2 m of the ground. Laboratory Rearing All sampled petioles were immediately brought into the laboratory for rearing. Petioles from three of the trees were placed in sealed plastic bags in an ambient laboratory. Petioles from the remaining 13 trees were placed into rearing chambers consisting of PVC tubes with mesh covers on one end and an emergence trap head on the opposite end. The rearing tubes were placed on racks in a laboratory under ambient conditions. Beetles were allowed two months to develop and emerge. After this period all emerged beetles were collected and a subset of the petioles was completely dissected to ensure that the emergence traps were capturing most of the beetles developing in the petioles. Beetles were then preserved in 70% EtOH and identified to species. Rearing in bags and tubes 73 produces similar results (J. Calabrese unpubl. data), but the tubes require much less maintenance during the rearing period. Analysis Aggregation To quantify the degree of intraspecific aggregation for species in this assemblage, we employed the mean crowding index (Ives 1988), which is defined B x = ? x 2 ? x 2 ? 1 ? x Eq.3.1 where ? and ? 2 are the mean and variance, respectively, of the number of individuals of the focal species over the resource units. The B index has range (-1, +?) and measures the degree to which individuals are aggregated over petioles relative to a Poisson (random) distribution. A Poisson distribution has B = 0, while B > 0 indicates intraspecific aggregation and B < 0 signifies regular dispersion. A B value of 1 is a 100% increase in intraspecific aggregation compared to a Poisson distribution of individuals over the petioles. For each focal species, we calculate two B statistics. We follow the approach of Krijger and Sevenster (2001) and define B x to be the B value for focal species x, while B y is the B value for the rest of the community excluding focal species x, treated as one ?super-species? (Shorrocks and Rosewell 1986, 1987, Sevenster and VanAlphen 1996, Wertheim et al. 2001, Krijger and Sevenster 2001). To quantify interspecific aggregation, we used the C index (Ives 1988), defined C xy = Cov xy ? x ? y Eq.3.2 74 where Cov xy is the covariance in the numbers of individuals of focal species x and ?super-species? y over the resource units. C xy measures the proportional increase in the mean number of heterospecifics (?super-species? y) that occur with each individual of focal species x in a petiole relative to the number expected if both x and y had independent Poisson distributions of individuals over petioles. Values of C are interpreted in the same way as those for B, with respect to heterospecifics. Sevenster (1996) has shown that B and C can be combined into T xy , a criterion that defines the necessary and sufficient conditions for the persistence of species x in competition with species y under the aggregation model. We follow the assemblage-level interpretation of T xy (Krijger and Sevenster 2001) for species x in competition with the rest of the assemblage, ?super-species? y, which is defined T xy = C xy +1 B y +1 <1 Eq.3.3 Values less than one predict long-term persistence of focal species x, while values greater than one indicate that species x will not be able to persist in competition with the rest of the assemblage. We present values of the B, C and T indices only for species represented by at least 15 individuals, because estimates for extremely rare species will not be robust. Parameter Estimation Volkov et al. (2003) derived an analytical solution for the distribution of abundance in a local community of Hubbell?s (2001) neutral model. The 75 analytical solution predicts the equilibrium number of species ? n with abundance n in a local community of size J (number of individuals, summed across species) under limited dispersal from the metacommunity for n = 1, 2, ?, J. The equilibrium local species abundance distribution is ? n =? J! n!(J ? n) ?(? ) ?(J +? ) ?(? ) ?(1+ y) ?(J ? n +? ? y) ?(? ? y) exp(?y? /? )dy 0 ? ? Eq.3.4 where ? = 2J M ? is the biodiversity parameter, J M is the size of the metacommunity, ? is the speciation rate, m is the per death probability of immigration from the metacommunity, ? (z) is the Gamma function and ? = m(J ?1) 1? m . We fit Eq. 3.4 to the species abundance data from all 16 trees pooled together to estimate values of the two free parameters,? and m, via maximum likelihood methods. Defining ? = (?*,m*) as a set of candidate values for the two free parameters of Eq. 3.4, the log likelihood function is L ? = q n ln E n (?) q n ? ? ? ? ? ? ? E n ?()? q n ? ? ? ? ? ? ? ? ? ?n=1 J ? Eq.3.5 where q n is the observed number of species with abundance n, E n ? ( ) is expected number of species with abundance n given parameter set ? (Hubbell 2001). We dealt with zeros in either q n or E n ? ( ) in the same way as McGill (2003), by ignoring the first term of the log likelihood function in these cases. We assessed the goodness of fit of the neutral model in several ways. First, we plotted the species abundance data as a Preston histogram (Preston 1948) and compared the goodness of fit of the neutral model predictions given the maximum likelihood 76 estimates via Chi-square. Next, we obtained 95% confidence intervals on the point estimates for both parameters from their likelihood profiles (Hillborn and Mangel 1997). Finally, we conducted a more rigorous assessment of goodness of fit by adopting a process/sampling error framework described below. Process and Sampling Error Neither of the first two methods for evaluating the fit of the neutral model to the data account for process and sampling error. Indeed, a shortcoming of the Volkov analytical solution is that it provides no estimate of process error. Given that the neutral model is inherently stochastic, an estimate of the process error is required to assess the fit of the model to the data. Similarly, sampling error will also affect the comparison between model and data, and therefore must be accounted for as well. The impacts of these sources of error can be assessed by considering how they enter into the data. Assuming the neutral model is true, process error can be represented as N act,i = E(N act,i |?, m, J)+W i Eq.3.6 where N act,i is the actual or ?true? abundance of the i th ranked species (i = 1, 2,?,S, rank one being the most common) in the local community, E(N act,i |?,m,J) is the expected abundance of the i th ranked species under the neutral model given ? , m and J, and W i is a random variable describing the process error. Because we lack an analytical expression for W, we simulated a neutral community using the observed community size and the maximum likelihood estimates of ? and m according to the methods described by Hubbell (2001) and Hubbell and Borda-de-Agua (2004). To obtain reliable estimates the 77 frequency distributions of the N act,i for each species i in the local community, we ran 10,000 independent simulations, each for 10,000,000 birth/death cycles (with 1 death per cycle). To incorporate sampling error, we first demonstrate that the negative binomial distribution (NBD) is a good descriptor for the spatial distributions of each species represented by at least 15 individuals in the sample. To do this, we fit the NBD to the frequency distribution of number of individuals per petiole using the method of moments (Bliss and Fisher 1953, Hilborn and Mangel, 1997) and assessed the goodness of fit via Chi-square. We then note that the abundance of each species in the sample represents the sum of h independent draws from its NBD distribution of individuals over petioles with mean ? and aggregation parameter k, where h, in this case, represents the number of petioles sampled. The distribution of the sum of h independent negative binomial random variables is also negative binomial with mean ? ? = h?? and aggregation parameter ? k = h?k (Anscombe 1950). Therefore the frequency distribution of observed abundances for a given species that would be produced by taking repeated, independent samples of h petioles would be negative binomial, as would the sampling error probability distribution. This can be represented as N obs,i = N act,i +U i Eq.3.7 where U ~ NBD( ? ? , ? k ) and N act,i is also a random variable, the distribution of which is resolved by simulating the neutral model as described above. Thus, the distribution of the N obs,i is, for each species, the frequency distribution of the abundances that would be observed by repeatedly sampling h petioles under both 78 process and sampling error, assuming the neutral model is true and that the species have negative binomial distributions of individuals over resource units. To estimate the shape of this distribution, we used for each species i, ? ? i = N act,i and ? k i =160?k obs,i , and generated one negative binomial random number for each of the 10,000 values of N act,i (generated by simulating the neutral model). To assess the fit of the neutral model to the data when both process and sampling error are accounted for, we estimated the 95% confidence intervals for each ranked species i from the frequency distribution of N obs,i . We then plotted the mean abundance predicted by the neutral model, E(N act,i |?,m,J), from the simulations with these confidence intervals against the empirical abundance data as a dominance-diversity curve (sensu Whittaker 1965). This method allows us to assess the fit of each species in the sample relative to the model predicted value taking into account the process error in the neutral model and sampling error caused by the unique spatial distribution of each species. Plotting the data in this manner also allows us to determine if the observed data deviate systematically from the neutral model predictions. Finally, to quantify and test for the significance of the overall deviation of the data from the dominance-diversity curve predicted by the neutral model, we adapt the neutral model deviation statistic originally described by Ewens (1972) and introduced to ecology by Caswell (1976), defined V = ?H ? E( ?H ) ?( ?H ) Eq. 3.8 79 where H? is the Shannon-Weiner diversity index, E( ? H ) is the expected value of H? under the neutral model and ?( ? H ) is the standard deviation of E(H?) under neutrality. Ewens (1972) derived this statistic for a different neutral model (the infinite alleles model from population genetics), and thus the analytical formulae he used to calculate E( ? H ) and ?( ? H ) cannot be used here. Instead, we calculated H? for each of the 10,000 communities (each community being a realization of N act,i , for i = 1, ?, S) used to develop the joint error distributions. From the 10,000 H? values thus calculated, we can estimate E( ? H ) and ?( ? H )and then substitute these values into Eq. 3.8. This is essentially recalibrating V for Hubbell?s neutral model and the sampling error caused by the observed spatial distributions of the beetle species in the sample. As in Caswell?s (1976) treatment, V > 0 signifies excess evenness in the distribution compared to the neutral model, V = 0 perfect agreement with the neutral model and V < 0 indicates excess dominance. As V is a normalized sum of random variables, by the Central Limit Theorem we expect the distribution of V to be approximately normal with mean 0 and standard deviation 1 (Ewens 1972, Caswell 1976, Hilborn and Mangel 1997). To estimate the shape of the distribution of the V statistic we bootstrap resampled the 10,000 H? values, 50,000 times, calculating the V statistic each time and plotted the frequency distribution of the V values (Hilborn and Mangel 1997). From this distribution, we assessed the statistical significance of the departure of the data from the model predicted abundance distribution. 80 RESULTS Aggregation The sample produced 2013 individual beetles distributed across 19 species (Table 3.1). Of the 160 sampled petioles, 79% produced at least one beetle. All species represented by at least 15 individuals (nine species) show strong intraspecific aggregation of individuals over petioles (B x , Table 3.1). All of these species satisfy the necessary condition for persistence, that intraspecific aggregation is stronger than interspecific aggregation when all heterospecific individuals are treated as one super species (B x >C xy , Table 3.1). Furthermore, all nine of these species satisfy both the necessary and sufficient conditions for long- term persistence, T xy < 1, which suggests that aggregation is likely an important mechanism in mediating competitive relationships and facilitating coexistence in this assemblage, though two species are marginal in this regard (Table 3.1). The fit of the NBD to the frequency distributions of the number of individuals of each of the nine most abundant species per petiole was excellent (Chi-square tests, Table 3.1). These results suggest that the NBD is a good descriptor of the spatial distributions of these species, justifying our use of the NBD to model sampling error. Neutral Model Parameter Estimation The maximum likelihood estimates for the parameters ? and m were 3.04 and 0.47, respectively. When binned as a Preston histogram (sensu Preston 1948), the species abundance data do not differ significantly from the analytical solution with the maximum likelihood estimates of ? and m, but the fit is clearly rough 81 (? 11df 2 = 5.953, p = .653, Figure 3.1). Likelihood profiles for each parameter reveal reasonable confidence for ?, but not for m (Figure 3.2a, b). The likelihood profile for m suggests that it is not significantly different than 1, indicating that there is little evidence of dispersal limitation in this beetle assemblage under the neutral theory (Figure 3.2b). Joint Error Distributions Table 3.1 presents the means, ? ? , and aggregation parameters, ? k , of the distributions of N obs,i for all 19 species. For the nine most abundant species, ? k , the clumping parameter of the negative binomial distribution was estimated as ? k i =160?k obs,i , where k obs,i was estimated from the data via the method of moments. For the remaining 10 species that were too rare to allow accurate estimation of k obs,i , we assigned the average value of ? k from the top nine species (Table 3.1). For species of higher abundance, we found that the distribution of the N obs tended to approximate a normal distribution, becoming increasingly right- skewed as abundance declines (Figure 3.3). The 95% confidence intervals for each species were estimated from these frequency distributions by identifying the values of N obs for which 2.5% (250) of the 10,000 values fell to the left (lower bound) and for which 2.5% of the values fell to the right (upper bound). Model Fit Under Process and Sampling Error For all species except the rank one species, Scolytodes blandfordi, the empirical dominance diversity curve falls within the 95% confidence intervals estimated from the joint process/sampling error distribution for each species (Figure 4). Scolytodes blanfordi is significantly more common (p = .0209) than 82 would be expected given the neutral model predictions and sampling error caused by its aggregated spatial distribution (Table 1). Furthermore, except for the four rarest species, the overall dominance diversity curve is steeper than that predicted by the neutral model, suggesting that there may be some effect of asymmetric species interactions on the species abundance distribution. Finally, when the entire species abundance distribution is compared to the neutral model predictions using the modified V statistic, the systematic deviation towards higher dominance becomes apparent. For the data, V = -3.554 and the probability that this value came from the bootstrapped distribution of the V statistic under the neutral model and negative binomial sampling error is < .0003 (Figure 5). Thus overall, we reject the neutral model for the Cecropia assemblage. DISCUSSION Here, we tested the ability of Hubbell?s neutral model to explain species abundance patterns in a DEP assemblage. This test was performed based on the rationale that strong intraspecific aggregation, which clearly exists in the Cecropia assemblage, might cause species to behave as functional equivalents, thus allowing ecological drift to drive species abundance patterns. We treated the neutral theory as a null hypothesis for species abundance patterns that can be tested, at least partially, by fitting a neutral model to the species abundance data. Furthermore, we sought to refine methods for testing the neutral theory in situations where assumptions of random sampling of individuals are clearly violated because of the spatial distributions of the focal species. We found that, 83 whereas a standard Chi-square test of the fit of the analytical solution of the neutral model to the data does not allow rejection of the model, a more comprehensive model evaluation framework that incorporates both process and sampling error suggests that the neutral model is not a good descriptor of the Cecropia assemblage. We used a combined approach that employed both analytical and simulation/randomization methods to assess the fit of the neutral model to the data. The analytical solution is preferable for parameter estimation because it circumvents issues of convergence inherent in a simulation approach (see Volkov et al. 2003). However, the analytical solution of Volkov et al. (2003) does not allow assessment of the variability inherent in the underlying stochastic model that generates the mean predictions. To assess this source of variability we used the simulation methods described by Hubbell (2001) and Hubbell and Borda-de- Agua (2004). Furthermore, because intraspecific aggregation is such a prominent feature of DEP systems, we sought to explicitly account for the sampling error resulting from randomly sampling resource units that contain non-randomly distributed individuals (i.e., species have negative binomial, and not Poisson, distributions of individuals over resource units). Combining these two sources of error into an overall estimate of error allowed us to put the deviations of the empirical data from the model predictions into the proper context. To our knowledge, this is the first test of neutral theory that has used an analytical solution to estimate parameters and has accounted for both process and sampling error when assessing the fit of the model to the data. 84 The distributions of N obs,i represent, for each species i, an estimate of the frequency distribution of abundance that would be observed by taking repeated, independent samples of 160 petioles under the joint action of process and sampling error. The shapes of these distributions (Figure 4), tending towards normality at high abundance and becoming increasingly right skewed at low abundance, are expected based on the shapes of the constituent process and sampling error distributions from which they were generated. The process error distributions show the same overall patterns as mean abundance decreases but to a lesser extent (Results not shown). For the sampling error distributions, this behavior is expected because they are negative binomial, and as the value of the aggregation parameter ( ? k ) increases, a negative binomial approximates a Poisson distribution (Anscombe 1950). Formally, the NBD becomes asymptotically equivalent to the Poisson in the limit where k ?? (Anscombe 1950), but for practical purposes, k > 10 is more or less Poisson-like. As the mean of a Poisson becomes large, it comes to approximate a normal, hence higher abundance species have sampling error distributions that are roughly symmetrical while the rarer species have Poisson-like right skewed sampling error distributions. These results also demonstrate that it is instructive to consider the fit of the neutral model on a species-by-species basis as well as for the trend in the entire abundance distribution. For example, Figure 4 demonstrates that the most common species is too common to be consistent with the neutral model. Furthermore, the rest of the assemblage tends to fall on the lower end of the 95% joint confidence intervals. That the four rarest species fall on or above the 85 dominance-diversity curve predicted by the neutral model may reflect nothing more than the fact that all species in the sample have aggregated spatial distributions and therefore occur in clumps. Figure 4 also clearly shows the trend in the empirical data towards excess dominance in the entire distribution, with the exception of the rarest species. However, from Figure 4 alone, we would not necessarily conclude that the empirical abundance distribution is inconsistent with the neutral theory because 18 of 19 species do fall within the 95% confidence intervals. To assess the fit of the model to the data across the entire dominance- diversity curve, we must quantify the total deviation between the two curves in their entireties, which is exactly what the V statistic does. The reason that we need to assess the fit for the entire abundance distribution is that, while the sampling error acts independently from species to species, process error does not. This is because the assumed constant community size in the neutral model causes negative correlations among species abundances (zero-sum rule, Hubbell 2001). Said another way, for one species to increase its abundance in a neutral community, another must decrease. Some fingerprint of these negative correlations will be retained in the joint error distributions, and thus we must consider how likely it is to observe a given set of species abundances together (and not just the individual abundances of each species) under the specified process and sampling error assumptions. Though Ewens (1972) originally derived the V statistic under the infinite alleles model and an assumption of perfect random sampling, it can be recalibrated to measure deviation under a different 86 neutral model and under different types of sampling. The recalibration of V presented here is specific to the details of this particular problem, but the approach is general and could easily be adapted to other situations. The advantage of such an overall deviation statistic is that it quantifies both the direction and magnitude of the deviation from neutrality and provides a formal statistical test for the entire species abundance distribution. The reason that V obs is so extreme in the present example is that, while 18 of 19 species fall within the 95% joint confidence intervals on a species-by-species basis, the probability of having a single assemblage that is so strongly skewed towards high dominance (large, negative V value) is extremely small. By rejecting the neutral model as a null hypothesis for the Cecropia assemblage, we accept the vaguely defined alternative hypothesis that some species, particularly S. blandfordi, retain some competitive advantage. The vagueness of the alternative hypothesis highlights the current state of the art in testing neutral theory and emphasizes the need for non-neutral models that incorporate key elements from neutral theory such as demographic stochasticity. Tilman?s (2004) stochastic niche model is a step in this direction, but is inappropriate for the present application because it focuses on nutrient competition in plant communities. Neutral theory clearly gets the species abundance predictions ?in the ballpark?, but explaining the excess dominance in the empirical abundance distribution will require the addition of non-neutral mechanisms. Such a range of neutral to progressively more non-neutral models will hopefully be developed in the future. Both Hubbell (2001) and Chave (2004) 87 have called for such models to reconcile the so-called ?dispersal-assembled? and ?niche-assembled? perspectives in community ecology, but the development of these models awaits further theoretical advances. Unfortunately, not much is currently known about S. blandfordi or the rest of the Cecropia petiole beetles (but see Jordal and Kirkendall 1998, Jordal 1998). From what is known, S. blandfordi does not appear to have a clear advantage over other species in terms of fecundity, tolerance to variable environmental conditions or other factors (Jordal and Kirkendall 1998). Recently, Fargione et al. (2003) and Wootton (2005) have highlighted the importance of testing predictions of the neutral theory experimentally whenever possible. A forthcoming paper will examine experimentally the abilities of species in the Cecropia assemblage to disperse and colonize new petioles. The colonization process is crucial in the Cecropia system because the petioles decay rapidly under field conditions, thus limiting the window of opportunity for successful brood development. Clear differences among species in colonization ability would violate the equivalence assumption and may help to explain the deviation from the neutral model predictions. Though we have rejected the neutral model based on the available evidence in this example, we stress that it should be tested on a wider range of DEP assemblages in the future to: 1) determine if it can be rejected in DEP systems in general, and 2) quantify the direction and magnitude of deviations from the neutral model predictions to assess the nature and strength of any non- neutral mechanisms that may be influencing community structure. Such tests 88 would likely lead to progress on a theory of abundance for DEP assemblages by elucidating when and to what extent non-neutral mechanisms are required to explain abundance patterns in DEP systems. 89 Table 3.1 Abundances and spatial distribution statistics for the 19 species present in the sample. The ? ? values are the abundances expected under the neutral model and negative binomial sampling error. The ? ? i are rounded to the nearest integer here, but are presented in decimal form in figure 3. The ? k values represent k obs values estimated from fitting a negative binomial to each species? distribution over sampled petioles multiplied the number of sampled petioles (i.e., ? k = k obs *160). Species with abundance < 15 were assigned the average ? k of the species with abundance ?15. The degrees of freedom (DF) of the Chi-square tests vary among species because different binnings were convenient for each species. The binnings do not strongly affect the goodness of fit. Species Abundance ? ? ? k ? 2 DF p B x C xy T xy Scolytodes blandfordi 1133 815 35.408 0.221 23 >.995 4.490 0.630 0.536 Xylosandrus morigerus 323 418 20.522 0.079 15 >.995 7.745 0.080 0.287 Coccotrypes cyperi 192 251 43.674 0.187 11 >.995 3.635 0.317 0.371 Lechriops LE132 103 162 9.444 0.038 11 >.995 16.826 2.120 0.998 Araptus costaricensis 97 107 34.656 0.074 10 >.995 4.578 0.350 0.403 Scolytodes atratus 42 76 19.644 0.010 5 >.995 8.070 0.658 0.514 Scolytodes caudatus 36 52 11.285 0.082 8 >.995 14.062 0.544 0.479 Eulechriops SL07 24 37 125.802 0.009 2 >.995 1.222 0.508 0.470 Araptus laevigatus 15 27 2.412 0.123 8 >.995 65.844 1.987 0.944 Coccotrypes advena 12 20 33.650 N/A N/A N/A N/A N/A N/A Hypothenemus sp006 8 14 33.650 N/A N/A N/A N/A N/A N/A Pseudolechriops PL05 8 11 33.650 N/A N/A N/A N/A N/A N/A Cryptorhynchus honestus 5 8 33.650 N/A N/A N/A N/A N/A N/A Araptus sp404 3 6 33.650 N/A N/A N/A N/A N/A N/A Cerambycidae sp01 3 4 33.650 N/A N/A N/A N/A N/A N/A Conoderinae sp07 3 3 33.650 N/A N/A N/A N/A N/A N/A Curculionidae sp01 2 2 33.650 N/A N/A N/A N/A N/A N/A Scolytodes parvulus 2 2 33.650 N/A N/A N/A N/A N/A N/A Eulechriops SL01 2 1 33.650 N/A N/A N/A N/A N/A N/A 90 Figure 3.1 The maximum likelihood fit of the analytical solution to the species abundance data. The curve represents the analytical solution and the bars are empirical data. Data are binned in Preston-type doubling classes of abundance. The maximum likelihood parameter estimates for ? and m are 3.04 and 0.47, respectively. The probability of obtaining this fit by chance is 0.347, suggesting the neutral model fits rather poorly, however the difference is not statistically significantly (? 11df 2 = 5.953, p = .653). 91 Figure 3.2 Likelihood profiles (dashed curves) and 95% confidence intervals for each parameter. Confidence limits occur at the intersections of the dashed likelihood profiles and the solid horizontal lines. Every value of the focal parameter whose likelihood falls on or below the solid lines is in the 95% confidence interval. The stars on the x-axis of each plot highlight the maximum likelihood value of the focal parameter. Panel A represents the likelihood profile for ? taken at the maximum likelihood estimate for m (0.47). The profile indicates reasonable support for ?. Panel B is the profile for m taken at the maximum likelihood estimate for ? (3.04). Based on the likelihood profile, m is not significantly different than one, indicating that there is little evidence of dispersal limitation in the beetle assemblage under the neutral theory. 92 Figure 3.3 Dominance-diversity curves for the simulated neutral model predictions (dashed line) and the beetle data (solid line). The error bars represent the 95% joint confidence intervals on the neutral model predicted abudances. All species except for the rank one species, Scolytodes blandfordi, fall within the joint confidence intervals. The overall trend in the curve is towards increased dominance relative to the neutral model predictions. 93 Figure 3.4 Frequency distributions of N obs,i under the joint action of process and sampling error for four species. Panels A, B, C and D represent the distributions for the Rank 1, 5, 10 and 15 species respectively. The distributions become progressively more right skewed as mean abundance decreases. 94 Figure 3.5 Bootstrapped frequency distribution of the modified V statistic. Results are based on 50,000 bootstrap resamplings of the 10,000 communities used to generate the joint error distributions. As expected, the distribution of V is approximately normal and has mean ? 0 and standard deviation?1. The vertical arrow indicates the position of V obs in the distribution. The probability of obtaining V obs = -3.554 from this distribution is < .0003. Thus we reject the neutral model for the beetle data. The negative value of V indicates strong deviation in the direction of excess dominance. 95 LITERATURE CITED Abrams, P.A. 2001. The unified neutral theory of biodiversity and biogeography. Nature 412: 858-859. Alcock, J. 1996. Male size and survival: The effects of male combat and bird predation in Dawson's burrowing bees, Amegilla dawsoni. Ecological Entomology 21(4): 309-316. Alcock, J. 1997. Small males merge earlier than large males in Dawson's burrowing bee (Amegilla dawsoni) (Hymenoptera: Anthophorini). Journal of Zoology 242: 453-462. Alcock, J. 1999. The nesting behavior of Dawson's burrowing bee, Amegilla dawsoni (Hymenoptera: Anthophorini), and the production of offspring of different sizes. Journal of Insect Ecology 12(3): 363-384. Allee, W.C., A.E. Emerson, O. Park, T. Park, and K.P. Schmidt. 1949. Principles of animal ecology. W.B. Saunders, Philadelphia, PA. Alonso, D., and A.J. McKane. 2004. Sampling Hubbell's neutral theory of biodiversity. Ecology Letters 7: 901-910. Andreassen, H.P., and R.A. Ims. 2001. Dispersal in patchy vole populations: Role of patch configuration, density dependence, and demography. Ecology 82: 2911-2926. Anscombe, F.J. 1950. Sampling theory of the negative binomial and logarithmic series distributions. Biometrika 37: 358-382. Anstett, M.C., G. Michaloud, and F. Kjellberg. 1995. Critical population size for fig/wasp mutualisms in a seasonal environment: effect and evolution of the duration of female receptivity. Oecologia 103: 453-461. Atkinson, W.D., and B. Shorrocks 1981. Competition on a divided and ephemeral resource: a simulation model. Journal of Animal Ecology 50: 461-471. Atkinson, W.D., and Shorrocks, B. 1984. Aggregation of larval Diptera over discrete and ephemeral breeding sites ? the implications for coexistence. American Naturalist 124: 336-351. Augspurger, C. K. 1981. Reproductive synchrony of a tropical shrub: experimental studies on effects of pollinators and seed predators on Hybanthus prunifolius (Violaceae). Ecology 62(3): 775-788. 96 Beaver, R.A. 1979. Non-equilibrium 'island' communities: a guild of tropical bark beetles. Journal of Animal Ecology 48: 987-1002. Bell, G. 2000. The distribution of abundance in neutral communities. American Naturalist 155: 606-617. Bell, G. 2001. Neutral macroecology. Science 293: 2413-2418. Bender, D.J, L. Tischendorf, and L. Fahrig 2003. Using patch isolation metrics to predict animal movement in binary landscapes. Landscape Ecology 18(1): 17-39. Berec, L., D.S. Boukal, and M. Berec. 2001. Linking the Allee effect, sexual reproduction, and temperature-dependent sex determination via spatial dynamics. American Naturalist 157: 217-230. Bliss, C.I., and R.A. Fisher. 1953. Fitting the negative binomial distribution to biological data ? note on the efficient fitting of the negative binomial. Biometrics 9: 176-200. Boukal, D.S., and L. Berec. 2002. Single-species models of the Allee effect: extinction boundaries, sex ratios, and mate encounters. Journal of Theoretical Biology 218: 375-394. Bullock, S.H., and K. S. Bawa. 1981. Sexual dimorphism and the annual flowering pattern in Jacaratia dolichaula (D. Smith) Woodson (Caricaceae) in a Costa Rican rain-forest. Ecology 62: 1494-1504. Bunn, A.G., D.L. Urban, and T.H. Keitt. 2000. Landscape connectivity: a conservation application of graph theory. Journal of Environmental Management 59: 265-278. Cabeza, M. 2003. Habitat loss and connectivity of reserve networks in probability approaches to reserve design. Ecology Letters 6: 665-672. Cantwell, M.D., and R.T.T. Forman. 1993. Landscape graphs: ecological modeling with graph theory to detect configurations common to diverse landscapes. Landscape Ecology 8: 239-255. Carvalho, M.C., P.C.D. Queiroz, and A. Ruszczyk. 1998. Protandry and female size-fecundity variation in the tropical butterfly Brassolis sophorae. Oecologia 116: 98-102. Case, T.J. 2000. An illustrated guide to theoretical ecology. University Press, Oxford. 97 Caswell, H. 1976. Community structure: a neutral model analysis. Ecological Monographs 46: 327-354. Chave, J., H.C. Muller-Landau, and S.A. Levin 2002. Comparing classical community models: theoretical consequences for patterns of diversity. American Naturalist 159: 1-23. Chave, J. 2004. Neutral theory and community ecology. Ecology Letters 7: 241- 253. Chesson, P. 1994. Multispecies competition in variable environments. Theoretical Population Biology 45: 227-276. Chisholm, R.A., and M.A. Burgman. 2004. The unified neutral theory of biodiversity and biogeography: Comment. Ecology 85: 3172-3174. Clark, J.S., M. Lewis, and L. Horvath. 2001. Invasion by extremes: Population spread with variation in dispersal and reproduction. American Naturalist 157: 537-554. Cooper, W.S., and R.H. Kaplan. 1982. Adaptive "coin-flipping": a decision- theoretic examination of natural selection for random individual variation. Journal of Theoretical Biology 94: 135-151. Courchamp, F., T. Clutton-Brock, and B. Grenfell. 1999. Inverse density dependence and the Allee effect. Trends in Ecology and Evolution 14(10): 405-410. Cresswell, J.E., J.L. Osbourne, and D. Goulson. 2000. An economic model of the limits to foraging range in central place foragers with numerical solutions for bumblebees. Ecological Entomology 25: 249-255. Cushman, J.H., C.L. Boggs, S.B. Weiss, D.D. Murphy, A.W. Harvey, and P.R. Ehrlich. 1994. Estimating female reproductive success of a threatened butterfly: influence of emergence time and host-plant phenology. Oecologia 99: 194-200. Davis, G.A.N., J.F.D. Frazer, and A.M. Tynan. 1958. Population numbers in a colony of Lysandra Bellargus Rott. (Lepidoptera: Lycaenidae) during 1956. Proceedings of the Royal Entomological Society of London 33: 31- 36. del Castillo, R.C., and J. Nunez-Farfan. 1999. Sexual selection on maturation time and body size in Sphenarium purpurascens (Orthoptera : Pyrgomorphidae): Correlated response to selection. Evolution 53: 209- 215. 98 Dennis, B. 1989. Allee effects: population growth, critical density, and the chance of extinction. Natural Resource Modeling 3(4): 481-538. Dennis, B. 2002. Allee effects in stochastic populations. Oikos 96: 389-401. D'Eon, R.G., S.M. Glenn, I. Parfitt, and M. Fortin. 2002. Landscape connectivity as a function of scale and organism vagility in a real forested landscape. Conservation Ecology 6: Article 10. Doherty Jr., P.F., G. Sorci, J.A. Royle, J.E. Hines, J.D. Nichols, and T. Boulinier. 2003. Sexual selection affects local extinction and turnover in bird communities. Proceedings of the National Academy of Science 100: 5858-5862. Dowdeswell, W.H., R.A. Fisher, and E.B. Ford. 1940. The quantitative study of populations in the Lepidoptera. I. Polyommatus icarus. Annals of Eugenics 10: 123-136. Dunham, J.B., and B.E. Rieman. 1999. Metapopulation structure of bull trout: Influences of physical, biotic, and geometrical landscape characteristics. Ecological Applications 9: 624-655. Engen S., R. Lande, and B.E. Saether. 2003. Demographic stochasticity and Allee effects in populations with two sexes. Ecology 84(9): 2378-2386. Etienne, R.S., and H. Olff. 2004. A novel genealogical approach to neutral biodiversity theory. Ecology Letters 7: 170-175. Etienne, R.S. 2005. A new sampling formula for neutral biodiversity. Ecology Letters 8: 253-260. Evelyn, O.B., and R. Camirand. 2003. Forest cover and deforestation in Jamaica: an analysis of forest cover estimates over time. International Forestry Review 5: 354-363. Ewens, W.J. 1972. The sampling theory of selectively neutral alleles. Theoretical Population Biology 3: 87-112. Fagan, W.F. 2002. Connectivity, fragmentation, and extinction risk in dendritic metapopulations. Ecology 83: 3243-3249. Fagan, W. F., E. Meir, J. Prendergast, A. Folarin, and P. M. Kareiva. 2001. Characterizing vulnerability to extinction for 758 species. Ecology Letters 4: 132-138. 99 Fagan, W.F., P.J. Unmack, C. Burgess and W.L. Minckley. 2002. Rarity, fragmentation, and extinction risk in desert fishes. Ecology 83: 3250-3256. Fagan, W.F., and J.M. Calabrese. 2005. Quantifying connectivity: Balancing metric performance with data requirements. In K. Crooks and M.A. Sanjayan, eds. Connectivity conservation. Cambridge University Press. Fargione, J., C.S. Brown, and D. Tilman. 2003. Community assembly and invasion: An experimental test of neutral versus niche processes. Proceedings of the National Academy of Sciences (USA) 101: 414-414. Fortin, M.J., B. Boots, F. Csillag and T.K. Remmel. 2003. On the role of spatial stochastic models in understanding landscape indices in ecology. Oikos 102(1): 203-212. Gillis, E.A., and C.J. Krebs. 1999. Natal dispersal of snowshoe hares during a cyclic population increase. Journal of Mammalogy 80: 933-939. Gillis, E.A., and C.J. Krebs. 2000. Survival of dispersing versus philopatric juvenile snowshoe hares: do dispersers die? Oikos 90: 343-346. Godfray, H.C.J., M.P. Hassell, and R.D. Holt. 1994. The population dynamic consequences of phenological asynchrony between parasitoids and their hosts. Journal of Animal Ecology 63: 1-10. Green, D.M. 2003. The ecology of extinction: population fluctuation and decline in amphibians. Biological Conservation 111: 331-343. Groom, M.J. 1998. Allee effects limit population viability of an annual plant. American Naturalist 151 (6): 487-496. Haddad, N.M. 1999. Corridor use predicted from behaviors at habitat boundaries. American Naturalist 152: 215-227. Haines-Young, R., and M. Chopping. 1996. Quantifying landscape structure: a review of landscape indices and their application to forested landscapes. Progress in Physical Geography. 20(4): 418-445. Hanski, I. 1981. Coexistence of competitors in patchy environments with and without predation. Oikos 37: 306-312. Hanski, I. 1991. Metapopulation dynamics: brief history and conceptual domain. Biological Journal of the Linnean Society 42: 3-16. 100 Hanski, I. 1994. A practical model of metapopulation dynamics. Journal of Animal Ecology 63: 151-162. Hanski, I. 2001. Spatially realistic theory of metapopulation ecology. Naturwissenschaften 88: 372-381. Hanski, I., and Y. Camberfort 1991. Competition in dung beetles, In Hanski, I., and Y. Camberfort, eds. Dung Beetle Ecology. Princeton University Press, Princeton. Pp. 305-329. Hanski, I., A. Moilanen, T. Pakkala and M. Kuussaari. 1996. The quantitative incidence function model and persistence of an endangered butterfly metapopulation. Conservation Biology 10: 578-590. Hanski, I. and O. Ovaskainen. 2000. The metapopulation capacity of a fragmented landscape. Nature 404: 755-758. Hanski, I., and O. Ovaskainen. 2003. Metapopulation theory for fragmented landscapes. Theoretical Population Biology 64: 119-127. Harte, J. 2003. Tail of death and resurrection. Nature 424: 1006-1007. Hartley, S., and B. Shorrocks 2002. A general framework for the aggregation model of coexistence. Journal of Animal Ecology 71: 651-662. Hastings, N.A.J., and J.B. Peacock. 1975. Statistical distributions. Butterworths, London, UK. Hassell, M.P., and H.N. Comins. 1976. Discrete-Time Models for 2-Species Competition. Theoretical Population Biology 9: 202-221. Havel, J.E., J.B. Shurin, and J.R. Jones. 2002. Estimating dispersal from patterns of spread: Spatial and local control of lake invasions. Ecology 83: 3306- 3318. Hilborn, R. and M. Mangel 1997. The ecological detective: Confronting models with data. Princeton University Press, Princeton. Holzapfel, C.M., and W.E. Bradshaw. 2002. Protandry: The relationship between emergence time and male fitness in the pitcher-plant mosquito. Ecology 83: 607-611. Horn, H.S. and R.H. MacArthur 1972. Competition among Fugitive Species in a Harlequin Environment. Ecology 53: 749-752. 101 Hubbell, S.P. 1997. A unified theory of biogegraphy and relative species abundance and its application to tropical rain forests and coral reefs. Coral Reefs 16 Supplement: S9-S21. Hubbell, S.P. 2001. The unified neutral theory of biodiversity and biogeography. Princeton University Press, Princeton. Hubbell, S.P. and L. Borda-De-Agua. 2004. The unified neutral theory of biodiversity and biogeography: Reply. Ecology 85: 3175-3178. Hutchinson, G.E. 1959. Homage to Santa-Rosalia or why are there so many kinds of animals. American Naturalist 93:145-159. Ims, R.A., and N.G. Yoccoz. 1997. Studying transfer processes in metapopulation: emigration, migration, and colonization. In Hanski, I.A., and M.E. Gilpin eds. Metapopulation biology: ecology, genetics, and evolution. Academic Press. San Diego, CA. Pp 247-265. Ives, A.R. and R.M. May 1985. Competition within and between species in a patchy environment: relations between microscopic and macroscopic models. Journal of Theoretical Biology 115: 65-92. Ives, A.R. 1988. Aggregation and the coexistence of competitors. Annales Zoologici Fennici 25: 75-88. Ives, A.R. 1991. Aggregation and coexistence in a carrion fly community. Ecological Monographs 61: 75-94. Iwasa, Y., F.J. Odendaal, D.D. Murphy, P.R. Ehrlich, and A.E. Launer. 1983. Emergence patterns in male butterflies: a hypothesis and a test. Theoretical Population Biology 23: 363-379. Iwasa, Y., and S.A. Levin. 1995. The timing of life-history events. Journal of Theoretical Biology 172: 33-42. Johnson, C.N. 2002. Determinants of loss of mammal species during the Late Quaternary ?megafauna? extinctions: life history and ecology, but not body size. Proceedings of the Royal Society of London, Series B 269: 2221- 2227. Johnson, M.L., and M.S. Gaines. 1985. Selective basis for emigration of the prairie vole, Microtus ochrogaster: open field experiment. Journal of Animal Ecology 54:399-410. 102 Jordal, B.H. 1998. A review of Scolytodes Ferrari (Coleoptera : Scolytidae) associated with Cecropia (Cecropiaceae) in the northern Neotropics. Journal of Natural History 32: 31-84. Jordal, B.H. and L.R. Kirkendall 1998. Ecological relationships of a guild of tropical beetles breeding in Cecropia petioles in Costa Rica. Journal of Tropical Ecology 14: 153-176. Kelly, M.G. and D.A. Levin 1997. Fitness consequences and heritability aspects of emergence date in Phlox drummondii. Journal of Ecology 85(6): 755- 766. Keitt, T.H., D.L. Urban, and B.T. Milne. 1997. Detecting critical scales in fragmented landscapes. Conservation Ecology 11: Article 4. Kokko, H., and R. Brooks. 2003. Sexy to die for? Sexual selection and the risk of extinction. Annales Zoologici Fennici 40: 207-219. Kot, M., M.A. Lewis, and P. van den Dreissche. 1996. Dispersal data and the spread of invading organisms. Ecology 77: 2027-2042. Krijger, C.L., and J.G. Sevenster. 2001. Higher species diversity explained by stronger spatial aggregation across six neotropical Drosophila communities. Ecology Letters 4: 106-115. Kunin, W.E. 1998. Extrapolating species abundance across spatial scales. Science 281: 1513-1515. Lederhouse, R.C. 1983. Population structure, residency and weather related mortality in the black swallowtail butterfly, Papilio polyxenes. Oecologia 59: 307-311. Lepsch-Cunha, N., and S.A. Mori. 1999. Reproductive phenology and mating potential in a low density tree population of Couratari multiflora (Lecythidaceae) in central Amazonia. Journal of Tropical Ecology 15: 97- 121. Levins, R. 1969. Some demographic and genetic consequences of environmental heterogeneity for biological control. Bulletin of the Entomological Society of America 15: 237-240. Levins, R., and D. Culver. 1971. Regional coexistence of species and competition between rare species. Proceedings of the National Academy of Sciences (USA) 69: 1246-1248. 103 MacArthur, R. and R. Levins 1967. The limiting similarity, convergence, and divergence of coexisting species. American Naturalist 101: 377-385. Marzluff, M., and K. P. Dial. 1991. Life history correlates of taxonomic diversity. Ecology 72: 428-439. Matziris, D.I. 1994. Genetic-variation in the phenology of flowering in black pine. Silvae Genetica 43(5-6): 321-328. McCarthy, M.A. 1997. The Allee effect, finding mates and theoretical models. Ecological Modelling 103: 99-102. McGarigal, K., S.A. Cushman, M.C. Neel, and E. Ene. 2002. FRAGSTATS: Spatial pattern analysis program for categorical maps. Computer software program produced by the authors at the University of Massachusetts, Amherst. Available at the following web site: www.umass.edu/landeco/research/fragstats/fragstats.html Meegan, R.P., and D.S. Maehr. 2002. Landscape conservation and regional planning for the Florida panther. Southeastern Naturalist 1: 217-232. Merriam, G. 1984. Connectivity: a fundamental ecological characteristic of landscape pattern. In Brandt, J., and P. Agger, eds. Proceedings of the first international seminar on methodology in landscape ecological research and planning. Roskilde University. Roskilde, Denmark. Pp. 5-15. Moilanen, A. 2002. Implications of empirical data quality to metapopulation model parameter estimation and application. Oikos 96: 516-530. Moilanen, A., and I. Hanski. 2001. On the use of connectivity measures in spatial ecology. Oikos 95: 147-151. Moilanen, A., and M. Nieminen. 2002. Simple connectivity measures in spatial ecology. Ecology 83: 1131-1145. M?ller A. P. 2003. Sexual selection and extinction: why sex matters and why asexual models are insufficient. Annales Zoologici Fennici 40: 221-230. Morris, W., and D. Doak. 2002. Quantitative conservation biology: theory and practice of population viability analysis. Sinauer Press, Sunderland, MA, USA. Nicholson, A.J., and V.P. Bailey. 1935. The balance of animal populations. Proceedings of the Royal Society of London. 3: 551-598. 104 Nikkanen, T. 2001. Reproductive phenology in a Norway spruce seed orchard. Silva Fennica 35(1): 39-53. Oh, R.J. 1979. Repeated copulation in the brown planthopper, Nilaparvata lugens Stal. (Homoptera: Delphacidae). Ecological Entomology 4: 345-353. Ollerton, J., and A. Diaz. 1999. Evidence for stabilising selection acting on flowering time in Arum maculatum (Araceae): the influence of phylogeny on adaptation. Oecologia 119: 340-348. Ovaskainen, O., and I. Hanski. 2001. Spatially structured metapopulation models: Global and local assessment of metapopulation capacity. Theoretical Population Biology 60: 281-302. Peterson, M.A. 1995. Phenological isolation, gene flow, and developmental differences among low- and high-elevation populations of Euphilotes enoptes (Lepidoptera: Lycaenidae). Evolution 49: 446-455. Philippi, T., and J. Seger. 1989. Hedging one?s evolutionary bets, revisited. Trends in Ecology and Evolution 4: 41-44. Pimm, S.L., H.L. Jones, and J. Diamond. 1988. On the risk of extinction. American Naturalist 132: 757-785. Pollard, E. 1981. Aspects of the ecology of the meadow brown butterfly, Maniola jurtina (L.) (Lepidoptera: Satyridae). Entomologist's Gazette 32: 67-74. Pollard, E., and T. J. Yates. 1993. Monitoring butterflies for ecology and conservation: Conservation biology series. Chapman and Hall, London. Porter, W.P., S. Budaraju, W.E. Stewart, and N. Ramankutty. 2000. Physiology on a landscape scale: Applications in ecological theory and conservation practice. American Zoologist 40: 1175-1176. Post, E., S.A. Levin, Y. Iwasa, and N.C. Stenseth. 2001. Reproductive asynchrony increases with environmental disturbance. Evolution 55(4): 830-834. Preston, F.W. 1948. The commoness, and rarity, of species. Ecology 29: 254-283. Primack, R.B. 1980. Variation in the phenology of natural populations of montane shrubs in New Zealand. Journal of Ecology 68: 849-862. Purrington, C. B. 1993. Parental effects on progeny sex-ratio, emergence, and flowering in Silene latifolia (Caryophyllaceae). Journal of Ecology 81: 807-811. 105 Purrington, C.B., and J. Schmitt. 1998. Consequences of sexually dimorphic timing of emergence and flowering in Silene latifolia. Journal of Ecology 86: 397-404. Ricklefs, R.E. 2003. A comment on Hubbell's zero-sum ecological drift model. Oikos 100: 185-192. Rossiter, M. C. 1991. Maternal effects generate variation in life history: consequences of egg weight plasticity in the gypsy moth. Functional Ecology 5: 386-393. Ruckelshaus, M., C. Hartway, and P. Kareiva. 1997. Assessing the data requirements of spatially explicit dispersal models. Conservation Biology 11: 1298-1306. Ruckelshaus, M., C. Hartway, and P. Kareiva. 1999. Dispersal and landscape errors in spatially explicit population models: a reply. Conservation Biology 13: 1223-1224. Saether, B.E. 1997. Life history variation, population processes, and priorities in species conservation: towards a reunion of research paradigms. Oikos 77: 217-226. Satake, A., A. Sasaki, and Y. Iwasa. 2001. Variable timing of reproduction in unpredictable environments: adaptation of flood plain plants. Theoretical Population Biology 60: 1-15. Schtickzelle, N., E. Le Bouleng?, and M. Baguette. 2002. Metapopulation dynamics of the bog fritillary butterfly: demographic processes in a patchy population. Oikos 97: 349-360. Schumaker, N.H. 1996. Using landscape indices to predict habitat connectivity. Ecology 77: 1210-1225. Seger, J., and H.J. Brockmann. 1987. What is bet-hedging? Oxford Survey of Evolutionary Biology 4: 182-211. Sevenster, J.G. 1996. Aggregation and coexistence. 1. Theory and analysis. Journal of Animal Ecology 65: 297-307. Sevenster J.G., and J.J.M. VanAlphen 1996. Aggregation and coexistence. 2. A neotropical Drosophila community. Journal of Animal Ecology 65: 308- 324. Shorrocks, B. 1991. Competition on a divided and ephemeral resource: a cage experiment. Biological Journal of the Linnean Society 43: 211-220. 106 Shorrocks B., J. Rosewell, K. Edwards, and W. Atkinson 1984. Interspecific competition is not a major organizing force in many insect communities. Nature 310: 310-312. Shorrocks, B., and J. Rosewell 1986. Guild size in drosophilids: a simulation model. Journal of Animal Ecology 55: 527-541. Shorrocks, B., and J. Rosewell 1987. Spatial patchiness and community structure: coexistence and guild size of drosophilids on ephemeral resources. In Gee, J.H.R., and P.S. Giller, eds. Organisation of communities: Past and present. Blackwell Scientifc Publications, Oxford. Pp. 29-51. Shorrocks, B., and M. Bingley 1994. Priority effects and species coexistence ? experiments with fungal-breeding Drosophila. Journal of Animal Ecology 63: 799-806. Siitonen, P., A. Tanskanen, and A. Lehtinen. 2002. Method for selection of old- forest reserves. Conservation Biology 16: 1398-1408. Siitonen, P., A. Tanskanen, and A. Lehtinen. 2003. Selecting forest reserves with a multiobjective spatial algorithm. Environmental Science and Policy 6: 301-309. Simmons, A.M., and M.O. Johnston. 1997. Developmental Instability as a bet- hedging strategy. Oikos 80(2): 401-406. Singer, M.C., and P.R. Ehrlich. 1979. Population dynamics of the checkerspot butterfly Euphydryas editha. Fortschritte der Zoologie 25: 53-60. Singleton, P.H., W.L. Gaines, and J.F. Lehmkuhl. 2002. Landscape permeability for large carnivores in Washington: a geographic information system weighted-distance and least-cost corridor assessment. USDA Forest Service Pacific Northwest Research Station Research Paper 549: 1-+. Southwood, T.R.E. 1978. Ecological methods. 2 nd edn. Chapman and Hall, London. Stephenson, A.G., and R.I. Bertin. 1983. Male competition, female choice, and sexual selection in plants. In L. Leal, ed. Pollination biology. Academic Press, New York, USA. Pp. 109-149. Sowter, F.A. 1949. Arum maculatum L. Journal of Ecology 37: 207-219. 107 Sutherland, W.J. 1996. Predicting the consequences of habitat loss for migratory populations. Proceedings of the Royal Society of London Series B ? Biological Sciences 263: 1325-1327. Tammaru, T., K. Ruohom?ki, and I. Saloniemi. 1999. Within-season variability of pupal period in the autumnal moth: a bet-hedging strategy. Ecology 80(5): 1666-1677. Taylor, B.W., C.R. Anderson, and B.L. Peckarsky. 1998. Effects of size at metamorphosis on stonefly fecundity, longevity, and reproductive success. Oecologia 114: 494-502. Taylor, P.D., L. Fahrig, K. Henein, and G. Merriam. 1993. Connectivity is a vital element of landscape structure. Oikos 68: 571-573. Tikkanen, O.P. and P. Lyytikainen-Saarenmaa 2002. Adaptation of a generalist moth, Operophtera brumata, to variable budburst phenology of host plants. Entomologia Experimentalis Et Applicata 103(2): 123-133. Tilman, D. 1994. Competition and biodiversity in spatially structured habitats. Ecology 75: 2-16. Tilman, D. 2004. Niche tradeoffs, neutrality, and community structure: A stochastic theory of resource competition, invasion, and community assembly. Proceedings of the National Academy of Sciences (USA) 101: 10854-10861. Tischendorf, L. 2001a. On the use of connectivity measures in spatial ecology: a reply. Oikos 95: 152-155. Tishcendorf, L. 2001b. Can landscape indices predict ecological processes consistently? Landscape Ecology. 16(3): 235-254. Tischendorf, L., and L. Fahrig. 2000. On the usage and measurement of landscape connectivity. Oikos 90: 7-19. Tischendorf, L., D.J. Bender, and L. Fahrig. 2003. Evaluation of patch isolation metrics in mosaic landscapes for specialist vs. generalist dispersers. Landscape Ecology 18: 41-50. Turchin, P. 1998, Quantitative analysis of movement: measuring and modeling population redistribution in animals and plants, Sinauer Associates, Inc. Urban, D., and T. Keitt. 2001. Landscape connectivity: a graph-theoretic perspective. Ecology 82: 1205-1218. 108 van Langevelde, F. 2000. Scale of habitat connectivity and colonization in fragmented nuthatch populations. Ecography 23: 614-622. Veech, J.A., T.O. Crist and K.S. Summerville 2003. Intraspecific aggregation decreases local species diversity of arthropods. Ecology 84: 3376-3383. Volkov, I., J.R. Banavar, S.P. Hubbell, and A. Maritan 2003. Neutral theory and relative species abundance in ecology. Nature 424: 1035-1037. Wahlberg, N., T. Klemetti, V. Selonen, and I. Hanski. 2002. Metapopulation structure and movements in five species of checkerspot butterflies. Oecologia 130: 33-43. Waldbauer, G.P. 1978. Phenological adaptation and polymodal emergence patterns of insects. In H. Dingle, ed. Evolution of insect migration and diapause. Springer-Verlag, New York, NY, USA. Pp. 127-144 Warren, M.S. 1987a. The ecology and conservation of the heath fritillary butterfly, Mellicta athalia. I. Host selection and phenology. Journal of Applied Ecology 24: 467-482. Warren, M.S. 1987b. The ecology and conservation of the heath fritillary butterfly, Mellicta athalia. II. Adult population structure and mobility. Journal of Applied Ecology 24: 483-498. Warren, M.S., E. Pollard, and T.J. Bibby. 1986. Annual and long-term changes in a population of the wood white butterfly Leptidea sinapis. Journal of Animal Ecology 55: 707-719. Wells, H., E.G. Strauss, M.A. Rutter, and P.H. Wells. 1998. Mate location, population growth and species extinction. Biological Conservation 86: 317-324. Wiens, J.A. 1997. The emerging role of patchiness in conservation biology. In Pickett, S.T.A., R.S. Ostfeld, M.Shachak, and G.E. Likens, eds. The ecological basis of conservation: heterogeneity, ecosystems, and biodiversity. Chapman and Hall. New York, NY. Pp. 93-107. Wertheim, B., J.G. Sevenster, I.E.M Eijs, and J.J.M. Van Alphen 2000. Species diversity in a mycophagous insect community: the case of spatial aggregation vs. resource partitioning. Journal of Animal Ecology 69: 335- 351. Whittaker, R.H. 1965. Dominance and diversity in land plant communities. Science 147: 250-260. 109 Wiklund, C., and T. Fagerstr?m. 1977. Why do males emerge before females? A hypothesis to explain the incidence of protandry in butterflies. Oecologia 31: 153?158. Wiklund, C., and C. Solbreck. 1982. Adaptive versus incidental explanations for the occurrence of protandry in a butterfly, Leptidea sinapsis L. Evolution 36: 56?62. Wootton, J.T. 2005. Field parameterization and experimental test of the neutral theory of biodiversity. Nature 433: 309-312.