RESEARCH ARTICLE | ECOLOGY Social copying drives a tipping point for nonlinear population collapse Daniel Oroa,b,1 , Lluís Alsedàc,d, Alan Hastingsb,e , Meritxell Genovarta , and Josep Sardanyésd Edited by Jorge X. Velasco-Hernandez, Universidad Nacional Autonoma de Mexico, Juriquilla, Mexico; received August 25, 2022; accepted January 21, 2023 by Editorial Board Member Pablo A. Marquet Sudden changes in populations are ubiquitous in ecological systems, especially under perturbations. The agents of global change may increase the frequency and severity of Significance anthropogenic perturbations, but complex populations’ responses hamper our under- standing of their dynamics and resilience. Furthermore, the long-term environmental Among the complex dynamics and demographic data required to study those sudden changes are rare. Fitting dynamical arising in all living systems, models with an artificial intelligence algorithm to population fluctuations over 40 y in sudden population collapses are a social bird reveals that feedback in dispersal after a cumulative perturbation drives a one of the most fascinating. population collapse. The collapse is well described by a nonlinear function mimicking Understanding the mechanisms social copying, whereby dispersal made by a few individuals induces others to leave that may cause these collapses is the patch in a behavioral cascade for decision-making to disperse. Once a threshold fundamental to the conceptual for deterioration of the quality of the patch is crossed, there is a tipping point for a study of population dynamics. social response of runaway dispersal corresponding to social copying feedback. Finally, We fit dynamical models to dispersal decreases at low population densities, which is likely due to the unwillingness of the more philopatric individuals to disperse. In providing the evidence of copying population fluctuations over 40 y for the emergence of feedback in dispersal in a social organism, our results suggest a in a social bird that showed an broader impact of self-organized collective dispersal in complex population dynamics. unexpected collapse after a This has implications for the theoretical study of population and metapopulation non- perturbation press that linear dynamics, including population extinction, and managing of endangered and progressively eroded harvested populations of social animals subjected to behavioral feedback loops. environmental conditions at the Downloaded from https://www.pnas.org by 158.109.55.145 on September 7, 2023 from IP address 158.109.55.145. world’s most suitable breeding tipping points | runaway dispersal | nonlinear population dynamics | social behavior | feedback patch. We demonstrated that the Understanding abrupt declines in a populations’ response to environmental perturbations density-dependence feedback is crucial for the theoretical study of population extinction and for managing harvested related to the simple behavior of and endangered species, especially under the impacts of global change (1–3). Under the social copying for dispersal to conceptual framework of the logistic model, the transition from a population level near other patches explained the carrying capacity to collapse in populations subjected to perturbations should occur collapse. The significance of our through a negative exponential decay (4–6). However, the logistic model has several study lies in showing that limitations because it assumes both a linear association between density and growth rate stochastic perturbations may and a lack of time lag in the response of individuals to changes in density (5). Time-lagged trigger a tipping point by responses may generate transient phenomena, which can explain abrupt regime shifts that runaway dispersal, driving are not directly associated with environmental changes (7, 8). populations to a state of A crucial driver of the dynamics in animal populations is individual behavior, which quasi-extinction. influences relevant demographic parameters such as dispersal and the evolution of life histories (9). Information theory frames how populations may persist in stochastic envi- ronments when behavior is targeted for gathering information and reducing environmental uncertainty for making optimal decisions (e.g., to be philopatric or to disperse) (10). Individuals living in social groups have two sources of information for making decisions that may have fitness consequences: The first is private information that each individual gathers and stores on its own, and the second is public information in the form of decisions Author contributions: D.O. and M.G. designed research; and performances others make. Examples of public information are flight from a predator D.O. performed research; L.A. and J.S. contributed new reagents/analytic tools; D.O., L.A., A.H., M.G., and J.S. approaching the group or average fertility in a particular patch indicating its suitability analyzed data; and D.O. wrote the paper. to a prospecting individual (11–13). The use of public information may generate densi- The authors declare no competing interest. ty-dependent behavioral responses, such as runaway dispersal from a patch to escape from This article is a PNAS Direct Submission. J.X.V.-H. is a a perturbation. This is why dramatic sudden collapses in populations are likely to occur guest editor invited by the Editorial Board. in social organisms. Social groups are complex systems in which the number of interactions Copyright © 2023 the Author(s). Published by PNAS. This article is distributed under Creative Commons within a group is not additive but grows in a factorial manner resulting in behavioral Attribution-NonCommercial-NoDerivatives License 4.0 feedback loops, such as those arising for information gathering, social copying, and group (CC BY-NC-ND). cohesion. These feedback loops are autocatalytic, meaning that as more individuals perform 1 To whom correspondence may be addressed. Email: d.oro@csic.es. a particular behavior, more other individuals will decide to follow suit. The density of This article contains supporting information online at conspecifics drives individual decision-making for staying or dispersing in a patch when https://www.pnas.org/lookup/suppl/doi:10.1073/pnas.​ the trade-off between the benefits of foraging and the risks of predation are considered 2214055120/-/DCSupplemental. (14). This is linked with collective behavioral loops that trigger informational cascades Published March 6, 2023. PNAS 2023 Vol. 120 No. 11 e2214055120 https://doi.org/10.1073/pnas.2214055120 1 of 8 used in the locating resources and threats (15). Collective behavior a density-dependent, nonlinear process (20, 21) (Fig. 1). may influence complex population responses such as tipping Determining how populations of social organisms may collapse points for regime shifts, boom-bust dynamics, and chaotic dynam- and diagnosing what behavioral mechanisms are involved is ics, all of them having potential consequences for population hampered by the challenges of collecting long-term data on both extinction (16–19). population fluctuations and environmental changes, including A key process linking individual and collective behavior to a perturbations (22). Perturbations may have different temporal demographic trait with consequences for populations is deci- features, such as stochastic pulses and presses, but these dynamics sion-making for dispersal (19). How dispersal is influenced by are crucial to understanding the impact and responses they may copying behavior in social animals, including humans, remains generate on populations. unexplored, but empirical and theoretical studies show that cop- An unprecedented opportunity to investigate the effects of ying is favored under perturbations (4, 15). Less is even known social feedback on population nonlinear collapses is afforded by about the consequences of that copying in dispersal and other 40 y of population and environmental data of the colonial long- behavioral feedbacks for the occurrence of tipping points and lived Audouin’s gull (Ichthyaetus audouinii) breeding seasonally at the collapse of populations. The main reason for this is that Punta de la Banya (Ebro Delta, western Mediterranean) (Fig. 2) theoretical population dynamic models including a mathemat- (23). Despite being a species that has evolved a slow life history, ical expression of dispersal, such as coupled logistic equations, the patch held almost 50% of the total world population after mostly consider dispersal as simple passive diffusion and not as only 6 y since colonization, due to its habitat suitability allowing A Downloaded from https://www.pnas.org by 158.109.55.145 on September 7, 2023 from IP address 158.109.55.145. B Fig. 1. (A) Schematic diagram of dispersal represented as a dynamic landscape. Depending on the patch quality (here represented as a perturbed patch with the presence of predators), individuals may undergo different dispersal dynamics (represented by δ) depending on population density (expressed as decreasing values over the x axis), e.g., density-independent (Left) or dispersal by social copying (Right). For social copying, a low population threshold may be found, accelerating such a dispersal. (B) Conceptual surface of dispersal strength as a function of population size and social copying feedback. In this study, we tested whether dispersal made by a few individuals induces others to leave the patch in a behavioral cascade once a tipping point is exceeded. This social copying for dispersal should be mimicked by a negative density-dependent function, where dispersal increases at a decreasing population density. 2 of 8 https://doi.org/10.1073/pnas.2214055120 pnas.org Fig. 2. (A) Population dynamics of gulls breeding at La Banya since colonization (1981 to 2021), expressed as density of reproducing females. The orange areas show the duration of the perturbed regime caused by the arrival of carnivores. In the final years, the population showed a slow recovery after predators were removed for conservation purposes. Numbers over the line are percentages of the total world population breeding here. (B) Spatial distribution of breeding patches, where 90% of the total world population breeds. Red star shows La Banya, red dots are the patches occupied before the perturbed regime, and green dots are the patches colonized since then. (Inset) Population dynamics for the same period 1981 to 2021 of the whole metapopulation, which holds ca. 45 breeding patches. Fluctuations in metapopulation density were greatly influenced by the dynamics of La Banya, but they did not show an abrupt decline because dispersed individuals were redistributed over the colonized patches. for extremely high rates of survival, fertility, and immigration from in explaining the population dynamics observed after the arrival other patches (24, 25). Low interference competition with other of predators. One dispersal form was positively density-dependent, sympatric species of the same ecological guild and the absence of whereas the other was negatively density-dependent mimicking terrestrial predators enhanced these rates (26, 27). The population the behavior of social copying (4). To assess how population den- continued to increase through a rapid transition, ending in a seem- sity influenced dispersal, we used an Elliot sigmoid function, Downloaded from https://www.pnas.org by 158.109.55.145 on September 7, 2023 from IP address 158.109.55.145. ingly dynamic equilibrium state, which occurred when predators which allowed us to test for different shapes of negative densi- invaded the patch (Fig.2A). Predator densities were very low and ty-dependence, including the extreme case of density-independent did not vary over the years, but a perturbation was pressing over dispersal. This function typically increases at decreasing popula- the breeding seasons eroding habitat quality and causing not only tion values (i.e., dispersal would increase with decreasing popula- a slight decrease in adult survival but also a high dispersal to other tion density), but it allows some other forms, such as sigmoidal occupied patches, especially to empty patches for colonization shapes and sudden nonlinear changes (SI Appendix, sections S2 (28–30) (Fig. 2B). Because we accurately monitored the metap- and S4, Figs. S8–S10). For instance, nonlinear shapes with a dras- opulation at a large spatial scale, we quantified dispersal processes tic increase in dispersal beyond a threshold value of population of emigration and immigration among the different patches, density would fit with a tipping point for runaway dispersal by including extinction and colonization events over the years social copying (4, 31). On the contrary, an Elliot function fitting (Fig. 2B). Just after peaking at a maximum density due to a sto- density-independent dispersal would suggest that dispersal chastic increase in food availability during 2005 to 2006, holding occurred only after using private information. Finally, we devel- 73% of the total world population, the local population under- oped an optimization technique to fit the population models to went a sustained decline leading to collapse, where it held only our field data using artificial intelligence algorithms. 3% of that total world population (Fig. 2A). In summary, we tested the hypothesis that the arrival of predators triggered a tip- Results and Discussion ping point by runaway dispersal driven by social copying behavior that caused the population to collapse. The geometric mean of the population growth rate over the whole To test our hypothesis, we first built a population model that study (expressed as xt+1/xt, being xt the size of the population at described the local population dynamics as a single-patch system time t since colonization, t0 = 1,981) was high, and it showed large and considered immigration and dispersal of individuals within fluctuations (mean = 1.106; SD = 0.835), likely due to the influ- the metapopulation, which was not explicitly considered. The ence of large rates of both immigration and dispersal (SI Appendix, model also took into account other ecological processes, such as Fig. S1). The model for the initial phase before the arrival of pred- intraspecific competition for resources and density-independent ators (1981 to 1997) given by Eq. 1 with ρ = λ = 0 (i.e., with no death rates. In previous studies, we accurately estimated death dispersal) clearly showed a logistic behavior with an initial expo- rates by capture-recapture modeling of +30,000 banding birds nential growth, an average growth rate = 0.349 y−1 and carrying and +28,000 resights over 32 y encompassing 69 local patches of capacity = 18,828 birds (least-squares LS = 2,593.05, coefficient the metapopulation (+90% of the total world population) (27, 29). of determination R2 = 0.976, N = 17) (SI Appendix, section S3.2.2) Metapopulation dynamics and most dispersal processes, in terms (Fig. 3). During this initial phase, the population was in a transient of the density of individuals, were governed by the dynamics of state, and it did not reach the predicted equilibrium value even the study patch, because it held a large percentage of the total though the patch already held a large proportion of the total world population (Fig. 2) (29). Most importantly, for modeling world population. This was due to the arrival of predators the dispersal out of the study patch in our model, we included (Fig. 3B and SI Appendix, Fig. S6). We confirmed this after we different forms of dispersal to assess which one played a major role fitted field data with the logistic model from 1981 to 2004 to test PNAS 2023 Vol. 120 No. 11 e2214055120 https://doi.org/10.1073/pnas.2214055120 3 of 8 A B Fig. 3. (A) Population dynamics from patch colonization in 1981 to the arrival of predators in 1997 and the fit of the logistic model displayed with a blue line, which shows an exponential increase in density (SI Appendix, section S3). The fitted parameters suggest that the population was approaching the expected equilibrium value when predators were established at the patch. The estimated carrying capacity is shown with a dashed line, in conjunction with the location of the equilibrium points (stable in orange, unstable in white), predicted by the estimated parameters (SI Appendix, section S3.2). (B) Potential function of the logistic model showing the predicted equilibrium (orange marble), which was not achieved due to the perturbation starting in 1997 (red marble) (SI Appendix, Eq. S3.2 and Fig. S5). whether the model performed well without considering the impact During the first years of the collapse period (initial slow phase, of predators, but the fitting was worse (LS = 5,639.34, R2 = 0.928, Fig. 4B), dispersal increased slightly with a decrease in population N = 24) (SI Appendix, section S4.3 lines 669 to 688, SI Appendix, density. Interestingly, dispersal accelerated once a tipping point Figs. S16 and S17). of population density was crossed (acceleration phase, Fig. 4B). After the onset of perturbation phase in 1997, population fluc- After attaining its maximum rate, dispersal decreased with lower tuations likely increased due to an increase in demographic noise, population densities (slowing down phase, Fig. 4B). These results caused by the negative density-dependent dispersal found in our contrast with the simple linear dispersal process assumed in sim- model. The population also showed a slight decrease between 1998 pler models. Downloaded from https://www.pnas.org by 158.109.55.145 on September 7, 2023 from IP address 158.109.55.145. and 2004 (Fig. 2A and SI Appendix, Figs. S1 and S18). This We used a mathematical technique for fitting a population decrease could not be attributed to a deterioration in environmen- model to real data to explore which type of dispersal could explain tal conditions, because none of the environmental parameters the collapse of a local population at a metapopulation level. Our influencing resource availability and predators (densities and kill- capacity to study how perturbations may cause population col- ing rates) changed significantly (SI Appendix, Fig. S1.2, sections lapses has increased in recent decades, thanks to the greater avail- S1.2, S1.3, and S4.3) (4, 29). Despite the low amount of data for ability of long-term data on both populations and environmental this phase, the parameters estimated in our model using the results variability (33–35). Furthermore, more detailed population mod- obtained for the full collapse phase (see below) indicated that the els are considering the complex nature of population dynamics by dominant dispersal was social copying, whereas positive densi- incorporating simultaneously density-dependent mechanisms and ty-dependent dispersal was negligible (SI Appendix, section S4.3, the different types of stochasticity (mainly environmental and lines 689 to 716 SI Appendix, Fig. S18). demographic) (36–39). Nevertheless, disentangling the role of For the final phase of population collapse (2006 to 2017), the endogenous (e.g., density-dependence) and exogenous (density-in- full model given by Eq. 1 fitted the population decline with high dependent, environmental perturbations) drivers on population accuracy (LS = 2,566.99, R2 = 0.975, N = 12) (SI Appendix lines declines remains a challenge, especially when environmental 646 to 668, SI Appendix, Table S4). The model considered a neg- changes are not accurately monitored (40). Further challenges ative density-dependent dispersal mimicking social copying for include modeling dispersal processes explicitly to understand local dispersal (4) (Fig. 4). Contrarily to what is predicted by simple and metapopulation dynamics and to assess the role of density-­ theory (32), models considering only density-independent or pos- dependence on those processes, which may be particularly crucial itive density-dependent dispersal had a worse fit for the full col- during phases of colonization and declining population (i.e., lapse phase (LS = 3,391.96 and R2 = 0.957 for the former model, immigration and emigration, respectively) (24, 41). Here, we built SI Appendix, section S5). For the positive density-dependent dis- a mechanistic population model to explain the population dynam- persal model, no optimum exists, but rather a sequence of pro- ics from colonization through collapse in a social bird subpopu- gressively better positive density-dependent functions with lation that went from holding 73% of the total world population quadratic errors converging to LS = 3,300. These functions were to only 3% in just 10 y. The results show that the collapse was essentially indistinguishable from the density-independent case caused by a nonlinear negative density-dependent dispersal, sug- (SI Appendix, Fig. S21), showing that the positive density-depend- gesting the behavior of social copying. The generality of copying ent dispersal played no role in reproducing the field data for the in social animals is supported by the bulk of theoretical studies full collapse phase (SI Appendix, section S5.2). Therefore, models and empirical evidence showing the importance of social infor- considering standard dispersal modes did not improve the fit of mation exchange when making decisions (15, 28, 42, 43). the model with dispersal by social copying. The set of parameters Copying also occurs in human behavior when decision-making that better fitted field data for describing the negative density-­ is at play and it is studied under the framework of the theory of dependence dispersal function showed that dispersal varied non- conformity (44). Other than primates, social copying has been linearly with population density over three different phases (Fig. 4). found in a broad range of animals, even in simple organisms, and 4 of 8 https://doi.org/10.1073/pnas.2214055120 pnas.org Fig. 4. (A) The entire population dynamics since colonization to collapse obtained from the field data (solid dots) and the fit of the model without (brown line) and with (blue lines) dispersal. The initial phase (from 1981 to 1997) was fitted using Eq. 1 with ρ = λ = 0 (i.e., with no dispersal). Here, the population showed a logistic growth before the arrival of predators in 1997. Later, the field data showed a first slight decline after this arrival, and a very abrupt decline from 2006 to 2017, when predators were removed. These two declining phases were fitted with the full model given by Eq. 1 including asocial and social dispersal terms. (B) Dispersal of individuals in the full collapse phase varied nonlinearly with population density, showing three phases: an initial slow phase (gray area), an acceleration of dispersal with a social tipping point (blue area), and a slowing down phase found at a low population size (yellow area). it may have important consequences for the individual, e.g., for collapse phase lasted nearly a generation. The time elapsed under Downloaded from https://www.pnas.org by 158.109.55.145 on September 7, 2023 from IP address 158.109.55.145. mate choice and thus for its fitness (45, 46). In an ecological the perturbation press regime was consistent with predictions of context, theoretical models exploring how individuals weigh private the transient theory for fast–slow dynamics: While population and social information to make decisions showed that the tendency growth after colonization of new patches can be extremely fast to follow the majority option is enhanced by deteriorating environ- (4), it only occurs once individuals gather enough information mental conditions (e.g., the presence of a predator) (15). about the suitability of alternative patches for what is termed Beyond the behavioral consequences that social copying may “informed dispersal” (12, 42). This behavioral process of pros- have for individual fitness, very little attention has been devoted pecting to gather information may be slow in empty patches due to the potential effects that copying may have for generating non- to the lack of public information, which is used as a cue for linear density-dependence dispersal and population dynamics assessing the quality of occupied patches (SI Appendix, Fig. S2) (4, 47). In colonial, long-lived species with spatially structured (4, 28, 43, 52). When predators first invaded La Banya, only a populations, individuals dispersing may use inadvertent social few small other patches in the metapopulation were occupied information available at each patch (in the form of conspecific and their local population growth was constrained by competi- density and the fertility of conspecifics) to assess its suitability and tion for breeding sites. Therefore the dispersal to these patches, make decisions about where to settle each reproductive sea- where social information was available, was challenging (Fig. 2B) son (11, 48, 49). Thus, social copying would influence the deci- (28, 29). However, little is still known about how individuals sion-making in the trade-off between staying in the patch and copy others in the process of dispersal and whether they are able dispersing to other patches, which could also affect dispersal at to track a decreasing density at the actual patch (53). the level of the group and the dynamics of local populations and Interindividual behavioral synchronization for collective disper- metapopulations (e.g., periodic dynamics, boom-bust, rescue sal has been seldom described (54). Group dispersal may occur effects, source-sink dynamics). Furthermore, some individuals of at different temporal scales and using several alternative, nonex- the population (e.g., specific ages, sexes) may be more prone to clusive, behavioral mechanisms. First, previous studies show that disperse and generate demographic heterogeneity among patches, individuals use public information by tracking density at some thus increasing stochasticity and spatial variability in extinction patches at a metapopulation scale, and density is a clue for patch risk (30, 37). suitability and therefore decision-making on where to settle for The first nonlinear response of the studied population was breeding (24, 29, 48, 55). Patches with larger population den- related to social processes that occurred after the onset of a per- sities generate higher conspecific attraction due to the advantages turbation. Population fluctuations increased mainly due to an of grouping (e.g., risk dilution, social cohesion, and collective unexpected sharp and short increase in the population leading vigilance). Individuals inhabiting patches subjected to press per- to a maximum density, and this corresponded to what resilience turbations and decreasing densities may disperse to avoid the theory calls demographic compensation (50). Compensation was costs of living in small groups (e.g., incurring in Allee effects) likely caused by high local recruitment and a decrease in inter- (56–58). We also know that birds and mammals may follow specific competition favored by a pulse of higher food availability others while foraging and commuting between feeding patches recorded at the patch (Fig. 2A and SI Appendix, Fig. S2) (27, 50, (59, 60). This process might also operate while individuals are 51). The transient arrival of predators in the beginning of the prospecting different breeding patches for gathering information PNAS 2023 Vol. 120 No. 11 e2214055120 https://doi.org/10.1073/pnas.2214055120 5 of 8 about where to disperse (4, 52). We have some lines of evidence population and metapopulation nonlinear dynamics, including of Audouin’s gulls copying each other when leaving to go on population extinction and resilience, predator-prey dynamics, and patch prospecting trips (52). the management of endangered and harvested populations of social It is clear that population dynamics of social organisms in dete- animals subjected to behavioral feedback loops. riorating environments may be more influenced by transient phe- nomena and density-dependent mechanisms than solitary, Materials and Methods territorial species (7). Empirical data on declining populations of Study System and Species. Demographic population monitoring of nonsocial organisms suggest their dynamics should follow a den- Audouin's gulls (Ichthyaetus audouinii) at Punta de la Banya (Ebro River Delta: sity-independent decrease toward extinction (61, 62). The exist- 40°34010:8900N, 0°39034:2800E) started in 1981 when the patch was colo- ence of social behavioral feedback, from very simple (e.g., social nized, and it has been monitored constantly for the last four decades (1981 to copying) to more complex (e.g., cooperation, hierarchies), may 2021). La Banya is a 2,500-ha sandy peninsula formed by sandy dunes covered increase stochasticity in local population processes like abrupt by halophilous vegetation. The site, once effectively protected in the early 80s, extinctions and the appearance of tipping points for collapses held extremely suitable environmental conditions for breeding ground-nesting (4, 14, 17, 63, 64). Nevertheless, disentangling the processes that, waterbirds (i.e., large availability of food and reproductive habitat free of preda- alone or in interaction, may drive population nonlinear declines, tors). Audouin's gulls are long-lived social birds with a bet-hedging life history. They and the mechanisms behind those processes that may cause the have evolved to cope with ephemeral habitats typical of Mediterranean marshes tipping point to appear, remain a challenge (14, 65, 66). Outside (71, 72). Consequently, Audouin's gulls have nomadic behavior between breeding of mesocosmos experimental studies, interpreting complex local seasons. When patch conditions change or worsen, then individuals are more population dynamics in empirical studies needs accurate moni- prone to disperse to other occupied sites, but in recent years, and after the arrival toring of both metapopulation dynamics and environmental sto- of invasive carnivores in the late 90s, colonization rate of new patches has largely chasticity. The exceptional fitting of our population model during increased. Ground-nesting gulls have not developed evolutionary defenses to cope the phase of collapse is an example of how an environmental with terrestrial predators like carnivores. This is why they select breeding patches stochastic perturbation (i.e., the arrival of predators) triggered a isolated and protected against carnivore invasions. The long-term monitoring at tipping point by runaway dispersal, i.e., a negative density-de- the metapopulation scale allowed reliable knowledge of the system in terms of pendent accelerating dispersal caused by positive behavioral feed- population fluctuations (Fig. 2B and SI Appendix, Fig. S1), interference competition back that drove the population to a state of quasi-extinction. This with species of the same ecological guild, demographic parameters (e.g., survival, type of tipping, which occurs by the progressive loss of the resil- fertility, recruitment, and dispersal), and environmental stochasticity, including ience of the current population state, can potentially be anticipated the invasion of predators, the occurrence of other perturbations, and fluctuations in food availability over the years (SI Appendix, lines 4 to 26) (24, 26, 27, 29). when using the tools of early warning signals (67). Fieldwork Methods. Results regarding the methods used to monitor demo- Downloaded from https://www.pnas.org by 158.109.55.145 on September 7, 2023 from IP address 158.109.55.145. Conclusions graphic parameters and population size were published previously (27), and exhaustive details concerning the fieldwork protocols and sites can be found in The critical significance of our results is that social copying, a simple the previous literature and in SI Appendix (lines 27 to 74). Briefly, population den- behavioral mechanism operating in all social organisms, may gen- sity, which corresponded to the number of breeding females, has been estimated erate nonlinear population collapses by the occurrence of tipping since colonization in 1981 by counting the nests containing eggs just before points in dispersal processes. These dynamics fit well with what has hatching, when biases due to the variance in laying dates were minimal. Errors at been observed in populations of social animals subjected to press counting nests were quantified and were small and constant over the years (<5%). perturbations, such as superabundant species subjected to culling The same procedure was followed to estimate the population density of the main and humans during warfare (4). Under these circumstances, there competitor of Audouin’s gull, the larger yellow-legged gull (Larus michahellis). is an initial reluctance to leave the patch due to a large availability The density of carnivores (mostly foxes, Vulpes vulpes) and the number of killings were estimated every year since their arrival in 1997 thanks to the tracks in the of information, the force of social cohesion, and the evolutionary sand and the intense monitoring of roosting and breeding dens. Other impacts advantages of being philopatric as opposed to being disperser, affecting demographic parameters, such as mortality by fishery bycatch, were also which is riskier in terms of fitness prospects (28, 68, 69). In duly monitored. Finally, we used the statistics of landings on trawling boats in the humans, this is reinforced by a sunk-cost effect, which prevents harbors close to the study site. Landings are highly correlated with the amounts people from abandoning their ways of living, cultures, and beliefs of fish discarded, which represent up to 70% of biomass ingested by the gulls (70). However, once a threshold for deterioration of the quality of and greatly influence their fertility. We showed that the population collapse was the patch has trespassed, there is a tipping point for a social response not caused by a deterioration of environmental conditions but by the stochastic of runaway dispersal corresponding to social copying feedback. fluctuations of carnivore densities (SI Appendix, Fig. S2, lines 75 to 97). Finally, dispersal decreases at low population densities, which is likely due to the unwillingness of some individuals to disperse after Mathematical Model and Data Fitting. To analyze population dynamics from occupying a suitable patch for extended periods. At an ecoevolu- colonization through collapse, we used an ordinary differential equation logistic tionary scale, this corresponds to the individual heterogeneity for population model and density-dependent migration terms (SI Appendix, lines dispersal within populations, with some individuals being extremely 98 to 139). Due to the very low number of predators and their negligible effect philopatric. This phenomenon may generate long queues of qua- on predation (28), we did not consider predator–prey dynamics explicitly as an ecological process causing population decline by direct mortality. The model, si-extinction states observed in some populations of social animals, which describes the population dynamics in a single-patch system, but considers and that would boost immigration by conspecific attraction once immigration and dispersal of individuals, is given by environmental conditions at the patch are recovered (4). Our study emphasizes the value of combining long-term demographic and d ( x(t) ) [ ] environmental data with mathematical modeling to uncover the x(t) = (𝜗 + 𝜔)x(t) 1 − − 𝜀x(t) − 𝜌x(t) + 𝜆D(x(t), 𝜇, 𝜎, 𝛿) , dt K behavioral mechanisms driving nonlinear responses of populations [1] of social organisms under perturbations. Our results suggest a where x(t) denotes population density at time t, 𝜗 is the intrinsic reproduction broader impact of self-organized collective dispersal in complex rate, 𝜔 is the immigration rate from other patches, and 𝜀 is the annual death population dynamics, such as tipping points and transient phe- rate estimated from field data [0.11 y−1, estimated from long-term monitoring nomena. This may have implications for the theoretical study of (1988 to 2015) using capture–recapture modeling (27)]. The model introduces 6 of 8 https://doi.org/10.1073/pnas.2214055120 pnas.org intraspecific competition for resources constrained by a logistic function with correspond to different nonlinear types of negative density-dependent dispersal, carrying capacity K, and importantly, it includes two distinct dispersal terms. The all indicating that dispersal occurred by social copying (SI Appendix, section S4.1, first dispersal term 𝜌 is the linear (exponential) dispersal to other patches, which Figs. S7–S9). Interestingly, we note that for some special conditions (e.g., 𝜎 close corresponds to a positive density-dependent dispersal. The second dispersal term to 0), the Elliot function represents a density-independent dispersal, which would is represented by the function D(x(t), 𝜇, 𝜎, 𝛿 ), which is an Elliot sigmoid function mimic individuals dispersing only using private information (SI Appendix, Remark (see below and SI Appendix, lines 250 to 355). This function introduces a negative S4.2, Figs. S7, S9, and S10). We searched for the model parameters that matched density-dependent term of dispersal (i.e., the fewer individuals at the patch, the our population data using artificial intelligence methods (SI Appendix, lines faster the dispersal, see Fig. 1), which mimics social copying of individuals leaving 520 to 645). Because the number of years for the full collapse phase was larger the patch, with dispersal rate 𝜆. We applied the model in Eq. 1 to fit empirical pop- than for the onset of the perturbation phase, we used the longer period to fit ulation fluctuations separately for three phases in which we divided the dynamics the parameters of the Elliot sigmoid function (SI Appendix, section S4.2). Then of the population (SI Appendix, lines 123 to 139, SI Appendix, Fig. S3). The first we used those parameters to fit the model to the onset of perturbation phase phase goes from the patch colonization to the arrival of predators (initial phase: (SI Appendix, lines 689 to 716, SI Appendix, Fig. S18). 1981 to 1997). In the absence of predators, we expected this initial phase to be characterized by an early exponential growth and a logistic behavior when Data, Materials, and Software Availability. Population time series data have approaching the carrying capacity. Here, we did not consider dispersal triggered been deposited in DIGITAL.CSIC (https://digital.csic.es/handle/10261/211284). by the perturbation because predators were not present in the patch during that ACKNOWLEDGMENTS. The authors thank the many fieldworkers who have phase, and previous studies showed that population dynamics were dominated contributed to data collection over the years, particularly Albert Bertolero and by high immigration rates from outside (24). Thus, we used Eq. 1 and removed the technical staff of the Ebro Delta Natural Park. We also want to thank Clara the two dispersal terms (i.e., 𝜌 = 𝜆 = 0). For fitting the population model to this Alsedà for her help in designing figures. Data collection and research have been phase, we used a standard trust region method and the Levenberg–Marquardt funded by the Spanish Ministry of Science, the Spanish State Research Agency algorithm (SI Appendix, lines 232 to 249, SI Appendix, Fig. S6). (AEI), The Fulbright Commission, and European Regional Development Fund Second, we analyzed the two phases in which predators were present at the (FEDER) [last grants: CGL2017-85210-P, PID2021-122893NB-C21, and Salvador patch. One phase (onset of perturbation phase) encompassed from the estab- de Madariaga-Fulbright (all to D.O.)]. This work is also supported through the lishment of predators in 1997 to 2004, before population density suddenly Severo Ochoa and María de Maeztu Program for Centers and Units of Excellence increased due to a stochastic pulse of food availability in the years 2005 to 2006 in R&D (CEX2020-001084-M) funded by (MCIN/AEI/10.13039/501100011033). (SI Appendix, Fig. S3). The last phase encompassed the maximum density attained We thank Centres de Recerca de Catalunya (CERCA) Programme (Generalitat to the full collapse phase between 2006 and 2017 (Fig. 2 and SI Appendix, de Catalunya) for institutional support. J.S. has been funded by the Ramón y Fig. S3, sections S4.2 and S4.3). For these two phases, we used the full Eq. 1, Cajal grant (RYC-2017-22243 funded by MCIN/AEI/10.13039/501100011033 which included the two distinct dispersal terms (SI Appendix, Eq. S4.1). We "European Social Fund (FSE) invests in your future"). Ll.A. has been also funded by tested whether there was a negative density-dependent term of dispersal by the Spanish Ministry of Science and Innovation (grant PID2020-118281GB-C31 social copying using different values of the parameters contained in function and MCIN/AEI/10.13039/501100011033). Two anonymous reviewers and an D(x(t), 𝜇, 𝜎, 𝛿). The parameters of this function drove the sharpness of the jump Downloaded from https://www.pnas.org by 158.109.55.145 on September 7, 2023 from IP address 158.109.55.145. Associate Editor provided very helpful comments for improving the manuscript. and its smoothness in the relationship between density and dispersal: 𝜇 reflects the tendency of dispersal for small population sizes, 𝜎 defines the sharpness and smoothness of the dispersal function, and 𝛿 models the transition between small Author affiliations: aTheoretical and Computational Ecology Laboratory, Centre d’Estudis and large population sizes (SI Appendix, lines 266 to 305). The Elliot sigmoid func- Avançats de Blanes, Consejo Superior de Investigaciones Científicas, 17300 Girona, tion has two features that make it suitable for testing what type of dispersal would Spain; bDepartment of Environmental Science and Policy, University of California, Davis, CA 95616; cDepartament de Matemàtiques, Universitat Autònoma de Barcelona, 08193 explain the population collapse. First, the function typically increases at decreas- Bellaterra, Spain; dCentre de Recerca Matemàtica, 08193 Bellaterra, Spain; and eSanta Fe ing population values, and second, different combinations of those parameters Institute, Santa Fe, NM 87501 1. A. J. Underwood, The analysis of stress in natural populations. Biol. J. Linn. Soc. 37, 51–78 (1989). 18. V. Dakos, A. Hastings, Editorial: Special issue on regime shifts and tipping points in ecology. Theor. 2. M. S. Boyce, C. V. Haridas, C. T. Lee, N. S. D. W. Group, Demography in an increasingly variable world. Ecol. 6, 253–254 (2013). Trends Ecol. Evol. 21, 141–148 (2006). 19. S. Arganda, A. Pérez-Escudero, G. G. de Polavieja, A common rule for decision making in animal 3. M. G. Burgess et al., Range contraction enables harvesting to extinction. Proc. Natl. Acad. Sci. U.S.A. collectives across species. Proc. Natl. Acad. Sci. U.S.A. 109, 20508–20513 (2012). 114, 3945–3950 (2017). 20. A. Hastings, Complex interactions between dispersal and dynamics: Lessons from coupled logistic 4. D. Oro, Perturbation, Behavioural Feedbacks, and Population Dynamics in Social Animals (Oxford equations. Ecology 74, 1362–1372 (1993). University Press, UK, 2020). 21. J. M. Travis, D. J. Murrell, C. Dytham, The evolution of density-dependent dispersal. Proc. R. Soc. 5. F. B. Hanson, H. C. Tuckwell, Logistic growth with random density independent disasters. Theor. Lond. SeriesB 266, 1837–1842 (1999). Popul. Biol. 19, 1–18 (1981). 22. Z. Ratajczak et al., Abrupt change in ecological systems: Inference and diagnosis. Trends Ecol. Evol. 6. B.-E. Sæther, S. Engen, Pattern of variation in avian population growth rates. Philos. Trans. R. Soc. 33, 513–526 (2018). Lond. B Biol. Sci. 357, 1185–1195, (2002). 23. D. Oro, LTER bird census data Ebro Delta. (2018) Available at https://digital.csic.es/handle/​ 7. A. Hastings et al., Transient phenomena in ecology. Science 361, eaat6412 (2018). 10261/211284. 8. R. Lande, Risks of population extinction from demographic and environmental stochasticity and 24. D. Oro, G. D. Ruxton, The formation and growth of seabird colonies: Audouin’s gull as a case study. J. random catastrophes. Am. Nat. 142, 911–927 (1993). Anim. Ecol. 70, 527–535 (2001). 9. W. J. Sutherland, From Individual Behaviour to Population Ecology (Oxford University Press, 1996). 25. G. Tavecchia, R. Pradel, M. Genovart, D. Oro, Density-dependent parameters and demographic 10. M. I. O’Connor et al., Principles of ecology revisited: Integrating information and ecological theories equilibrium in open populations. Oikos 116, 1481–1492 (2007). for a more unified science. Front. Ecol. Evol. 7, 219 (2019). 26. D. Oro et al., Interference competition in a threatened seabird community: A paradox for a successful 11. E. Danchin, T. Boulinier, M. Massot, Conspecific reproductive success and breeding habitat selection: conservation. Biol. Conserv. 142, 1830–1835 (2009). Implications for the study of coloniality. Ecology 79, 2415–2428 (1998). 27. M. Genovart, D. Oro, S. Tenan, Immature survival, fertility, and density dependence drive global 12. J. Clobert, J.-F.L. Galliard, J. Cote, S. Meylan, M. Massot, Informed dispersal, heterogeneity in animal population dynamics in a long-lived species. Ecology 99, 2823–2832 (2018). dispersal syndromes and the dynamics of spatially structured populations. Ecol. Lett. 12, 197–209 28. A. Payo-Payo et al., Colonisation in social species: The importance of breeding experience for (2009). dispersal in overcoming information barriers. Sci. Rep. 7, 42866 (2017). 13. S. R. X. Dall, L.-A. Giraldeau, O. Olsson, J. M. McNamara, D. W. Stephens, Information and its use by 29. A. Fernández-Chacón et al., When to stay, when to disperse and where to go: Survival and dispersal animals in evolutionary ecology. Trends Ecol. Evol. 20, 187–193 (2005). patterns in a spatially structured seabird population. Ecography 36, 1117–1126 (2013). 14. M. A. Gil, M. L. Baskett, S. B. Munch, A. M. Hein, Fast behavioral feedbacks make ecosystems 30. A. Payo-Payo et al., Predator arrival elicits differential dispersal, change in age structure and sensitive to pace and not just magnitude of anthropogenic environmental change. Proc. Natl. Acad. reproductive performance in a prey population. Sci. Rep. 8, 1971 (2018). Sci. U.S.A. 117, 25580–25589 (2020). 31. E. H. van Nes et al., What do you mean, ‘Tipping point’? Trends Ecol. Evol. 31, 902–904 (2016). 15. A. Pérez-Escudero, G. G. de Polavieja, Adversity magnifies the importance of social information in 32. O. Ovaskainen, B. Meerson, Stochastic models of population extinction. Trends Ecol. Evol. 25, decision-making. J. R. Soc. Interface 14, 20170748 (2017). 643–652 (2010). 16. L. Avilés, P. Abbot, A. D. Cutter, Population ecology, nonlinear dynamics, and social evolution. I. 33. N. Reid et al., Regime shift tipping point in hare population collapse associated with climatic and Associations among Nonrelatives. Am. Nat. 159, 115–127 (2002). agricultural change during the very early 20th century. Glob. Change Biol. 27, 3732–3740 (2021). 17. L. Avilés, Cooperation and non-linear dynamics: An ecological perspective on the evolution of 34. J. A. Hutchings, J. D. Reynolds, Marine fish population collapses: Consequences for recovery and sociality. Evol. Ecol. Res. 1, 459–477 (1999). extinction risk. BioScience 54, 297–309 (2004). PNAS 2023 Vol. 120 No. 11 e2214055120 https://doi.org/10.1073/pnas.2214055120 7 of 8 35. J. F. McLaughlin, J. J. Hellmann, C. L. Boggs, P. R. Ehrlich, The route to extinction: Population 54. J. Cote et al., Behavioural synchronization of large-scale animal movements – disperse alone, but dynamics of a threatened butterfly. Oecologia 132, 538–548 (2002). migrate together? Biol. Rev. 92, 1275–1296 (2017). 36. P. Turchin, A. D. Taylor, Complex dynamics in ecological time series. Ecology 73, 289–305 (1992). 55. D. Oro et al., Lessons from a failed translocation program with a seabird species: Determinants of 37. B. A. Melbourne, A. Hastings, Extinction risk depends strongly on factors contributing to success and conservation value. Biol. Conserv. 144, 851–858 (2011). stochasticity. Nature 454, 100–103 (2008). 56. R. Etienne, B. Wertheim, L. Hemerik, P. Schneider, J. Powell, The interaction between dispersal, the 38. W. W. Murdoch, Population regulation in theory and practice. Ecology 75, 271–287 (1994). Allee effect and scramble competition affects population dynamics. Ecol. Model. 148, 153–168 39. Z. Zeng, R. M. Nowierski, M. L. Taper, B. Dennis, W. P. Kemp, Complex population dynamics in the (2002). real world: Modeling the influence of the time-varying parameters and time lags. Ecology 79, 57. D. Serrano, D. Oro, E. Ursúa, J. L. Tella, Colony size selection determines adult survival and dispersal 2193–2209 (1998). preferences: Allee effects in a colonial bird. Am. Nat. 166, E22–E31 (2005). 40. T. J. Hefley, M. B. Hooten, J. M. Drake, R. E. Russell, D. P. Walsh, When can the cause of a population 58. A. Martínez-Abraín, P. S. Tomillo, J. Mouriño, S. Tenan, D. Oro, Delayed predator-prey collapses: decline be determined? Ecol. Lett. 19, 1353–1362 (2016). The case of black-legged kittiwakes and Iberian sardines. Mar. Ecol. Prog. Ser. 631, 201–207 41. S. D. Bona et al., Spatio-temporal dynamics of density-dependent dispersal during a population (2019). colonisation. Ecol. Lett. 22, 634–644 (2019). 59. S. S. Urmy, Visual trail following in colonial seabirds: Theory, simulation, and remote observations. 42. M. M. Delgado, K. A. Bartoń, D. Bonte, J. M. J. Travis, Prospecting and dispersal: Their eco- Ecol. Monogr. 91, e01429 (2021). evolutionary dynamics and implications for population patterns. Proc. R. Soc. B Biol. Sci. 281, 60. T. Kashetsky, T. Avgar, R. Dukas, The cognitive ecology of animal movement: Evidence from birds and 20132851 (2014). mammals. Front. Ecol. Evol. 9 (2021). 43. M. A. Gil, A. M. Hein, O. Spiegel, M. L. Baskett, A. Sih, Social information links individual behavior to 61. B.-E. Sæther et al., Time to extinction of bird populations. Ecology 86, 693–700 (2005). population and community dynamics. Trends Ecol. Evol. 33, 535–548 (2018). 62. B. C. Lister, A. Garcia, Climate-driven declines in arthropod abundance restructure a rainforest food 44. B. D. Bernheim, A theory of conformity. J. Polit. Econ. 102, 841–877 (1994). web. Proc. Natl. Acad. Sci. U.S.A. 115, E10397–E10406 (2018). 45. E. Danchin et al., Cultural flies: Conformist social learning in fruitflies predicts long-lasting mate- 63. B. C. Lister, Information, behaviour and population dynamics. Oikos 123, 1431–1438 (2014). choice traditions. Science 362, 1025–1030 (2018). 64. A. R. Bhowmick, B. Saha, J. Chattopadhyay, S. Ray, S. Bhattacharya, Cooperation in species: Interplay 46. A. D. Davies, Z. Lewis, L. R. Dougherty, A meta-analysis of factors influencing the strength of mate- of population regulation and extinction through global population dynamics database. Ecol. Model. choice copying in animals. Behav. Ecol. 31, 1279–1290 (2020). 312, 150–165 (2015). 47. M. Sapage, S. A. M. Varela, H. Kokko, Social learning by mate-choice copying increases dispersal and 65. C. N. Anderson et al., Why fishing magnifies fluctuations in fish abundance. Nature 452, 835 (2008). reduces local adaptation. Funct. Ecol. 35, 705–716 (2021). 66. V. Dakos et al., Ecosystem tipping points in an evolving world. Nat. Ecol. Evol. 3, 355 (2019). 48. D. Parejo, D. Oro, E. Danchin, Testing habitat copying in breeding habitat selection in a species 67. L. Dai, D. Vorselen, K. S. Korolev, J. Gore, Generic indicators for loss of resilience before a tipping adapted to variable environments. Ibis 148, 146–154 (2006). point leading to population collapse. Science 336, 1175–1177 (2012). 49. D. Serrano, J. L. Tella, Dispersal within a spatially structured population of lesser kestrels: The role of 68. M. Genovart et al., Decrease in social cohesion in a colonial seabird under a perturbation regime. spatial isolation and conspecific attraction. J. Anim. Ecol. 72, 400–410 (2003). Sci. Rep. 10, 18720 (2020). 50. P. Capdevila, I. Stott, M. Beger, R. Salguero-Gómez, Towards a comparative framework of 69. M. del M. Delgado, I. I. Ratikainen, H. Kokko, Inertia: The discrepancy between individual and demographic resilience. Trends Ecol. Evol. 35, 776–786 (2020). common good in dispersal and prospecting behaviour. Biol. Rev. 86, 717–732 (2011). 51. P. A. Abrams, When does greater mortality increase population size? The long history and diverse 70. M. Scheffer, Anticipating societal collapse; Hints from the Stone Age. Proc. Natl. Acad. Sci. U.S.A. mechanisms underlying the hydra effect. Ecol. Lett. 12, 462–474 (2009). 113, 10733–10735 (2016). 52. D. Oro, J. Bécares, F. Bartumeus, J. M. Arcos, High frequency of prospecting for informed dispersal 71. A. Martínez-Abraín, D. Oro, M. G. Forero, D. Conesa, Modeling temporal and spatial colony-site and colonisation in a social species at large spatial scale. Oecologia 197, 395–409 (2021). dynamics in a long-lived seabird. Popul. Ecol. 45, 133–139 (2003). 53. A. M. Berdahl et al., Collective animal navigation and migratory culture: From theoretical models to 72. D. Oro, G. Tavecchia, M. Genovart, Comparing demographic parameters for philopatric and immigrant empirical evidence. Philos. Trans. R. Soc. B Biol. Sci. 373, 20170009 (2018). individuals in a long-lived bird adapted to unstable habitats. Oecologia 165, 935–945 (2011). Downloaded from https://www.pnas.org by 158.109.55.145 on September 7, 2023 from IP address 158.109.55.145. 8 of 8 https://doi.org/10.1073/pnas.2214055120 pnas.org