Ecological opportunity and upward prey-predator radiation cascades

A general goal in community ecology and evolutionary biology is to understand how diversity has arisen. In our attempts to reach such goals we become increasingly aware of interacting ecological and evolutionary processes shaping biodiversity. Ecological opportunity and adaptive radiations can, for example, drive diversification in competitive communities but little is known about how such processes propagate through trophic levels in adaptive radiation cascades. I use an eco-evolutionary model of trait-based ecological interactions and micro-evolutionary processes to investigate the macro-evolutionary aspects of predator diversification in such cascades. Prey diversification facilitates predator radiation through predator feeding opportunity and disruptive selection. Predator radiation, however, often disconnects from the prey radiation as the diversification progresses. Only when predators have an intermediate niche width, high predatory efficiency, and high evolutionary potential can radiation cascades be maintained over macro-evolutionary time scales. These results provide expectations for predator response to prey divergence and insight into eco-evolutionary feedbacks between trophic levels. Such expectations are crucial for future studies that aim for a better understanding of how diversity is generated and maintained in complex communities.

www.nature.com/scientificreports www.nature.com/scientificreports/ the abundance of all other populations to whom it may interact see also 15,16,33 . The fitness of a predator is a function of its traits, the traits, and abundance of its prey and the traits, and abundance of other predators to which the focal population may compete with.
I implement the ecological model in an eco-evolutionary context with connected predator-prey adaptive radiations as emergent model outputs 20,34 . Prey-resource, prey-prey and predator-prey trait matchings dictate resource utilization, competition, and trophic interactions respectively. Predators will be selected to match the trait of their prey and prey will be selected to mismatch predators while trait matching within trophic levels will be selected against due to competition for resources. A mutant (the source of phenotypic variation in the model) that occurs in trait space where many similar and abundant populations already exist or where resources are low thus tend to have low invasion fitness. Contrary, mutants in trait space where resources are available and competition is low have high invasion fitness. I utilize these properties of the model to simulate emerging adaptive radiations and I focus on the properties of predators (e.g. niche width, efficiency, and evolvability) that can lead to upward radiation cascades while at the same time considering the prey properties that may facilitate downward effects from predators on prey. More specifically, I investigate the fitness landscape of predators throughout macro-evolutionary history to assess if predators are filling niche space that is constituted of already diversified prey or if predator diversification is driven by synchronized predator-prey diversification. I also match the detailed mechanistic drivers of radiations across trophic levels to a more empirically tangible measure of congruence between prey and predator phylogenetic trees.

Results
Irrespective of other parameters (e.g. predator efficiency (b max ), predator mutation rate (μ pred ) or predator niche width (σ a )) most of the predator diversity occurs when prey niche width (σ α ) is low (Figs. 2a,b, and S1 and S2). This is expected as a low σ α facilitates the co-existence of multiple prey populations and prey branching 20 . Prey diversity is also a prerequisite for predator diversification as the distribution of prey constitutes predator niche space and there needs to exist more than one phenotypically distinct prey for disruptive selection on predators to occur (Eq. 9 in the methods section). The distribution of prey tends to be wide with many phenotypically similar prey species when σ α is low (see prey radiation in Fig. 3a) which explains the positive relationship between predator diversity and σ α . Another striking result is a hump-shaped relationship between predator diversity and predator niche width (σ a ) (Figs. 2c and S2). Such patterns make sense, as specialized predators are sensitive to prey that may out-evolve them. If the prey evolves to peripheral parts or outside of the predator's niche the probability of extinction increases. In contrast, if predators are generalists the risk for prey evolving out of the predator niche width is low but multiple predator co-existence is reduced. Furthermore, high predator mutation rates (μ pred ) tend to reduce diversity, compared to low and intermediate μ pred (Fig. 2c). This is in line with Pontarp and Petchey 20 who show that high predator mutation rates tend to disrupt diversification both in prey and predators. Top consumers (a) with some trait z (e.g. birds of prey with body size z) and competitive prey (b) with trait u (e.g. granivorous birds with beak size u) that interact in a local habitat (c) defined by some implicit resource distribution with peak abundance as uopt and width σK. The three trophic levels are distributed on the same trait dimension (e.g. size) here illustrated by color. Competition between species is dictated by their niche width (black and gray Gaussian kernels), and I assume that populations with similar traits interact more than less similar ones. The invasion fitness of a mutant is thus a function of its traitmatching to its resources, the traits of its competitors on the same trophic level and their niche widths. Image created in Adobe Illustrator CS6 Version 16.0.0 (64-bit).
www.nature.com/scientificreports www.nature.com/scientificreports/ The results presented above provide insights into the eco-evolutionary processes of predator diversification (see for example [12][13][14] ) but they remain silent on whether such diversification is synchronized among trophic levels due to co-evolution in adaptive radiation cascades.
Focusing on radiation cascades, results show that when predator efficiency (b max ) and evolvability (μ pred ) are low (Fig. 3a), radiation cascades across macroevolutionary time scales do not exist. Indeed, it follows from the model assumptions (Eqs. 2 and 9) that prey diversity is required for disruptive selection on predators to exist. All predator radiations can thus be viewed as initiated by upward radiation cascades but predators radiate largely independently, filling up niche space without co-evolution with any particular prey species. This disconnection between radiations is confirmed by large regions of high invasion fitness in the predator fitness landscapes (Fig. 3a) and low congruence among phylogenetic trees throughout evolutionary time (light gray line in Fig. 4). www.nature.com/scientificreports www.nature.com/scientificreports/ Elucidating the effect of predator niche width (σ a ), efficiency (b max ) and evolvability (μ pred ) on co-radiation patterns and congruence I find that a change in σ a does not affect prey radiations when b max and μ pred are low. The disconnect in predator-prey radiations thus remains, irrespective of σ a . Congruence is, however, positively related to σ a but this is seemingly not due to increased co-radiation but rather due to a positive relationship between predator diversity and congruence. Fewer leaves in trees also lead to fewer leaves removed to reach isomorphism. However, when predator efficiency (b max ) is increased (Fig. 3b) predators start to affect prey radiations, disrupting prey branching and leaving gaps in prey niche space (see also 20 ). As a result, the radiations reveal a more distinct pattern of predators co-evolving with prey. With this being said, the prey is still radiating outside the niche range of predators, and large parts of the prey and predator radiations are thus disconnected. This is revealed by the predator fitness landscape which shows large regions of positive invasion fitness (Fig. 3). Congruence is dramatically increased ( Fig. 4) but in this case, such an increase seems to be due to the reduced prey diversity rather than actual co-radiation. When b max is even larger (b max = 0.0005-0.0007) predator radiations are disrupted altogether ( Fig. S3a,b) as the predator pushes prey radiations away from niche space. Only when b max is large in combination with a large σ a can the predator seemingly co-radiate with prey but still with some prey radiating outside the range of the predators (Fig. S3c,d). These results show an interaction effect between b max and σ a . An efficient specialist predator does not radiate as easy as an efficient generalist (results shown by the shift in diversity peak along the σ a axis in Fig. S2).
Finally, the effect of an increased μ pred affects prey radiation by increasing diversification rate but still when b max is low such radiation remains largely disconnected from the prey radiation (Fig. 3c). Congruence increases even further (Fig. 4), plausibly due to a combination of increased co-radiation and decreased overall species richness. Ultimately, a combined increase in μ pred and b max provides what looks most like large co-radiations that persist throughout macro-evolutionary history (Fig. 3d). The predators radiate and they spread out in traits space such that the proportion of independently radiating prey is low. Congruence is also high, although diversity is quite low because predators reduce the radiation of prey 20 .

Discussion
Brodersen, et al. 5 list three potential predator responses to the diversification of prey. Predators may remain adapted to and specialized on one prey, predators may evolve to be a generalist, or the predator may diversify. In this paper, I expand on the mechanisms for each of these outcomes by explicitly focusing on predator properties. I focus on a combination of predator niche width, evolvability, and efficiency. Each of these properties affects the probability of predator survival and diversification and they have been shown to affect prey diversification in a downward radiation cascade 20 . In concordance with classic adaptive dynamics and radiation theory 9 , I conclude that an intermediate predator niche width facilitates predator radiation. Predator efficiency and evolvability also dictate how connected the prey and predator radiations are. Predator radiations are initiated by prey diversification but the connection between prey and predator radiations can break down soon after the first predator branching. Prey radiation is thus a prerequisite for predator diversification as it creates niche space in which www.nature.com/scientificreports www.nature.com/scientificreports/ predators can diversify. Whether the radiations represent radiation cascades depends on the definition of a cascade. Brodersen, et al. 5 use a wide definition of cascades including seemingly disconnected cichlid radiations in Lake Victoria as an example. In this sense, my modeled radiation also falls into the category of cascades, although it can be argued that such predator radiations fall in the category of "classical" radiations as it is mainly driven by adaptation to empty niche space and competition for prey even though the niche space may be constructed by radiating prey.
In my results, radiations become more and more connected with a combined increase in predator evolutionary potential (related to mutation rate μ pred ) and efficiency (b max ). As predator efficiency increases, predator abundance tends to increase and predator extinction probability decreases. Increased predator population sizes also increase evolutionary potential, as do increased mutation rate. Predator evolvability thus also facilitates radiation cascades (Figs. 3 and 4). With this being said, a too efficient or a too fast-evolving predator can exclude prey from niche space, disrupting predator diversification altogether (Fig. S3). Interestingly, the diversity of the prey is also dramatically decreased in scenarios with efficient and fast-evolving predators, indicating a clear two-way process where predators affect the prey radiation and vice versa. Pontarp and Petchey 20 go into depth about the mechanistic underpinning of how predator diversification can restrict prey radiations through decreased population sizes, decreased disruptive selection and the exclusion of prey from parts of niche space. Such a downward effect has been shown before and it questions the concept of one-way (upward or downward) radiation cascades. Does a one-way radiation cascade ever exist, and if so how can we quantify it? Such questions also open for a discussion about the possibility to determine and quantify radiation cascades from data. In a simple model with known processes, like the one analyzed here, it is relatively easy to determine if cascades occur through combined visual inspection of the radiations, fitness landscape analysis and congruence analysis. Such a wealth of information is however not available in most empirical systems which make quantification difficult. Congruence may be the only measure that is available and my analysis point to some difficulties associated with such measures. Congruence analysis is, for example, highly dependent on community diversity and null model approaches may be needed to account for such effects. www.nature.com/scientificreports www.nature.com/scientificreports/ Despite the open questions about the definition and quantification of radiation cascades, the above provides sought after theory and testable expectations for when upward adaptive radiation cascades occur. For example, aquatic systems 35,36 have been suggested to include upward radiation cascades and data are thus available for such tests. With this being said, results should be viewed in the context of the model design and specific assumptions. Although explicitly designed to model adaptive radiations in predator-prey systems the model is quite general and simplistic as it models asexual organisms, linear functional responses, and I assume a constant environment and resource availability. The model also assumes a one-dimensional trait space and it includes explicit assumptions about competition and predation through trait matching. It is thus important to note that the results presented here are relevant for when a prey trait (e.g. size or some other ecological trait) dominates and drives the response in predator traits or vice versa. The signal of radiation cascades may be less prevalent when multiple correlated functional traits 37,38 affect eco-evolutionary dynamics within and among trophic levels. As noted in the introduction it may also seem unrealistic to compare communities that contain species with non-evolving narrow or broad niche widths, respectively. It has been shown that niche filling can occur through either evolved generalization or diversification and the outcome depends on resource diversity and resource-acquisition trade-off 39,40 . Given the scope of this paper where I explicitly focus on diversification, I exclude the possibility of generalization on the expense of diversification. I model wide resource diversity, I assume non-evolving niche widths, I do not model an acquisition trade-off, and I set parameters such that diversification can occur. The expectations about radiation cascades presented here are thus most relevant in the context of wide resource diversity and high cost of obtaining a generalist niche as this is a prerequisite for diversification 39 .
Despite the assumptions and simplifications listed above, the results provide an understanding of the fundamental link between eco-evolutionary processes and diversification and it reveals several of the mechanisms that are operating in diversifying complex communities. Such a detailed investigation is virtually impossible in experimental or field studies. By excluding some of the complications listed above the model reveals several interesting points associated with predator properties that may promote radiation cascades, the likely two-way process between co-radiating predator-prey radiations, and issues related to the quantification or causality of such cascades from data. This said I am sure this will not be the last words written about the complex and intriguing topic of diversification in complex communities.

Methods ecological model. Following Pontarp and Petchey 20 I base my trait-based model of adaptive predator-prey
radiations on the assumption that the diverging traits (e.g. body size) also have a direct influence on ecological interactions. I follow established approaches 21,22,24,32,33 by modeling resources along some generally defined trait dimension and I assume that both predator and prey populations are defined by some resource utilization trait, distributed along the same trait dimension (Fig. 1). The per capita growth (fitness) of a focal prey individual associated with a given population is thus a function of its resource utilization trait, the abundance of the individual's population, the local resource distribution and the abundance of all other populations to whom it may interact see also 15,16,33 . The fitness of a predator is a function of its trait, the traits, and abundance of its prey and the traits, and abundance of other predators to which the focal population may compete with. Mathematically the above formulate as: where N i and P k denote prey and predator population size respectively ultimately modeling the ecological dynamics, in per capita form, of n prey populations and p predator populations for i =1 to n, k =1 to p. The model is essentially built on the generalized Lotka-Volterra (GLV) model 41 but extended to include trait-based interactions as carrying capacity, the prey interactions and predator-prey interactions are formulated as trait dependent functions: Here, as in Pontarp and Petchey 20 , K i (u i , u opt ) represents the carrying capacity for a monomorphic population of prey individuals with trait value u i in a habitat characterized by a resource distribution with its peak resource availability at the point u opt . The parameter K 0 denotes the maximal carrying capacity (at u = u opt ) and it follows from Eq. 3 that the resource availability declines symmetrically as u deviates from u opt according to the width www.nature.com/scientificreports www.nature.com/scientificreports/ of the resource distribution (σ K ). The interaction coefficient, α ij (u i ,u j ), is a function that models the ecological interaction between the focal prey population (defined by its trait u i ) and its competitors (defined by their traits u j ). The competition coefficients is standardized such that α ii =1 and 0 <α ij < 1 (u i ≠u j ). The parameter σ α determines the degree of competition between individuals given certain utilization traits and can thus be viewed as the niche width of the prey. Similarly, a ik (u i ,z k ) models the interaction between a focal predator population k with trait value z and a prey population i with trait value u. The parameter b max denotes the maximum attack rate obtained when u i =z k and this rate then falls of symmetrically as u i deviates from z k according to a Gaussian function with variance σ a. Similar to the σ α parameter, σ a can be viewed as the niche width of the predator. eco-evolutionary simulation implementation. I implement the model in an eco-evolutionary context 20,34 . The prey-resource, prey-prey and predator-prey trait matchings (Eqs. 3-5) dictate resource utilization, competition, and trophic interactions respectively. Trait matching will be selected for between trophic levels while trait matching within trophic levels will be selected against due to competition for resources. A mutant that occurs in trait space where many similar and abundant populations already exist or where resources are low thus tend to have low invasion fitness. Contrary, mutants in trait space where resources are available and competition is low have high invasion fitness. I utilize these properties of the model to simulate adaptive radiations through alternating phases of determining equilibrium population abundances and the introduction of mutants as a source of phenotypic variation 34,42 .
I start the simulations by setting up the resource distribution with parameters u opt = 0, σ K = 1 and K 0 =10000. I implement the assumption of ecological opportunity by seeding the system with two monomorphic and perfectly matching populations (one prey and one predator) with trait values u = z = u opt , and I compute the equilibrium population sizes by integrating over Eqs. 1 and 2 until equilibrium or a steady state is reached (integrating from time 0 to 1000 is proven sufficient 34 ). I then allow one individual of a selected population to mutate and change its traits to u' or z', depending on what population was picked. More specifically, a single mutant population is drawn at random according to the product of the population sizes and mutation probability (µ prey = 0.01 and µ pred = 0.005 -0.1). Note that I am assuming a constant µ prey throughout this study while evolvability is one of the predator properties that I aim to evaluate in the context of radiation cascades, hence the range of µ pred analyzed during different model realizations. Mutation size for both predators and prey are assumed to be a random trait value drawn from a normal distribution with mean equal to the trait value of the mutating population and a variance σ mut = 0.02.
Following the adaptive dynamics framework 9,43 I compute trait dependent invasion fitness for an arbitrary prey and predator mutant with trait value u' or z' according to: where u' denotes the trait value of the mutant prey and z' denotes the trait value of a mutant predator. The vectors u, z, N, and P are defined as above containing the resident community trait distributions and abundances. If the mutant has positive invasion fitness, I continue the analysis with a mutual invasibility test, evaluating if a single individual with the mother trait can invade a population with the mutant trait at equilibrium population size. If the mutant invasion fitness is positive but mutual invasibility does not exist, the mutant will replace the resident. However, if mutual invasibility does exist, the resident and the mutant can co-exist. After the mutant is either introduced alongside the resident or is allowed to replace the resident, I recalculate the equilibrium, remove populations that may have gone extinct due to the introduction of the mutant population. Finally, before iterating into the next cycle of computation that constitutes an evolutionary time step I assign each population to a species id using a trait-based speciation definition where species are defined as populations having common descent and a continuous distribution of traits (no gaps in the trait distribution> 3* σ µ ) 15,16 . If a gap> 3* σ µ is detected in the trait distribution within an existing species, a speciation event is registered (i.e. one species branching into two). Thereafter I progress into the next evolutionary step, repeating the whole procedure for 5000 iterations.
Model analysis. It is known from previous studies with similar models that diversity is contingent on niche width 9,10,13,20 . I build on such insights and investigate prey and predator diversity as a function of evolutionary time and prey niche width (σ α ) ranging from 0.1-0.7 with increments of 0.2. I also analyze predator diversity as a function of predator niche width (σ a ) ranging from 0.1-0.7 with increments of 0.2, predator efficiency (b max ) ranging from 0.0001-0.0007 with increments of 0.0002, and predator mutation probability (μ pred ) of 0.005, 0.01 and 0.1. Such analyses provide a general understanding of the diversification process, which I use as a basis for the investigation of radiation cascades. I initiate simulations with constant resource distribution, constant mutation model for prey and constant prey niche widths (σ α = 0.1). This initiation ensures diversification of prey and I use it as a basis for my investigations of predator properties in the context of upward adaptive radiation cascades. The parameter space presented above is relevant for the diversification of predators as well as downward radiation cascades 20 . To quantify radiation cascades I analyze the fitness landscape (Eq. 9) and congruence between prey and predator phylogenetic trees measured as the maximum number of leaves removed to reach isomorphism among trees. Phylogenetic trees were constructed using data on the distance to the nearest common ancestor between the species and the MATLAB function seqlinkage.
www.nature.com/scientificreports www.nature.com/scientificreports/ Other constant model parameters for the simulations were r = 1, d = 0.2 and c = 0.3. These, as well as the parameters already presented above, produce diverse enough communities to analyze adaptive radiations within reasonable computational time. Given the stochastic components in the implementation of the deterministic ecological model as an eco-evolutionary simulation model, some variation between model realizations can occur. Each parameter combination throughout parameter space was thus replicated 10 times, the variation in predator and prey diversification was quantified (Fig. S4) and a robustness check across replicates was done on congruence (Fig. S5).