text:
Multi-species interactions and the evolution of biological systems by Leithen K. M'Gonigle B.Sc. Mathematics (Honours), University of Victoria, 2005 a thesis submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in the faculty of graduate studies (Zoology) The University Of British Columbia (Vancouver) December 2011 c Leithen K. M'Gonigle, 2011 Abstract In this thesis I develop several models examining how genetic evolution can affect evolutionary processes at a broader scale. First, I ask how evolution would proceed at a locus that governs the mutation rate between alleles mediating interactions between hosts and parasites. By relaxing several simplifying assumptions I am able to explore the affects of sex and recombination. I find that, when the modifier locus is completely linked, the mutation rate evolves toward the optimum rate. With looser linkage, however, lower mutation rates evolved. This work can potentially explain the high rates of antigenic switching observed in many asexual taxa. Second, I investigate how ploidy levels and the genetic model underlying species interactions affect how evolution proceeds from a free-living to a par- asitic life-history. I find that the transition to parasitism occurs over a broader range of parameters when the parasite is haploid. The role of host ploidy is more complicated, depending on the model governing host-parasite interac- tions. These results provide a first characterization of how genetic architecture affects selection on life-history in antagonistic species interactions. Third, I develop a model of sexual selection in an environment with spatial variation in the carrying capacity, but no variation in resource type. I show that, when searching for a mate is costly, this variation can stabilize demographic fluctuations, facilitating long-term coexistence of species differing only in sex- ual traits. This is the first study to demonstrate the existence of conditions under which sexual selection alone can promote the long-term coexistence of ecologically equivalent species in sympatry. Finally, I develop a model characterizing the effects of mating preferences on species interactions in hybrid zones. I find that the spatial distribution of geno- types observed in many "mosaic" hybrid zones might be better explained by species-specific differences in mating than by differences in ecology (the com- ii mon explanation). In addition, I develop a statistical method that can be applied to empirical hybrid zone data to estimate how "mosaic" the hybrid zone is. I test this statistic on data from the Mytilus edulis and M. galloprovincialis hybrid zone. iii Preface I conceived and developed the model in chapter 2, with guidance and assistance from Sarah P. Otto. An undergraduate student, Jane Shen, assisted me on the project as part of her honours thesis. Sarah P. Otto helped with writing and editing the manuscript. This article has been published: M'Gonigle, L.K., Shen, J.J., and Otto, S.P. 2009. Mutating away from your en- emies: The Evolution of Mutation Rate in a Host-Parasite Model. Theor. Popul. Biol.: 75: 301–311. I conceived and developed the model in chapter 3, with guidance and assistance from Sarah P. Otto. She also helped with writing and editing the manuscript. The article has been published: M'Gonigle, L.K., and Otto, S.P. 2011. Ploidy and the evolution of parasitism. Proc. R. Soc. B: 278: 2814–2822. I conceived the model in chapter 4, in collaboration with Ulf Dieckmann when I was a participant in the 2008 Young Scientists Summer Program at the Interna- tional Institute for Applied Systems Analysis in Laxenburg, Austria. Sarah P. Otto, Ulf Dieckmann, and Rupert Mazzucco all provided helpful input in the development and analysis of the model and the editing of the manuscript. This article has been submitted. I conceived the model in chapter 5, and developed and analyzed it in collabo- ration with Richard FitzJohn. Sarah P. Otto also provided helpful input in the development and analysis of the model and both Sarah P. Otto and Richard FitzJohn helped in the editing of the manuscript. Nicolas Bierne kindly pro- iv vided data for analysis and inclusion in the manuscript. The article has been published. M'Gonigle, L.K. and FitzJohn, R.G. 2010. Assortative mating and spatial struc- ture in hybrid zones. Evolution: 64: 444–455. v Table of Contents Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . ii Preface . . . . . . . . . . . . . . . . . . . . . . . . . . . iv Table of Contents . . . . . . . . . . . . . . . . . . . . . . . vi List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . viii List of Figures . . . . . . . . . . . . . . . . . . . . . . . . ix Acknowledgments . . . . . . . . . . . . . . . . . . . . . . xii 1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.1 Genetics and host-parasite co-evolution . . . . . . . . . . . . . . . . 2 1.2 Sexual selection and co-existence . . . . . . . . . . . . . . . . . . . . 3 1.3 Hybrid zones . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 1.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2 Mutating away from your enemies: the evolution of mutation rate in a host-parasite system . . . . . . . . . . . . . . . . . . . . . . 6 2.1 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 2.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 2.3 Methods and results . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 2.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 3 Ploidy and the evolution of parasitism . . . . . . . . . . . . . . 26 3.1 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 3.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 3.3 Model summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 3.4 Analytical results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 vi 3.5 Simulation model summary . . . . . . . . . . . . . . . . . . . . . . . 33 3.6 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 4 Sexual selection enables long-term coexistence despite ecological equivalence . . . . . . . . . . . . . . . . . . . . . . . . . 45 4.1 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 4.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 4.3 Model summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 4.4 Results and discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 48 5 Assortative mating and spatial structure in hybrid zones . . . . . . . 55 5.1 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 5.2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 5.3 Model description . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57 5.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 5.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 6 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . 76 6.1 Future directions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77 Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . 81 A Mutating away from your enemies (Chapter 2) . . . . . . . . . . . 91 A.1 Cost of deleterious mutations . . . . . . . . . . . . . . . . . . . . . . 91 A.2 Solving the recursion equations for the disequilibrium . . . . . . . 93 A.3 Different generation times . . . . . . . . . . . . . . . . . . . . . . . . 94 A.4 One-species model . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96 B Ploidy and the evolution of parasitism (Chapter 3) . . . . . . . . . 99 B.1 Additional analyses . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 B.2 Simulations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 C Sexual selection enables co-existence (Chapter 4) . . . . . . . . . . 114 C.1 Model description . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114 C.2 Model extensions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119 D Hybrid zones (Chapter 5) . . . . . . . . . . . . . . . . . . . . 132 D.1 Hybrid zone structure likelihood method . . . . . . . . . . . . . . . 132 vii List of Tables Table 3.1 Description of parameters . . . . . . . . . . . . . . 39 Table 3.2 Invasion matrices . . . . . . . . . . . . . . . . . 39 Table 3.3 The fitness advantage of the modifier allele when genetic associations are weak . . . . . . . . . . . . . . . . . . . 40 Table 3.4 Invasion condition for a modifier that increases the level of parasitism in MAM and IMAM . . . . . . . . . . . . . . . 41 Table B.1 Equations for FA,H and FA,P when yH and yP are not as- sumed to be near 1 . . . . . . . . . . . . . . . . . . . . 102 Table B.2 Full equations for w̄diff, without assuming high levels of sex- ual reproduction in either species . . . . . . . . . . . . . . 103 Table C.1 Model parameters and model variables. . . . . . . . . 124 viii List of Figures Figure 2.1 Period of co-evolutionary cycles as a function of the strength of selection in parasites . . . . . . . . . . . . . . . . . . 21 Figure 2.2 Phase difference between host and parasite cycles as a func- tion of the difference in their mutation rates . . . . . . . . . . 22 Figure 2.3 Evolutionarily attracting mutation rate in hosts as a function of the strength of selection . . . . . . . . . . . . . . . . . 23 Figure 2.4 Evolutionarily attracting mutation rate in hosts, with multi- ple alleles . . . . . . . . . . . . . . . . . . . . . . . . 24 Figure 2.5 Co-evolutionarily trajectories in simulations with differing strengths of selection in hosts and parasites . . . . . . . . . . 25 Figure 3.1 Sample trajectories from simulations in the matching-alleles model with complete parasitism . . . . . . . . . . . . . . . 42 Figure 3.2 Invasion conditions in the matching-alleles and inverse-matching- alleles model . . . . . . . . . . . . . . . . . . . . . . . 43 Figure 3.3 Invasion conditions in the gene-for-gene model . . . . . 43 Figure 3.4 Evolutionary convergent level of parasitism as a function of the strength of selection in parasites under the matching-alleles model and the inverse-matching-alleles model . . . . . . . . . 44 Figure 4.1 Sexual selection can enable long-term coexistence of ecolog- ically equivalent species . . . . . . . . . . . . . . . . . . 51 Figure 4.2 Distributions of allele frequencies at the display locus through time . . . . . . . . . . . . . . . . . . . . . . . . . . 52 Figure 4.3 Time until loss of polymorphism at the display locus when females are choosy . . . . . . . . . . . . . . . . . . . . 53 Figure 4.4 Four representative model runs in a patchy landscape. . . 54 ix Figure 5.1 An example model fitted to hypothetical allele frequency data along a transect through a hybrid zone . . . . . . . . . . 70 Figure 5.2 Example simulated hybrid zones with varying strengths of assortative mating. . . . . . . . . . . . . . . . . . . . . 70 Figure 5.3 Mosaicity as a function of long distance dispersal distance in the one-locus linear-preference model. . . . . . . . . . . . 71 Figure 5.4 Mosaicity, linkage disequilibrium, and bimodality plotted over time . . . . . . . . . . . . . . . . . . . . . . . . 72 Figure 5.5 Mosaicity as a function of mate preference under different preference models. . . . . . . . . . . . . . . . . . . . . 73 Figure 5.6 Mosaicity at each locus as a function of the strength of se- lection against hybrids . . . . . . . . . . . . . . . . . . . 74 Figure 5.7 Best fit models for three diagnostic loci in theMytilus edulis/M. galloprovincialis hybrid zone . . . . . . . . . . . . . . . . . 75 Figure A.1 Co-evolutionarily trajectories in simulations with differing costs associated with mutation rate modifiers. . . . . . . . . . 98 Figure B.1 Evolutionary convergent level of parasitism in MAM and IMAM . . . . . . . . . . . . . . . . . . . . . . . . . 104 Figure B.2 Evolutionary convergent level of parasitism in MAM and IMAM with aH = 0.01 . . . . . . . . . . . . . . . . . . . 105 Figure B.3 Evolutionary convergent level of parasitism in MAM and IMAM with stronger selection in hosts . . . . . . . . . . . . 106 Figure B.4 Evolutionary convergent level of parasitism in MAM and IMAM with smaller population sizes . . . . . . . . . . . . . 107 Figure B.5 Evolutionary convergent level of parasitism in MAM and IMAM with some hosts reproducing asexually . . . . . . . . . 108 Figure B.6 Evolutionary convergent level of parasitism in MAM and IMAM with different generation times in hosts and parasites . . . 109 Figure B.7 Evolutionary convergent level of parasitism in MAM and IMAM with three alleles at the A-locus in each species . . . . . . 110 Figure B.8 Evolutionary convergent level of parasitism in the GFG with conditional costs to virulence . . . . . . . . . . . . . . . . 111 Figure B.9 Evolutionary convergent level of parasitism in the GFG with unconditional costs to virulence and differing initial levels of parasitism 112 x Figure B.10 Evolutionary convergent level of parasitism in the GFG with unconditional costs to virulence, differing initial levels of parasitism, and small populations . . . . . . . . . . . . . . . . . . . 113 Figure C.1 Variation in three components of fitness as a function of the local carrying capacity . . . . . . . . . . . . . . . . . . . 125 Figure C.2 Search costs experienced by females for the model run in figure 4.1D . . . . . . . . . . . . . . . . . . . . . . . 126 Figure C.3 Distributions of allele frequencies at the display locus through time with and without mating-dependent dispersal . . . . . . . 126 Figure C.4 Minimum level of variation in local carrying capacity needed to ensure that coexistence is maintained for, on average, at least 5N generations . . . . . . . . . . . . . . . . . . . . . . . 127 Figure C.5 Minimum level of variation in local carrying capacity needed to ensure that coexistence is maintained for, on average, at least 5N generations . . . . . . . . . . . . . . . . . . . . . . . 128 Figure C.6 Time until loss of polymorphism at the display locus when females are choosy in the fecundity model . . . . . . . . . . . 128 Figure C.7 Time until loss of polymorphism at the display locus when females are choosy . . . . . . . . . . . . . . . . . . . . 129 Figure C.8 Effects of changes in genetic architecture in a two-dimensional bimodal landscape . . . . . . . . . . . . . . . . . . . . 130 Figure C.9 Effects of asymmetric fitness costs of display traits in a two- dimensional bimodal landscape . . . . . . . . . . . . . . . 131 xi Acknowledgments This thesis would not have been possible without the support of numerous friends and family. First and foremost, I would like to thank Sally Otto. Words cannot express how grateful I am for the support, both academic and personal, that you have given me over these past six years. Thank you for making gradu- ate school the wonderful experience it has been. I would like to thank my committee members, Michael Doebeli, Loren Riese- berg, and Jeannette Whitton for providing me with helpful advice of all sorts. Several other faculty at both UBC and SFU have also provided help in some form or another. In no particular order: Colin Brauner, Elizabeth Elle, Darren Irwin, Bernie Roitberg, Dolph Schluter, and Mike Whitlock. I would also like to thank, for technical and administrative assistance, Alistair Blachford, Andy LeBlanc, Sanja LeBlanc, Richard Sullivan, and Alice Liou. I have been fortunate to be a member of a truly fantastic lab. Dilara Ally, Rich FitzJohn, Aleeza Gerstein, Jessica Hill, Kay Hodgins, Crispin Jordan, Karen Magnuson-Ford, Itay Mayrose, Jasmine Ono, Kate Ostevik all contributed to making my time at UBC both fun and productive. In particular, I would like to thank Rich for spending so much time helping me with the many computer- related problems I have encountered along the way, from the most basic pro- gramming to the formatting of this thesis. Ulf Dieckmann, thank you for enabling me to spend three summers in one of my favourite cities, and for welcoming me into your research group at IIASA. My wonderful memories of summers at IIASA were made possible by the good friends and colleagues I met through working there: Jan Ohlberger, Josh Payne, Davnah Urbach, Barbara Fischer, Rupert Mazzucco, Rebecca Whitlock, Daniel Falster and Jacob Johansson. A number of friends and colleagues have also contributed invaluably to the completion of this thesis, whether through scientific discussions, or providing xii company while cycling to conferences: Aneil Agrawal, Arianne Albert, Lisa An- derson, Rowan Barrett, Mirjam Barreuto, Alistair Blachford, Ella Bowles, Evan Braswell, Christina Burch, Anne Dalziel, Julie Granka, Philip Greenspoon, Luke Harmon, Emma Hodsgon, Travis Ingram, Nathan Kraft, Julie Lee-Yaw, Milica Mandic, Kerry Marchinko, Blake Matthews, Jon Mee, J.S. Moore, Scott Nuismer, Daniel Ortiz-Barrentos, Antoine Paccard, Seth Rudman, Kieren Samuk, Alana Schick, Maria Servedio, Jane Shen, Matt Siegle, Laura Southcott, Ben Speers- Roesch, Andrea Stephens, Dave Toews, Jabus Tyerman, Tim Vines, John Willis, Sam Yeaman. Finally I would like to thank the friends and family that provided me with the support outside of academic life that made it possible for me to complete this degree. Eric Stoehr, thank you for repeatedly reminding me that I do not have a real job. Deanna Bicego, thank you for helping me make the transition to Vancouver and for remaining a good friend over these past few years. Carla Crossman, thank you for making my last year in Vancouver the best one. I think you deserve most of the credit for my making me not want to finish this degree earlier. Lastly and perhaps most importantly, I would like to thank my parents, Michael and Wendy, for providing me with endless encouragement and support and also for helping to remind me that there are greater problems in the world than the ones I have encountered. xiii Chapter 1 Introduction To fully understand the complexities of the biological world it is necessary that we have both empirical documentation of a particular process, as well as a con- ceptual understanding of the mechanisms causing that process to occur. For many simple biological processes it is possible to easily intuit the outcome of a particular interaction. However, in many cases, such intuition is not possi- ble without some aid. By providing us with additional computational power and/or methods of analysis, theoretical models can provide such an aid, and thus facilitate our understanding in these more difficult cases. In addition to helping us understand documented biological phenomena, theoretical models can change how we understand existing data and, in some cases, lead to empirical discoveries that would not have occurred otherwise. For example, theoretical investigations of interactions between hybridizing species have led to a number of predictions about the expected spatial genetic patterns in the region where the two types overlap and form hybrid offspring (i.e., in their "hybrid zone"; Barton and Hewitt, 1985, 1989). One such prediction is that a hybrid zone should often move in space, until finally settling in a region of low population density. Several empirical studies have since found convinc- ing evidence that such a process really does occur in nature. For example, Saitoh and Katakura (1996) found that the boundaries between three parapatric species of flightless leaf beetles mostly occur along streams or cliffs where pop- ulation densities are low to non-existent. Similarly Barrowclough et al. (2005) found that the range of overlap between the northern spotted owl and the Cal- ifornia spotted owl is restricted to a narrow band within the region of lowest documented owl density. Such empirical findings as these serve to illustrate the critical role theoretical models can play in biology. A central focus in evolutionary biology is on understanding how processes occurring at one level of biological organization affect evolutionary dynamics 1 Chapter 1 at another level (e.g., how cellular processes affect organismal behaviour or how the behaviour of an individual may affect the dynamics of a population or group; Maynard Smith and Szathmary, 1997; Okasha, 2007). Because inter- actions between multiple levels of organization are complex and can easily be obscured by more direct interactions, many insights on this topic derive from theoretical work. At the most general level, the chapters of my thesis can be summarized as a series of theoretical models focusing on these interactions be- tween different evolutionary levels. In chapters 2 and 3 I have examined the relationship between species-level interactions and the genetic architecture of the traits regulating those interactions, both in the context of host-parasite in- teractions, and in chapters 4 and 5 I have examined the relationship between population-level interactions and individual behaviour (namely female mating behaviour). I will now briefly summarize each of my chapters. 1.1 Genetics and host-parasite co-evolution Historically host-parasite interactions have received substantial attention be- cause of the implications any findings may have for disease control and/or pest management. Host-parasite interactions are of particular interest to evo- lutionary biologists, however, because they represent one of the most complex and highly-evolved relationships between separate species. In many cases, the two interacting species are locked in a perpetual evolutionary "arms race", with selection constantly favouring novel mechanisms allowing hosts to evade their parasites and parasites to invade these same hosts. Understanding how these interactions can affect selection at the most basic genetic level is critical if we hope to understand how parasite genomes evolve in response to their biotic environment. Mutation is the fundamental source of genetic variation, without which adaptation could not occur. Models of the evolution of mutation rates in a con- stant environment have shown that asexual species evolve an optimal mutation rate that balances the costs and benefits of producing adaptive and deleterious mutations. In a host-parasite system, however, the fitness effect of a mutation is not constant but instead depends on the current composition of the other species. In chapter 2 I use modifier theory to ask how evolution would proceed at a locus that governs the mutation rate between alleles that mediate host- parasite interactions. In my analysis I am able to relax a number of common simplifying assumptions, and also explore the effects of sex and recombination. I find that when the modifier locus is completely linked to the locus mediat- ing the host-parasite interaction, the mutation rate evolves toward the optimum 2 Chapter 1 rate. With looser linkage, however, lower mutation rates evolve. Complementary to studies that focus on how host-parasite interactions se- lect for changes in genetic parameters (like that above) are those that investi- gate how genetic parameters affect evolution of host-parasite interactions. Such studies are necessary if we hope to fully understand the origins of broad-scale macro-evolutionary patterns. In the context of hosts and parasites, such an un- derstanding is essential as we try to characterize the traits that make a particular species more or less prone to adopting a parasitic life-style. One assumption that is common in models of host-parasite co-evolution is that species are either completely parasitic or completely non-parasitic. While many species do fit this assumption, there are numerous examples for which this assumption is not appropriate. For example, a number of species from a range of taxonomic groups have been shown to be facultatively parasitic: ciliates (Reynolds, 1936; Thompson and Moewus, 1964), flatworms (Hooge and Tyler, 1999), fungi (Morin et al., 1993), nematodes (Benham, 1974). These species are parasitic if the opportunity arises, but are otherwise free-living and capable of reproduction without the aid of a host species. In chapter 3 I use modifier the- ory to investigate how evolution proceeds from a free-living to a parasitic life- history and, in particular, I investigate how the ploidy levels of the interacting species and the genetic model underlying species interactions affects whether selection favours a more or less parasitic life-cycle. In general, I find that the transition to parasitism occurs over a broader range of parameters when the parasite is haploid. The role of host ploidy is more complicated, depending on the model governing host-parasite interactions. 1.2 Sexual selection and co-existence Just as genetic architecture can affect higher-level species interactions, so can the behaviour of an individual affect higher-level population dynamics. For example, the mating behaviour of females can sometimes select for particu- lar traits in males which, under the right conditions, can split a population into multiple reproductively isolated groups (Lande, 1981; Kirkpatrick, 1982; Seger, 1985). In some situations, therefore, sexual selection (an individual-level process) seems to promote the generation of diversity (a population-level pro- cess). While appearing to be capable of generating diversity, sexual selection has not been thought to be capable of promoting the long-term co-existence of distinct types, and thus not capable of maintaining biodiversity. This is because changes in mating behaviour or preferences do not lead to niche divergence, which is generally considered to be a necessary precondition for co-existence 3 Chapter 1 when species' ranges overlap (Weissing et al., 2011). In chapter 4 I develop a model of sexual selection in an environment with variation in resource quantity but no heterogeneity in resource type across space. I show that the resultant variation in the number of individuals can stabilize demographic fluctuations, which then allows for the stable long term co-existence of populations differing only in their mating preferences. While costs associated with searching for a mate do not enable co-existence on their own, they turn out to be critical in facilitating co-existence in this spatially ex- plicit model. 1.3 Hybrid zones In addition to maintaining biodiversity amongst or within populations, female mating preferences are also likely to have an impact on the spatial distribu- tion of diverse types. The most natural place to investigate spatial interactions between distinct genotypes is in the region where their ranges meet. If repro- ductive isolation is not complete, then hybridization can occur. These hybrid zones are, therefore, an ideal setting for observing the forces that maintain and/or break down species boundaries (Barton and Hewitt, 1989, 1985). The majority of theoretical investigations of hybrid zones have focused on ecologi- cal differences between species, which are often assumed to exist despite mixed empirical support (Harrison and Rand, 1989; Cain et al., 1999; Bridle et al., 2001). In chapter 5 I develop and analyze a model characterizing the poten- tial effects of mating preferences on spatial interactions between distinct mating types in hybrid zones. Interestingly, I find that the spatial distribution of geno- types observed in many "mosaic" hybrid zones (those characterized by patchy species distributions) is consistent with patterns generated by species-specific differences in mating preferences. In cases where there is little evidence for eco- logical differences between the interacting species, divergent mating preferences may, therefore, provide a better explanation for this mosaic pattern. In addition to the primary model, I develop a new statistical method that can be applied to empirical hybrid zone data to estimate how "mosaic" the hybrid zone is. I test this statistic on data from three loci from the hybrid zone of the marine musselsMytilus edulis andM. galloprovincialis. The estimated "mosaicity" is sig- nificantly higher for one of these three loci, suggesting that this particular locus or a linked region may, at least partially, underlie assortative mating and/or local adaptation.. 4 Chapter 1 1.4 Conclusions In addition to providing insight into several interesting topics in biology, the chapters in this thesis also illustrate the general importance of theoretical re- search. Theory has enabled us to fully understand many complex processes and, additionally, revealed many interesting avenues for future theoretical and empirical work. In each of the following chapters I will discuss the empiri- cal implications of my work, as well as any possible future directions for both theoretical and empirical research. In doing so I hope to demonstrate the impor- tance of theoretical work to empirical biology, and vice versa. Without theory, many potentially interesting components of the natural would remain unknown to us and without empirical tests, the relevance of theoretical research (and its potential contradictions) would be greatly under-appreciated. 5 Chapter 2 Mutating away from your enemies: the evolution of mutation rate in a host-parasite system 2.1 Summary The rate at which mutations occur in nature is itself under natural selection. While a general reduction of mutation rates is advantageous for species inhabit- ing constant environments, higher mutation rates can be advantageous for those inhabiting fluctuating environments that impose on-going directional selection. Analogously, species involved in antagonistic co-evolutionary arms-races, such as hosts and parasites, can also benefit from higher mutation rates. We use mod- ifier theory, combined with simulations, to investigate the evolution of mutation rate in such a host-parasite system. We derive an expression for the evolution- ary stable mutation rate between two alleles, each of whose fitness depends on the current genetic composition of the other species. Recombination has been shown to weaken the strength of selection acting on mutation modifiers, and ac- cordingly, we find that the evolutionarily attracting mutation rate is lower when recombination between the selected and the modifier locus is high. Cyclical dynamics are potentially commonplace for loci governing antagonistic species interactions. We characterize the parameter space where such cyclical dynam- ics occur and show that the evolution of large mutation rates tends to inhibit cycling and thus eliminates further selection on modifiers of mutation rate. We then find using computer simulations that stochastic fluctuations in finite pop- ulations can increase the size of the region where cycles occur, creating selection for higher mutation rates. We finally use simulations to investigate the model behaviour when there are more than two alleles, finding that the region where 6 Chapter 2 cycling occurs becomes smaller and the evolutionarily attracting mutation rate lower when there are more alleles. 2.2 Introduction Mutation is the fundamental source of genetic variation, without which adap- tation could not occur. The rate at which mutations occur has, therefore, been of long-standing interest to evolutionary biologists (Muller, 1928; Mukai, 1964; Drake et al., 1998). The notion that mutation rates themselves are subject to nat- ural selection has been well documented empirically (Sniegowski et al., 1997, 2000; Baer et al., 2007). It would seem at first that there should be a strong selec- tive pressure to eliminate mutation altogether, due to the high probability that any particular mutation will have deleterious effects. This general reduction principle has been shown to hold in theoretical models pioneered by Sam Kar- lin and colleagues (Karlin and McGregor, 1974; Liberman and Feldman, 1986) for populations at equilibrium in non-changing environments, where there is no cost to increased replication fidelity. However, when higher replication fidelity entails a cost (due to, for example, increased energy allocation during tran- scription), then selection favours an intermediate mutation rate that balances the fitness cost associated with deleterious mutations and the energy savings associated with a higher mutation rate (André and Godelle, 2006). While reducing mutation rates is evolutionarily favoured in a constant en- vironment, were it not for costs, the same is not true in a novel environment. In a novel environment a mutator lineage (one with a higher mutation rate) still carries the burden of an increased deleterious mutation load. However, it now also has a higher probability of experiencing a beneficial mutation. When such a beneficial mutation occurs in a mutator lineage it can pull the mutator allele to a higher frequency (a mechanism referred to as genetic hitch-hiking; Maynard Smith and Haigh, 1974). The above ideas have been confirmed in a number of theoretical models (Kimura, 1967; Leigh, 1970; Johnson, 1999). Because modifiers of mutation rate are primarily subject to indirect selection via their effects at other loci (the cost of replication fidelity being the excep- tion to this rule), modifier dynamics are highly sensitive to their rate of recom- bination with the target loci that are the subject of both mutation and direct selection. Recombination has been shown to weaken the strength of indirect selection acting on mutation modifiers (Kimura, 1967; Leigh, 1970; Sniegowski et al., 1997). This can be understood as follows. Suppose a beneficial mutant arises in a mutator lineage. With a high enough recombination rate, it is proba- ble that the beneficial allele will recombine into the non-mutator lineage during 7 Chapter 2 its spread. In this way the mutator lineage "shares" the benefits of the benefi- cial mutants it creates. Meanwhile deleterious mutations will accumulate on the mutator background more rapidly than on the non-mutator background. While recombination allows some of these mutator lineages to escape the deleterious effects of their mutations, a large number of mutator lineages will be eliminated by purifying selection before any recombination occurs. Consequently, because beneficial mutations have a longer sojourn time within a population, on aver- age, than deleterious mutations, the benefits of producing beneficial mutations are dissipated over time by recombination, while the costs of producing dele- terious mutations are more immediately felt. This intuitive argument has been considered explicitly in models of directional selection (Johnson, 1999). A novel environment is not the only scenario that has been shown capable of selecting for higher mutation rates. For example, a number of theoretical mod- els have shown that increased mutation rates can be advantageous in fluctuating environments, where the direction of selection periodically changes (see Leigh, 1970, 1973; Ishii et al., 1989). In this situation individuals are repeatedly un- der pressure to adapt. Antagonistic co-evolutionary interactions, such as those that occur between hosts and parasites, can also create selection for increased mutation rates. Because adaptive changes in one species often have detrimental effects on the other, hosts and parasites repeatedly create "novel" environments for one another (this is commonly referred to as the "Red Queen" hypothe- sis; Van Valen, 1973). Within host-parasite interactions, it has been shown that mean fitness is optimized by high or non-zero mutation rates, in a manner simi- lar to that found with fluctuating abiotic environments (Nee, 1989; Sasaki, 1994; Haraguchi and Sasaki, 1996). There is also empirical evidence that host-parasite co-evolution can favour increased mutation rates. For example, Pal et al. (2007) recently found that co-evolution with viruses drove up the mutation rates in the bacterium Pseudomonas fluorescens. Host-parasite interactions are commonly mediated through antigen molecules, which are expressed on the surface of parasite cells. For hosts, a pathogen's antigen molecule provides a useful target to aid in its detection and thus elim- ination. Consequently, hosts have evolved sophisticated mechanisms enabling them to detect and respond to a wide array of possible antigen types, while parasites have evolved complex mechanisms allowing them to regularly pro- duce offspring with antigens differing from their own. This process of "antigen switching" in parasites has been well documented empirically (Brannan et al., 1994; Turner, 1997; Frank, 2002), and consequently many of the mechanisms that have evolved in order to make the process more efficient are now well understood. A good summary of many of these is available in Frank (2002). 8 Chapter 2 Frank describes, for example, the process of "gene conversion," whereby some parasites can copy one of several non-expressed archival antigen alleles into a single expressed site. In this way they are able to store many variant antigen alleles in their genome, while only ever expressing a single one. Because each mutant is drawn from a pool of alleles that have been historically exposed to selection, switching in this way is more likely than random mutation to pro- duce viable mutants. Such a strategy can be found in Trypanosoma brucei, which causes "sleeping sickness" in several African mammals. This particular parasite has been shown to carry hundreds of alternative loci in its genome (Pays and Nolan, 1998). For discussions of similar mechanisms see also Donelson (1995) (in Borrelia hermsii); Svard et al. (1998) (in Giardia lamblia); Kusch and Schmidt (2001) (in free-living protozoa). In a theoretical model of host-parasite co-evolution using mean fitness argu- ments, Nee (1989) found that antagonistic interactions can lead to an indefinite escalation of mutation rates in both species, provided that the selection induced by these interactions is strong relative to mutation rates. Haraguchi and Sasaki (1996) later extended this model to consider the effects of deleterious mutations and found that only a small amount of unconditionally deleterious mutation was sufficient to prevent the above reported indefinite escalation of mutation rates. Like Haraguchi and Sasaki (1996), we consider a modifier extension to the model introduced in Nee (1989). Unlike Haraguchi and Sasaki, however, we allow recombination to occur between the mutation modifier and the selected locus. We find that under certain conditions, this addition results in a qualita- tively different outcome, with recombination acting to reduce the evolutionarily attracting mutation rate by distributing the benefits of advantageous mutations. We closely follow the methods of Gandon and Otto (2007), who analyzed a modifier of the rate of recombination between two loci mediating host-parasite interactions. We supplement our analytical model with computer simulations, which allow us to investigate extensions to finite populations, as well as to more alleles. 2.3 Methods and results 2.3.1 Model description We consider two co-evolving haploid species: a host and a parasite. We follow haplotype frequencies through a life cycle consisting of a census, selection, and reproduction. Mating is random within each species, and (except where noted) 9 Chapter 2 population sizes are assumed constant and large enough that drift can be ig- nored. Antagonistic interactions between species are mediated through a single locus in hosts and a single locus in parasites, each with two alleles (Ah/ah in hosts and Ap/ap in parasites). A second locus in each species determines the rate of mutation between Ai and ai, where forward and backward mutation rates are assumed equal. Throughout i denotes the species type: h for hosts, p for parasites. We denote the two alleles at this mutator locus by Mi and mi. The four haplotypes within a species are thus {aimi, aiMi, Aimi, AiMi}. We let xi,j denote the frequency of the jth genotype in species i (labelled in the order given above). Fitness is determined by a matching-alleles model, introduced by Hamilton (1980). Specifically, when the host genotype matches that of the parasite (e.g., Ap parasite and Ah host, or ap parasite and ah host), the parasite experiences a fitness increase of ap > 0 and the host a fitness decrease of ah < 0. The fitness of the jth genotype in species i is then given by wi,j = 1+ ai(z j(xī,1 + xī,2) + (1 z j)(xī,3 + xī,4)) (2.1) where z j is equal to 1 for j 2 {1, 2} and 0 for j 2 {3, 4} and where an overbar denotes the other species type (h̄ = p, p̄ = h). If µi,M (µi,m) is the probability that mutation occurs at the A-locus in an individual of species i carrying the Mi (mi) allele at the M-locus (with µi,m and µi,M both 1/2), then the frequencies after selection and mutation, but before sex and recombination are given by x0i,1 = (1 µi,m)(wi,1/wi)xi,1 + µi,m(wi,3/wi)xi,3 x0i,2 = (1 µi,M)(wi,2/wi)xi,2 + µi,M(wi,4/wi)xi,4 x0i,3 = (1 µi,m)(wi,3/wi)xi,3 + µi,m(wi,1/wi)xi,1 x0i,4 = (1 µi,M)(wi,4/wi)xi,4 + µi,M(wi,2/wi)xi,2 (2.2) where wi (the mean fitness of species i) is given by wi = 4  j=1 wi,jxi,j. In the text, we assume that the modifier allele only alters the mutation rate at the A-locus, but in appendix A.1 we also include a cost to the modifier of producing unconditionally deleterious alleles. Letting yi denote the product of the probability that a haploid individual engages in sexual reproduction and the probability of recombination between the A and the M loci in species i, 10 Chapter 2 genotype frequencies after sex and recombination are given by x00i,1 = x0i,1 yiD0i x00i,2 = x0i,2 + yiD0i x00i,3 = x0i,3 + yiD0i x00i,4 = x0i,4 yiD0i (2.3) where D0i denotes the linkage disequilibrium in species i after mutation and selection. It can be computed as D0i = x0i,1x0i,4 x0i,2x0i,3. Denoting the frequency of A by pi,A, and the frequency of M by pi,M, we have pi,A = xi,3 + xi,4 pi,M = xi,2 + xi,4 (2.4) In the matching-alleles model, allele frequencies typically fluctuate around 1/2 over time. The key simplification that we make is that these fluctuations are relatively small, so that the allele frequencies remain near 1/2. If we define di,A as the departure of the frequency of the A allele from 1/2, so that di,A = (xi,3 + xi,4) 1/2, we may then convert the recursion equations described in equation (2.3) into a new system of recursion equations involving di,A, pi,M, and Di. Assuming that di,A is small, specifically that it is of the order of a small term e, and that the modifier also has a small effect on the mutation rate, such that µi,M µi,m is also of order e, we find that disequilibrium rapidly evolves to a level that is of order e2 (see appendix A.2 for details). With the above assumptions it is possible to simplify the recursion for di,A. Up to order e we get dh,A[t+ 1] dp,A[t+ 1] ! =M dh,A[t] dp,A[t] ! (2.5) where M = (1 2µh,M) ah(2+ah) (1 2µh,M) ap (2+ap) (1 2µp,M) (1 2µp,M) ! . While this approach requires that disequilibrium is small (order e2), it is not necessary that it remains at the steady state values predicted by the current state of the population, as is the case in a quasi-linkage equilibrium analysis (for more information on quasi-linkage equilibrium analyses see Barton, 1995). Thus equation (2.5) can be applied even when selection is large relative to the 11 Chapter 2 rates of sex and recombination. The point dh,A = dp,A = 0 is an equilibrium of equation (2.5), as well as of the full recursion equations. A local stability analysis shows that this point is unstable with complex eigenvalues as long as R = q det(M) = s (1 2µp,M)(1 2µh,M) (1+ ap/2+ ah/2)(1+ ap/2)(1+ ah/2) (2.6) is greater than one, where det(M) is the determinant of the matrixM and R is the magnitude of the leading eigenvalue ofM. When R from is greater than one, the frequencies in a species cycle sinu- soidally outwards with time. However, when R is near one the cycles remain small for extended periods of time. Without loss of generality we let t = 0 denote the time at which the allele frequency in the host first passes 1/2 (equiv- alently dh,A first passes zero). This allows us to express the general solution to (2.5), after simplifying, as dh,A[t] = Rt 1(1 2µh,M) ah(ah + 2) sin[ft] sin[f] dp,A[0] (2.7a) dp,A[t] = Rt 1 q R2 (1 2µh,M)(1 2µp,M)sin[ft+ s]sin[f] dp,A[0] (2.7b) where f denotes the speed of evolutionary cycles (the period is 2p/f), and s is the phase difference between the host and parasite cycles: f = cos 1 (1 µp,M µh,M) R (2.8a) s = cos 1 24 µp,M µh,Mq R2 (1 2µp,M)(1 2µh,M) 35 (2.8b) The asymmetry between (2.7a) and (2.7b) is a consequence of our starting time when dh,A[0] = 0. Recalling that ah < 0 and ap > 0 it can be shown that R increases with the strength of selection (magnitude of ai). Thus, we can see from equation (2.8a) that cycle period decreases as selection becomes stronger in either species. This is illustrated in figure 2.1. Similarly the period decreases as µi,M (the mutation rate) increases in either species. Figure 2.2 shows how the phase difference between host and parasite cycles (s) depends on the mutation rates in the two species. Our results qualitatively confirm those of Nee (1989). When the mutation rates are equal the two species 12 Chapter 2 cycle 90 out of phase, regardless of the strength of selection in either species (see equation 2.8b). When the hosts have a higher mutation rate the cycles are > 90 out of phase, and when the parasites have a higher mutation rate they are < 90 out of phase. If all else is equal, the phase difference s is larger when selection is weak. This is because changes in one species take longer to induce a response in the other species. Interestingly, in this model only small deviations from 90 are possible in the parameter region where cycling occurs (R > 1), as indicated by the small range of the y-axis in figure 2.2. Nee (1989) argued that large phase shifts cannot occur, because they require that at some point in each cycle one of the two species evolves away from the currently favoured allele. For the parameters investigated here (e.g., those in figure 2.1) the phase difference from a 90 shift is always less than what would result from a single generation of evolution in the wrong direction. Thus, the phase differences from 90 that occur in this model result from slight overshooting based on allele frequencies from the previous generation. Having described the dynamics of the selected locus, we turn now to the dynamics of the modifier locus. After selection and sexual reproduction, the change in frequency of the modifier allele is exactly given by Dpi,M[t] = pi,M[t+ 1] pi,M[t] = 4Di[t]aidī,A[t] 2+ ai(1+ 4di,A[t]dī,A[t]) (2.9) This is analogous to equation (29a) in Gandon and Otto (2007). In appendix A.2, we solve for the disequilibrium in species i. We show that when di,A is assumed small (of order e), the effect of the mutation rate modifier is assumed weak (of order e), and the disequilibrium in the previous generation (Di) has reached order e2, the general solution for the disequilibrium is Di[t] = 2(1 yi)(µi,m µi,M)(1 pi,M)pi,M t  t=1 Xt 1i ✓ di,A[t t] + ai(ai + 2)dī,A[t t] ◆ (2.10) where Xi = (1 yi)(1 2µi,M) and the influence of the initial conditions is assumed to have dissipated. Substituting (2.10) into (2.9) yields Dpi,M = 8ai(1 yi)(µi,m µi,M)(1 pi,M)pi,Mdī,A[t] 2+ ai(1+ 4di,A[t]dī,A[t]) ⇥ t  t=1 Xt 1i ✓ di,A[t t] + ai(ai + 2)dī,A[t t] ◆ (2.11) 13 Chapter 2 which is the analog for a mutation modifier of equation (33) for the recombi- nation modifier in Gandon and Otto (2007). Because Xi < 1 whenever there is some sex/recombination or mutation, the sum may be evaluated explicitly us- ing equation (2.7). We can then average the change in modifier allele frequency over one co-evolutionary cycle, ignoring transient dynamics due to the initial conditions. We also average over all possible starting points in the host-parasite cycle. Doing so, we find that the modifier allele M that increases the mutation rate in species i spreads whenever the following is positive sign[ Dpi,M ] ' sign[c1Xi(1 µī,M µi,M)/2 (µī,M µi,M)(Xi (1 2µī,M))] (2.12) where c1 is the positive constant ahap/(2+ ah + ap). The first term in equa- tion (2.12) is always positive and thus always favours higher mutation rates. It is largest when sex/recombination is rare (small yi and large Xi) and when selection is strong (large c1). We may re-write equation (2.12) as sign[ Dpi,M ] ' sign[(R2 (1 2µi,M)(1 µi,M µī,M))Xi + R2(µī,M µi,M)] (2.13) in order to facilitate interpretation. Because R > 1 when there are cycles, the first term must be positive. It tends to dominate when linkage is tight (Xi large) and when the mutation rates are similar (second term small). Conditioning upon the existence of cycles (R > 1), equation (2.13) can become negative and favour lower mutation rates only when µi,M is sufficiently greater than µī,M. It follows that selection will tend to increase the mutation rate in a species, unless that mutation rate is substantially higher than the mutation rate in the other species. This creates a "ratcheting-up" effect, where an increase in the mutation rate in one species creates selection for a corresponding increase in the other species. 2.3.2 Evolutionarily attracting mutation rate Assuming mutation rates start small we expect, based on our reasoning above, that they will increase until the second term in equation (2.13) becomes nega- tive. Thus setting equation (2.13) (or equivalently equation 2.12) equal to zero allows us to find the evolutionarily attracting level of mutation (µ⇤i ), for a given mutation rate in the other species. Defining the positive constants c2,i = c1 2 c1 ✓ 1+ 1 2(1 yi) ◆ + yi 2(1 yi) 14 Chapter 2 allows us to simplify µ⇤i , giving µ⇤i = µī,M + (1/2 µī,M) c2,i + s c22,i + 2 ✓ c1 2 c1 ◆! (2.14) As the terms in parentheses are positive, equation (2.14) shows that the evolu- tionarily attracting mutation rate is higher than that of the antagonistic species. Figure 2.3 plots the evolutionarily attracting mutation rate in hosts (solid line) as a function of the strength of selection in parasites (left panels) and in hosts (right panels). The shaded region indicates where cycling ceases (R < 1) and thus where selection on mutation modifiers becomes neutral. We expect pop- ulations with low initial mutation rates to evolve until they reach µ⇤i , or until mutation rates become large enough that cycles at the selected locus disappear (R 1). In figure 2.3, therefore, we would expect population mutation rates starting at or near zero to proceed upwards until they reach the minimum of the black line or the border of the grey region. 2.3.3 Different generation times A main assumption in our model is that host and parasite generation times are equal. We relaxed this assumption using two methods (described in more detail in appendix A.3). In both cases we varied host generation time, while assuming parasite dynamics to be governed by the same equations as in the main text. We first assumed that all hosts were subject to selection and mutation at each time step, but that only a subset of the individuals reproduced sexually. We then changed the assumption that all individuals were subject to mutation at each time step, and assumed instead that mutations only occurred when individuals reproduced sexually. In both cases we found that, when hosts lived longer, the evolutionarily attracting mutation rate in hosts was higher than it was when generation times were equal. This is because individuals who live longer undergo less recombination per unit time, and thus remain linked to the beneficial mutations they produce for more time steps (see Introduction). 2.3.4 One-species model To help understand the effects of the speed versus the amplitude of evolutionary cycles we considered a simplified one-species model where selection fluctuated deterministically over time (see appendix A.4). We found that the speed of cy- cles was important in deciding whether a modifier of mutation rate was selected for or against (see equation A.4), with faster cycles selecting for higher mutation 15 Chapter 2 rates. In contrast to the results of our two-species model, however, a one-species model revealed that the strength of selection at the selected locus only affected the strength of selection at the modifier locus, not the direction of selection. This qualitative difference occurs as a result of the dependence of cycle speed (f) on selection strength (ah and ap) in the two-species model, whereas in the one-species model cycle speed is fixed. As in the two-species model we also find that higher recombination rates reduce the evolutionary attracting mutation rate in the one-species model. 2.3.5 Numerical simulations We ran computer simulations to investigate separately the effects of two ad- ditional factors: drift in finite populations, and more than two alleles at the selected locus. We considered a Wright-Fisher model with constant population size. Each generation consisted of selection followed by recombination and then mutation. Host-parasite interactions were governed by the same fitness matrix as was used in the analytical model. Mutation rate was first allowed to evolve in only one species at a time. While it is theoretically clearer this way, it may also be biologically realistic in the case where one of the species is more con- strained in its ability to modify mutation rates. All individuals in the evolving population were initialized with mutation rates equal to zero. A single novel mutator allele was introduced at the M-locus in a randomly drawn individual whenever fixation at that locus occurred. The mutation rate of the novel muta- tor was drawn from a Gaussian distribution centered at the current population mean and with a standard deviation of 0.01 (negative mutation rates were set to 0). All simulations were run for 107 generations. For very large populations (109 individuals in each species), where drift is negligible, simulations matched up perfectly with analytical predictions (see solid black points in figure 2.3). In smaller populations (1000 individuals in each species), however, stochastic fluctuations in allele frequency created additional opportunities that favoured the evolution of higher mutation rates (for example, by creating fluctuations when cycles would not be expected based on analytical predictions). Thus, for small populations higher mutation rates evolved (see open circles in figure 2.3) We also used simulations to investigate the effects of higher numbers of al- leles at the selected locus. With k possible alleles at the selected locus, mutation was set so that an allele of type i had a uniform probability of mutating to any of the other k 1 types. There was a significant reduction in the size of the region in which cycling occurred with more alleles (see figure 2.4), and where cycling 16 Chapter 2 did occur, the period was longer. The combination of these two changes led to an overall reduction in the evolutionarily attracting mutation rate. However, the dependence on the strength of selection and recombination exhibited qualita- tively similar results to those for the two-allele case. For example, in figure 2.4a high recombination causes a reduction in the evolutionarily attracting mutation rate into the region where cycles persist. Finally, we investigated co-evolution between hosts and parasites by allow- ing mutation rates to evolve simultaneously in both species. As before, a single novel mutator allele was introduced into either host or parasite populations when either became monomorphic at the mutation rate locus. We focused on the two-allele case at the locus mediating species interactions. When popula- tions were large, escalation of mutation rates in both species led to the eventual cessation of cycles, at which point evolution stopped (figure 2.5a). As before, however, when populations were small stochastic fluctuations drove mutation rates to higher values than expected, based on our analytical model (figure 2.5b). When the strength of selection differed between the species, mutation rates rose to higher levels in the species with the larger selection coefficient. In essence, the species with more at stake kept ahead in the co-evolutionary arms race (figure 2.5). 2.4 Discussion We have used analytical and simulation models to investigate the evolution of mutation rate in a co-evolving host-parasite system. In order to maintain con- sistency with other similar articles we have based our models on the matching- alleles model of host-parasite interactions, where each host allele can be matched by one parasite allele (whether matching causes infection or resistance is imma- terial in the haploid two-allele model because we are free to define the alleles such that Ap parasites are able to infect Ah hosts; Otto and Michalakis (1998)). We followed host and parasite populations over multiple generations, where interactions between the two species occurred at random at each time step. We derived an expression for the evolutionary attracting mutation rate be- tween two alleles, each of whose fitness was dependent on the current genetic composition of the other co-evolving species (see equation 2.14). We found, in accordance with previous literature for directional selection (Kimura, 1967; Leigh, 1970; Sniegowski et al., 1997), that lower mutation rates are expected to evolve with higher levels of recombination (figure 2.3). In order to account for the higher number of deleterious mutations that often accompany increased mutation rates, several authors have imposed additional 17 Chapter 2 fitness costs on individuals with higher mutation rates (Sasaki, 1994; Haraguchi and Sasaki, 1996). In our model individuals with a higher mutation rate suffer a fitness cost at certain points in the cycle through their elevated probability of mutating from the more fit to the less fit allele. However, to maintain consis- tency with other authors, we also considered a general cost to increased mu- tation rates (see appendix A.1). We found, not surprisingly, that higher costs reduced the evolutionarily attracting mutation rate. Because an imbalance in allele frequencies in one species creates selection for an imbalance in the opposite direction in the other species, cyclical fluctua- tions can occur at the selected locus. When these fluctuations are absent allele frequencies at the selected locus converge to a stable polymorphism, at which point selection at the mutation rate locus disappears. We characterized the pa- rameter space in which cycling does or does not occur (see non-shaded and shaded regions in figures 2.3 and 2.4) and showed that large mutation rates can inhibit cycling (a result previously shown by both Seger 1988 and Nee 1989). We, therefore, expect that evolution will lead populations with initially small mutation rates to our predicted evolutionarily attracting mutation rate (µ⇤i from equation 2.14) or until the point where cycles disappear, should this occur first. This prediction was tested with computer simulations, which confirmed our expectations when population sizes were very large (see solid circles in figure 2.3). However, when population sizes were small, higher than expected muta- tion rates evolved. Drift in small populations can create stochastic fluctuations in allele frequencies (Seger and Hamilton, 1988). Such fluctuations extend the region where selection on a modifier occurs, accounting for an increased evolu- tionarily attracting mutation rate. From figures 2.3 and 2.4, it is apparent that the evolutionarily attracting mu- tation rate predicted in our model, when interpreted as a "per site" mutation rate is extraordinarily high. When viewed as an antigenic "switching rate", or a "per trait" mutation rate, however, it is not as unrealistic and is consistent with theoretically predicted switching rates for parasites based on within host dynam- ics (Sasaki, 1994). By not including an explicit cost to high mutation rates in the text, we have implicitly assumed that the modifier has a very localized effect on the mutation rate at the locus mediating the host-parasite interaction. While mutation rates can be highly site specific (Frank, 2002), it is likely that mutation rate modification is an intricate process, and therefore probable that modifiers of mutation rate will have pleiotropic deleterious effects. Accounting for this in our model by including additional costs reduced the predicted mutation rate toward which the system evolved (see appendix A.1), as was the case in both Sasaki (1994) and Johnson (1999). It is also possible that very high mutation 18 Chapter 2 rates would be evolutionarily favoured at loci that mediate host-parasite inter- actions, but that such high rates take a long time to evolve (or may not evolve at all because of genetic constraints), so that real populations may not have reached such high levels. A major assumption in our analytical model is that there are only two possi- ble alleles at the selected locus. Seger (1988) has shown that dynamics can differ qualitatively with multiple alleles. We, therefore, used computer simulations to investigate more than two alleles. We found that the evolutionarily attracting mutation rate decreased, and in accordance with Seger (1988), the region where cycling occurred became smaller with more alleles (see figure 2.4). As in the case for two alleles, mutation rates did not always evolve simply to the point where cycles stopped. This can be seen by comparing points in figure 2.4a and observing that when selection is strong and recombination is high the evolution- arily attracting mutation rate lies in the region where cycling persists. Because hosts and parasites usually have a large number of potential beneficial alleles (Frank, 2002), we would predict, based on the above, that most real populations would exhibit much lower mutation rates than those predicted by a two-allele model. Factors that led to longer cycles decreased the evolutionarily attracting muta- tion rate. Cycles were longer when selection was weak in either species (see fig- ure 2.1), when the difference between host and parasite mutation rates was large (see equation (2.8a)), and when there were more alleles (results not shown). Be- cause cycles in one species create a fluctuating selective pressure in the other species, the observed decrease in the evolutionarily attracting mutation rate with each of the above factors is consistent with previous findings that in fluc- tuating environments the optimal mutation rate is proportional to the inverse of the length of selective episodes (Leigh, 1970, 1973; Ishii et al., 1989). When cycles persisted in our model, we observed that hosts and parasites were usually very close to 90 out of phase. This is consistent with the results of Nee (1989) who argued that large deviations from 90 cannot be sustained, because they would require one of the two species to regularly evolve in the "wrong" direction. While in our model, the phase shift s can theoretically range from to 0 to 180 , we found that cycles disappear well before s becomes too different from 90 and large deviations were, therefore, never observed (see figure 2.2). For large population sizes, cycles also never appeared drastically more or less than 90 out of phase in our simulations. However, in small pop- ulations, where cycles were largely driven by genetic drift, large phase shifts were observed (not shown), as was also noted by Gandon and Nuismer (2009) for a spatially structured population. 19 Chapter 2 For a wide range of parameters, cycles were ephemeral, disappearing once mutation rates had evolved to sufficiently high levels (see figure 2.3). Further- more, the range of parameters where cycles persist became smaller as the num- ber of alleles increased (see figure 2.4). These findings question the legitimacy of the pervasive assumption that loci governing host-parasite interactions exhibit cyclical dynamics. A lack of empirical evidence for recurrent cycles further questions this assumption's legitimacy. To our knowledge, there is only a single empirical study that presents data supporting cyclical dynamics (i.e. changes in allele frequency which eventually return to the same value; several addi- tional studies show change consistent with either cycles or directional selective sweeps). Stahl et al. (1999) used coalescence theory to argue that a 9.8 million year old polymorphism for disease resistance in Arabidopsis thaliana must have been maintained by frequency-dependent selection, and furthermore, showed signs of having historically differed in frequency. While consistent with cycling, their results are also consistent with an alternative explanation; a recent change in the environment could have altered the relative fitness cost of the resistance allele, resulting in a shift in allele frequencies. More long-term empirical data are, therefore, needed to determine whether true cycles occur in host-parasite systems. 20 Chapter 2 0 2 4 6 8 10 0 20 40 60 80 Strength of selection in parasites ( αp ) Pe rio d of e vo lut ion ar y c yc les ( 2 ππ / φφ ) αh −0.25 −0.5 −0.75 −1 Figure 2.1: Period of co-evolutionary cycles (in generations) as a function of the strength of selection in parasites. Other parameters were µh,M = µp,M = 1⇥ 10 3 per generation. 21 Chapter 2 0.00 0.05 0.10 0.15 0.20 0.25 90 95 100 105 110 Mutation rate difference ( µh − µp ) Ph as e dif fe re nc e (1 80 σσ / ππ ) αh −0.25 −0.5 −0.75 −1 Figure 2.2: Phase difference (in degrees) between host and parasite cycles as a function of the difference in their per generation mutation rates. Other parame- ters were µp,M = 1⇥ 10 3 and ap = (1+ ah) 1 1 (for each curve). Curves are only plotted where R > 1. 22 Chapter 2 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● 0 0.1 0.2 0.3 a ψh = ψp = 0 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● b ψh = ψp = 0 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● 0 0.1 0.2 0.3 c ψh = ψp = 0.35 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● d ψh = ψp = 0.35 0 5 10 15 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● 0 0.1 0.2 0.3 e ψh = ψp = 0.5 0.0 0.2 0.4 0.6 0.8 1.0 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● f ψh = ψp = 0.5 Strength of selection in parasites ( αp ) Strength of selection in hosts (− αh ) M ut at ion ra te in h os ts ( µµ h ) Figure 2.3: Evolutionarily attracting per generation mutation rate in hosts as a function of the strength of selection in parasites (a, c, e) and in hosts (b, d, f). The rate at which sex and recombination (yi) occurs is given in each panel. Solid black lines denote the predicted mutation rate corresponding to equation (2.14) and shaded regions indicate where cycling does not occur (R < 1). ah = 1 in panels a, c, e, and ap = 100 in panels b, d, f. Filled and open circles correspond to simulations with population sizes of 109 and 103 (in both hosts and parasites), respectively. µp was 0.002 in all panels (results are very similar unless µp 0). Evoloved mutation rates in hosts are the mean of 50 replicate runs, each of which was averaged over the last 5⇥ 106 of 107 generations. 23 Chapter 2 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● a 0 5 10 15 0 0.1 0.2 0.3 3 alleles ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● b 0 5 10 15 10 alleles ● ● ψh = ψp = 0.0 ψh = ψp = 0.5 Strength of selection in parasites ( αp ) M ut at ion ra te in h os ts ( µµ h ) Figure 2.4: Evolutionarily attracting per generation mutation rate in hosts with multiple alleles at the selected locus. The shaded light grey regions indicate where cycling did not occur in the three (panel a) and ten (panel b) allele model, and the shaded dark grey region indicates where cycling did not occur in the two allele model and is included for reference. Shaded regions were generated using simulations; cycling was considered absent whenever the allele frequency range over the last 500 generations of the simulations fell below 0.005. With two alleles, this region corresponded well to the region where R < 1. For the large population sizes used here (109 hosts and 109 parasites), a clear transition from cycling to no cycling occurs. Recombination rates were set to 0.5 (filled circles) and 0 (open circles), ah = 1, and µp = 0.002. Each point is the mean of 10 replicate runs, each of which was averaged over the last 5⇥ 106 of 107 generations. 24 Chapter 2 0.0 0.1 0.2 0.3 0.0 0.1 0.2 0.3 Index NU LL ● ● ● ● ● ● ● ● ● a Selection stronger in parasites Equal strength of selection Selection stronger in hosts 0.0 0.1 0.2 0.3 Index NU LL ● ● ● ● ● ● ● ● ● b Host Mutation Rate Pa ra sit e M ut at ion R at e Figure 2.5: Co-evolutionarily trajectories in simulations with differing strengths of selection in hosts and parasites. Each curve represents the mean of 100 repli- cate runs. Initial mutation rates were either: 0 in both hosts and parasites, 0 in hosts and 0.15 in parasites, or 0.15 in hosts and 0 in parasites. Points repre- sent final values after 107 generations. The light shaded regions indicate where cycles are not expected in the analytical model (R < 1) for ah = 0.9, ap = 9 and the dark shaded regions indicate where cycles are not expected for both ah = 0.9, ap = 99 and ah = 0.99, ap = 9. Population sizes in both hosts and parasites were 109 in (a) and 104 in (b). When the strength of selection was equal in both species we used ah = 0.9 and ap = 9; when selection was stronger in the parasites we used ah = 0.9 and ap = 99; finally, when selection was stronger in the hosts we used ah = 0.99 and ap = 9. The dotted line has a slope of one and is included for reference. Recombination rates were set to 0. 25 Chapter 3 Ploidy and the evolution of parasitism 3.1 Summary Levels of parasitism are continuously distributed in nature. Models of host- parasite co-evolution, however, typically assume that species can be easily char- acterized as either parasitic or non-parasitic. Consequently, it is poorly under- stood which factors influence the evolution of parasitism itself. We investigate how ploidy level and the genetic mechanisms underlying infection influence evolution along the continuum of parasitism levels. In order for parasitism to evolve, the selective benefits to successful invasion of hosts must outweigh the losses when encountering resistant hosts. However, we find that exactly where this threshold occurs depends not only on the strength of selection, but also on the genetic model of interaction, the ploidy level in each species, and the nature of the costs to virulence and resistance. With computer simulations we are able to incorporate more realistic dynamics at the loci underlying species interactions and to extend our analyses in a number of directions, including finite population sizes, multiple alleles, and different generation times. 3.2 Introduction Understanding the complex ecological and evolutionary interactions between parasites and their hosts has long been a central focus in the biological sciences. This is largely due to the important consequences that advances in this field have had on the development of new strategies for disease and pest manage- ment. The continued need for progress has led to high levels of communication between theoreticians and empiricists, which has helped propel research in both 26 Chapter 3 fields (e.g., Frank, 2002; Galvani, 2003). Consequently, there are numerous the- oretical models covering a wide range of topics, including: the evolution of virulence (e.g., May and Anderson, 1983; Frank, 2002; Galvani, 2003), sex (e.g., Hamilton et al., 1990; Peters and Lively, 2007), recombination and mutation rates (e.g., Gandon and Otto, 2007; M'Gonigle et al., 2009), the evolution of host resistance (e.g., Boots and Bowers, 1999; Miller et al., 2005), and local adaptation (e.g., Gandon, 2002). One typical assumption of theoretical host-parasite models is their treatment of species as either parasitic or non-parasitic (e.g., Sasaki, 1994; Haraguchi and Sasaki, 1996; Gandon, 2002; Day and Proulx, 2004; Nuismer and Otto, 2004; Gandon and Otto, 2007). In other words, models typically operate under the assumption that a species lives strictly as a parasite. While many species do fit this assumption (e.g., those for whom the very completion of their life-cycle depends on the successful infection of a host, such as the plasmodium species that cause malaria), there are numerous examples of species for whom this as- sumption is not appropriate. For example, a number of species from a range of taxonomic groups have been shown to be "facultatively parasitic"; e.g., ciliates (Reynolds, 1936; Thompson and Moewus, 1964), flatworms (Hooge and Tyler, 1999), fungi (Morin et al., 1993), nematodes (Benham, 1974). These species are parasitic if the opportunity arises but are otherwise free-living and capable of reproduction without the aid of a host species. Levels of parasitism should thus be seen as distributed along a continuum in which "completely parasitic" and "completely non-parasitic" define the extreme cases. One question that then arises is, how does evolution occur along this continuum, and what are the main factors that determine whether evolution occurs toward higher or lower levels of parasitism? Empirical work on a number of different host-parasite systems has uncov- ered a variety of genetic mechanisms employed by hosts and parasites to gener- ate the phenotypic variation needed to defend against and invade one another (Frank, 2002). For example, a single allele in flax (Linum usitatissimum) causes resistance to the fungal pathogen Melampsora lini, and a single virulent allele in the pathogen allows infection of both non-resistant and resistant strains of flax (a "gene-for-gene" interaction (Flor, 1942, 1955, 1956)). Host-parasite inter- actions have also been shown to exert strong selection on the underlying genes that modulate species interactions (e.g., favouring changes in expression level (Nuismer and Otto, 2005) or ploidy level (Nuismer and Otto, 2004)). That there are many ways species can interact on a genetic level, and that these interac- tions have been shown to be under selection, suggests that the nature of the genetic interactions between species also exerts a selective force on the degree 27 Chapter 3 of parasitism. Here we ask how ploidy level, an important component of such genetic interactions, influences evolutionary transitions along the continuum from free-living to parasitic life-histories. Using a combination of analytical models and simulations, we examine evo- lution at a locus that modifies the amount of time a facultatively parasitic species spends parasitizing its host species. This is done in the context of each of the three models of host-parasite interactions that are thought to describe a large number of host-parasite systems (Nuismer and Otto, 2004). 3.3 Model summary We consider two interacting species, denoted H and P for hosts and parasites, respectively. The term "parasite" is used loosely here, as the species in question can spend anywhere from 0 to 100% of its time as a parasite. Species interactions are governed by a single locus (from here on referred to as the A-locus) with two alleles in each species (AH and aH in hosts, and AP and aP in parasites). We suppose that parasites spend a proportion of their life cycle parasitizing hosts, and the remaining proportion as free-living organisms. A second locus (from here on referred to as the M-locus or the "modifier" locus) determines how a parasitic individual partitions its time between these two strategies; individuals of genotype i spend a proportion fi of their life cycle as parasites (see table C.1 for a complete list of parameters and their descriptions). We consider here three models of host-parasite interactions. The matching- alleles model (abbreviated MAM) is based on a system of self/non-self recogni- tion (Hamilton, 1980; Frank, 1994; Peters and Lively, 1999; Grosberg and Hart, 2000), as typically occurs in immune systems that develop via the elimination of self-compatible MHC molecules. In this model hosts are susceptible to parasites carrying only alleles that mimic or "match" their own cell signals and are re- sistant to parasites possessing any non-matching alleles. The inverse-matching- alleles model (abbreviated IMAM) is essentially the opposite of the MAM; hosts can defend against parasites carrying any matching alleles and are susceptible to parasites carrying only non-matching alleles (Frank, 1994). This model describes components of the vertebrate MHC system, where host alleles influence the ar- ray of antigen molecules that can be detected. Hosts can only defend against parasites whose antigens they can detect. In the gene-for-gene model (abbrevi- ated GFG), avirulent parasite alleles produce signal molecules that bind to cell surface receptors on resistant host cells, triggering an immune response and thus unsuccessful invasion (Albersheim and Anderson-Prouty, 1975; Gabriel and Rolfe, 1990). Virulent pathogens, however, are able to suppress the pro- 28 Chapter 3 duction of these elicitors and are, therefore, able to invade both resistant and non-resistant hosts. These systems are typically characterized by dominant re- sistance alleles and recessive virulence alleles (Frank, 2002); we shall assume these dominance interactions throughout. Gene-for-gene interactions were first discovered and have since been shown to be quite common in plant-pathogen interactions (Gabriel and Rolfe, 1990). Because we are interested in the effects of ploidy on species interactions, we will consider all combinations of haploid and diploid hosts and parasites. We let {xH,1[t], . . ., xH,k[t]} and {xP,1[t], . . ., xP,l [t]} denote the frequencies of the k host and l parasite genotypes at time t. As a free-living organism, each individual has some basal fitness, which we arbitrarily set to 1. Selection coefficients for other life stages are then measured relative to this fitness. Parasitic individuals that successfully infect hosts experience a fitness gain of aP, while those that encounter resistant hosts experience a fitness loss of bP. If a parasite attempts to find a host, but fails, and if it can no longer reproduce as a free-living organism, then fitness would be lower in both cases. Infection by a parasite is assumed to lower host fitness by aH. We define the indicator variable hi,j to equal 1 if parasites of genotype i can infect hosts of genotype j, and 0 otherwise. Table 3.2 summarizes the infection patterns for each of the models considered here. The fitness of a genotype i parasite at time t is then given by wP,i[t] = (1 fi) + fi k  j=1 (1+ ahi,jP · ( bP)(1 hi,j))xH,j[t], (3.1) and the fitness of a genotype i host is given by wH,i[t] = 1 aH l  j=1 hj,i f jxP,j[t]. (3.2) The above ignores demographic fluctuations and assumes that each individual engages in at most one host-parasite interaction per time-step. Costs of resistance and virulence have been demonstrated in some GFG sys- tems (Tian et al., 2003; Thrall and Burdon, 2003). Without such costs, we would expect the resistant host alleles and/or virulent parasite alleles to spread to fixation. We, therefore, assume that both the virulent parasite allele and the resistant host allele are costly. In hosts we assume that the resistant allele (AH) reduces the fitness of its carriers by an amount cH. In parasites, we consider two types of costs: a conditional cost (cP,c) that only impacts individuals involved in host-parasite interactions (e.g., reduces growth within a host), and an un- 29 Chapter 3 conditional cost (cP,u) that impacts all virulent individuals (e.g., reduces growth within and outside of hosts). The effects of these costs act additively, so that the fitness of a virulent individual of genotype i is reduced by ( ficP,c + cP,u). The frequency of genotype i in species j (j = H or j = P) after selection may then be computed as x0j,i = xj,i[t]wj,i[t] w̄j[t] , (3.3) where w̄j[t] = Âi xj,i[t]wj,i[t] is the mean fitness of species j (the sum is taken over all genotypes). While we largely focus on the effects of ploidy and the model of genetic interaction, it is worth mentioning that the above model also captures possi- ble ecological changes in the opportunity for parasitism; if the environment clearly favours one life-history strategy over another (as may occur, for exam- ple, when a new host species becomes available), then parasitism would be expected to evolve, regardless of the genetic architecture underlying species interactions. This possibility would be captured by high values of aP (large advantages of successful invasion) and potentially low values of bP (weak host defenses against the parasite). In cases where the environmental forces favour- ing parasitism are not absolute, however, our analysis will help predict how the underlying genetics shapes the course of evolution. We let yH and yP denote the proportion of hosts and parasites, respectively, that undergo sexual reproduction at each time step, and we assume that the remaining individuals consist of surviving parents or asexual offspring. We let x00H,i and x00P,i denote the frequency of genotype i individuals in hosts and para- sites, respectively, formed through random mating within the parental genera- tion after selection. In both hosts and parasites all sexual individuals contribute their gametes to a general gamete pool, out of which offspring are selected at random. Recombination between the modifier and the A-locus occurs dur- ing meiosis in parasites at rate r. After reproduction, genotype frequencies in species j are then given by xj,i[t+ 1] = (1 yj)x0j,i + yjx00j,i . (3.4) 3.4 Analytical results We make the assumption that selection is weak (aH, aP, and bP are all on the same order as some small term e), and that most individuals are sexual in both species (yH and yP are on the order of 1-e; this assumption is relaxed in the Appendix). We also assume that the modifier has a small effect (i.e. we set 30 Chapter 3 fM = fm + DM or fMm = fmm + DMm and fMM = fmm + DMM and then assume that the D's are also of order e). Performing a change of variables allows us to describe the system in terms of the departure from a frequency of 0.5 at the A locus in each species (dH [t] in hosts and dP[t] in parasites), the frequency of the modifier in parasites (pM[t]), and several higher order association measures, such as the departure from Hardy-Weinberg equilibrium and linkage disequi- librium (as defined in Barton and Turelli, 1991). A quasi-linkage equilibrium analysis (Barton and Turelli, 1991) showed that all genetic associations are of order e2 or higher, and that changes in allele fre- quency at the M-locus are governed largely by terms of order e, which describe differences in fitnesses of the different genotypes. Specifically, a modifier M of parasitism level will spread only if the difference between the marginal fit- nesses of alleles M and m, which we denote by w̄diff = w̄M w̄m, is positive. The expressions for w̄diff for an allele that increases parasitism for the four com- binations of host/parasite ploidy levels are given in table 3.3. As is typical in models of host-parasite co-evolution, dynamics at the A- locus are characterized by cyclical fluctuations (figure 3.1). In the matching- alleles and the inverse-matching-alleles models, these cycles are symmetric about an allele frequency of 0.5 (equivalently, about dH = 0 and dP = 0). By assuming (for these cases) that these cycles are small (e.g., that both dH and dP are also on the order of the small term e), we are able to find simple conditions under which selection favours increased levels of parasitism (table 3.4). It is clear from the expressions in table 3.4 that the fitness effects of match- ing versus not matching the genotype of the host, aP/bP, must be sufficiently beneficial for parasites to adopt a less free-living life-cycle in both MAM and IMAM. However, where this threshold occurs depends on both the model of genetic interactions and the ploidy level of each species (figure 3.2). In gen- eral, the matching-alleles model tends to favour parasitism more strongly than the inverse-matching-alleles model (compare figure 3.2A to figure 3.2B), mainly because it is easier for a parasite to mimic hosts that are heterozygous diploid (MAM) than to evade detection by them (IMAM). In both MAM and IMAM, the transition to parasitism occurs over a broader range of parameters when the par- asite is haploid, because such parasites express only one antigen allele (compare solid to dashed lines in figure 3.2). The role of host ploidy is more complicated, however. Diploidy allows for the appearance of heterozygous hosts that are in- fected by any type of parasite in MAM, but resistant to every type of parasite in IMAM. Thus host diploidy favours (disfavors) the evolution of parasitism in MAM (IMAM) (compare thick to thin lines in figure 3.2). Because cycles in the gene-for-gene model are not typically centered around 31 Chapter 3 0.5 we take a slightly different approach in this case. We first solve for the equilibrium dH and dP, and then substitute these into the expressions for w̄diff. Assuming weak selection, we are again able to simplify the expressions for w̄diff. For all ploidy combinations, we find w̄diff = (aP cP,c) f 2aH cHcP,u f 2aH , (3.5) where f denotes the resident parasitism level ( f = fm for haploid parasites and fmm for diploid parasites). Unlike for MAM and IMAM models, there are no major effects of ploidy on the evolution of parasitism in the GFG. While quali- tative dynamics at the A-locus differ between cases, parasitism is favoured for the same combinations of selection and cost parameters across all ploidy levels. This contrasting result for GFG is a consequence of our empirically motivated assumption of complete dominance. With both the resistant allele in hosts and the virulent allele in parasites completely dominant, the two species are essen- tially composed of only two types, and thus effectively interact as haploids. In contrast to ploidy, the nature of costs of virulence are critically important to the evolution of parasitism in the GFG (figure 3.3). Consider setting cP,u equal to zero. With just conditional costs (cP,c), we find, w̄diff = aP cP,c , (3.6) and thus parasitism should evolve whenever the benefits to successful invasion, aP, are greater than the conditional cost of the virulent allele cP,c. In contrast, when cP,c equals zero, w̄diff = ✓ aP cP,u cHf 2aH ◆ . (3.7) Here the cost term is weighted by 1/ f 2. With lower resident parasitism levels (smaller f values), a larger selective benefit to parasitism (aP) is required in or- der for selection to favour further increases in parasitism. This makes it exceed- ingly difficult for parasitism to evolve from initially low levels when costs are unconditional. Intuitively, because the unconditional cost is paid by all virulent individuals, it is unlikely that any fitness gains acquired through parasitism will sufficiently compensate for the costs of virulence when the chance of infecting a host is low. 32 Chapter 3 3.5 Simulation model summary We ran computer simulations to investigate the robustness of our model to vio- lations of its assumptions, such as small cycles (figures B.1, B.2), weak selection (figure B.3), infinite population sizes (figure B.4), and high rates of sexual repro- duction (figure B.5). We also investigated the effect of differences in generation times (figure B.6) and of multiple alleles (figure B.7) on the evolution of par- asitism. In each case we employed a Wright-Fisher model with constant and finite population size. Each time step consisted of selection followed by sex and recombination (with r=0.5). Because population sizes were finite, mutation between alleles at the A-locus was necessary to ensure that allelic variation at this locus was not permanently lost. Mutation between alternative alleles at the A-locus occurred in both species at rate µ per generation. In order to investigate the evolution of parasitism we initialized populations with low levels of parasitism ( f = 0.1), and tracked evolution at the modifier lo- cus. We also ran simulations initialized with high levels of parasitism ( f = 0.9), but because the final level of parasitism attained was typically similar, these results are not presented. Where this change did affect the final outcome, we provide a more detailed discussion. All individuals were initially identical at the modifier locus, and whenever fixation occurred, a novel modifier allele was introduced at low frequency (we used 0.01), and in linkage equilibrium with the A-locus. The parasitism level ( f ) corresponding to the novel modifier was drawn from a Gaussian distribution centered on the current level of parasitism with a standard deviation of 0.1 (parasitism levels were re-drawn if they fell out- side the range [0, 1]). While introducing a mutant allele into the population at linkage equilibrium is not biologically realistic, it eliminates unwanted artifacts that may result from biased initial associations between the modifier and the A-locus. It is also worth mentioning that the initial frequency of the new modi- fier and the standard deviation used to draw new mutants did not qualitatively affect results, but they did affect the speed of the simulations. We examined a number of extensions to our model (see Appendix). Most extensions had little effect on our results. Here we focus on only the simplest and most informative extensions. Unless specified otherwise, simulations were run for 106 generations and initial frequencies at the A-locus in each species were drawn independently from a uniform distribution. 33 Chapter 3 3.5.1 MAM/IMAM simulation results Small cycles We begin by examining the simpler case when cycles at the A-locus are small in amplitude. To constrain cycle size in the simulations we increased the muta- tion rate, which pushes allele frequencies toward intermediate values and thus dampens cycles (exposure to multiple parasites per time step has also been shown to dampen cycles (Lively, 2010)). With intermediate to high mutation rates, cycles were characterized by smooth sinusoidal curves (figure 3.1C, 3.1D ; note that for very high mutation rates cycles were absent altogether, as in figure 3.1A, 3.1B). When cycles were absent or small, increased levels of parasitism evolved as predicted in table 3.4 (columns 1 and 2 in figure B.1). Large cycles Large amplitude cycles are observed with lower mutation rates (figure 3.1E, 3.1F). Mutation to other serotypes in Borrelia hermsii has been estimated to oc- cur at a rate somewhere between 10 3 and 10 4 per generation (Stoenner et al., 1982), and new variant surface glycoproteins arise in trypanosomes at a rate somewhere between 10 2 and 10 6, per cell doubling time (Turner and Barry, 1989). We thus set µ = 10 5 to investigate realistic mutation rates. With large amplitude cycles a few cases did not match the small-cycle analytical approxi- mations in table 3.4 (figure 3.4), although the large-cycle conditions in table 3.3 continued to hold, given the observed dynamics for dH and dP (results available upon request). We will describe these cases in turn. When parasites were diploid and hosts were haploid, large cycles led to a reduction in the size of the region where parasitism evolved (figure 3.4C and 3.4D). Because in MAM and IMAM heterozygous parasites could not invade either haploid host, parasites responded slowly to changes in allele frequency in the hosts. Consequently, the proportion of time parasites spent "losing" the host-parasite arms race grew as cycle size increased (left column in figure 3.1), and thus the region where parasitism was favoured shrunk. When parasites were haploid and hosts were diploid, the opposite scenario occurred. Here it was the diploid hosts that were slow to respond to allele frequency changes in the haploid parasites. Furthermore, because heterozy- gous hosts in IMAM are more resistant than homozygous hosts, cycles tended to dampen (remaining near di = 0), whereas a slow coevolutionary response in hosts was observed in MAM (right column, figure 3.1). The region where parasitism was favoured thus grew slightly with MAM (figure 3.4E). 34 Chapter 3 When both species were diploid, general conclusions could not be drawn about which species would lag behind in the arms race. Unlike the previous comparisons, whether the region where parasitism was favoured slightly grew or shrunk depended more sensitively on the strength of selection in hosts (aH, see figures B.2 and B.3). 3.5.2 GFG simulation results Because high mutation rates drive allele frequencies to 1/2, which is not gener- ally the equilibrium in the gene-for-gene model, we only consider low mutation rates and thus large cycles in this case. With conditional costs our simulations exactly matched our predictions, and parasitism evolved whenever the fitness benefit of successfully invading a host, aP, was greater than the cost of the virulent allele, cP,c (figures 3.3A, B.8). As predicted for unconditional costs, the initial level of parasitism present in the population strongly affected which parameter combinations favoured further evolution of parasitism (figures 3.3B, B.9). In large populations (N = 106 individuals) and initially low levels of par- asitism ( f = 0.1), increased parasitism never evolved, as expected. With high initial levels of parasitism, however, evolution of a more parasitic life history was possible. In other words, the system exhibited bistability. Interestingly, in regions where the evolution of a free-living life cycle was expected, the GFG system would often converge to an M/m polymorphism fixed for allele a. That is, the initial mA/ma polymorphism involving a costly virulent allele and a sensitive allele was replaced with an Ma/ma polymorphism involving sensitive alleles with higher and lower levels of parasitism (explaining why the regions in figures B.8 and B.9 were grey rather than white). With unconditional costs the effect of initial conditions described above disappeared altogether in small populations (figure B.10). Stochastic fluctuations in allele frequency at the in- teraction locus, combined with drift at the modifier locus, allowed occasional excursions into the parameter space in which further evolution of parasitism became advantageous. 3.6 Discussion We have used analytical and simulation methods to investigate the evolution of parasitism in a pair of co-evolving species. Our results provide an initial characterization of how genetic architecture affects selection on life-history in antagonistic species interactions. By and large, the evolution of parasitism depends on the mean fitness of 35 Chapter 3 allelic variants at a locus governing how much time a species spends as a para- site and is not strongly influenced by genetic associations. By comparing mean fitness of these allelic variants we were able to characterize the conditions un- der which high levels of parasitism were expected to evolve. While the fitness effects of matching or not matching the genotype of the host had to be suffi- ciently beneficial in order for parasites to adopt a more parasitic life-cycle, the exact threshold depended on both the model of genetic interactions and, in most situations, the ploidy level of each species. In situations where hosts are only able to defend against parasites for which they have the correct allele, as with IMAM, hosts that carry a larger suite of alleles (diploids) or parasites that carry few alleles (haploids) tend to thrive. Thus lower ploidy levels in either species tend to increase the benefits to parasitism. In contrast, in situations where par- asites must match host genotypes in order to invade (e.g., MAM), diploid hosts can be infected by a greater number of parasite types, and thus diploidy in hosts tends to favour parasitism, while haploidy in parasites is again most con- ducive to further evolution of parasitism. With GFG interactions ploidy had little impact on the evolution of parasitism because of the complete dominance assumed. The above predictions were derived under a number of assumptions, the most significant being intermediate allele frequencies at the locus governing host-parasite interactions (i.e., small cycles). Using simulations we were able to investigate our model's behaviour when no such constraints were imposed on allele frequencies. The predictions based on small cycles were altered slightly under some conditions (panels C-E, G, H in figure 3.4), although the main qual- itative results continued to hold under all conditions. The differences from our predictions occurred mostly when host and parasite fluctuations were not 90 out of phase with one another (figure 3.1). Typically this occurred when het- erozygotes of one species had low fitness (e.g., hosts in MAM and parasites in IMAM). These low fitness heterozygotes reduced the efficacy of selection in this species, as beneficial alleles, when rare, were found almost exclusively in the heterozygous form. As a result this species responded slowly to changes in allele frequency in the other species. This meant that more time was spent in a population configuration that favoured the faster responding species, and thus the region of parameter space where parasitism evolved was shifted in favour of that species. Violations of our other main assumptions (infinite population sizes, weak selection, and primarily sexual populations) were also tested using simulations and shown to have only minor effects (see Appendix). In nature, parasites typically have much shorter generation times than their hosts, and furthermore, many host-parasite interactions are governed by more 36 Chapter 3 than two alleles (e.g., trypanosomes are known to possess hundreds of allelic antigen variants (Van Der Ploeg et al., 1982)). Using simulations we investigated how these extensions changed our general conclusions. While neither led to any qualitative changes across ploidy combinations, more alleles at the interaction locus had significant and opposite effects in the matching-alleles and inverse- matching-alleles models. Because higher genetic diversity among hosts with more alleles makes them resistant to a larger number of parasites in MAM, more alleles were less conducive to the evolution of parasitism. Similarly, with MAM, high diversity in parasites tends to help hosts recognize their parasites as genetically distinct. The opposite held true in the IMAM, where greater genetic diversity in hosts allows parasites to invade a greater proportion of host genotypes and greater genetic diversity in parasites allows them to remain undetected by more host genotypes. Thus, the more alleles segregating at the genes mediating host-parasite interactions, the more conducive IMAM systems are to the evolution of parasitism. Another factor found to have a large influence on the evolution of para- sitism was the nature of the costs to virulence in the GFG models. Interestingly, conditional costs were much more conducive to the evolution of parasitism. When parasitism is rare, unconditional costs of virulence typically outweigh the benefits of being parasitic and result in the spread of sensitive parasites and resistant hosts, which prevents the evolution of further parasitism. Had uncon- ditional costs been weak enough to allow parasitism to increase when low, then they would have been too weak to prevent the fixation of virulent alleles once parasitism levels were high. In the absence of factors such as strong genetic drift, which may stochastically shift parasitism levels upward, a predominantly free-living life history is thus expected with substantial unconditional costs of virulence. Previous theoretical work has shown that transitions between haploidy and diploidy are expected as a consequence of host-parasite interactions (Nuismer and Otto, 2004). In particular, haploidy is most favoured in parasites because of the advantage of reducing antigenic expression to a single allele, while diploidy is more often favoured in hosts because of the advantage (in many cases) of het- erozygous hosts being able to recognize multiple parasites. In accordance with the above theoretical predictions, a survey of empirical data revealed an associa- tion between ploidy and life history (Nuismer and Otto, 2004); parasitic protists are three to four times as likely as non-parasitic protists to be haploid. This pattern would, however, be consistent with either parasites evolving more hap- loid life-cycles (Nuismer and Otto, 2004), or haploids evolving more parasitic life cycles (herein). Indeed, if transitions in parasitism occur more frequently 37 Chapter 3 than transitions in ploidy, transitions in parasitism may be more important in explaining the association between parasitism and haploidy. Some groups of species today are almost wholly parasitic (e.g., Apicom- plexa), while others contain a mixture of both free-living and parasitic individ- uals (e.g., dinoflagellates) (Als et al., 2004; Moran and Wernegreen, 2000), and many are wholly non-parasitic. In groups such as dinoflagellates, the ability to photosynthesize and thus produce ones own food may make the switch between parasitic and free-living life-styles relatively easy, whereas in other groups it appears that the ability to regain a free-living lifestyle has been altogether lost (e.g., no Borrelia sp. proliferating in an environment outside of a vertebrate or invertebrate host has been observed (Barbour and Hayes, 1986)). A comparative phylogenetic analysis of closely related groups of species, which differ in their proportions of parasitic species, would provide additional insight into exactly what sorts of traits facilitate acquisition or loss of parasitism, and furthermore, just how common such transitions have been. There are a number of worthwhile extensions to our model. Ample em- pirical evidence suggests that many, if not most, host-parasite interactions are governed by more than a single locus (May and Anderson, 1983). For exam- ple, the brown planthopper (Nilaparvata lugens), a pest on rice in South East Asia, was originally assumed to be engaged in a GFG interaction, but it has since been shown to contain several biotypes, each determined by different co- adapted gene complexes (Thompson and Burdon, 1992). Extending our model to include multiple interaction genes would allow us to consider the buildup of the co-adapted gene complexes that facilitate life-history transitions. Fur- thermore, the model presented here assumes that some level of parasitism is already present or that at least the genetic architecture is already in place for proper parasitic invasion of hosts. De novo evolution of parasitism realistically requires more than a single mutational event, perhaps mediated by intermedi- ate stages involving mutualistic or trophic interactions. More detailed models on these early stages could provide insight into how parasitic life styles have evolved out of non-parasitic ones. 38 Chapter 3 Table 3.1: Description of parameters. Parameter Description fi proportion of time spent as a parasite by genotype i xi,j frequency of genotype j in species of type i (i = H or P) aP fitness gained by a parasite that successfully infects a host bP fitness lost by a parasite that attempts but fails to infect a resistant host aH fitness lost by hosts when they are infected hi,j indicator variable defined to equal 1 if parasites of genotype i can infect hosts of genotype j, and 0 otherwise wi,j fitness of genotype j in species of type i cH cost of the resistant allele in hosts (GFG only) cP,c conditional cost of the virulent allele in parasites (GFG only) cP,u unconditional cost of the virulent allele in parasites (GFG only) yi proportion of species of type i that reproduce sexually r recombination rate in parasites DM effect size of the modifier (haploid parasites) DMm effect size of the modifier when present in heterozygotes (diploid parasites) DMM effect size of the modifier when present in homozygotes (diploid parasites) di deviation from a frequency of 0.5 at the A-locus in species i pM frequency of the modifier in parasites w̄M, w̄m marginal fitnesses of alleles M and m w̄diff difference between marginal fitnesses (i.e., w̄M w̄m) µ mutation rate at the A-locus in both species (simulations only) Table 3.2: Invasion matrices. Each entry represents the outcome of interactions in the three models in the following order {matching-alleles, inverse-matching- alleles, gene-for-gene}. I is used to denote infection (hi,j = 1 in equations (3.1) and (3.2)) and R resistance (hi,j = 0). Both haploid and diploid hosts and para- sites are included in the table. Host AH or AHAH AHaH aH or aHaH Parasite AP or APAP {I,R,I} {I,R,I} {R,I,I} APaP {R,R,R} {I,R,R} {R,R,I} aP or aPaP {R,I,R} {I,R,R} {I,R,I} 39 C h apter 3 Table 3.3: The fitness advantage of the modifier allele, w̄diff = w̄M w̄m, when genetic associations are weak. We have dropped a factor DM from the haploid parasite cases and (pM(DMM DMm) + (1 pM)DMm) from the diploid parasite cases; these terms can be interpreted as the average effect size of the modifier. Model Host Parasite w̄diff ploidy ploidy MAM 1 1 (aP bP)/2+ (aP + bP)(2dHdP) 1 2 (aP 3bP)/4+ (aP + bP)(2dHdP + d2P) 2 1 (3aP bP)/4+ (aP + bP)(2dHdP d2H) 2 2 (5aP 3bP)/8+ (aP + bP)(4dHdP(1+ dHdP) 3d2H + d2P)/2 IMAM 1 1 (aP bP)/2 (aP + bP)(2dHdP) 1 2 (aP 3bP)/4 (aP + bP)(2dHdP d2P) 2 1 (aP 3bP)/4 (aP + bP)(2dHdP d2H) 2 2 (aP 7bP)/8 (aP + bP)(4dHdP(1 dHdP) d2H d2P)/2 GFG 1 1 (3aP bP)/4+ (aP + bP)(dP dH(1 2dP))/2 cP,c(1+ 2dP)/2 1 2 (7aP bP)/8+ (aP + bP)(2(1+ 2dH)(1 dP)dP dH)/4 cP,c(3/4+ (1 dP)dP) 2 1 (5aP 3bP)/8+ (aP + bP)(3dP + 2(1 2dP)(d2H dH))/4 cP,c(1+ 2dP)/2 2 2 (13aP 3bP)/16+(aP+ bP)(d2H(1 2dP)2 dH(1 2dP)2+ 3(1 dP)dP)/4 cP,c(3/4+ (1 dP)dP) 40 Chapter 3 Table 3.4: Invasion condition for a modifier that increases the level of parasitism in MAM and IMAM, assuming small cycles around allele frequencies of 1/2 (or dH = dP = 0). Model Host Parasite Invasion ploidy ploidy condition MAM 1 1 aP > bP 1 2 aP > 3bP 2 1 aP > bP/3 2 2 aP > 3bP/5 IMAM 1 1 aP > bP 1 2 aP > 3bP 2 1 aP > 3bP 2 2 aP > 7bP 41 Chapter 3 0 10 20 30 40 50 0.48 0.50 0.52 AHost Parasite 0 10 20 30 40 50 B 0 10 20 30 40 50 0.0 0.5 1.0 C 0 10 20 30 40 D 0 50 100 150 200 250 300 350 0.0 0.5 1.0 E 0 50 100 150 200 250 300 F Generation Fr eq ue nc y o f A Haploid Hosts, Diploid Parasites Diploid Hosts, Haploid Parasites µ = 10−1 µ = 10−3 µ = 10−5 Figure 3.1: Sample trajectories from simulations (after a burn-in period) in the matching-alleles model with complete parasitism ( f = 1). Column labels indi- cate ploidy levels and row labels mutation rates. The background is coloured grey to indicate when the parasite is "losing" the arms race with the host, and white when it is "winning". Other parameters were aP = 0.05, bP = 0.05, aH = 0.05,yH = yP = 1, r = 0.5, and population sizes were 106 in both species. 42 Chapter 3 Index NU LL 0 0.05 0.1 0 0.05 0.1 A Index NU LL 0 0.05 0.1 B αP βP MAM IMAM Figure 3.2: Invasion conditions in the matching-alleles (A) and inverse- matching-alleles model (B). Solid (dashed) lines correspond to haploid (diploid) parasites, and thick (thin) lines to haploid (diploid) hosts. For a given case, par- asitism is expected to evolve when selection is such that the point (aP, bP) lies below the corresponding line. The slopes of these lines can be inferred from the invasion conditions in table 3.4. Note that there are two lines with the same slope in panel (B). Index NU LL 0 0.05 0.1 0 0.05 0.1 A cP,c 0.01 0.05 Index NU LL 0 0.05 0.1 B cP,u 0.01 0.05 αP βP GFG Conditional Costs GFG Unconditional Costs Figure 3.3: Invasion conditions in the gene-for-gene model with conditional (A) and unconditional costs (B). For a given case, parasitism is expected to evolve for all ploidy combinations when aP lies to the right of the plotted line. The three curves for each unconditional cost in panel (B), from left to right, correspond to different initial parasitism levels of f = 0.9, 0.7, and 0.5. aH = 0.05 and cH = 0.01. 43 Chapter 3 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 A B C D E F G H _P `P MAM IMAM Haploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Haploid Parasites Diploid Hosts / Diploid Parasites 0.0 0.5 1.0 Figure 3.4: Evolutionary convergent level of parasitism ( f ) as a function of the strength of selection in parasites under the matching-alleles model (column 1) and the inverse-matching-alleles model (column 2). The mutation rate was (µ = 10 5), and thus large cycles occurred in all cases except in panel F. Dashed lines denote the analytical invasion condition assuming small cycles (table 3.4). Cells are shaded based on the mean level of parasitism present in the population after 106 generations of evolution in a single simulation (darker = higher, see grayscale in panel H). Initial frequencies at the A-locus were randomly drawn for each cell. Different ploidy combinations are indicated on the right hand side. Other parameters were as in figure 3.1. 44 Chapter 4 Sexual selection enables long-term coexistence despite ecological equivalence 4.1 Summary Empirical data indicate that sexual preferences are critical for maintaining species boundaries (Eberhard, 1985; Seehausen and Van Alphen, 1999; Gray and Cade, 2000; Wilson et al., 2000; Irwin et al., 2001; Huber, 2003), yet theoretical work has suggested they can play only a minimal role in maintaining biodiversity on their own (Turner and Burrows, 1995; Panhuis et al., 2001; Van Doorn et al., 2004; Johansson and Ripa, 2006; Weissing et al., 2011). This is because long-term coexistence within overlapping ranges is thought to be unlikely in the absence of ecological differentiation Weissing et al. (2011). Here we challenge this widely held view by generalizing a standard model of sexual selection to include two ubiquitous features of populations with sexual selection: spatial variation in local carrying capacity and mate-search costs in females. We show that, when these two features are combined, sexual preferences can single-handedly main- tain coexistence. Remarkably, coexistence can occur for spatial variation in lo- cal carrying capacity that is so slight that it might go unnoticed empirically. This is the first study to demonstrate that sexual selection alone can promote the long-term coexistence of ecologically equivalent species with overlapping ranges, thus providing a novel explanation for the maintenance of biodiversity. 45 Chapter 4 4.2 Introduction A central objective of evolutionary ecology is to understand the mechanisms that allow species to coexist. One such mechanism is ecological differentiation. By occupying different niches, species in overlapping ranges are able to reduce direct competition among one another (Schluter, 2000). While there are numer- ous examples of closely related species occupying different ecological niches, many recently diverged and coexisting taxa are known to differ most dramati- cally in their secondary sexual characters, exhibiting few, if any, ecological dif- ferences (Eberhard, 1985; Seehausen and Van Alphen, 1999; Gray and Cade, 2000; Wilson et al., 2000; Irwin et al., 2001; Huber, 2003). It seems, therefore, that sexual selection is an important mechanism for maintaining coexistence. Indeed, models of sexual selection have shown that populations of choosy fe- males and their preferred males can arise and, under various conditions, form reproductively isolated mating groups (Fisher, 1930; Lande, 1981, 1982; Kirk- patrick, 1982; Seger, 1985). However, because sexual selection does not lead to ecological differentiation, species differing only in their mating preferences compete for the same ecological niche. This has traditionally led to the conclu- sion that, if their ranges overlap, one of these species will eventually displace the other (Turner and Burrows, 1995; Panhuis et al., 2001; Van Doorn et al., 2004; Johansson and Ripa, 2006; Weissing et al., 2011). Coexistence can become less difficult when species are able to reduce their range overlap. Sexual selection provides a natural mechanism whereby they may accomplish this. Any process that creates spatial variation in female pref- erences indirectly also creates selection on male display traits, locally favouring those males that are most preferred by the local females. As a consequence, spa- tially segregated mating domains, characterized by the co-occurrence of match- ing display and preference traits, can emerge from populations with an initially random spatial distribution. Once segregated, interactions between different mating types are limited to individuals at the peripheries of these domains. In finite populations, however, the mating domains may shrink or grow, and the in- terface between them may drift randomly in space. Such fluctuations eventually lead to one mating domain replacing all others (Fig. 4.1A, C). In a pioneering study, Payne and Krakauer (Payne and Krakauer, 1997) argued that lower dis- persal in males with better mating prospects facilitates spatial segregation and maintains coexistence. In finite populations, however, such mating-dependent dispersal fails to stabilize coexistence (Fig. C.3). Given these difficulties associ- ated with sexual selection, a recent review concluded that sexually divergent, but ecologically equivalent, species cannot coexist for significant lengths of time 46 Chapter 4 (Weissing et al., 2011). 4.3 Model summary Here we report model results that suggest the contrary and demonstrate that sexual selection can promote long-term coexistence, even without any ecological differentiation. Building on a standard model of sexual selection (Kirkpatrick, 1982), we develop an individual-based model, which naturally allows us to con- sider finite populations and to examine the long-term fate of species differing only in their secondary sexual characters in an ecologically neutral context. We assume a simple genetic structure with individuals characterized by two un- linked haploid loci: the first locus (with alleles Q and q) governs a display trait that is expressed only in males, while the second (with alleles P and p) gov- erns a preference trait that is expressed only in females (we later generalize this to more than two alleles and quantitative mating traits; see Supplemen- tary Information (SI) and Fig. 4.4D–E). Because we are interested in coexistence, and not speciation, we assume that the genetic variation at both loci is already present, for example, due to recent migration from allopatric ranges. All else being equal, females bearing a P (p) allele prefer (Kirkpatrick, 1982; Seger, 1985; Payne and Krakauer, 1997) to mate with males carrying a Q (q) allele by a factor a, and a female's preference for a given male attenuates with increasing distance between them. Likewise, competition decreases as the distance between indi- viduals increases. Competition is assumed to reduce an individual's probability of surviving until it reaches reproductive maturity. Importantly, hybrids suffer no intrinsic fitness costs, other than potentially carrying mismatched preference and trait alleles. Further model details are provided in the SI. Mating domains can be lost either through movement of the interface be- tween them (as described above), or when individuals of one mating type colo- nize the domain of another mating type. In particular, because selection at the preference locus disappears when there is no variation at the display locus, for- eign preference alleles may drift into regions with low variation in male display alleles, eventually leading one mating type to displace the other. Loss of mating domains can, however, be prevented by the inclusion of two features ubiquitous in populations experiencing sexual selection: spatial variation in local carrying capacity and mate-search costs in females. Spatial variation in carrying capacity is present in most, if not all, biological systems (see Figs. 4.1 and 4.4 and the SI for model details). Mate-search costs occur if a female spends time and energy looking for a suitable mate and rejecting non-preferred males, thereby reducing her ability to invest in offspring. To account for such costs we assume that the 47 Chapter 4 fecundity of a particular female increases from 0 to a maximum level with the local density of available males, weighted according to her preference (SI). 4.4 Results and discussion Our model confirms the longstanding view that sexual selection in homoge- neous spatial models, without mate-search costs, does not facilitate coexistence and can, in fact, hasten the demise of species (compare Fig. 4.2A to 4.2B). Spatial variation in local carrying capacity, on its own, also has little, if any, effect in stabilizing populations (compare Fig. 4.2B to 4.2C). Sexual selection with mate- search costs slightly prolongs coexistence in a spatially homogeneous environ- ment by helping to prevent mixing of the mating domains, but this effect is weak (Fig. 4.2D). However, in an environment with spatial variation in local carrying capacity, sexual selection with mate-search costs dramatically increases coexis- tence times (compare Fig. 4.2E to Fig. 4.2B and also Fig. 4.1A, C to Fig. 4.1B, D). In this case, mate-search costs curb the neutral drift of preference alleles, thus preventing the dilution of mating domains, while areas of high local carrying capacity provide spatial "anchors", stabilizing the mating domains at a rela- tively constant size and location (Fig. 4.1B, D). This anchoring is a consequence of the net flow of migrants that occurs each generation from the high-density regions into the low-density regions. While neither spatial variation in local carrying capacity nor mate-search costs suffice on their own to stabilize populations, surprisingly little of both can be enough to ensure the long-term persistence of divergent mating types (Fig. 4.3). When mate-search costs in females are high, coexistence can be main- tained with less than 20% spatial variation in local carrying capacity. When mate-search costs are low, 50% spatial variation in local carrying capacity is sufficient to stabilize mating domains. Throughout our analyses, we have kept population sizes relatively small, so as to exacerbate the challenge of coexis- tence owing to drift in finite populations. When population sizes are larger, we find that less than 10% variation in local carrying capacity is needed to stabilize mating domains (Fig. C.4D). Because movement and demographic stochastic- ity obscure spatial variation in local carrying capacity across space, levels of variation in this range would be virtually undetectable in nature. The stabilizing effect of spatial variation in local carrying capacity and mate- search costs readily extends to more realistic and natural landscapes (Fig. 4.4) and also to multiple genotypes (Fig. 4.4D–E) . While the spatial distribution of mating types becomes less predictable in complex environments, coexistence is again greatly facilitated. As long as spatial variation in local carrying capac- 48 Chapter 4 ity does not become so insignificant that it hardly affects the landscape, or so asymmetric that a single local population dominates, different mating domains are maintained in mosaic sympatry Mallet (2008); Mallet et al. (2009) (Fig. C.7). Our findings are also robust to changes in female preference strength, mate- search distance, movement distance, and competition distance (Figs. C.4A, C.5), to changes in the relative importance of ecological competition versus sexual selection (Fig. C.4B–C), to changes in the genetic architecture of the display and preference traits (Fig. C.8), and to the inclusion of fitness differences caused by the male display trait (Fig. C.9). Generally, coexistence will be maintained if fe- male preferences are sufficiently strong to prevent extensive interbreeding, and if individuals move and interact on a spatial scale such that they are affected by spatial variation in local carrying capacity. This phenomenon can be inter- preted in terms of a more general mechanism: whenever positive frequency dependence creates multiple possible stable states, global coexistence of these states can become possible in a spatially structured environment if this structure allows the states to become anchored in space. Our results in Fig. 4.4 also ex- tend a previous finding from theoretical work on hybrid zones, predicting that the spatial interface between species moves in space until settling in a region of low population density (Barton and Hewitt, 1985, 1989). Similarly, earlier theo- retical work (Dieckmann, 2004) using reflective boundaries for such anchoring, has shown that ecologically equivalent types can coexist when fecundity drops, or mortality or mobility rise, in the company of heterospecifics. Because both spatial variation in local carrying capacity and costs associated with mate search are ubiquitous in nature, our model may provide an explana- tion for the coexistence of many species whose reproductive barriers are primar- ily sexual. For example, local habitat availability and quality vary around the shoreline of Lake Victoria. It is, therefore, possible that the mechanism reported here could explain how ecologically similar cichlid species can coexist in such vast diversity. That sexual differences have been the primary force maintaining species boundaries in this group is supported by the increasing frequency of hybridization that is occurring as a consequence of high turbidity levels, which reduce a female's ability to discern different male phenotypes (Seehausen et al., 1997). Similar explanations could plausibly be applied to other species pairs that seem to be largely maintained by sexual selection (e.g., species of fruit flies (Hollocher et al., 1997), weakly electric fish (Feulner et al., 2008), frogs (Ryan and Wilczynski, 1988), crickets (Gray and Cade, 2000), and grasshoppers (Tregenza et al., 2000), among others). To test this hypothesis, one could analyse spatial associations between mating domains and local carrying capacity: Fig. 4.4 sug- gests that boundaries of mating domains should often align with troughs of 49 Chapter 4 low local carrying capacity (occasionally mating domains span such troughs, and very rarely they abut where ridges of high local carrying capacity are nar- row). Our work demonstrates that under reasonable conditions, with variation in local carrying capacity over space and costs to females that encounter few preferred mates, sexual selection can maintain species that are not ecologically differentiated. This is in stark contrast to the widespread opinion that sexual selection, on its own, is unable to maintain ecologically equivalent species that overlap in space. Throughout, we have deliberately avoided making any claims about the emergence of diversity and speciation, choosing instead to focus on the coexistence of mating types. Further theoretical work is, therefore, needed to determine which conditions are most conducive to the initial appearance of variation in mating preferences, as well as to the maintenance of these variants. 50 Chapter 4 Figure 4.1: Sexual selection enables long-term coexistence of ecologically equiv- alent species. We consider a population distributed across a continuous habitat in one dimension (columns A, B) or two dimensions (columns C, D) with a lo- cal carrying capacity that is either spatially uniform (A, C: top panels) or that exhibits two peaks (B, D: top panels). Each peak is of Gaussian shape with stan- dard deviation sk. The level v of spatial variation may be altered by changing the height of these peaks relative to the base. A value of v = 0.25, as shown in B and D, means that local carrying capacity is elevated by 25% at the peaks. The four lower rows show model runs through time. Each generation, individuals survive after a round of local competition and reproduce after a round of local mating, followed by offspring movement and the death of all parents. Competi- tion between individuals decreases with their distance according to a Gaussian function with standard deviation ss. Coloured curves in A and B show the effective local density of competitors of each type (weighted by their competi- tive effect, SI, Eq. C.4), across the one-dimensional arena, while dots in C and D show surviving adults. Individuals are coloured according to their geno- type at the display locus (similar patterns are observed at the preference locus; Fig. C.2). Females are a times more likely to mate with a preferred male, when encountered. Males are encountered with a probability that decreases with the distance between male and female according to a Gaussian function with stan- dard deviation sf. Female fecundity declines with the strength of mate-search costs m (SI). Movement distances are drawn from a Gaussian function with standard deviation sm, centered at 0, with wrap-around boundaries. The total carrying capacity is K = 500, supporting the survival of approximately half of the N = 1000 offspring produced each generation; other parameters: sk = 0.1, ss = 0.05, a = 5, sf = 0.05, sm = 0.05, and m/K = 1 (roughly halving fecundity, Fig. C.1). 51 Chapter 4 0 0.5 1 A α = 1 m K = 0 v = 0 0 0.5 1 B α = 5 m K = 0 v = 0 0 0.5 1 C α = 5 m K = 0 v = 0.25 D α = 5 m K = 1 v = 0 0 0.5 1 0 500 1000 E α = 5 m K = 1 v = 0.25 , , Generation, t Fr eq ue nc y o f Q Figure 4.2: Conditions for long-term coexistence. Panels show distributions of allele frequencies at the display locus through time across 1000 model runs in a two-dimensional landscape; coexistence occurs only while these frequencies remain intermediate. Inset panels depict the spatial variation in local carrying capacity as viewed along transects at y = 0.25. A. Homogeneous environment with no sexual selection (a = 1). B. Same as A, except that females are choosy (a = 5). C. Same as B, except with variation in local carrying capacity (v = 0.25). D. Same as B, except with mate-search costs in females (m/K = 1). E. Same as B, except with spatial variation in local carrying capacity (v = 0.25) and mate- search costs in females (m/K = 1); only when both features are combined is long-term coexistence observed. To focus on the maintenance of coexistence, we begin with two equally sized and spatially segregated populations of PQ and pq genotypes (all individuals on the left half of the arena initially have the PQ genotype, while all individuals on the right initially have the pq genotype). This mimics a scenario in which types that previously arose in allopatry come back into contact, revealing the conditions under which they can persist in sympatry. All other parameters are as in Fig. 4.1. 52 Chapter 4 0 1 2 0 0.25 0.5 B C D E > 5000 4000 − 5000 3000 − 4000 2000 − 3000 1000 − 2000 0 − 1000 Strength of mate−search costs in females, m K Sp at ial va ria tio n in ca rry ing ca pa cit y, v Figure 4.3: Conditions for long-term coexistence. Shading indicates the number of generations that polymorphism at the display locus persists when females are choosy (a = 5) in a two-dimensional landscape (darker = longer). Each cell represents the mean time to loss of polymorphism for 10 replicate model runs. Letters indicate parameter combinations used to generate the lower four panels in Fig. 4.2. Side panels illustrate the extent of spatial variation in local carrying capacity for the three parameter values shown along the vertical axis. Model runs are initialized as in Fig. 4.2. All other parameters are as in Fig. 4.1. 53 Chapter 4 Figure 4.4: Mosaic sympatry. Four representative model runs in a patchy two- dimensional landscape with random variation in local carrying capacity. Panel A depicts the underlying spatial variation in local carrying capacity, while pan- els B–E show results from independent model runs after 10, 000 generations overlaid on the local carrying capacity. Panels B and C are initialized with two types, whereas panels D and E are initialized with ten display alleles and ten corresponding preference alleles, all at equal frequencies and distributed randomly across the arena (SI, Section C.2.2). Some of these alleles are then lost during the colonization phase. As in Fig. 4.1, individuals are coloured ac- cording to their genotype at the display locus. The spatial arena is eight times larger than in Fig. 4.1 and the total carrying capacity is K = 4000, supporting the survival of approximately half of the N = 8000 offspring produced each generation. All other parameters are as in Fig. 4.1. 54 Chapter 5 Assortative mating and spatial structure in hybrid zones 5.1 Summary The spatial genetic composition of hybrid zones exhibits a range of possible patterns, with many characterized by patchy distributions. While several hypo- thetical explanations exist for the maintenance of these "mosaic" hybrid zones, they remain virtually unexplored theoretically. Using computer simulations we investigate the roles of dispersal and assortative mating in the formation and persistence of hybrid zone structure. To quantify mosaic structure we de- velop a likelihood method, which we apply to simulation and empirical data. We find that long distance dispersal can lead to a patchy distribution that as- sortative mating can then reinforce, ultimately producing a mosaic capable of persisting over evolutionarily significant periods of time. By reducing the mat- ing success of rare males, assortative mating creates a positive within-patch frequency-dependent selective pressure. Selection against heterozygotes can similarly create a rare-type disadvantage and we show that it can also preserve structure. We find that mosaic structure is maintained across a range of as- sumptions regarding the form and strength of assortative mating. Interestingly, we find that higher levels of mosaic structure are sometimes observed for inter- mediate assortment strengths. The high incidence of assortment documented in hybrid zones suggests that it may play a key role in stabilizing their form and structure. 55 Chapter 5 5.2 Introduction Hybrid zones provide a natural setting in which to study the effects of selection and gene flow on alleles or combinations of alleles (Barton and Hewitt, 1985, 1989). Theoretical studies have provided methods to infer important quantities such as the strength of selection against hybrids and dispersal distance from em- pirical measurements of the spatial distribution of genotypes within a hybrid zone. Consequently, genetic spatial structure can be an informative property of a hybrid zone (Barton, 1979; Barton and Hewitt, 1985, 1989). To date the major- ity of theory has assumed that a monotonic change in genetic composition will be observed along a transect through a hybrid zone; that is, the hybrid zone will be "clinal" (Bazykin, 1969; Barton and Hewitt, 1985, 1989). However, a number of hybrid zones exhibit significant departures from a cline, with patches alter- nating in species composition (e.g., Howard, 1986; Harrison and Rand, 1989; Bierne et al., 2003). These "mosaic" hybrid zones are a spatial patchwork of populations, each fixed (or nearly fixed) for only one of the parental species' types. Of the factors capable of generating mosaic patterns, habitat heterogeneity has received the most attention (Harrison and Rand, 1989; Cain et al., 1999; Bri- dle et al., 2001; Bridle and Butlin, 2002). This hypothesis assumes that the hybrid zone consists of alternating patches of different habitats and that the individuals preferentially occupy their own parent's habitat. The patchy species distribution then reflects the patchy environmental distribution. This hypothesis has been tested in a number of empirical studies (e.g., Howard and Harrison, 1984b,a; Howard et al., 1993; Harrison and Bogdanowicz, 1997; Bridle et al., 2001; Bridle and Butlin, 2002; Vines et al., 2003). However, where ecological patterns reflect- ing the patchiness of the hybrid zone have been found, they often explain only a small proportion of the deviations from a clinal distribution (e.g., Bridle et al., 2001; Bridle and Butlin, 2002; Bierne et al., 2002b). Long distance dispersal during colonization has also been shown to be ca- pable of creating patchy population structures (see Nichols and Hewitt, 1994; Ibrahim et al., 1996; Le Corre et al., 1997; Bialozyt et al., 2006). However, with- out any reinforcing mechanism, the spatial structure created by this process will tend to be transient, reverting to a cline over many generations. The pattern of dispersal on its own is not, therefore, sufficient to explain the persistence of structure observed in many hybrid zones. Several authors (e.g., Cruzan and Arnold, 1993; Jiggins and Mallet, 2000; Bridle and Butlin, 2002; Bailey et al., 2004) have suggested that assortative mating may also contribute to stabilizing the observed mosaic distribution in hybrid zones, but this idea has yet to be 56 Chapter 5 theoretically tested. In this paper we develop a simulation model to investigate the role of as- sortative mating and colonization patterns in the formation and maintenance of mosaic spatial structure within a hybrid zone without habitat heterogene- ity. We also develop a likelihood method to quantify the level of "mosaicity" within a hybrid zone. We find that when dispersal during colonization contains sufficiently many long distance dispersers, assortative mating can reinforce the initial spatial distribution of colonists to enhance the mosaic pattern. This mo- saic structure is then able to persist for prolonged periods of time, even in the face of high levels of dispersal. 5.3 Model description We consider the population dynamics of a one-dimensional hybrid zone us- ing a computer simulation to investigate the role of migration and assortative mating in generating mosaic structure. Our simulation proceeds with discrete generations, tracking genotypes at two loci for diploid individuals along 30 ecologically identical patches, flanked on the left and right by pure patches of two different species. Individuals at the left boundary were assumed fixed for allele a at the first locus and allele b at the second, while those at the right were assumed fixed for alleles A and B. The A-locus is always assumed to af- fect assortative mating, while the role of the B-locus changes depending on our analyses. Recombination occurs between the loci at rate r. Each generation, individuals reproduce sexually within their patch to create the next generation and then die. The offspring then may migrate to a different patch, and new individuals migrate in from the pure edge populations. Hybrids are assumed to be completely viable and fertile for the majority of our analyses, but we also consider viability selection against hybrids through under-dominance at the B- locus. All interior patches were initially empty, and initial colonization occurred via dispersal from the edge populations. 5.3.1 Reproduction and population growth We assume that population growth within each patch is logistic, with growth rate r and carrying capacity K (we used r = 1.25 and K = 1000 for the major- ity of our simulations). Each individual has a Poisson distributed number of offspring with mean 1+ r(K n)/K 57 Chapter 5 where n is the current patch population size. When viability selection acts, indi- viduals heterozygous at the B-locus survive to reproduce with probability 1+ s (s being the strength of selection for or against hybrids). Assortative mating occurs through female preference. Females mate assortatively using the "best of N" scheme presented in Seger (1985). In this model a female surveys N males and then chooses a mate. For brevity we will often use the term "lek" to refer to the group of N males sampled by a female, even though this is not the precise meaning of the term. We chose the best of N mating scheme for its generality. Where N = 1, this model reduces to random mating for all fe- male preference strengths, whereas at the other extreme where N is the whole patch, the model becomes equivalent to the "fixed relative-preference" scheme of Kirkpatrick (1982). We will use N = frp to denote this case. Assortment usu- ally occurs with respect to the A-locus, however, when we consider two-locus assortment (see below), then both loci are assumed to affect female preferences. We let yi,j denote the probability that a female of genotype i mates with a male of genotype j from among the males in her lek, and c denote the strength of a female's preference, ranging from one (random mating) to infinity (mating only with preferred males, if present). yi,j can be interpreted as a combination of all pre-zygotic isolating factors, including both pre-mating (e.g., female choosiness) and post-mating, pre-zygotic factors (e.g., conspecific sperm precedence). We consider the following three models of assortative mating. Linear-preference (k loci) While we investigate assortative mating based on at most two loci, we will describe a general version of this preference model for k loci. In this model preferences of parental species towards hybrids are intermediate between their preferences towards the pure types. If a female of type i and a male of type j share a total of t alleles across the k bi-allelic loci that underlie assortative mating (that is, both have t alleles from the same source population), then the probability that female i mates with male j in her lek is yi,j = ✓ 1+ t 2k (c 1) ◆ f j (5.1) where f j denotes the frequency of genotype jmales in the female's sample of N, and the yi,j's are standardized across all males. For example, in the one-locus case, if c = 2, then, if there were equal genotype frequencies, a homozygous female is twice as likely to mate with a conspecific than a heterospecific and 1.5 times as likely to mate with a heterozygous individual. A heterozygous female 58 Chapter 5 in this case would be 1.5 times as likely to mate with another heterozygote than with either con-specific. Self-preference In this model a female prefers males of her own genotype by an equal factor (c) over all other males in her lek, discriminating equally against any other genotype. The probability that a female of genotype i mates with a male of genotype j is then given by yi,j = cdi,j f j (5.2) where fj again denotes the frequency of genotype jmales in the female's sample of N, and di,j is Kronecker's Delta, which is equal to one when i = j and zero otherwise. Again the yi,j's are standardized across all males. Dominant-preference In this model we assume that heterozygous males and females are indistinguish- able from homozygous individuals characteristic of the source populations on the right side of the hybrid zone. Equation (5.2) then applies where di,j is one if the phenotypes of i and j are the same. 5.3.2 Dispersal There are two types of dispersal in our model: internal dispersal among the 30 patches and dispersal of new immigrants from the pure edge patches. In the internal patches, after reproduction each individual disperses with probability m. Individual dispersal distances were drawn randomly from a mixture of two exponential distributions, generating a leptokurtic (fat-tailed) dispersal kernel (Clark, 1998): pE(µs) + (1 p)E(µl) (5.3) where E(µ) is an exponential distribution with mean µ. We assumed that µs < µl , so that p represents the proportion of short-distance dispersers. Dispersal distance was measured as the number of patches away from the focal patch, with an equal probability of being in either direction. Individuals that dispersed beyond the edge patches and into the pure source populations were assumed to have a negligible impact, and thus were removed from the system. In addition to local migration between patches, a fixed number I of individuals arrived each generation from each of the pure patches. These migrated according to the same dispersal kernel as above, relative to the edge of the pure patches. 59 Chapter 5 5.3.3 Analyzing simulation results To analyze the results of our simulation we developed a likelihood method to quantify the degree of "mosaicity" of a hybrid zone (i.e., deviation from a monotonic cline). We fit a series of horizontal steps to allele frequency data along a transect through the hybrid zone (figure 5.1). In order to compute the allele frequency within a patch we sampled every individual. A model with k steps is defined by k step locations s = {s1, s2, · · · , sk} and k 1 step heights h = {h1, h2, · · · , hk 1}, where hi is the height between steps si and si+1. Because we assumed that the edge patches were fixed for each parental type, we added an initial step of height 0 and a final step of height 1 (figure 5.1). In appendix D.1 we show how the likelihood of the observed allele frequency data given this model can be computed. For a set of data and given number of steps, we found the set of step locations and heights that maximize this likelihood and then used likelihood ratio tests to find the number of statistically significant steps required to best explain the data (see appendix D.1 for details). A similar approach to ours was used by Macholán et al. (2008) in order to fit step models through data collected from theMus musculus musculus/M. m. domesticus hybrid zone. However, Macholán et al. constrained their step heights to change monotonically through the hybrid zone, and thus their approach cannot be used to make inferences as to the level of mosaic structure a hybrid zone displays. Once the step-wise model is fit to the allele frequencies, any measure of mosaicity may be investigated. Here, we focus on the sum of the magnitudes of the downward step sizes as a measure of the "mosaicity" of a data set, given by: M = k 2  i=1 max(0, hi hi+1) (5.4) (figure 5.1). This quantity has a minimum of zero for a clinal model (monotoni- cally increasing steps), regardless of the steepness of the cline, and grows as the number and size of reversals in step height increases, attaining its maximum possible value when the hybrid zone consists of patches alternating between fixation on one or the other allele. For such cases (where each patch is fixed for either the a or the A allele), M is equal to the number of times the frequency of the A allele changes in frequency from 1 to 0 and back to 1 again. Because wider hybrid zones (e.g., those with more patches) have a greater number of possible locations where reversals in allele frequency may occur, M may be inflated in hybrid zones with finer levels of sampling. Accounting for the number of samples when computing M, however, would lead to the same hybrid zone having different mosaicity scores when sampled at different scales 60 Chapter 5 (an undesireable outcome). Thus we have assumed that the scale of sampling has been chosen so that it appropriately captures the spatial distribution of alle- les across the hybrid zone, that is, no patches are missed during sampling. We have not, therefore, included the number of patches in our measure of mosaic- ity. Four sample best fit models with associated M values are shown in figure 5.2. Other measures of structure Two additional measures that are commonly applied to quantify other aspects of genetic structure within a hybrid zone are linkage disequilibrium and bi- modality. Linkage disequilibrium measures the association between alleles at different loci, and bimodality measures the lack of heterozygotes at a single lo- cus, compared with the expectation under Hardy-Weinberg equilibrium (Jiggins and Mallet, 2000). We also assess these measures in our simulations. To compute linkage disequilibrium we use Lewontin's D0, which is link- age disequilibrium standardized by its maximum possible value given the ob- served allele frequencies (see Lewontin, 1988, for details). This statistic varies from -1 to 1, with -1 representing a population composed entirely of A/b and a/B genotypes, 1 representing a population composed entirely of A/B and a/b genotypes, and 0 representing a population with no associations between loci. To measure bimodality we use FIS, as presented in Jiggins and Mallet (2000), which is computed as FIS = 1 pAa2pApa (5.5) where pAa is the frequency of heterozygotes and pA and pa denote the frequen- cies of the two alleles. This statistic lies between -1 and 1, with 0 implying that a population is at Hardy-Weinberg equilibrium, and 1 (-1) representing a complete lack of heterozygotes (homozygotes). To compute both FIS and D0 for a given population we pooled individuals from all patches. Alternatively, we could have computed a patch average of each statistic. Because many patches were often fixed or nearly fixed for a single genotype, however, we found that such an approach did not provide meaningful results. In particular a mosaic hybrid zone consisting of patches alternating in state, from fixation on one allele to fixation on the other, would have a maximum mosaicity score, but in each patch both FIS and D0 would be undefined. At the other extreme, if hybrids never survive, FIS and D0 would be maximal within each patch, regardless of whether the mosaicity score were low or high. Thus, we see that the mosaicity score provides a distinct measure of the "patchiness" of the hybrid zone. 61 Chapter 5 5.4 Results Because our findings are largely insensitive to which model of assortative mat- ing was used, we primarily present results from the one-locus linear-preference model, and will only explicitly state when this is not the case. We compute mo- saicity for the one-locus linear model at the A-locus, which governs assortative mating (figure 5.2). We found that a combination of long distance dispersal and assortative mat- ing allow mosaic population structure to form (figures 5.2, 5.3) and persist for long periods of time (figure 5.4A). The initial mosaic structure is generated by individuals leap-frogging over heterospecific populations during colonization to found new populations. These processes have been discussed in detail else- where (e.g., Nichols and Hewitt, 1994; Ibrahim et al., 1996). In contrast with these previous studies, which assumed random mating, we found that once a mosaic population structure establishes, female preference acts to preserve the mosaic patterns (figure 5.4A). This is a consequence of a within-patch mating advantage to males of the more abundant female type, leading to within-patch fixation (or near fixation) on whichever mating type attains a higher frequency during the early stages of colonization. Once patches have reached carrying capacity, invasion by the other mating type becomes unlikely, and thus the fi- nal mosaic pattern persists, even with relatively high levels of between-patch dispersal (m = 0.01,Nm ⇡ 10). The strength of assortative mating strongly affected the final level of mosaic structure observed in any particular hybrid zone (figure 5.2). Where mating was random, the constant arrival of foreign migrants biases patches on the left toward one allele and patches on the right toward the other, creating a grad- ual, roughly monotonic cline (figure 5.2A). With weak female preference mo- saic structure was not stable over time. While sexual selection was a sufficiently strong force to reduce the extent of within patch co-existence of the two mat- ing types, it was not strong enough to prevent the occasional invasion by the rare mating type. The constant arrival of different migrant types from each end of the zone, combined with migration among patches eventually removes any traces of the founding population structure, producing a cline that was both steeper and more monotonic than when mating was random (figure 5.2B). Weak female preference thus resulted in slightly lower final mosaicity values than those obtained in the neutral case (figure 5.4A). For stronger assortative mating the hybrid zone became a mosaic of alternating patches of each geno- type (figure 5.2C, 5.2D). Due to the near fixation of each patch on one allele type or the other, the mosaicity scores in these cases can be interpreted as approxi- 62 Chapter 5 mately counting the number of complete reversals in patch genotype frequency (e.g., figure 5.2D has three nearly complete down-steps and a mosaicity score of 2.91). The effects of assortative mating extended to other loci in the genome. Link- age disequilibrium between the assortative mating locus and an unlinked neu- tral locus persisted at higher levels with stronger female preference (figure 5.4B). Tighter linkage further increases the level of disequilibrium (figure 5.4B). Bi- modality at both loci was also higher with stronger female preferences (figure 5.4C), achieving nearly the maximum possible value at the assortative mating locus when c = 10. When migration distances were short, dispersal did not create a sufficient level of genetic patchiness for assortative mating of any level to reinforce, and thus no mosaics were observed. Increasing the mean of the long distance dis- persal component (µl) led to a roughly linear increase in the observed level of mosaicity, when assortative mating was sufficiently strong enough to maintain a mosaic pattern (figure 5.3). To account for changes in mean dispersal dis- tance that occur when increasing µl we ran a separate set of simulations where we decreased µs whenever we increased µl , in order to maintain a constant mean dispersal distance. Results remained qualitatively identical and thus have not been included. Similarly, increasing the proportion of individuals within a patch that disperse each generation (m), or the number of foreign migrants arriving each generation (I), sped up the rate at which colonization occurred, but did not qualitatively change results. Varying either population growth rate (r) or carrying capacity (K) also had little effect, and thus further exploration of these parameters has not been included. The level of mosaic structure not only depended on the strength of female preference, but also on her lek size (figure 5.5A). Whenever mating was non- random (e.g., for lek sizes greater than 1), the level of mosaicity increased with female mate preferences from c ⇡ 1.5. For small leks mosaicity scores remained high throughout the entire range of strong female preferences (e.g., N = 2, 5 in figure 5.5A). Surprisingly, however, large leks displayed a decline in mosaicity with very strong female preferences (e.g., N = 20, frp in figure 5.5A). This implies that when a female samples only a small number of individuals from the population before mating, mosaic structure is more stable than when she has access to many males. This occurs because, when females only sample a few males, there is a high chance that the lek of a rare female will not contain any males of her preferred type and she will, therefore, mate with a more common male, which further reduces the expected mating success of the rare males. Conversely, when females sample many males, the lek of a rare female will often 63 Chapter 5 contain at least one male of her preferred type. A few males and females may, in this case, be sufficient to colonize a patch occupied by heterospecifics. At very strong mate preferences and large lek sizes the population structure reduces to a noisy cline, with different alleles often co-existing in the same patch. However, the mosaicity level is approximately 2-3 times that for random mating (figure 5.5A). The higher mosaicity score is due to a complete lack of heterozygotes, whose presence in the random mating case acts to reduce differences in allele frequencies between adjacent patches and thus the average step size. A drop in mosaicity values below those of the random mating case can be observed for intermediate mating preferences in the limiting case where a female samples the entire patch of males before mating (see figure 5.5 and the c = 1000 curve in figure 5.3). With intermediate c values, rare types are selected against, but not so strongly that they do not occasionally fix. The constant arrival of migrants from the boundaries thus eventually overwhelms any traces of the founding population structure, eventually creating a steep monotonic cline. Changing the model of female preference did not qualitatively affect re- sults (see figure 5.5B). However, with both the dominant and self-preference models the strength of female preferences required to overwhelm the mosaic structure was weaker than in the one-locus linear model (compare the c ⇡ 102 region in the dominant and self-preference models to the c ⇡ 103 region in the one-locus linear model). This occurs because the weak discrimination against hybrids combined with the reduced preference strength of hybrids in the linear model effectively reduces the efficacy of assortative mating for a given prefer- ence strength, compared to the other two models. Mean mosaicity scores were lower with the dominant-preference model, due to there being a slight bias to- ward patches fixing the dominant allele. In the two-locus linear model results were qualitatively unchanged when female preferences were weak. When pref- erences were stronger, however, curves appeared qualitatively more similar to those from the one-locus model with a smaller lek size than used in the sim- ulations. This shift is likely a consequence of recombination breaking down genotypes of rare individuals, making co-existence and thus invasion more dif- ficult for that species. So far we have treated hybrids as having equal fitness as parentals. In a similar manner to assortative mating, an intrinsic reduction in hybrid viability or fertility can create a frequency-dependent selective pressure against the rarer of the two parental species within a patch (figure 5.6). This is because each hy- bridization event represents a larger fraction of the rare species' matings, and thus each hybrid death reduces the fitness of the rare species by a proportionally larger amount. Barton and Whitlock (1997) showed that stabilizing selection on 64 Chapter 5 polygenic traits can similarly induce frequency-dependent selection against rare genotypes and thus, through an analagous process to that described above, can also lead to the maintenance of different allelic combinations between popula- tions, provided migration rates are low enough. As expected, we found that as the strength of selection against heterozygotes increased a tighter correla- tion emerged between mosaicity scores at the assortative mating locus and the viability locus. When hybrids were lethal, recombination never occured and the two loci were in essence completely linked. With strong assortative mating there was higher concordance between mosaicity scores at both loci, even when viability selection was relatively weak (figure 5.6C). The marginal decline in mo- saicity observed for strong viability selection arises due to a slight reduction in the initial establishment of highly mosaic populations and does not represent a reduction in the stability of mosaics, once established. In the simulation results presented, we assumed that allele frequencies were determined from the entire patch (K = 1000 individuals). The majority of em- pirical hybrid zone data most certainly contain many fewer samples. In order to investigate this sampling effect we sampled 100 individuals from each patch (without replacement) and then fit the best model to the sampled data. Results remained largely unchanged, with sampled data usually producing a model missing a few small steps but having a nearly identical mosaicity score. Although we have only presented models consisting of 30 patches here, it is worth mentioning that the width of the cline or mosaic for any particular run often fluctuated well within the limits of these 30 patches, typically settling down to a much narrower final width than is possible with 30 patches. The resultant clines or mosaics were thus flanked on either side by large regions consisting of pure parental genotypes (e.g. figure 5.2). To assess whether the width of the hybrid zones would differ noticeably with more patches, we ran simulations with twice the number of patches (60). Indeed the hybrid zones in these cases settled into a comparable width as with 30 patches, but with larger flanking regions of pure parental species on either side. Thus our results are not likely to be affected by the number of patches in the hybrid zone. 5.4.1 Application to Mytilus edulis and M. galloprovincialis hybrid zone data The smooth-shelled mussels Mytilus edulis and M. galloprovincialis form a mo- saic hybrid zone that stretches around the coast of western Europe. While dif- ferences in temperature, salinity, and wave exposure affect species composition within the hybrid zone (Gardner, 1994; Bierne et al., 2002b), Bierne et al. (2002b) argued that local adaptation was not sufficient on its own to explain the ob- 65 Chapter 5 served pattern in this hybrid zone. In a separate study Bierne et al. (2002a) suggested that the presence of assortative fertilization has likely contributed in maintaining the current population structure. To test our methods we fit models to three separate loci from this hybrid zone, using data presented in Bierne et al. (2003) (figure 5.7). Mosaicity scores and 95% bootstrap confidence intervals were 3.66 (3.08, 3.85) for Glu-5', 2.43 (1.88, 2.70) for mac-1', and 2.72 (2.14, 2.97) for Efbis. The mosaicity score for Glu-5' is significantly higher than that for both mac-1' and Efbis, despite similar sample sizes, suggesting that Glu- 5' is possibly more closely linked to a locus influencing assortative mating or experiencing under-dominant selection. Furthermore, this difference between mosaicity scores demonstrates our method's ability to detect differences among loci from the same species with realistic levels of empirical sampling. 5.5 Discussion Our results show that assortative mating, when coupled with long-distance dis- persal during colonization, can lead to the stable persistence of mosaic patterns, even in the absence of ecological differences between incipient species. Fur- thermore, these results are robust to changes in a variety of assumptions about migration, female mating behaviour, and the fitness of hybrids. It has been previously demonstrated that founder effects caused by long distance dispersal into vacant habitats during colonization can create spatially mosaic populations (Nichols and Hewitt, 1994; Ibrahim et al., 1996), and it is this process that drives the initial mosaic patterns observed in our simulations. However, this process alone can not explain the long-term persistence of spatial structure in these populations (figure 5.4). Rather, assortative mating and/or hybrid inviability are essential to stabilize the mosaic structure in the absence of environmental heterogeneity. We have further shown that the effects of assortative mating will likely carry over to other regions of the genome, even with high levels of recombination. The combination of continued immigration of pure AB and ab genotypes into the hybrid zone, and the slow decay of linkage disequilibrium which has been shown to characterize stepping-stone models (De and Durret, 2007), led to a non-zero final value of D0, even when mating was random (figure 5.4). These equilibrium values of linkage disequilibrium were noticeably higher when fe- male preference was stronger and/or when linkage to the assortative mating locus was tighter (figure 5.4B). Furthermore, higher levels of bimodality can be maintained at neutral loci with stronger female preferences (figure 5.4C). Sam- pling at a neutral locus will, therefore, not tend to reveal just how mosaic the 66 Chapter 5 hybrid zone may be at loci directly involved in assortative mating. However, such data can still potentially provide insight into whether a hybrid zone does have an underlying mosaic structure. While our results indicate that assortative mating can stabilize mosaic hy- brid zones for extensive periods of time, it is likely that other factors work in concert with this process. For example, habitat heterogeneity could strengthen the mosaic effect observed here. The combination of ecological differences be- tween incipient species and a patchy environment could create a small degree of initial spatial segregation within a population. Assortative mating could then help push sub-populations towards fixation on one or the other type, whichever is locally more abundant. Similarly, in some cases the combined effects of via- bility selection against hybrids (which we have shown can also preserve mosaic structure) and assortative mating may allow for the preservation of a highly mosaic structure, where each force in isolation would not be sufficient to do so. We have also presented a method that can be used to fit step-wise models through one-dimensional empirical hybrid zone data in order to objectively es- timate their level of mosaicity. Our mosaicity statistic measures the number of reversals in allele frequency. It is not proposed as an alternative to the types of cline fitting that have been traditionally used in mosaic hybrid zones (e.g., see Bridle et al., 2001), but instead provides a complementary measure. To test our method we applied it to data from the Mytilus edulis/M. galloprovincialis hybrid zone (Bierne et al., 2003). Our best fit models exhibited high mosaicity scores at all three loci, demonstrating that it is informative when applied to empirical data sets. Interestingly, we found a statistically higher mosaicity score at Glu- 5', suggesting linkage to a locus involved in assortment, hybrid inviability, or ecological adaptation. How common is assortment likely to be in hybrid zones? Assortative mat- ing is thought to evolve as a consequence of divergence, with recently diverged species being more likely to hybridize (Felsenstein, 1981; Coyne and Orr, 1997). Bailey et al. (2004) documented strong assortative mating between the field grasshoppers Chorthippus brunneus and C. jacobsi, which form a mosaic hybrid zone in northern Spain. The reported "isolation index" of I = 0.59 between these species corresponds to a c of 3.9 in our model. Bridle et al. (2006) have also recently argued that this assortative mating, through preferences for male song, likely plays an important role in maintaining the observed structure in this hybrid zone. Howard and Gregory (1993) conducted sperm competition experiments between the ground crickets Allonemobius fasciatus and A. socius, which form a mosaic hybrid zone in northeastern United States (Britch et al., 2001). When mated to both types of males A. socius females exhibited a con- 67 Chapter 5 specific sperm precedence of at least 95%, and A. fasciatus females of at least 98%. These translate into c values of approximately 42.5 and 49 respectively. The strength of assortative mating documented in both of these hybrid zones is sufficiently strong to lead to a high level of mosaicity in our model, provided that patchiness was initially present. In some systems, however, assortative mating will likely play only a small role, if one at all. In the Bombina bombina – B. variegata (fire-bellied toad) hybrid zone, MacCallum et al. (1998) found that B. bombina-like hybrids were most often associated with pond habitats, whereas B. variegata-like hybrids were most often found in puddles. This strong habitat specialization explains most of the observed spatial variation, leaving little need for additional processes. Mosaic hybrid zones may be more common than has been reported. The pattern observed in any particular hybrid zone may reflect the scale at which individuals are sampled (Schilthuizen, 2000; Ross and Harrison, 2002). When sampling is too coarse a mosaic pattern can appear clinal. Harrison and Bog- danowicz (1997) provided a simple characterization of hybrid zones, based on the shape of the genotypic distribution at the cline center. Mosaic hybrid zones tend to exhibit an overabundance of parental types, relative to hybrids and thus have a bimodal genotypic distribution. Resampling our mosaic hybrid zones on a coarser spatial scale would, on average, lead to a clinal pattern with a bimodal distribution of genotypes. This suggests that many bimodal hybrid zones may be mosaic at a finer spatial scale than measured. Jiggins and Mal- let (2000) found a strong positive correlation between the measured strength of assortative mating and the level of bimodality in a survey of several empir- ical hybrid zones. While bimodality does not necessarily imply that a hybrid zones is mosaic (Cruzan and Arnold, 1993; Emms and Arnold, 1997), the study by Jiggins and Mallet (2000) suggests that an underlying correlation between mosaicity and assortativity may exist. Understanding the mechanisms by which assortative mating occurs in hy- brid zones would provide insight into the applicability of our model, as well as possibly reveal interesting theoretical extensions to other models of sexual selection. Our findings depend on sexual selection inducing positive frequency- dependent selection; rare males are always at a disadvantage in the best of N mating scheme. This is not necessarily the case with other models of assorta- tive mating. In a "grouping-based" model, females mate within some group of individuals with a fixed probability and otherwise mate with a male drawn ran- domly from the population (Felsenstein, 1981; Otto et al., 2008). Because group membership may be frequency independent (e.g., groups could be chosen based on spatial or temporal proximity), this model does not necessarily induce a rare- 68 Chapter 5 type disadvantage, and thus we would expect it to yield a qualitatively different outcome. Our model assumes an underlying demic system. Consideration of a sim- ilar model in continuous space may provide an interesting avenue for future research. In a homogeneous environment, where there are no patch boundaries stabilizing the sizes of pure (or nearly pure) populations, we would expect fluc- tuations in population size to lead to the eventual loss of mosaic structure. This may not be the case, however, in heterogenous environments. Ecologically het- erogeneous environments that favour some degree of local adaptation, or envi- ronments where the carrying capacity and thus the density of individuals varies in space, may stabilize population sizes of the different mating types, and thus allow for the preservation of mosaic structure. Temporal fluctuations in popula- tions size, or regularly occuring local extinctions followed by re-colonization via long distance dispersal may also allow for the long term persistence of mosaic structure, although, in this case it would vary spatially in time. Despite having received significant empirical attention, mosaic hybrid zones have remained largely unexplored in the theoretical literature, and the majority of models so far have assumed an underlying clinal structure. Both theoreti- cal and empirical work have demonstrated, however, that clinal models do not always make accurate predictions about evolutionary processes occurring in a mosaic hybrid zone. For example, Cain et al. (1999) found that reinforcement evolved under a much wider set of circumstances for a mosaic hybrid zone than for a clinal one, and the empirical hybrid zones described in Bridle and Butlin (2002) and Cruzan and Arnold (1993) did not conform well to the clinal model expectations. Given that our current estimate of the prevalence of mosaic hybrid zones in nature is probably underestimated (because their detection is sensitive to the scale of sampling) and that predictions based on clinal models may not apply to mosaic hybrid zones, it is important that theory be developed to help understand the forces creating and maintaining mosaic hybrid zones and the effects that these mosaic hybrid zones have on evolutionary processes. 69 Chapter 5 ● ● ● ● ● ● ● ● ● ● ● ● 0.0 0.5 1.0 Al lel e fre qu en cy s1 s2 s3 s4 s5 Spatial position h1 h2 h3 h4M = h2 − h3 Figure 5.1: An example model fitted to hypothetical allele frequency data along a transect through a hybrid zone. This model has five fitted steps. The height between step si and step si+1 is hi. Because there is only a single down-step the mosaicity score, M, is simply the magnitude of this step (i.e., h2 h3, see equation 5.4). ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● x p 0.0 0.5 1.0 c =1 M =0.66 A ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● x p c =1.2 M =0.13 B ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● x p 0 10 20 30 0.0 0.5 1.0 c =2 M =0.86 C ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● x p 0 10 20 30 c =10 M =2.91 D Spatial position Al lel e fre qu en cy Figure 5.2: Example simulated hybrid zones with varying strengths of assorta- tive mating (increasing across panels A to D), after 1000 generations. The most likely stepwise model is indicated by the line through the data. Assortative mating was based on the one-locus linear-preference model. Parameters were µs = 1, µl = 15, m = 0.01, p = 0.75, r = 1.25, K = 1000, I = 50, and N = 5. 70 Chapter 5 5 10 15 0 1 2 3 Mean of long distance dispersal ( µl ) M os aic ity c = 1 c = 2 c = 10 c = 1000 Figure 5.3: Mosaicity as a function of long distance dispersal distance (µl) in the one-locus linear-preference model. Each point is the mean of 100 replicate sim- ulations. Different curves correspond to different female preference strengths (indicated by labels). Error bars denote ± one standard error. Other parameters were as in figure 5.2 with N = frp. Curves were plotted after 1000 generations. 71 Chapter 5 0 1 2 M os aic ity ( M ) A 0.0 0.2 0.4 0.6 0.8 1.0 Lin ka ge D ise qu ilib riu m ( D' ) B ρ =0.1 ρ =0.5 0.0 0.2 0.4 0.6 0.8 1.0 Bi m od ali ty ( F IS ) Time (generations) 50 250 2500 25000 C c 1 1.2 1.5 2 10 Assortative mating locus Neutral locus Figure 5.4: Mosaicity at the assortative mating locus (panel A), linkage disequi- librium between the assortative mating locus and a freely recombining neutral second locus (panel B) and bimodality at each locus (panel C) plotted over time. Each curve is the mean of 100 replicate simulations. Different curves correspond to different female preference strengths in the one-locus linear-preference model (the key in C applies to all panels). Recombination rates were set to r = 0.5 (free recombination) for all curves except the topmost curve in panel B, which corre- sponds to r = 0.1. In panel C, under random mating (c = 1) both loci display nearly identical trajectories, and thus the A-locus has been omitted for clarity. Error bars denote ± one standard error and are only shown at generation 25000 for clarity. Other parameters were as in figure 5.2 with N = frp. 72 Chapter 5 0 1 2 M os aic ity 1 101 102 103 104 105 frp N = 1 N = 2 N = 5 N = 20 A 1 101 102 103 104 105 Linear−preference (1 locus) Linear−preference (2 locus) Self−preference Dominant−preference B Conspecific mate preference ( c ) Figure 5.5: Mosaicity as a function of mate preference for different female lek sizes in the linear-preference model (A) and under different preference models with N = frp (B). Each point is the mean of 100 replicate simulations. Error bars denote the standard error and are only plotted for the rightmost points. Other parameters were as in figure 5.2 with N = frp and r = 0.5. Simulations were run for 1000 generations. 73 Chapter 5 0 1 2 A c = 0 0 1 2 B c = 1.5 0 1 2 C c = 10 Assortative mating locusHybrid viability locus 0 −0.2 −0.4 −0.6 −0.8 −1 Strength of selection against hybrids ( s ) M os aic ity Figure 5.6: Mosaicity at each locus as a function of the strength of selection against hybrids in the one-locus preference model. Each curve is the mean of 250 replicate simulations. Different panels correspond to different female preference strengths (indicated by labels). Error bars denote the standard error. Other parameters were as in figure 5.2 with N = frp and free recombination (r = 0.5). 74 Chapter 5 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●0.0 0.2 0.4 0.6 0.8 1.0 Al lel e fre qu en cy Glu−5' M =3.66 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●0.0 0.2 0.4 0.6 0.8 1.0 Al lel e fre qu en cy mac−1' M =2.43 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●0.0 0.2 0.4 0.6 0.8 1.0 Al lel e fre qu en cy 1 6 11 16 21 26 31 Efbis M =2.72 Location Figure 5.7: Best fit models for three diagnostic loci in the Mytilus edulis/M. galloprovincialis hybrid zone (data from Bierne et al., 2003). Mosaicity scores and 95% bootstrap confidence intervals were 3.66 (3.08, 3.85) for Glu-5', 2.43 (1.88, 2.70) for mac-1', and 2.72 (2.14, 2.97) for Efbis. 75 Chapter 6 Conclusions Overall, the above models illustrate how mathematics can be a useful tool in guiding our understanding when interactions between populations or species lead to complex dynamics. In chapter 2 I used modifier theory to investigate the evolution of mutation rate at a locus regulating host-parasite interactions. I found that lower mu- tation rates evolved when recombination occurred between those loci and the loci regulating mutation rate. This finding can potentially help to explain the high rates of antigenic switching that have evolved in many asexual taxa, where linkage is complete. In this model, higher mutation rates tended to dampen the cycles in allele frequency characteristic of host-parasite interactions. When mutation rates evolved to be sufficiently high, cycles disappeared altogether, effectively eliminating selection at the modifier locus. In small populations, however, stochastic fluctuations in allele frequency still occurred and thus led to higher mutation rates than expected from the deterministic theory. In chapter 3 I developed a model to investigate how ploidy and the archi- tecture of the genetic interactions between hosts and their parasites might affect the evolution of parasitism itself. Because parasites who possess only a single allele are often able to better evade detection by hosts (i.e., homozygotes or hap- loids often outperform heterozygotes), I found that the transition to parasitism occurred over a broader range of parameters when the parasite was haploid. The role of host ploidy was more complicated and depended on the model gov- erning host-parasite interactions. These results provide a first characterization of how genetic architecture affects selection on life-history strategies in antago- nistic species interactions. In chapter 4 I investigated the implications of mating interactions for long- term co-existence of biological diversity. I found that, when two criteria were met (spatial variation in the carrying capacity and search costs associated with 76 Chapter 6 rarity in females), sexual selection dramatically prolonged co-existence times. Essentially, these criteria together ensure that rare males and rare females arriv- ing into a high density area of another type have low fitness and that such areas are relatively fixed in space and time. This is the first study to demonstrate the existence of conditions under which sexual selection alone can promote the long-term co-existence of ecologically equivalent populations with overlapping ranges even in the face of drift. This work thus marks an important contribution to our understanding of the role played by sexual selection in the maintenance of biodiversity. In chapter 5 I considered the effects of mating preferences on the spatial distribution of genotypes observed in hybrid zones. I found that the distri- bution of genotypes observed in many so-called "mosaic" hybrid zones might be better explained by species-specific differences in mating preferences rather than by differences in ecology, the more commonly invoked cause. In analyzing the model I also developed a statistic to quantify the level of "mosaicity" of a particular hybrid zone. I tested this statistic on empirical data from the coastal hybrid zone between the mussel species Mytilus edulis and M. galloprovincialis (Bierne et al., 2003). I found that "mosaicity" was significantly higher for one of the three loci analyzed. This suggests that that locus or a linked locus, at least partially, underlies assortative mating and/or local adaptation. Interest- ingly, previous observations indicate that assortative mating does occur in this hybrid zone, consistent with the potential role of mating preferences in shaping the mosaic nature of the zone (Bierne et al., 2002a). 6.1 Future directions Because the natural world is vastly complex, models are necessarily simplifi- cations of the actual processes occurring in real ecosystems. For this reason, mismatches between models and the systems they are meant to represent are to be expected. Mismatches should not, however, be viewed as "failures" of the theoretical models, but instead, as opportunities for future insights (e.g., under- standing why there was a discrepancy). Only through direct empirical tests can we identify when and where such mismatches occur. Such empirical tests are, therefore, indispensable as we seek to develop models that better map onto the world around us. Many theoretical models address questions that are currently beyond the realm of empirically testable hypotheses. For example, there is little, if any, em- pirical support for or against the ubiquitous theoretical assumption that allele frequencies at loci regulating host-parasite interactions display cyclical dynam- 77 Chapter 6 ics. Thus, there is still a need for experimental tests of many of the most basic assumptions and/or findings of some of the earliest and simplest theoretical models. With this caveat, however, I will discuss some of the most promising empirically testable hypotheses that emerge from the models discussed above. The primary finding in chapter 2 that lower mutation rates evolve with a higher recombination rate between the loci regulating host-parasite interactions and the loci regulating mutation rate could possibly be tested experimentally. For example, comparing mutation assays for antibiotic resistance in closely re- lated species of microbes that differ in their level of sexuality would allow one to evaluate whether a correlation between sexuality and mutation rate exists in nature. The secondary finding that higher mutation rates evolved in small populations could also be experimentally tested by, for example, controlling population size in a microbial host-parasite system, and then measuring the rate of accumulation of mutations allowing hosts or parasites to invade and/or evade one another. In chapter 4, the primary finding that sexually divergent, but ecologically equivalent, species can co-exist depended critically on two assumptions. The first of these (that there exists spatial variation in the carrying capacity) is surely met in most, if not all habitats. However, the same is not necessarily true for the second assumption, namely that rare mating types suffer a fitness cost. It is well known that female preferences can create frequency-dependent selection against males possessing non-preferred traits. While an analogous cost to rarity in females is often assumed to exist, to my knowledge, there is little empir- ical support for or against such an assumption. Such a cost seems plausible (e.g., females who rarely encounter their preferred male phenotypes may have to spend longer searching for a mate), however, an empirical justification for such an assumption would be valuable. Once the more basic assumptions re- garding costs to rarity have been experimentally verified, it seems realistic that an experiment could directly test the main result presented in this paper (e.g., using a spatially structured grid of connected vials with varying resources and Drosophila strains with divergent preferences). In addition to mismatches with empirical work, the simplicity of theoreti- cal models should not be viewed as necessarily limiting their potential power. Often, it is this simplicity that can lead to general wide-ranging hypotheses and findings. For example, because the models I have developed above are not based on any particular species (or pairs of species in the case of hosts and par- asites), the main findings should apply to any species (or pairs of species) that satisfy the central assumptions. This creates opportunities for potentially inter- esting comparative analyses testing whether species that fit the assumptions of 78 Chapter 6 a model confirm its predictions better than less relevant species. Furthermore, those cases where the models do not seem to apply provide avenues for future theoretical exploration. There are a number of comparative analyses and the- oretical studies that would provide interesting avenues for future research as follow-ups to the models presented here. In chapter 3 I investigated how the genetic architecture regulating host- parasite interactions might select for increases or decreases in the level of para- sitism of the parasitic species. I found that, in general, parasitism should more easily evolve in haploids than in diploids. This is similar to a finding by Nuis- mer and Otto (2004), who used a theoretical model to show that, in order to better evade detection by hosts, parasites should more often evolve toward hap- loid genomes than diploid ones. A survey of protists confirmed their prediction (i.e., there was a correlation between parasitism and haploidy). Such a finding is, however, also consistent with my model (i.e., haploids more readily evolve parasitism). Comparative analyses, such as those presented in Pagel (1994), could be used to evaluate whether parasites more often evolve haploidy or hap- loids more often evolve parasitism. In addition to the above comparative study, there is a challenging, but po- tentially interesting theoretical follow-up. In the model, as presented, evolution of parasitism from a non-parasitic life-style can occur in a single evolutionary step. This is because we assumed from the outset that all of the machinery was in place for parasitism to evolve (e.g., only a mutation at a modifier locus was necessary in order for parasitism to initially appear; all the alleles necessary for successful invasion and resistance were already present in both species). As a result, this model does not fully address the initial stages in the evolution of parasitism. Such a model is needed in order for us to appreciate why some species might be more likely to undergo large shifts in life-history strategies. It is possible that parasitism could evolve through a number of complex evolu- tionary pathways (e.g., degradation of a mutualistic relationships or a gradual exploitation of a newly encountered species are two such possibilities). In addi- tion to helping understand the origins of parasitism and life-history shifts, such a model may also provide insight into how adaptation proceeds when many complex genetic changes are required. In chapter 4 I investigated how sexual selection might enable long-term co-existence of ecologically equivalent species. In doing so, however, I as- sumed that genetic variation in mating preferences and trait displays was al- ready present in the population. How the most basic assumptions of my model inhibit or facilitate the process of speciation and the generation of this variation is, therefore, necessary if we hope to gain a full appreciation of the model's 79 Chapter 6 generality. If it turns out that costs to rarity in females and variation in the local carrying capacity prevent sympatric speciation by sexual selection, then my model is most applicable in cases of secondary contact between sexually divergent mating types, but it cannot help explain origins of diversity. If, on the other hand, the parameter space in which long-term co-existence occurs in my model overlaps with that where speciation can occur by sexual selection, then my model may provide a general mechanism whereby sympatric specia- tion can produce species capable of stable co-existence, without any disruptive ecological selection or allopatry. Additionally, generalizing the model (e.g., to quantitative or multi-dimensional preference and display traits) would provide insight into the broader applicability and robustness of these results. In chapter 5 I suggested that variation in mating preferences between species might provide a better explanation for the spatial structure observed in some mosaic hybrid zones. In order to assess whether such mating preferences really do provide a better explanation it would be worthwhile to conduct a meta- analysis of hybrid zone data, quantifying the level of "mosaicity" within each hybrid zone, and testing for a correlation between that trait and the docu- mented level of assortative mating, versus ecologically relevant attributes. In some species, it might also be possible to conduct reciprocal transplants among patches in a mosaic hybrid zone. If individuals survive to reproduction equally well in the different patches, but suffer reduced mating success, that would indicate that the explanation provided here may be the more appropriate one. In most complex biological systems, intuitively understanding the outcome of a particular process can be difficult, if not impossible. In my thesis I have developed models demonstrating how the use of additional theoretical tools can aid in our pursuit of a fuller understanding. Often these models can provide unexpected insights and, additionally, exciting directions for future research, both empirical and theoretical. 80 Bibliography Albersheim, P. and A. J. Anderson-Prouty, 1975. Carbohydrates, proteins, cell-surfaces, and the biochemistry of pathogenesis. Annu. Rev. Plant Physiol. 26:31–52. Als, T. D., R. Vila, N. P. Kandul, D. R. Nash, S. Yen, Y. Hsu, A. A. Mignault, J. J. Boomsma, and N. E. Pierce, 2004. The evolution of alternative parasitic life histories in large blue butterflies. Nature 432:386–390. André, J. and B. Godelle, 2006. The evolution of mutation rate in finite asexual populations. Genetics 172:611–626. Baer, C. F., M. M. Miyamoto, and D. R. Denver, 2007. Mutation rate variation in multicellular eukaryotes: causes and consequences. Nature Reviews: Genetics 8:619–631. Bailey, R. I., C. D. Thomas, and R. K. Butlin, 2004. Premating barriers to gene exchange and their implications for the structure of a mosaic hybrid zone between Chorthippus brunneus and C. jacobsi (Orthoptera: Acrididae). J. Evol. Biol. 17:108–119. Barbour, A. G. and S. F. Hayes, 1986. Biology of Borrelia species. Microbiol. Rev. 50:381–400. Barrowclough, G. F., J. G. Groth, L. A. Mertz, and R. J. Gutierrez, 2005. Genetic structure between northern and Californian spotted owls (Strix occidentalis). Mol. Ecol. 14:1109–1120. Barton, N. H., 1979. The dynamics of hybrid zones. Heredity 43:341–359. ———, 1995. A general model for the evolution of recombination. Genetical Research 65:123–144. Barton, N. H. and G. M. Hewitt, 1985. Analysis of hybrid zones. Annu. Rev. Ecol. Syst. 16:113–148. ———, 1989. Adaptation, speciation and hybrid zones. Nature 341:497–503. Barton, N. H. and M. Turelli, 1991. Natural and sexual selection on many loci. Genetics 127:229–255. 81 Bibliography Barton, N. H. and M. C. Whitlock, 1997. The evolution of metapopulations. Pp. 183–210, in I. A. Hanski and M. E. Gilpin, eds. Metapopulation biology: ecology, genetics, and evolution. Academic Press, San Diego. Bazykin, A. D., 1969. Hypothetical mechanism of speciation. Evolution 23:685–687. Benham, G. S., 1974. A synopsis of the obligate and facultative insect parasitic nematodes. J. Invertebr. Pathol. 24:263–270. Bialozyt, R., B. Ziegenhagen, and R. J. Petit, 2006. Contrasting effects of long distance seed dispersal on genetic diversity during range expansion. J. Evol. Biol. 19:12–20. Bierne, N., P. Borsa, C. Daguin, D. Jollivet, F. Viard, F. Bonhomme, and P. David, 2003. Introgression patterns in the mosaic hybrid zone between Mytilus edulis and M. galloprovincialis. Mol. Ecol. 12:447–461. Bierne, N., P. David, P. Boudry, and F. Bonhomme, 2002a. Assortative fertilization and selection at larval stage in the mussels Mytilus edulis and M. galloprovincialis. Evolution 56:292–298. Bierne, N., P. David, A. Langlade, and F. Bonhomme, 2002b. Can habitat specialisation maintain a mosaic hybrid zone in marine bivalves? Mar. Ecol-Prog. Ser. 245:157–170. Boots, M. and R. G. Bowers, 1999. Three mechanisms of host resistance to microparasites – avoidance, recovery and tolerance – show different evolutionary dynamics. J. Theor. Biol. 201:13–23. Brannan, L. R., C. M. R. Turner, and R. S. Phillips, 1994. Malaria parasites undergo antigenic variation at high rates in vivo. Proc. R. Soc. Lond. B 256:71–75. Bridle, J. R., S. J. E. Baird, and R. K. Butlin, 2001. Spatial structure and habitat variation in a grasshopper hybrid zone. Evolution 55:1832–1843. Bridle, J. R. and R. K. Butlin, 2002. Mating signal variation and bimodality in a mosaic hybrid zone between Chorthippus grasshopper species. Evolution 56:1184–1198. Bridle, J. R., C. I. Saldamando, W. Koning, and R. K. Butlin, 2006. Assortative preferences and discrimination by females against hybrid male song in the grasshoppers Chorthippus brunneus and Chorthippus jacobsi (Orthoptera: Acrididae). J. Evol. Biol. 19:1248–1256. Britch, M. L., M. L. Cain, and D. J. Howard, 2001. Spatio-temporal dynamics of the Allonemobius fasciatus—A. socius mosaic hybrid zone: a 14-year perspective. Mol. Ecol. 10:627–638. 82 Bibliography Cain, M. L., V. Andreasen, and D. J. Howard, 1999. Reinforcing selection is effective under a relatively broad set of conditions in a mosaic hybrid zone. Evolution 53:1343–1353. Clark, J. S., 1998. Why trees migrate so fast: Confronting theory with dispersal biology and the paleorecord. Am. Nat. 152:204–224. Coulson, T., D. R. Macnulty, D. R. Stahler, B. Vonholdt, R. K. Wayne, and D. W. Smith, 2011. Modeling effects of environmental change on wolf population dynamics, trait evolution, and life history. Science 334:1275–1278. Coyne, J. A. and H. A. Orr, 1997. "Patterns of speciation in Drosophila" revisited. Evolution 51:295–303. Cruzan, M. B. and M. L. Arnold, 1993. Ecological and genetic associations in an iris hybrid zone. Evolution 47:1432–1445. Day, T., 2000. Sexual selection and the evolution of costly female preferences: spatial effects. Evolution 54:715–730. Day, T. and S. R. Proulx, 2004. A general theory for the evolutionary dynamics of virulence. Am. Nat. 163:E40–E63. De, A. and R. Durret, 2007. Stepping-stone spatial structure causes slow decay of linkage disequilibrium and shifts the site frequency spectrum. Genetics 176:969–981. Dieckmann, U., 2004. Assortative mating and spatial coexistence. Pp. 306–307, in U. Dieckmann, J. A. J. Metz, M. Doebeli, and D. Tautz, eds. Adaptive Speciation. Cambridge University Press. Doebeli, M. and U. Dieckmann, 2003. Speciation along environmental gradients. Nature 421:259–264. Donelson, J. E., 1995. Mechanisms of antigenic variation in Borrelia hermsii and African trypanosomes. The Journal of Biological Chemistry 270:7783–7786. Drake, J. W., B. Charlesworth, D. Charlesworth, and J. F. Crow, 1998. Rates of spontaneous mutation. Genetics 148:1667–1686. Eberhard, W. G., 1985. Sexual Selection and Animal Genitalia. Harvard University Press, Cambridge, Mass. Emms, S. K. and M. L. Arnold, 1997. The effect of habitat on parental and hybrid fitness: Transplant experiments with Louisiana irises. Evolution 51:1112–1119. Felsenstein, J., 1981. Skepticism towards Santa Rosalı́a, or why are there so few kinds of animals? Evolution 35:124–138. 83 Bibliography Feulner, P. G. D., F. Kirschbaum, and R. Tiedemann, 2008. Adaptive radiation in the Congo River: an ecological speciation scenario for African weakly electric fish (Teleostei; Mormyridae; Campylomormyrus). J. Physiol. Paris 102:340–346. Fisher, R. A., 1930. The Genetical Theory of Natural Selection. 1st ed. Oxford: Clarendon Press. Flor, H. H., 1942. Inheritance of pathogenicity in Melamspora lini. Phytopathology 32:653–669. ———, 1955. Host-parasite interaction in flax rust - its genetics and other implications. Phytopathology 45:680–685. ———, 1956. The complementary genic systems in flax and flax rust. Advances in Genetics 8:29–54. Frank, S., 1994. Recognition and polymorphism in host-parasite genetics. Philos. Trans. R. Soc. Lond., Ser. B 346:283–293. Frank, S. A., 2002. Immunology and the evolution of infectious diseases. Princeton University Press. Gabriel, D. W. and B. G. Rolfe, 1990. Working models of specific recognition in plant-microbe interactions. Annu. Rev. Phytopathol. 28:365–391. Galvani, A. P., 2003. Epidemiology meets evolutionary ecology. Trends Ecol. Evol. 18:132–139. Gandon, S., 2002. Local adaptation and the geometry of host-parasite coevolution. Ecol. Lett. 5:246–256. Gandon, S. and S. L. Nuismer, 2009. Interactions between genetic drift, gene flow, and selection mosaics drive parasite local adaptation. The American Naturalist 173:212–224. Gandon, S. and S. P. Otto, 2007. The evolution of sex and recombination in response to abiotic or coevolutionary fluctuations in epistasis. Genetics 175:1835–1853. Gardner, J. P. A., 1994. The structure and dynamics of naturally occurring hybrid Mytilus edulis and M. galloprovincialis (Bivalvia, Mollusca) populations: review and interpretation. Arch. Hydrobiol. 99 (supp):37–71. Gray, D. A. and W. H. Cade, 2000. Sexual selection and speciation in field crickets. Proc. Natl. Acad. Sci. USA 97:14449–14454. Grosberg, R. K. and M. W. Hart, 2000. Mate selection and the evolution of highly polymorphic self/nonself recognition genes. Science 289:2111–2114. Hamilton, W. D., 1980. Sex versus non-sex versus parasite. Oikos 35:282–290. 84 Bibliography Hamilton, W. D., R. Axelrod, and R. Tanese, 1990. Sexual reproduction as an adaptation to resist parasites (a review). Proc. Natl. Acad. Sci. USA 87:3566–3573. Haraguchi, Y. and A. Sasaki, 1996. Host-parasite arms race in mutation modifications: Indefinite escalation despite a heavy load? Journal of Theoretical Biology 183:121–137. Harrison, R. G. and S. M. Bogdanowicz, 1997. Patterns of variation and linkage disequilibrium in a field cricket hybrid zone. Evolution 51:493–505. Harrison, R. G. and D. M. Rand, 1989. Mosaic hybrid zones and the nature of species boundaries. Pp. 111–133, in J. A. Otte, D. Endler, ed. Speciation and its Consequences. Sinauer. Hollocher, H., C. Ting, F. Pollack, and C. Wu, 1997. Incipient speciation by sexual isolation in Drosophila melanogaster: variation in mating preference and correlation between sexes. Evolution 51:1175–1181. Hooge, M. D. and S. Tyler, 1999. Musculature of the facultative parasite Urastoma cyprinae (platyhelminthes). J. Morphol. 241:207–216. Howard, D. J., 1986. A zone of overlap and hybridization between two ground cricket species. Evolution 40:34–43. Howard, D. J. and P. G. Gregory, 1993. Post-insemination signalling systems and reinforcement. Philos. Trans. Roy. Soc. B 340:231–236. Howard, D. J. and R. G. Harrison, 1984a. Habitat segregation in ground crickets – experimental studies of adult survival, reproductive success, and oviposition preference. Ecology 65:61–68. ———, 1984b. Habitat segregation in ground crickets – the role of interspecific competition and habitat selection. Ecology 65:69–76. Howard, D. J., G. L. Waring, C. A. Tibbets, and P. G. Gregory, 1993. Survival of hybrids in a mosaic hybrid zone. Evolution 47:789–800. Huber, B. A., 2003. Rapid evolution and species-specificity of arthropod genitalia: fact or artifact? Org. Divers. Evol. 3:63–71. Ibrahim, K. M., R. A. Nichols, and G. M. Hewitt, 1996. Spatial patterns of genetic variation generated by different forms of dispersal during range expansion. Heredity 77:282–291. Irwin, D. E., S. Bensch, and T. D. Price, 2001. Speciation in a ring. Nature 409:333–337. Ishii, K., H. Matsuda, Y. Iwasa, and A. Sasaki, 1989. Evolutionarily stable mutation rate in a periodically changing environment. Genetics 121:163–174. Jiggins, C. D. and J. Mallet, 2000. Bimodal hybrid zones and speciation. Trends Ecol. Evol. 15:250–255. 85 Bibliography Johansson, J. and J. Ripa, 2006. Will sympatric speciation fail due to stochastic competitive exclusion? Am. Nat. 168:572–578. Johnson, T., 1999. Beneficial mutations, hitchhiking and the evolution of mutation rates in sexual populations. Genetics 151:1621–1631. Karlin, S. and J. McGregor, 1974. Towards a theory of the evolution of modifier genes. Theor. Popul. Biol. 5:59–103. Kimura, M., 1967. On the evolutionary adjustment of spontaneous mutation rates. Genetical Research 9:23–34. Kirkpatrick, M., 1982. Sexual selection and the evolution of female choice. Evolution 36:1–12. Kusch, J. and H. J. Schmidt, 2001. Genetically controlled expression of surface variant antigens in free-living protozoa. Journal of Membrane Biology 180:101–109. Lande, R., 1981. Models of speciation by sexual selection on polygenic traits. Proc. Natl. Acad. Sci. USA 78:3721–3725. ———, 1982. Rapid origin of sexual isolation and character divergence in a cline. Evolution 36:213–223. Le Corre, V., N. Machon, R. J. Petit, and A. Kremer, 1997. Colonization with long-distance seed dispersal and genetic structure of maternally inherited genes in forest trees: a simulation study. Genet. Res. 69:117–125. Leigh, E. G., 1970. Natural selection and mutability. The American Naturalist 104:301–305. ———, 1973. The evolution of mutation rates. Genetics 73 (supp):1–18. Lewontin, R. C., 1988. On measures of gametic disequilibrium. Genetics 120:849–852. Liberman, U. and M. W. Feldman, 1986. Modifiers of mutation rate: a general reduction principle. Theoretical Population Biology 30:125–42. Lively, C. M., 2010. An epidemiological model of host-parasite coevolution and sex. J. Evol. Biol. 23:1490–1497. MacCallum, C. J., B. Nurnberger, N. H. Barton, and J. M. Szymura, 1998. Habitat preference in the Bombina hybrid zone in Croatia. Evolution 52:227–239. Macholán, M., S. J. E. Baird, P. Munclinger, P. Dufková, B. Bı́mová, and J. Piálek, 2008. Genetic conflict outweighs heterogametic incompatibility in the mouse hybrid zone? BMC Evolutionary Biology 8:271. 86 Bibliography Mallet, J., 2008. Hybridization, ecological races, and the nature of species: empirical evidence for the ease of speciation. Phil. Trans. R. Soc. B 363:2971–2986. Mallet, J., A. Meyer, P. Nosil, and J. L. Feder, 2009. Space, sympatry and speciation. J. Evol. Bio. 22:2332–2341. May, R. M. and R. M. Anderson, 1983. Epidemiology and genetics in the coevolution of parasites and hosts. Proc. R. Soc. Lond. B 219:281–313. Maynard Smith, J. and J. Haigh, 1974. The hitchhiking effect of a favourable gene. Genetical Research 23:23–25. Maynard Smith, J. and E. Szathmary, 1997. The major transitions in evolution. Oxford University Press. M'Gonigle, L. K., J. J. Shen, and S. P. Otto, 2009. Mutating away from your enemies: the evolution of mutation rate in a host-parasite system. Theor. Popul. Biol. 75:301–311. Miller, M. R., A. White, and M. Boots, 2005. The evolution of host resistance: tolerance and control as distinct strategies. J. Theor. Biol. 236:198–207. Moran, N. A. and J. J. Wernegreen, 2000. Lifestyle evolution in symbiotic bacteria: insights from genomics. Trends Ecol. Evol. 15:321–326. Morin, L., B. A. Auld, and J. F. Brown, 1993. Interaction between Puccinia xanthii and facultative parasitic fungi on Xanthium occidentale. Biol. Control. 3:288–295. Mukai, T., 1964. The genetic structure of natural populations of Drosophila melanogaster. i. spontaneous mutation rate of polygenes controlling viability. Genetics 50:1–19. Muller, H. J., 1928. The measurement of gene mutation rate in Drosophila, its high variability, and its dependence upon temperature. Genetics 13:279–357. Nee, S., 1989. Antagonistic co-evolution and the evolution of genotypic randomization. Journal of Theoretical Biology 140:499–518. Nichols, R. A. and G. M. Hewitt, 1994. The genetic consequences of long-distance dispersal during colonization. Heredity 72:312–317. Nuismer, S. L. and S. P. Otto, 2004. Host-parasite interactions and the evolution of ploidy. Proc. Natl. Acad. Sci. USA 101:11036–11039. ———, 2005. Host-parasite interactions and the evolution of gene expression. PLoS Biol. 3:1283–1288. Okasha, S., 2007. Evolution and the levels of selection. Oxford University Press. Otto, S. P. and Y. Michalakis, 1998. The evolution of recombination in changing environments. Trends Ecol. Evol. 13:145–151. 87 Bibliography Otto, S. P., M. R. Servedio, and S. L. Nuismer, 2008. Frequency-dependent selection and the evolution of assortative mating. Genetics 179:2091–2112. Pagel, M., 1994. Detecting correlated evolution on phylogenies: a general method for the comparative analysis of discrete characters. Proc. R. Soc. Lond. B 255:37–45. Pal, C., M. D. Maciá, A. Oliver, I. Schachar, and A. Buckling, 2007. Coevolution with viruses drives the evolution of bacterial mutation rates. Nature 450:1079–1081. Panhuis, T. M., R. Butlin, M. Zuk, and T. Tregenza, 2001. Sexual selection and speciation. Trends Ecol. Evol. 16:364–371. Payne, R. J. H. and D. C. Krakauer, 1997. Sexual selection, space, and speciation. Evolution 51:1–9. Pays, E. and D. Nolan, 1998. Expression and function of surface proteins in Trypanosoma brucei. Molecular and Biochemical Parasitology 91:3–36. Peters, A. D. and C. M. Lively, 1999. The red queen and fluctuating epistasis: a population genetic analysis of antagonistic coevolution. Am. Nat. 154:393–405. ———, 2007. Short- and long-term benefits and detriments to recombination under antagonistic coevolution. J. Evol. Biol. 20:1206–1217. R Development Core Team, 2008. R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria. URL http://www.R-project.org. ISBN 3-900051-07-0. Reynolds, B. D., 1936. Colpoda steini, a facultative parasite of the land slug Agriolimax agrestis. J. Parasitol. 22:48–53. Ross, C. L. and R. G. Harrison, 2002. A fine-scale spatial analysis of the mosaic hybrid zone between Gryllus firmus and Gryllus pennsylvanicus. Evolution 56:2296–2312. Ryan, M. J. and W. Wilczynski, 1988. Coevolution of sender and receiver: effects of local mate preference in cricket frogs. Science 240:1786–1788. Saitoh, S. and H. Katakura, 1996. Strictly parapatric distribution of flightless leaf beetles of the Chrysolina angusticollis species complex (Coleoptera, Chrysomelidae) in the vicinity of Sapporo, northern Japan. Biol. J. Linn. Soc. 57:371–384. Sasaki, A., 1994. Evolution of antigen drift/switching: continuously evading pathogens. J. Theor. Biol. 168:291–308. Schilthuizen, M., 2000. Bimodal hybrid zones and the scale of a snail. Trends Ecol. Evol. 15:469. 88 Bibliography Schluter, D., 2000. The Ecology of Adaptive Radiation. Oxford University Press. Seehausen, O. and J. J. M. Van Alphen, 1999. Can sympatric speciation by disruptive sexual selection explain rapid evolution of cichlid diversity in lake victoria? Ecol. Lett. 2:262–271. Seehausen, O., J. J. M. Van Alphen, and F. Witte, 1997. Cichlid fish diversity threatened by eutrophication that curbs sexual selection. Science 277:1808–1809. Seger, J., 1985. Unifying genetic models for the evolution of female choice. Evolution 39:1185–1193. ———, 1988. Dynamics of some simple host-parasite models with more than two genotypes in each species. Philosophical Transactions of the Royal Society of London 319:541–555. Seger, J. and W. D. Hamilton, 1988. Parasites and sex. Pp. 176–193, in R. Michod and B. R. Levin, eds. The Evolution of Sex. Sinauer. Sniegowski, P. D., P. J. Gerrish, T. Johnson, and A. Shaver, 2000. The evolution of mutation rates: separating causes from consequences. BioEssays 22:1057–1066. Sniegowski, P. D., P. J. Gerrish, and R. E. Lenski, 1997. Evolution of high mutation rates in experimental populations of E. coli. Nature 387:703–705. Stahl, E. A., G. Dwyer, R. Mauricio, M. Kreitman, and J. Bergelson, 1999. Dynamics of disease resistance polymorphism at the Rpm1 locus of Arabidopsis. Nature 400:667–671. Stoenner, H. G., T. Dodd, and C. Larsen, 1982. Antigenc variation of Borrelia hermsii. J. Exp. Med. 156:1297–1311. Svard, S. G., T. Meng, J. M. Hetsko, M. L. McCaffery, and F. D. Gillin, 1998. Differentiation-associated surface antigen variation in the ancient eukaryote Giardia lamblia. Molecular Microbiology 30:979–989. Thompson, J. C. and L. Moewus, 1964. Miamiensis avidus n. g., n. sp., a marine facultative parasite in the ciliate order hymenostomatida. J. Protozool. 11:378–381. Thompson, J. N. and J. J. Burdon, 1992. Gene-for-gene coevolution between plants and parasites. Nature 360:121–125. Thrall, P. H. and J. J. Burdon, 2003. Evolution of virulence in a plant host-pathogen metapopulation. Science 299:1735–1737. Tian, D., M. B. Traw, J. Q. Chen, M. Kreitman, and J. Bergelson, 2003. Fitness costs of R-gene-mediated resistance in Arabidopsis thaliana. Nature 423:74–77. 89 Bibliography Tregenza, T., V. L. Pritchard, and R. K. Butlin, 2000. The origins of premating reproductive isolation: testing hypotheses in the grasshopper Chorthippus parallelus. Evolution 54:1687–1698. Turner, C. M. and J. D. Barry, 1989. High frequency of antigenic variation in Trypanosoma brucei rhodesiense infections. Parasitology 99:67–75. Turner, C. M. R., 1997. The rate of antigenic variation in fly-transmitted and syringe-passaged infections of Trypanosoma brucei. FEMS Microbiology Letters 153:227–231. Turner, G. F. and M. T. Burrows, 1995. A model of sympatric speciation by sexual selection. Proc. R. Soc. Lond. B 260:287–292. Van Der Ploeg, L. H. T., D. Valerio, T. De Lange, A. Bernards, P. Borst, and F. G. Grosveld, 1982. An analysis of cosmid clones of nuclear dna from Trypanosoma brucei shows that the genes for variant surface glycoproteins are clustered in the genome. Nucleic Acids Res. 10:5905–5923. Van Doorn, G. S., U. Dieckmann, and F. J. Weissing, 2004. Sympatric speciation by sexual selection: a critical reevaluation. Am. Nat. 163:709–725. Van Valen, L., 1973. A new evolutionary law. Evolutionary Theory 1:1–30. Vines, T. H., S. C. Kohler, A. Thiel, I. Ghira, T. R. Sands, C. J. MacCallum, N. H. Barton, and B. Nurnberger, 2003. The maintenance of reproductive isolation in a mosaic hybrid zone between the fire-bellied toads Bombina bombina and B. variegata. Evolution 57:1876–1888. Weissing, F. J., P. Edelaar, and G. S. Van Doorn, 2011. Adaptive speciation theory: a conceptual review. Behav. Ecol. Sociobiol. 65:461–480. Wilson, A. B., K. Noack-Kunnmann, and A. Meyer, 2000. Incipient speciation in sympatric Nicaraguan crater lake fishes: sexual selection versus ecological diversification. Proc. R. Soc. Lond. B 267:2133–2141. 90 Appendix A Mutating away from your enemies (Chapter 2) A.1 Cost of deleterious mutations Here we modify the model by including an explicit cost to modifiers that in- crease the mutation rate. This cost reflects the assumption that mutators will suffer from a higher load of unconditionally deleterious mutations at sites other than the A-locus. Rather than model these other sites explicitly, we incorpo- rate a cost Ciµi of having a mutation rate of µi. Technically this is equivalent if the unconditionally deleterious mutations are lethal. To incorporate costs of producing unconditionally deleterious mutations we add the following step be- tween mutation/selection (equation 2.2) and sex/recombination (equation 2.3). x00i,1 = (1 Ciµi,m)/C̄i x0i,1 x00i,2 = (1 Ciµi,M)/C̄i x0i,2 x00i,3 = (1 Ciµi,m)/C̄i x0i,3 x00i,4 = (1 Ciµi,M)/C̄i x0i,4 C̄i = 1 Ciµi,m(x0i,1 + x0i,3) Ciµi,M(x0i,2 + x0i,4) (A.1) A mutation modifier M that increases the mutation rate by (µi,M µi,m) will ex- perience a reduction in frequency each time step that is proportional to Ci(µi,M µi,m). In order that costs do not dominate dynamics at the modifier locus (i.e., they are not so large as to drive all mutation rates to 0) we must assume that the Ci are of order e2, where e is a small term, as discussed prior to equation (2.10). An analysis similar to that in the main text shows that a modifier that increases 91 Appendix A the mutation rate invades whenever the following is positive Dpi,M = c3[c1Xi(1 µī,M µi,M)/2 (µī,M µi,M)(Xi (1 2µī,M))] Ci(1 pi,M)pi,M(µi,M µi,m) (A.2) where c3 is the positive constant c3 = (R 4p f 1)f csc2[f]a2h(1 yi)(1 pi,M)pi,M(µi,M µi,m)(dp,A[0])2 p log[R](R2 2 cos[f]XiR+ X2i )(ai + 2)2(1 2µī,M) (A.3) Costs always select against higher mutation rates, but as long as the costs are sufficiently weak that equation (A.2) is positive when the mutation rate in a species is zero, then evolution will lead the system toward a non-zero level of mutation. The evolutionarily attracting mutation rate can be determined by setting Dpi,M equal to 0 using equation (A.2) and numerically solving for µi. We used simulations to investigate model behaviour when costs differed be- tween hosts and parasites (see main text for simulation methods). In accordance with Haraguchi and Sasaki (1996), we found that when mutation rates were al- lowed to co-evolve from low initial values (zero in each species), populations experienced an initial phase of selection for high mutation rates followed by a subsequent phase of selection for decreased rates in the species bearing the higher costs (see trajectories beginning at the origin in figure A.1). In contrast to Haraguchi and Sasaki (1996), however, mutation rates in our model remained at non-zero values in both species (convergence toward equilibrium values oc- curred from both directions; figure A.1). Furthermore, unlike Haraguchi and Sasaki (1996), we did not observe any asymmetries between hosts and parasites. This is because parasite populations did not require hosts for their survival in our model, whereas in the model of Haraguchi and Sasaki (1996), parasite sur- vival was dependent on their ability to find a host. Consequently, Haraguchi and Sasaki (1996) find that it is more often the hosts that tend to retreat towards a zero mutation rate, due to their having less at stake. Interestingly, we did not observe a large effect of population size when costs were included (compare left and right panels of figure A.1). Because sufficiently large costs can outweigh any indirect benefits to mutation rate modifiers, selec- tion for higher mutation rates is likely to disappear while rates are still low enough to permit co-evolutionary cycling (as was the case for the parameters presented in figure A.1). When this happens, the added indirect benefits created by stochastic allele frequency fluctuations in small populations have only a small effect and thus do not significantly change the final outcome of co-evolution. 92 Appendix A A.2 Solving the recursion equations for the disequilibrium Here we derive recursion equations for the disequilibrium in species i. By defi- nition D00i = x00i,1 x00i,4 x00i,2 x00i,3 We can then use equations (2.2) and (2.3) to write D00i in terms of genotype frequencies from the previous generation. Making the following substitutions xi,1 = 1/2 (1 pi,M)(1 2di,A) + Di xi,2 = 1/2 pi,M(1 2di,A) Di xi,3 = 1/2 (1 pi,M)(1+ 2di,A) Di xi,4 = 1/2 pi,M(1+ 2di,A) + Di gives us a recursion equation for Di. We next assume that di,A is small (specifically on the order of a small term e) and that the mutation rate modifier has a small effect (i.e. µi,M µi,m is also of order e). A simple analysis of the first three terms in in the Taylor series of Di (taken with respect to e) shows that the O(e2) term is the first that could possibly grow in a single generation (both the O(1) and O(e) terms only decay by a factor proportional to Xi = (1 yi)(1 2µi,M)). It follows that any initial disequilibrium in the system will rapidly evolve to a level that is at most of order e2. We, therefore, make the assumption that disequilibrium is already of this order in order to avoid any transient effects of initial conditions. Our recursion equation for Di can then be simplified to Di[t] = (1 yi) (1 2µi,M)Di[t 1]+ 2pi,M(1 pi,M)(µi,m µi,M) ✓ di,A[t 1] + ai(ai + 2)dī,A[t 1] ◆ +O(e3). We assume that enough time has passed that the initial disequilibrium Di[0] ex- erts negligible influence. With this assumption we can solve the above recursion equation to get Di[t] = 2(1 yi)(µi,m µi,M)(1 pi,M)pi,M t  t=1 ((1 yi)(1 2µi,M))t 1 ✓ di,A[t t] + ai(ai + 2)dī,A[t t] ◆ +O(e3). 93 Appendix A A.3 Different generation times Here we relax the assumption made in the main text that host and parasite generation times are equal. We assume that parasite dynamics remain exactly as they are in the main text, and then vary the generation times in hosts in two ways. In the first model we assume that all hosts are subject to mutation and selec- tion at every time step, but that only a proportion g take part in sexual repro- duction. This is incorporated by replacing equation (2.3) in the main text (for hosts only), with x00h,1 = (1 g)x0h,1 + g(x0h,1 yhD0h) = x0h,1 gyhD0h x00h,2 = (1 g)x0h,2 + g(x0h,2 + yhD0h) = x0h,2 + gyhD0h x00h,3 = (1 g)x0h,3 + g(x0h,3 + yhD0h) = x0h,3 + gyhD0h x00h,4 = (1 g)x0h,4 + g(x0h,4 yhD0h) = x0h,4 gyhD0h (A.1) It is clear from equation (A.1) that when g = 1 the model reduces to that pre- sented in the main text. It is also clear that this modification effectively scales recombination rates in hosts by the factor g. This has the net effect of reducing the recombination rate in hosts (by a factor g), which will ultimately increase the evolutionarily attracting mutation rate, provided it is predicted to lie within a region where cycles occur (see figure 2.3). Thus, increasing host generation time (smaller g) results in hosts evolving a higher mutation rate. In our second model we assume that all hosts are subject to selection every time step, but only a proportion g take part in sexual reproduction, at which point mutation and recombination occur. This can be incorporated by replacing our recurrence equation (2.3) in the main text (for hosts only) with x000h,1 = (1 g)x0h,1 + g(x00h,1 yhD00h ) x000h,2 = (1 g)x0h,2 + g(x00h,2 + yhD00h ) x000h,3 = (1 g)x0h,3 + g(x00h,3 + yhD00h ) x000h,4 = (1 g)x0h,4 + g(x00h,4 yhD00h ) (A.2) where x0h,j denotes the frequency of genotype j after selection and x 00 h,j denotes the frequency of genotype j after selection and mutation. Again it is clear that when g = 1 our model reduces to that in the main text. A similar analysis to 94 Appendix A that in the main text then yields dh,A[t] = Rt 1(1 2gµh,M) ah(ah + 2) sin[ft] sin[f] dp,A[0] dp,A[t] = Rt 1 q R2 (1 2gµh,M)(1 2µp,M)sin[ft+ s]sin[f] dp,A[0] where R, f and s are: R = s (1 2µp,M)(1 2gµh,M) (1+ ap/2+ ah/2)(1+ ap/2)(1+ ah/2) f = cos 1 (1 µp,M gµh,M) R s = cos 1 24 µp,M gµh,Mq R2 (1 2µp,M)(1 2gµh,M) 35 The recursion equation for Dh[t] (equation 2.10) also changes becoming: Dh[t] = 2(1 yh)g(µh,m µh,M)(1 ph,M)ph,M t  t=1 Xt 1i ✓ dh,A[t t] + ah(ah + 2)dp,A[t t] ◆ where Xh = 1 gyh 2gµh,M + 2gµh,Myh. An analysis following that in the main text shows that a modifier of mutation rate in hosts can invade whenever the following is positive sign[ Dph,M ] ' sign ⇥ c1Xh(1 µp,M gµh,M)/2+ (gµh,M µp,M)(Xh (1 2µp,M)) ⇤ (A.3) and in parasites sign[ Dpp,M ] ' sign ⇥ c1Xp(1 gµh,M µp,M)/2+ (µp,M gµh,M)(Xp (1 2gµh,M)) ⇤ (A.4) From the first of these equations, it can be seen that lengthening the generation time in hosts (reducing g) increases the relative importance of the first term in equation (A.3), which will tend to increase the evolutionarily attracting muta- tion rate in the host. The reverse is true in the parasite, where all else being equal a longer host generation time effectively reduces the host's mutation rate per unit time, which reduces selection for increased mutation rates in the para- site. 95 Appendix A A.4 One-species model Here we consider a simplified version of the model above where instead of two interacting species there is only a single species under a fluctuating selection regime. As before we consider two loci, each with two alleles. The four possible haplotypes are then am, aM, Am and AM. We let xj denote the frequency of the jth genotype (labelled in the order just given). The A-locus is assumed to be under a fluctuating selection regime. That is, the fitness of genotype j in the tth generation is given by wj[t] = ( 1+ a sin[bt] if j = 1, 2 1 a sin[bt] if j = 3, 4 (A.1) As in the two-species model we construct recursion equations, beginning with selection, then mutation, and finally sex and recombination. Performing a change of variables allows us to describe the system with recursions for the linkage disequilibrium (D[t]), the departure from a frequency of 0.5 at the A locus (d[t]), and the frequency of the modifier (pM[t]). In addition to the order assumptions made in the two-species model (see the main text preceding equa- tion 2.5), we must assume here that selection is also weak (a is on the order of e; this assumption wasn't necessary in the two-species model, but there we as- sumed dh,A and dp,A were both small). It can then be shown in a similar manner to that in appendix A.2 that, transient initial effects aside, the disequilibrium will remain of order at most e2. With these assumptions we can find a general solution for d[t] which, after simplification, is given by d[t+ 1] = a(1 2µM)((1 2µM) sin[bt] sin(b[t 1])) 4(2µ2M + (1 cos[b](1 2µM))) . (A.2) As in the two species model we next find a general solution for the linkage disequilibrium (D[t]), again ignoring transient effects of initial conditions. This can be simplified to D[t] = (1 y)(1 pM)pM(µm µM) t  t=1 ((1 y)(1 2µM))t 1(2d[t t] a sin[b(t t)]) (A.3) Finally we consider the change in frequency of a mutation rate modifier, av- eraged over a complete selective cycle. The expected change per generation is 96 Appendix A then E[DpM[t]] ⇡ (1 y)(µm µM)pM(1 pM) a 2(cos[b](1+ Xq) X q) (1+ X2 2X cos[b])(1+ q2 2q cos[b]) (A.4) where q = (1 2µM). It is apparent from equation (A.4) that the direction of selection on a modifier depends on the speed of evolutionary cycles (b), with faster cycles (larger b) selecting for higher mutation rates. In contrast to the two-species model, where the direction of selection was also dependent on the strength of selection (because the speed of evolutionary cycles f contained both ah and ap), we find that the strength of selection, a, in the one-species model does not affect the sign of selection on a modifier, although it does affect the strength of indirect selection. 97 Appendix A 0.00 0.05 0.10 0.15 0.00 0.05 0.10 0.15 ● ● ● 0.00 0.05 0.10 0.15 ● ● ● Cp < Ch Cp = Ch Cp > Ch Host Mutation Rate Pa ra sit e M ut at ion R at e Figure A.1: Co-evolutionarily trajectories in simulations with differing costs associated with mutation rate modifiers. Each curve represents the mean of 100 replicate runs. Solid points represent final values after 107 generations for trajectories initialized with both host and parasite mutation rates set to 0 and + symbols represent final values for populations initialized with mutation rates of 0 in hosts and 0.165 in parasites, or 0.165 in hosts and 0 in parasites. The shaded region indicates where cycles are not expected in the analytical model (R < 1). Population sizes were 109 (left panel) and 103 (right panel) in both hosts and parasites. When costs were equal we used Ch = Cp = 0.5 and when they differed we used Ch = 0.25,Cp = 0.5 and Ch = 0.5,Cp = 0.25. The dotted line has a slope of one and is included for reference. Recombination rates were set to 0, ah to 0.9, and ap to 1/(1+ ah) 1 (so that the strength of selection was the same in both species). The parameters in this figure are identical to figure 2.5 with ah = 0.9 and ap = 9, except for the addition of costs of mutation. 98 Appendix B Ploidy and the evolution of parasitism (Chapter 3) B.1 Additional analyses Quasi-linkage equilibrium analyses were also performed relaxing the assump- tion that yH and yP were near one. The main difference is that the departure from Hardy-Weinberg at the A locus (denoted FA,H in hosts and FA,P in para- sites) then becomes substantial (see table B.1 for the full expressions for FA,H and FA,P). Full expressions for w̄diff are given in table B.2. Because FA,H and FA,P are still on the same order as the strength of selection, however, these terms again drop out when we assume selection is weak and focus on the leading order terms. B.2 Simulations We summarize here a number of extensions to our model that were examined using simulations. Simulations matched our analytical predictions when mu- tation rates were high and cycles were small in every case, except when we considered more alleles at the interaction locus (discussed below). Where dis- crepancies were observed, they could be explained by accounting for the allele frequency dynamics (i.e., calculating dH and dP in every generation and using these in table 3.3). Also note that in each case, any shifts that did occur did not affect our main conclusions that parasitism is more likely to evolve under MAM than IMAM (but see results with three alleles) and that haploid parasites are more likely to evolve higher parasitism levels than diploids. • Figure B.1: simulations of MAM/IMAM comparing low and high muta- tion rates (µ = 10 1 versus µ = 10 5). 99 Appendix B • Figure B.2: same as figure B.1, except with selection in hosts reduced (aH = 0.01). • Figure B.3: same as figure B.1, except with stronger selection in hosts (aH = 0.5) and in parasites (0 < aP < 1). • Figure B.4: same as figure B.1, except with population sizes of 103 in both species. • Figure B.5: same as figure B.1, except with some hosts reproducing asexu- ally (yH = 0.2). Mutations were introduced at the same rate during sexual and asexual reproduction. • Figure B.6: same as figure B.1, except with different generation times in hosts and parasites. Only 20% of hosts reproduced at each time step, and thus had, on average, a generation time five times that of parasites. Mutations were introduced only during reproduction, so that hosts had 20% the rate of mutations per unit times as parasites. • Figure B.7: same as figure B.1, except with three alleles at the A-locus in each species. See below for a more detailed description of this case. • Figure B.8: evolutionary convergent level of parasitism in the GFG with conditional costs to virulence. • Figure B.9: evolutionary convergent level of parasitism in the GFG with unconditional costs to virulence. • Figure B.10: same as figure B.9, except with population sizes of 103 in both species. B.2.1 Three alleles at the interaction locus Qualitative shifts occurred in all cases when there were three alleles at the A- locus. These shifts could be described analytically by developing the model explicitly for multiple alleles. With MAM, the region where parasitism evolved shrunk with more alleles, whereas with IMAM it grew (figure B.7). This is because with MAM the heterozygous parasites can infect a lower proportion of the genotypes when there are more alleles present, whereas the opposite is true in the IMAM. For example, with three or more alleles, a host homozygous for an allele not present in a heterozygous parasite cannot be invaded by that parasite in the matching-alleles model, whereas when there are only two alleles 100 Appendix B a heterozygous parasite can invade all possible host genotypes. In the inverse- matching-alleles model the opposite is true, with new host genotypes providing additional targets for heterozygous parasites (in the case of two alleles, het- erozygous parasites cannot invade any host genotypes, but with three they can invade hosts homozygous for the allele that the parasite does not carry). Large cycles had the same effect with three alleles as they did with two, reducing and enlarging the same regions (figure B.7). 101 A ppen dix B Table B.1: Equations for FA,H and FA,P when yH and yP are not assumed to be near 1. Model Host Parasite w̄diff ploidy ploidy MAM/IMAM 1 2 FA,P = (1/4 d2P)2 fmm(aP + bP)(1 yP)/yP 2 1 FA,H = (1/4 d2H)2 fmaH(1 yH)/yH 2 2 FA,H = (1/4 d2H)2(3/4 d2P)2 fmmaH(1 yH)/yH, FA,P = (1/4 d2P)2(1/4+ d2H + FA,H)2 fmm(aP + bP)(1 yP)/yP GFG 1 2 FA,P = (1/4 d2P)2(cP,u + fmmcP,c fmm(aP + bP)(1/2+ dH))(1 yP)/yP 2 1 FA,H = (1/4 d2H)2(cH fmaH(1/2 dP))(1 yH)/yH 2 2 FA,H = (1/4 d2H)2(cH fmmaH(1/2 dP)2)(1 yH)/yH, FA,P = (1/4 d2P)2(cP,u + fmmcP,c + fmm(aP + bP)(d2H dH 3/4))(1 yP)/yP 102 A ppen dix B Table B.2: Full equations for w̄diff = w̄M w̄m, without assuming high levels of sexual reproduction in either species (e.g., without assuming yH and yP are on the order of 1 e). Model Host Parasite w̄diff ploidy ploidy MAM 1 1 (aP bP)/2+ (aP + bP)(2dHdP) 1 2 (aP 3bP)/4+ (aP + bP)(2dHdP + d2P + FA,P) 2 1 (3aP bP)/4+ (aP + bP)(2dHdP d2H FA,H) 2 2 (5aP 3bP)/8+ (aP + bP)[4dHdP(1+ dHdP) 3d2H + d2P+ (1+ 4d2H)FA,P (3 4d2P)FA,H + 4FA,PFA,H ]/2 IMAM 1 1 (aP bP)/2 (aP + bP)(2dHdP) 1 2 (aP 3bP)/4 (aP + bP)(2dHdP d2P FA,P) 2 1 (aP 3bP)/4 (aP + bP)(2dHdP d2H FA,H) 2 2 (aP 7bP)/8 (aP + bP)[4dHdP(1 dHdP) d2H d2P (1+ 4d2H)FA,P (1+ 4d2P)FA,H 4FA,HFA,P/]2 GFG 1 1 (3aP bP)/4+ (aP + bP)(dP dH(1 2dP))/2 cP,c(1+ 2dP) 1 2 (7aP bP)/8+ (aP + bP)(2(1+ 2dH)(1 dP)dP dH 2(2dH + 1)FA,P)/4 cP,c(3/4+ (1 dP)dP FA,P) 2 1 (5aP 3bP)/8+ (aP + bP)(3dP + 2(1 2dP)(d2H dH + FA,H))/4 cP,c(1+ 2dP) 2 2 (13aP 3bP)/16+ (aP + bP)[d2H(1 2dP)2 dH(1 2dP)2 3(dP 1)dP+ (1 2dP)2FA,H + (4d2H 4dH 3)FA,P + 4FA,HFA,P]/4 cP,c(3/4+ (1 dP)dP + FA,P) 103 Appendix B 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.10 0.05 0.1 0 0.05 0.1 A B C D E F G H I J K L M N O P _P `P +=10<1 +=10<5 MAM IMAM MAM IMAM Haploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Haploid Parasites Diploid Hosts / Diploid Parasites 0.0 0.5 1.0 Figure B.1: Evolutionary convergent level of parasitism ( f ) in MAM and IMAM. Right two columns are identical to figure 3.4. Left two columns report simula- tions with µ = 10 1 for comparison. Dashed red lines denote the analytical invasion condition assuming small cycles (table 3.4). Cells are shaded based on the mean level of parasitism present in the population after 106 generations of evolution in a single simulation (darker = higher, see grayscale in panel P). Parameters were aH = 0.05,yH = yP = 1, r = 0.5 and population sizes of 106 in both species. 104 Appendix B 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.10 0.05 0.1 0 0.05 0.1 A B C D E F G H I J K L M N O P αP βP µ = 10−1 µ = 10−5 MAM IMAM MAM IMAM Haploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Haploid Parasites Diploid Hosts / Diploid Parasites Figure B.2: Same as figure B.1, except with aH = 0.01. 105 Appendix B 0 0.5 1 0 0.5 1 0 0.5 1 0 0.5 1 0 0.5 1 0 0.5 10 0.5 1 0 0.5 1 A B C D E F G H I J K L M N O P αP βP µ = 10−1 µ = 10−5 MAM IMAM MAM IMAM Haploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Haploid Parasites Diploid Hosts / Diploid Parasites Figure B.3: Same as figure B.1, except with stronger selection in hosts (aH = 0.5) and in parasites (note change in axes ranges). 106 Appendix B 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.10 0.05 0.1 0 0.05 0.1 A B C D E F G H I J K L M N O P αP βP µ = 10−1 µ = 10−5 MAM IMAM MAM IMAM Haploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Haploid Parasites Diploid Hosts / Diploid Parasites Figure B.4: Same as figure B.1, except with population sizes of 103 in both species. 107 Appendix B 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.10 0.05 0.1 0 0.05 0.1 A B C D E F G H I J K L M N O P αP βP µ = 10−1 µ = 10−5 MAM IMAM MAM IMAM Haploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Haploid Parasites Diploid Hosts / Diploid Parasites Figure B.5: Same as figure B.1, except with some hosts reproducing asexually (yH = 0.2). Mutations were introduced at the same rate during sexual and asexual reproduction. 108 Appendix B 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.10 0.05 0.1 0 0.05 0.1 A B C D E F G H I J K L M N O P αP βP µ = 10−1 µ = 10−5 MAM IMAM MAM IMAM Haploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Haploid Parasites Diploid Hosts / Diploid Parasites Figure B.6: Same as figure B.1, except with different generation times in hosts and parasites. Only 20% of hosts reproduced at each time step, and thus had, on average, a generation time five times that of parasites. Mutations were in- troduced only during reproduction, so that hosts had 20% the rate of mutations per unit times as parasites. 109 Appendix B 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.10 0.05 0.1 0 0.05 0.1 A B C D E F G H I J K L M N O P αP βP µ = 10−1 µ = 10−5 MAM IMAM MAM IMAM Haploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Haploid Parasites Diploid Hosts / Diploid Parasites Figure B.7: Same as figure B.1, except with three alleles at the A-locus in each species. Dashed red lines denote the analytical invasion condition for the two- allele case, as given in table 3.4, and are included for comparison. 110 Appendix B 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 A B C D E F G H αP βP cP,u = 0 cP,c = 0.01 cP,c = 0.05 Haploid Hosts / Haploid Parasites Diploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Diploid Parasites Figure B.8: Evolutionary convergent level of parasitism ( f ) in the GFG with conditional costs to virulence, as indicated by column headings. Dashed red lines denote the value of aP for which w̄diff in eq. (3.6) is zero. Other parameters were as in figure B.1, along with µ = 10 5 and cH = 0.01. 111 Appendix B 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 A B C D E F G H αP βP cP,c = 0, cP,u =0.05 f i = 0.1 f i = 0.5 f i = 0.9 Haploid Hosts / Haploid Parasites Diploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Diploid Parasites Figure B.9: Evolutionary convergent level of parasitism ( f ) in the GFG with unconditional costs to virulence and differing initial levels of parasitism, as in- dicated by column headings. Dashed red lines denote the value of aP for which w̄diff in eq. (3.7) is zero. The poor fit for small bP (grey triangular regions to the right of the dashed red lines) is a consequence of selection being insufficiently strong to maintain the costly virulent allele, thereby reducing the advantage of being parasitic. Thus cycles do not occur and parasitism does not evolve. Other parameters were as in figure B.1, along with µ = 10 5 and cH = 0.01. 112 Appendix B 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 0 0.05 0.1 A B C D E F G H αP βP cP,c = 0, cP,u =0.05 f i = 0.1 f i = 0.5 f i = 0.9 Haploid Hosts / Haploid Parasites Diploid Hosts / Haploid Parasites Haploid Hosts / Diploid Parasites Diploid Hosts / Diploid Parasites Figure B.10: Evolutionary convergent level of parasitism ( f ) in the GFG with unconditional costs to virulence, differing initial levels of parasitism, and small populations. All parameters are as in figure B.9, except population sizes were 103. 113 Appendix C Sexual selection enables co-existence (Chapter 4) C.1 Model description We consider an individual-based model with discrete non-overlapping genera- tions in one- or two-dimensional continuous space with wrap-around bound- aries. Below, we describe the two-dimensional model, from which the corre- sponding one-dimensional model is readily generated by removing the spatial y-dimension. Each individual has a spatial location and is characterized by a display trait (expressed only in males) and a preference trait (expressed only in females). In our main set of model runs, these traits are assumed to be gov- erned by separate unlinked haploid loci, each with two alleles (display alleles are denoted by Q/q and preference alleles by P/p). Each generation, N individ- uals are produced and compete for resources, with those experiencing stronger competition being more likely to die before reaching reproductive maturity. Re- sources in our model may be interpreted in the broadest possible sense, describ- ing the biotic and abiotic factors that are subject to local ecological competition. Among the individuals surviving ecological competition, females choose mates, with the probability of a specific male being chosen depending on her mating preference and the spatial distance separating them. Females produce offspring in proportion to their fecundities. Offspring then disperse from their natal loca- tion and the parents die. Below we detail these steps in the order in which they occur. The names and descriptions of all parameters and variables are listed in Table C.1. 114 Appendix C C.1.1 Competition for resources The habitat at each location (x, y) is characterized by the local density k(x, y) of available resources. The total amount of resources over the spatial arena is given by K = RR k(x, y) dx dy. The function relating resource gain to survival is chosen such that if every individual received an equal share of these resources, the expected number of survivors would be K. Consequently, we refer to k(x, y) as the local carrying capacity and to K as the total carrying capacity. Except for Figs. 4.4 and C.7, we investigate a local carrying capacity that is bimodal, with peaks located at (x, y) = (0.25, 0.25) and (0.75, 0.25). To do so, we combine two Gaussian functions according to k(x, y) = exp( (x 0.25) 2 + (y 0.25)2 2s2k ) + exp( (x 0.75) 2 + (y 0.25)2 2s2k ) + b ! k0 , (C.1) where sk denotes the width of the two Gaussian peaks. The parameters b and k0 allow us to adjust the average height and degree of variation in k(x, y). Specif- ically, the height is adjusted such that the total carrying capacity equals K, and the degree of variation is adjusted to give the desired relation between peaks and troughs. Specifically, we measure the degree of spatial variation in local carrying capacity as v = max k(x, y) min k(x, y) min k(x, y) . (C.2) A value of v = 0.25 therefore means that the local carrying capacity is 25% higher at the peaks than at the troughs. For Fig. C.7, landscapes are generated in a similar way, except that the heights and widths of the two peaks differ. For Fig. 4.4, the landscape is generated by adding white noise to the baseline level, filtered to have a reasonable amount of spatial autocorrelation, with the highest peak set to twice the height of the lowest trough. Through competition, each individual obtains a share of the local carrying capacity, which we refer to as its resource share, ri = k(xi, yi) Âj nij , (C.3) where nij is the contribution of individual j to the effective density of competi- tors at the location of individual i, and the sum extends over all N individu- als. The competitive impact of individual j on individual i decreases with the distance dij separating them, according to a Gaussian function with standard 115 Appendix C deviation ss, nij = exp( d2ij/(2s2s ))/(2ps2s ) ; (C.4) in the one-dimensional model, the divisor is p 2pss. Note that the effect nii of an individual i on itself declines as ss increases, because the individual then competes for resources over larger distances and thus has less of a negative impact on its available resources. As defined, the resource share of an individual i is approximately K/N. This can be seen by assuming that the N individuals in the population are distributed over the arena according to the local carrying capacity, so that their expected density is N k(x, y)/K. Replacing the sum over individuals in Eq. C.3 with an integral over space, we obtain ri = k(xi, yi)RR N k(x,y) K exp( d2ij/(2s2s )) 2ps2s dx dy ⇡ K/N +O(v) , (C.5) where the second line assumes that spatial variation in the local carrying ca- pacity is low. In our individual-based model runs, departures from the above occur due to clumping, fecundity variation over space (Section C.1.4), as well as discrepancies due to replacing the sum in Eq. C.3 with the integral in Eq. C.5 (especially when ss is very small or large relative to the arena). That said, the mean resource share is typically close to K/N in our model runs. In Fig. C.1 we show the effect of spatial variation in local carrying capac- ity k(xi, yi) on various components of fitness, including the resource share, ri. Interestingly, ecological competition is weaker (ri is higher) in regions of low carrying capacity (Fig. C.1A), increasing the survival probability of individuals in these regions (si, Fig. C.1B). This occurs because females are less likely to encounter preferred males wherever the carrying capacity is low, causing their fecundity to be lower due to increased mate-search costs (ci; Section C.1.4 and Fig. C.1C). Consequently, fewer offspring are produced than expected based on the low local carrying capacity, resulting in weaker competition among those offspring. The net result of lower ecological competition and higher search costs in regions with low local carrying capacity is that females have roughly equal fitness across space. C.1.2 Survival We assume that individuals that gain more resources are more likely to survive to reproductive maturity. The probability si of such survival is assumed to 116 Appendix C be zero when an individual fails to gain any resources, to rise approximately linearly with its resource share ri when that share is small, and to taper off at a maximal survival probability of smax (ranging between 0 and 1). Specifically, we use a hyperbolic (or Holling type-2) function (Coulson et al., 2011) to relate resource share to the probability of survival, si = smax 1+ r/ri , (C.6) where r is the resource share that must be obtained for an individual to sur- vive with a probability equal to half the maximal survival probability. Unless stated otherwise, we assume that the maximum probability smax of surviving to reproductive maturity equals 1. The value of r is chosen to ensure that, on average, K individuals survive to reproduce if all individuals obtain an equal share of resources (ri = K/N). By setting the expected survival probability si to K/N in Eq. C.6 and substituting ri = K/N, we obtain r = smax K/N. With this choice of r, approximately K individuals survive each generation (with a variance that is typically small). For example, in Fig. C.1, the average survival probability was 0.484, close to the ex- pected value of K/N = 1/2. While competition for resources causes substantial mortality, survival probabilities across the arena differ only slightly (Fig. C.1B). Importantly, the survival of an individual does not depend on whether or not it is a hybrid. C.1.3 Mating Of the individuals that survive to mate, the probability that female i chooses male j as a mate depends on whether his display trait matches her preference trait and on the spatial distance separating them. Females bearing a P (p) al- lele prefer males bearing a Q (q) allele by a factor a. We assume that females encounter males in the vicinity of their home location. Specifically, each female spends a proportion of time at distance dij from her home that is described by a Gaussian distribution with standard deviation sf, so that her encounter probability eij with a male at distance dij is proportional to eij = exp( d2ij/(2s2f ))/(2ps2f ) ; (C.7) in the one-dimensional model, the divisor is p 2psf. In our main model, we assume that females encounter resources and males over the same spatial scales (i.e., sf = ss); we relax this assumption in Fig. C.5. The probability that female i 117 Appendix C chooses male j as a mate is proportional to pij = adij 1eij , (C.8) where dij equals 1 when the display trait of male j matches the preference trait of female i, and 0 otherwise. Once a female chooses a mate, we assume that all her offspring are sired by that male (monogamy). C.1.4 Reproduction The fecundity of a female i is given by: fi = fmax(1 ci) , (C.9) where fmax is the maximum fecundity and ci (ranging from 0 to 1) measures the cost associated with finding a preferred mate for female i. The factor 1 ci is assumed to be zero when there are no preferred males locally, to rise approxi- mately linearly with the local density of preferred males, µi =  males j pij , (C.10) and to taper off at 1 when preferred mates are readily encountered, resulting in maximal fecundity. Specifically, we use a hyperbolic (or Holling type-2) function (Doebeli and Dieckmann, 2003), 1 ci = 11+m/µi , (C.11) where m is the value of µi at which a female's fecundity is halved by mate-search costs. Because µi is obtained by summing over the entire male population, its value can be large, on the order of the number of surviving males, so values of m on the order of the surviving population's size K are needed for costs to be appreciable. This is why we express m relative to K, specifying the ratio m/K in the figures. We refer to ci as the mate-search cost of female i and to m as the strength of mate-search costs. Unless noted otherwise, we use m = 500. In our main simulations (with m/K = 1), mate-search costs reduce female fecundity by about 50%, on average, from the maximum fecundity (Fig. C.1C), with relatively minor differences in fecundity among females over space. Other values for m are explored in Fig. 4.3. For m = 0, all females have equal and maximal fecundity. As m is raised, fecundity declines, on average, and becomes more variable, with females in low- 118 Appendix C density regions or surrounded by non-preferred males having lower fecundity (Fig. C.2). After mating, offspring are produced. Inheritance at both loci is Mendelian, and we assume no linkage between the display and preference loci, except where noted (Section C.2.6). To allow us to explore various parameters re- lating to competition and mate-search costs independently, we hold the total number of offspring constant at N. For each offspring, a mother is chosen in proportion to the females' fecundities. Consequently, the maximum fecundity fmax only matters insofar as it is high enough to result in at least N offspring being produced across the population. Similar patterns are observed when fmax is fixed and offspring numbers are given by a Poisson distribution with a mean of fi for each female (data not shown). We consider N to be the total number of offspring surviving the phase during which resources are largely provided by the parents, after which the offspring migrate and begin the next phase of competition for resources. C.1.5 Movement Each offspring moves from its mother's location according to a distance drawn from a Gaussian function with mean 0 and standard deviation sm. Movements occur in all directions with equal probability. C.2 Model extensions To assess the robustness of our results, we consider several extensions and/or modifications to our main model described above. C.2.1 Incorporating mating-dependent dispersal To compare our results with those of Payne and Krakauer (1997), we con- sider mating-dependent dispersal. In their model, male movement distances are lower for males with better mating prospects, and we thus assume that the movement distance of male j is drawn from a Gaussian function with mean 0 and standard deviation sm,j = sm exp ✓ l Âi pij Âik pik ◆ , (C.12) where l determines how quickly movement distances decrease with increasing mating prospects and pij is given by Eq. C.8 in Section C.1.3. For l = 0, the above reduces to our main model. We find that the addition of mating-dependent 119 Appendix C dispersal in males extends coexistence times only marginally, if at all (compare Fig. C.3A to C.3B). We also examine the related case in which males with low mating prospects move farther, but again, coexistence times are not appreciably prolonged in our individual-based model. C.2.2 Introducing multiple allelic types To examine whether coexistence of more than two types is possible, we extend our main model so that one of n alleles p1, . . . , pn can occur at the preference locus and one of n alleles q1, . . . , qn can occur at the display locus. Specifically, in Fig. 4.4, we consider n = 10 preference and display types. A female with preference allele pi prefers males with display allele qi to all other males by the factor a. All other components of mate choice remain the same as for our main model with n = 2 mating types. C.2.3 Allowing competition to impact fecundity In our main model, competitive interactions reduce the survival probability of an individual. Alternatively, individuals that gain fewer resources might sur- vive, but have lower fecundity. To explore this possibility, we allow all N off- spring to survive, while reducing their reproductive success according to the impact of competition, as measured by si. Specifically, for males, the probabil- ity of being chosen as a mate is set to pij = adij 1eijsi. Likewise for females, fecundity is set to fi = fmax(1 ci)si. Such competition-dependent fecundity generates less demographic stochasticity, because all individuals reach repro- ductive maturity and can mate, albeit with reduced probability when their re- source share ri is low. Indeed, all else being equal, incorporating competitive effects on fecundity, rather than survival, enables long-term coexistence over a wider range of parameters (compare Fig. C.6 to Fig. 4.3). C.2.4 Altering the strength of density-dependent competition To measure the strength of density dependence on survival, we may define l = r/(1 K/N), with r = smax K/N (Section C.1.2). In our main model, the maximum survival rate smax is set to 1 so that l = 1, indicating that survival is strongly density-dependent. At the other extreme, if smax is set to K/N, all individuals survive with probability smax = K/N, regardless of their resource share, so there is no density-dependent effect on survival (l = 0). As shown in Fig. C.4B, coexistence does not occur in the absence of density dependence (l = 0); spatial variation in local carrying capacity then becomes irrelevant and 120 Appendix C cannot stabilize mating domains in space. As the importance of competition increases (larger l, or equivalently, larger smax), coexistence can occur over a broader parameter space. Once about half of the mortality is due to density- dependent competition (l > 0.5), results become similar to those for l = 1 (Fig. C.4B). We also explored the effects of density-dependent competition by varying the total carrying capacity K (Fig. C.4C), while holding the total number of off- spring constant. Because we are interested in the effects of population size per se, we also keep constant the relative strength of mate-search costs (m/K = 1), so the ease with which females encounter preferred mates remains unaffected by variation in K. Again, when density-dependent competition is weak (K near N), coexistence requires much higher levels of spatial variation in local carrying capacity. Conversely, when K is very small, stochasticity in survival becomes so great that coexistence is not maintained without high levels of spatial variation in local carrying capacity. Thus, intermediate values of K, relative to N, best facilitate co-existence. The effects of demographic stochasticity can also be seen in Fig. C.4D, where the strength of density dependence and the expected survival probability K/N are held constant (l = 1 and K/N = 1/2), while the total number N of offspring is varied, as is the time point at which co-existence is evaluated (at generation 5N). Again, we also keep constant the relative strength of mate-search costs (m/K = 1). All else being equal, larger population sizes facilitate the mainte- nance of coexisting types, as expected given the reduced stochasticity. C.2.5 Altering the spatial scale of competition, mate-search, and migration In the main model, we specified the spatial scale of several phenomena, in- cluding the breadth of the competition function (ss = 0.05), the breadth of the mate-search function (sf = 0.05), and movement distances (sm = 0.05). Fig. C.5 illustrates the minimum level of variation, v, in local carrying capacity required for coexistence to occur in our simulations for at least 5N generations. Coex- istence is easier to maintain when females search in smaller regions for males (small sf) and when movement is localized (small sm ) because mating types that predominate in different spatial locations remain more isolated. By con- trast, coexistence is easier to maintain when competition occurs across a broad spatial range (large ss) because individuals near the resource peaks compete more strongly for resources available in the troughs, reducing the population size in the troughs and promoting isolation of the mating types near each peak. 121 Appendix C C.2.6 Incorporating alternative genetic architectures Our main model assumes free recombination between the trait and preference loci. Fig. C.8 explores the effect of linkage, finding no substantive differences between complete linkage and free recombination between the trait and prefer- ence loci. To test whether our findings are robust to changes in the number of loci, we consider a quantitative genetic model in which an individual's preference and display trait are determined by two quantitative characters. This model can be interpreted as assuming that a large (infinite) number of additive loci code for each of the two traits. Complementing our main model, which features a finite number of alleles, this extension allows for arbitrarily many mating types. In this quantitative genetic model, the probability that female i mates with male j is proportional to pij = exp( (pi qj)2/(2s2p))eij , (C.13) where pi qj is the difference between the preference trait of female i and the display trait of male j, sp denotes the strength of female preference (smaller sp means females are choosier), and eij is proportional to the encounter probability between female i and male j, as defined in Eq. C.7. Offspring trait values are drawn from a Gaussian function centred at the mean of the parental phenotypes for each trait, with a standard deviation so that measures the variation among offspring due to segregation, recombination, and mutation. All other details of the quantitative genetic model are the same as for our main model. Despite the different genetic assumptions, the behaviour of the quantitative genetic model closely resembles that of the allelic model (Fig. C.8). Coexistence of mating domains is again possible over a wide range of parameters, provided female preferences are sufficiently strong (small sp). As in the allelic model, loss of mating domains in the quantitative model, when it happens, tends to occur through the replacement of one type by the other. Compared with the al- lelic model, the quantitative genetic model exhibits two additional mechanisms through which mating domains may be lost. First, when female preference is weak (large sp), interbreeding between adjacent mating domains may become so common that the resultant offspring form their own mating domains, fa- cilitating the merging of the original domains. Second, the random drift of matched trait and preference values in one mating domain may cause them to coincide by chance with the values in an adjacent mating domain, so the two originally separate domains may merge due only to the random genetic drift of quantitative mating traits that results from segregation, recombination, and mutation in finite populations. 122 Appendix C We note that our results bear some connection to those of Day (2000), who analyzed a spatial model of costly sexual preferences and male traits using a quantitative genetic model. A key difference is that Day (2000) incorporated ecological differences across space, with natural selection favoring different op- timal male traits as a function of position (essentially incorporating niche differ- ences across space). In parallel to our results, costly female preferences evolved and could maintain more diversity in male traits than expected based on the variation in the optimum trait across space. C.2.7 Incorporating asymmetric display costs Display traits can incur fitness costs in males. Our main model assumes that such costs, if present, affect all individuals equally. It may often be the case, however, that display traits differ in their effects on fitness. We therefore exam- ine what happens when the Q allele causes males to have a reduced survival probability relative to those carrying the q allele (i.e., for Q-bearing individ- uals, the survival probability si is reduced by a factor 1 a, with a ranging between 0 and 1). Provided that the resultant cost is not so strong that the sta- bilizing effect of spatial variation in local carrying capacity is overwhelmed by selection against Q-bearing males, our main findings remain largely unchanged (Fig. C.9). 123 Appendix C Symbol Eq. Description Model parameters a Strength of selection against Q-bearing males (only C.2.7) k(x, y) C.1 Local carrying capacity at location (x, y) l C.12 Strength of mating-dependence in male dispersal (only C.2.1) m C.11 Strength of mate-search costs smax C.6 Maximum survival probability v C.2 Spatial variation in local carrying capacity K Total carrying capacity N Number of offspring a C.8 Strength of female preference fmax C.9 Maximum female fecundity l Strength of density-dependent competition sf C.8 Width of mate-search distribution sk C.1 Width of peaks in local carrying capacity sm Width of movement distribution so Width of offspring distribution (only C.2.6) sp C.13 Width of female preference (only C.2.6) ss C.4 Width of competition distribution Model variables ci C.11 Mate-search costs of female i dij C.4 Spatial distance between individuals i and j eij C.7 Propensity for female i to encounter male j fi C.9 Fecundity of female i nij C.4 Competitive effect of individual j on individual i pij C.8 Propensity for female i to choose male j as a mate si C.6 Survival probability of individual i µi C.10 Local density of preferred males as seen by female i ri C.3 Resource share of individual i Table C.1: Model parameters and model variables. 124 Appendix C Figure C.1: Variation in three components of fitness as a function of the local carrying capacity experienced by each individual at t = 1000 for the model run in Fig. 4.1D. Individuals are coloured according to their genotype at the display locus. A. Resource share ri in males and females. B. Survival probability si of males and females. C.Mate-search costs ci of females. Lines show least-squares regression lines. 125 Appendix C Figure C.2: Mate-search costs for the model run in Fig. 4.1D. Panels in column A are identical to those in Fig. 4.1D, except that only females are shown and they are coloured according to their preference allele. Panels in column B show the costs associated with searching for a mate and rejecting non-preferred males for each female (Eq. C.9), as a function of her location y. For m/K = 1, female fecundity is typically only halved by mate-search costs. 0 0.5 1 A α = 5 m = 0, 0 0.5 1 0 500 1000 B α = 5 m = 0, Generation, t Fr eq ue nc y o f Q Figure C.3: Effects of mating-dependent dispersal in males. Panels show dis- tributions of allele frequencies at the display locus through time across 1000 replicate model runs in a two-dimensional homogeneous landscape; coexistence occurs only while these frequencies remain intermediate. Darker shading indi- cates a higher probability of observing a given frequency of the Q allele. PanelA is identical to Fig. 4.2B. Panel B is the same as A, except with mating-dependent dispersal in males (l = 100). Results for other values of l are qualitatively iden- tical. Model runs are initialized as in Fig. 4.2. All other parameters are as in Fig. 4.1B. 126 Appendix C Strength of female preference, α A 2 4 6 8 10 0.0 0.1 0.2 0.3 0.4 0.5 Strength of density−dependent competition, λ B 0 0.2 0.4 0.6 0.8 1 0.0 0.1 0.2 0.3 0.4 0.5 Expected survival probability, K N C 0 0.2 0.4 0.6 0.8 1 0.0 0.1 0.2 0.3 0.4 0.5 Total number of offspring, N D 500 1000 1500 2500 0.0 0.1 0.2 0.3 0.4 0.5 Sp at ial va ria tio n in ca rry ing ca pa cit y, v Figure C.4: Minimum level of spatial variation v in local carrying capacity needed to ensure long-term coexistence (grey regions) in a two-dimensional bi- modal landscape. v was increased until the average persistence time of 20 repli- cate runs exceeded 5N generations (vertical lines indicate standard errors). A. Effect of the strength a of female preference. Coexistence becomes more likely as female preferences become stronger (larger a), although once preference ex- ceeds a ⇡ 5, its impact is small. B. Effect of the strength l of density-dependent competition (varying smax while holding K = 500 and N = 1000). The limit l = 0 corresponds to completely density-independent survival, while the limit l = 1 corresponds to completely density-dependent survival. C. Effect of the expected survival probability K/N (N = 1000, l = 1, and m/K = 1 are held constant). Values near K/N = 0 correspond to very small mating populations, while the limit K/N = 1 corresponds to the absence of ecological competi- tion. D. Effect of population size N (holding constant K/N = 0.5, l = 1, and m/K = 1). All other parameters are as in Fig. 4.1D. 127 Appendix C 0.03 0.04 0.05 0.06 0.07 0.08 0.0 0.1 0.2 0.3 0.4 0.5 Movement distance, σm , competitive effect distance, σs , Sp at ial va ria tio n in ca rry ing ca pa cit y, v σm σf σs or mate−search distance, σf Figure C.5: Minimum level of spatial variation v in local carrying capacity needed to ensure long-term coexistence (grey regions) in a two-dimensional bimodal landscape. v was increased until the average persistence time of 20 replicate runs exceeded 5N generations (vertical lines indicate standard errors). The three curves show the effects of the competition width ss (red), the width of the mate-search distribution sf (green), and the width of the movement dis- tribution sm (blue), while holding all other parameters constant at their values in Fig. 4.1D. 0 1 2 0 0.25 0.5 > 5000 4000 − 5000 3000 − 4000 2000 − 3000 1000 − 2000 0 − 1000 Strength of mate−search costs in females, m K Sp at ial va ria tio n in ca rry ing ca pa cit y, v Figure C.6: Conditions for long-term coexistence with competition-dependent fecundity (Section C.2.3) in a two-dimensional bimodal landscape. All parame- ters are as in Fig. 4.3. 128 Appendix C 0.01 0.2 0 1 > 5000 4000 − 5000 3000 − 4000 2000 − 3000 1000 − 2000 0 − 1000 ● ● ● ● ● ● ● ● Width of peaks in local carrying capacity, σK Di ffe re nc e in pe ak h eig ht s Figure C.7: Effects of altering the shape of the local carrying capacity (Eq. C.1) in a two-dimensional bimodal landscape. Shading indicates how long polymor- phism persists at the display locus (darker = longer). Each cell represents the mean time to loss of polymorphism for 10 replicate model runs. Side panels indicate the extent of spatial variation in local carrying capacity along transects at y = 0.25 for nine parameter combinations indicated by the closest black cir- cle. The inset at the bottom center corresponds to the parameter combination used in Fig. 4.3. Spatial variation in local carrying capacity is relatively weak throughout this figure, with v ranging from 0.28 for sk = 0.01 (far left) to 0.038 for sk = 0.2 (far right). All other parameters are as in Fig. 4.1D. 129 Appendix C 0 0.25 0.5 0.0 0.1 0.2 0.3 A t = 5,000 0 0.25 0.5 B t = 25,000 Allelic, free recombination Allelic, no recombination Quantitative, σp = 0.1 Quantitative, σp = 0.2 Quantitative, σp = 0.3 Spatial variation in local carrying capacity, v Va ria nc e in dis pla y t ra it Figure C.8: Effects of changes in genetic architecture in a two-dimensional bi- modal landscape. Variance in display trait after 5, 000 (A) and 25, 000 (B) gen- erations for a variety of genetic architectures, averaged over 20 replicate model runs (vertical lines indicate standard errors). The dashed line indicates the max- imum possible variance in the allelic model (0.25). For determining variances in the allelic model, alleles Q and q are assigned trait values 0 and 1, respectively. In the quantitative genetic model, the initial preference/display trait values are set to 0/0 or 1/1 (corresponding to P/Q or p/q in the allelic model) with equal probability, yielding an initial variance of 0.25. Over time, the variance of 0.25 can be exceeded due to random genetic drift. For comparison, the red curve shows results of our main model. Model runs are initialized as in Fig. 4.2. All other parameters are as in Fig. 4.1; in the quantitative model, so = 0.01. 130 Appendix C 0 0.25 0.5 0.0 0.1 0.2 0.3 A t = 5,000 0 0.25 0.5 B t = 25,000 a = 0 a = 0.01 a = 0.05 a = 0.1 Spatial variation in local carrying capacity, v Va ria nc e in dis pla y t ra it Figure C.9: Effects of asymmetric fitness costs of display traits in the allelic model in a two-dimensional bimodal landscape. Variance in display trait after 5, 000 (A) and 25, 000 (B) generations when males bearing the Q allele have their survival lowered by a factor 1 a relative to males bearing the q allele, averaged over 20 replicate model runs (vertical lines indicate standard errors). The dashed line indicates the maximum possible variance in this allelic model (0.25). For comparison, the red curve (identical to that in Fig. C.8) shows results of our main model, corresponding to the limit a = 0. Model runs are initialized as in Fig. 4.2. All other parameters are as in Fig. 4.1. 131 Appendix D Hybrid zones (Chapter 5) D.1 Hybrid zone structure likelihood method Here we describe a likelihood method for fitting a series of steps to allele fre- quency data, as in figure 5.1. We consider a one-dimensional transect through a hybrid zone, with m "patches". We assume an initial step height at zero, before the first patch, and a final step height at one, after the last patch. A stepwise model through the hybrid zone will consist of k step locations s = {s1, s2, · · · , sk} and k 1 step heights (in addition to the first and last step heights fixed at zero and one), h = {h1, h2, · · · , hk 1}, where hi is the height between steps si and si+1. The stepwise model "partitions" the patches, and we will refer to a particular partition as the set of patches between two adjacent model steps. The step heights correspond to the estimated allele-frequency for each partition. Our method aims to quantify both the number and placement of these steps using a maximum likelihood approach. Suppose that the expected genotype frequencies within the ith patch are Ei[AA], Ei[Aa], and Ei[aa], and the ith patch has observed genotype counts of xi,AA, xi,Aa, and xi,aa for the genotypes AA, Aa and aa. The likelihood of sampling (with replacement) from the underlying frequencies is given by the multinomial probability✓ ni xi,AA xi,Aa xi,aa ◆ Ei[AA]xi,AAEi[Aa]xi,AaEi[aa]xi,aa where ni = xi,AA + xi,Aa + xi,aa is the number of individuals in the ith patch and ( nx1 x2 x3) is the multinomial coefficient defined as n!/(x1!x2!x3!). The expected genotypic frequencies are a function of the allele frequencies given by the partition allele frequency and the inbreeding coefficient for the patch, fi. Since fi is not directly of interest we set it to its most likely value, 132 Appendix D given the local genotypic values; that is, we define fi as fi = 1 xi,Aa2pi(1 pi)ni where pi is the observed allele frequency in the ith patch. If the ith patch is within the jth partition, it shares that partition's step height, hj, and its expected frequencies are Ei[AA] =h2j + fihj(1 hj) Ei[Aa] =2hj(1 hj)(1 fi) Ei[aa] =(1 hj)2 + fihj(1 hj) In cases where an expected frequency was negative, that frequency was set to zero and the other expectations were standardized appropriately. The likelihood of observing the data across all m patches can then be calcu- lated as Pr(x|s, h) = m ' i=1 ✓ ni xi,AA, xi,Aa, xi,aa ◆ Ei[AA]xi,AAEi[Aa]xi,AaEi[aa]xi,aa (D.1) where the product is taken over all m patches. For a given set of step locations s we find the heights h that maximize equation (D.1) using univariate optimiza- tion in R (R Development Core Team, 2008). For a given number of steps, k, equation (D.1) must be maximized with respect to step locations. For large data sets it is not feasible to exhaustively search for the best model, as the number of possible models is on the order of m! for m patches. Instead, we used a genetic algorithm to identify the best model for a given number of steps. The algorithm begins with a randomly gen- erated initial pool of k-step models. It then runs through multiple generations of mutation, recombination and selection. Mutation randomly replaces one or more steps within a model with other possible steps, recombination switches steps between different models (while maintaining step number), and selection samples the best models, weighted by their log-likelihood, to initiate the next generation. To find the number of statistically significantly steps required to best ex- plain the data, we started with the best single step model (k = 1) and added steps until the difference between the best k step model and k + 1 step model was not statistically significant following a likelihood ratio test. Each step re- quires two additional parameters (a step location and a step height), so we compared the likelihood ratio with a Chi-squared distribution with two de- 133 Appendix D grees of freedom. Simulations revealed that our use of two degrees of freedom was in fact a conservative assumption, while one degree of freedom was not. This method to fit mosaic hybrid zone data is available as an R package at http://www.zoology.ubc.ca/prog/mosaic/ 134
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Multi-species interactions and the evolution of biological systems M'Gonigle, Leithen K. 2011
Item Metadata
- Media
- 24-ubc_2012_spring_mgonigle_leithen.pdf [ 6.22MB ]
- Metadata
- JSON: 24-1.0072468.json
- JSON-LD: 24-1.0072468-ld.json
- RDF/XML (Pretty): 24-1.0072468-rdf.xml
- RDF/JSON: 24-1.0072468-rdf.json
- Turtle: 24-1.0072468-turtle.txt
- N-Triples: 24-1.0072468-rdf-ntriples.txt
- Original Record: 24-1.0072468-source.json
- Full Text
- 24-1.0072468-fulltext.txt
- Citation
- 24-1.0072468.ris
Full Text
Citation Scheme:
Usage Statistics
Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
async >
Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url: