Introduction

Modern genomic methods enable estimation of a lineage’s long-term effective population sizes back to its origins (Li & Durbin, 2011; Mazet et al., 2015; Nadachowska-Brzyska et al., 2015). This ability allows unprecedented opportunities to study the populational effects of past evolutionary and climatic phenomena on organismal lineages, and this new area of paleodemographics is likely to grow rapidly (Nadachowska-Brzyska et al. 2021, Germain et al. 2023a,b). We use this approach in a novel way to test hypotheses about how adoption of a major life-history trait affected lineages’ populations relative to those without the trait.

Seasonal migration is a common life-history strategy among the world’s birds, but its effects on long-term population sizes relative to lineages that don’t have this trait are largely unknown (Newton, 2008; Rappole, 2013). Demographic differences between seasonally migratory and resident lineages have been studied for decades, but, thus far, ecological, behavioral, and genetic approaches have necessarily been focused on recent and microevolutionary phenomena (e.g., Greenberg, 1980; Pérez-Tris et al., 2004; Sandercock & Jaramillo, 2002). Here we look much deeper, contrasting changes in effective population sizes through millions of years between closely related migratory and resident lineages.

Thrushes in the genus Catharus and this group’s sister Hylocichla mustelina (Aves: Turdidae) are a model songbird system for studying the evolutionary effects of seasonal migration. This group has both migratory and resident lineages, and it has been important in the study of seasonal migration, divergence, and speciation (e.g., Delmore et al,. 2015, 2016; Outlaw et al., 2003; Ruegg et al., 2014; Termignoni-Garcia et al., 2022; Voelker et al., 2013; Winker, 2000, 2010; Winker and Pruett, 2006). These animals are relatively small, omnivorous or insectivorous, forest-related birds with lineages that are long-distance seasonal migrants in North America and others that are resident in tropical North and South America (Collar, 2005). The breeding grounds of the migratory lineages span Siberia, Alaska, northern Canada, and the United States, and their wintering ranges include Mexico, the Caribbean, Central America, and tropical South America. In contrast, the resident lineages occupy Mexico, Central America, and northwestern South America (Collar, 2005; Fig. S1).

We made three predictions about effective population sizes (N_e) in seasonally migratory versus resident lineages of these thrushes: 1) N_ewould on average be higher among seasonal migrants; 2) variation in population size would also be higher; and 3) early population growth would be higher as seasonal migration opens novel ecological and geographic space, achieving the first prediction.

These predictions stem from the following reasoning. We can think of North America roughly as an inverted triangle, with substantially more geographic space in the north than in the south. Among these thrushes, the migratory lineages breed in that northern, on average, larger space. The resident lineages instead currently occupy, on average, smaller geographic ranges, even in South America, where they are largely montane (Fig. S1). Thus, larger breeding ranges among migrants, on average, should yield higher population sizes. But climatic variability is also higher with increased latitude, often temporarily rendering vast swaths of northern North America unoccupiable to many long-distance migrant species (Pielou, 1991). Thus, higher-latitude breeding occupancy should cause larger fluctuations in population size. Finally, seasonal migration is likely to cause enhanced ecological (and in this case also geographic) release (i.e., in the occupancy of new niche space) in relation to resident relatives. Thus, as migratory lineages’ mobility enabled them to take advantage of seasonal resource blooms in ecosystems at higher latitudes and expand their breeding ranges, their populations should have grown proportionally more and at a faster rate than sedentary lineages. Here we test these predictions by reconstructing effective population sizes of these thrush lineages through time using whole-genome data and the pairwise sequential Markovian coalescent (PSMC; Li & Durbin, 2011).

Results

The focal variables of our hypotheses all lacked phylogenetic signal. Mean effective population sizes and their variation were both higher among migratory lineages, as predicted (Table 1; U = 6 and 9, p < 0.01 and 0.05, respectively; one-tailed Mann-Whitney U-test). The degree of early population growth was also proportionally higher among migrants, as predicted (Table 1; U = 9, p < 0.05; one-tailed Mann-Whitney U-test), but the rate of this early growth was not different between the two groups (Table 1; U = 21.5, p > 0.05).

Table 1.

To examine why the degree of initial migrant population growth was higher but the rate of growth was not, we considered that this might be because that growth extended over a longer period of time. To our surprise, we found that among migrants these growth periods were remarkably long, ranging from 1.67-4.26 My and averaging 2.97 My (Table 1, Fig. 1). It appeared that migratory lineages did have a longer initial period of growth (migrants mean deltaT = 2.97 My vs. resident mean deltaT = 1.49 My; Table 1). Unlike our other variables, initial growth period (deltaT) showed significant phylogenetic signal (p = 0.02; Fig. S3), requiring PGLS analyses, which supported this post-hoc hypothesis, both with actual branch lengths and with branch lengths set to 1.0 (Table S2). We also used this approach to examine the relationship between deltaT and breeding latitude (also significant) and to confirm the relationship between migration and breeding latitude (Table S2).

Figure 1.

In addition to testing our hypotheses, our results reveal other important attributes of these songbirds. First, they appear to show few of the general population declines at the beginning of the last glacial period (LGP, ∼115 Kya) that Nadachowska-Brzyska et al. (2016) found among a global sampling of birds (Figs. 2, S4). Secondly, and more importantly, there is a general lack of temporal synchrony in historic population size changes (Figs. 2, S4), despite considerable geographic proximity and overlap of ranges (especially during the nonbreeding season). Finally, a pattern widely shared in the group is that most of the recent populations are below each lineage’s historic peak, and this trend predates the arrival of humans in the New World (Figs. 2, S4). A striking exception is C. ustulatus swainsoni, which has had remarkable long-term growth (Fig. S4).

Figure 2.

Discussion

To our knowledge, this is the first study to contrast long-term effective population size changes among closely related species to test a priori hypotheses about the population effects of a major life-history trait; here we focused on seasonal migration. This perspective offers several new insights into migrant demographics and evolutionary ecology. Migratory thrush lineages showed higher effective population sizes, greater population size variation, and proportionally larger initial population growth than resident relatives. A migratory life-history strategy thus imbues these lineages with population size characteristics that are fundamentally different, on average, from those of resident relatives.

Although migratory thrush lineages showed higher proportional initial growth than resident relatives, this was not achieved through rapid adaptation to new environments as we predicted. Population growth in migrants and residents was positive early in the lineages’ histories (deltaT) in all cases but one (C. aurantiirostris), and it occurred over fairly long periods of time (Table 1, Fig. 1). That this process initially (on average) involves millions of years of population growth in this group, apparently extended by a migratory life-history strategy that provided access to higher latitudes (Tables 1, S2), is an insight from our study that bears examination from a broader taxonomic perspective. While actual, census population sizes would likely have shown considerable variation across this time, especially after the glacial cycles of the Pleistocene began ∼2.6 Mya, the long-term effective population sizes of migrants tended to show steady increases, achieving peaks on average almost 3 My after the lineages’ PSMC inception (Figs. 1, S4; Table 1).

With deltaT being both long and showing phylogenetic signal (Table 1, Fig. S3), one can envision in this group a sort of ur-thrush trajectory in which unique, lineage-specific processes occur as each evolves, influenced by ancestry, geography, ecological factors, and life-history strategy. These long periods of early growth are concordant with the initial expansion or adaptation phases of taxon cycles, a conceptual framework of lineage expansion and contraction over time through interactions between biogeography and evolutionary ecology (Pepke et al., 2019; Ricklefs & Bermingham, 2002; Ricklefs & Cox, 1972; Wilson, 1961).

Another noteworthy pattern widely shared in this group is that most recent populations are below each lineage’s historic peak (Figs. 2, S4). This, too, is reminiscent of taxon cycles (smaller populations in later stages), and this trend both predates the arrival of humans in the New World (Fig. S4) and reflects a broader pattern of historic avian declines (Germain et al. 2023b). Yet another pattern is a general lack of temporal synchrony in historic population size changes, despite considerable geographic proximity and overlap of ranges (Figs. 2, S4). In their review of taxon cycles, Ricklefs & Bermingham (2002) inferred from such a lack of similarity among lineages that each is in a unique population-environment relationship, responding idiosyncratically in a coevolutionary relationship with factors such as parasites, predators, pathogens, and competitors.

Population structure

PSMC estimates of population size can be affected by population structure: a change in structure can produce the signal of a change in size (Mazet et al., 2015, 2016; Wakeley, 1999). Some of the population size changes we observed might therefore reflect changes only in levels of population connectivity. There are well-understood differences between seasonal migrants and residents in population structure. Migrants generally have lower levels of geographic partitioning of both phenotypic (subspecies) and genetic variation (Belliure et al., 2000; Delmore et al., 2020; Mayr, 1963; Montgomery, 1896). The resident lineages we studied are more prone to this structure than the migrant lineages, with migrants having an average of 2.5 subspecies per lineage and residents having an average of 5.9 (Collar, 2005). Migratory lineages therefore come closer to meeting the assumption of PSMC analyses that populations are panmictic. Thus, in addition to our geographic scenario differing between the migrant and resident lineages we studied, we also have a pervasive effect of population structure affecting our results, likely increasingly complementing geographic expansion among seasonal migrants as lineage ages approach the recent.

Accounting for population structure in these models to distinguish between the two phenomena is a difficult computational problem, and the correlation is not simple, and it might not be tractable with data from a single individual (Mazet et al., 2016; Nadachowska-Brzyska et al., 2016, 2022). However, each individual has a genome amalgamated from the entire lineage metapopulation, including individuals and populations that have not been sampled (Mazet et al., 2015). When not focusing on relatively recent or rapid changes (e.g., in human or other bottlenecked systems, where this has been most studied), these issues are probably not as important (Mazet et al., 2016). For example, sister lineages usually show identical effective population size attributes at their origins, reflecting their shared histories (Delmore et al., 2020; Li & Durbin, 2011; Nadachowska-Brzyska et al., 2016). Inasmuch as many of our most interesting results are in deep time, where we expect current individuals of a lineage to have the most thorough mixing or amalgamation of lineage-specific metapolulation history, the effects of population structure are likely to become less important and the results reflect (as we have necessarily interpreted) more of a lineage-wide phenomenon than a variably structured, demically biased one. However, it is clear that the latter bias becomes stronger as lineages progress toward the present, and for more recent times population structure will be an important aspect of comparing effective population size estimates and their changes between migrant and resident lineages (in this respect we were conservative in ignoring the most recent 50 Kyr).

At present, there is no way to fully disentangle the effects of population structure and geographic space on our results, but we make two observations leading us to infer that general differences in population structure are not the main factor driving our results. First, structure will have smaller effects on the PSMC signal in the deep past than in more recent history. Second, climatic changes are known to cause coordinated shifts among codistributed taxa, e.g., in promoting the partitioning of genetic variation among glacial refugia (Hewitt, 2000). A lack of coordinated shifts in our PSMC results suggests that climate-driven shifts in population structure are not a major feature (Fig. S4).

Geography and evolutionary ecology

Early population growth among migrant lineages is likely affected by some degree of ecological release as these birds engage with the novel environments their increased movements expose them to and occupy new niche space. However, in this system we cannot readily decouple the expanded geographic and ecological opportunities that seasonal migration provides. The relationship between avian range size and migration is complex, but correlation with latitude is notably strong in North America (Pegan & Winger, 2020). Two aspects of our results suggest that geography operating alone is unlikely to explain the larger migrant population sizes observed (Table 1, Fig. S4). The first is that early population growth is positive in all but one of the lineages we studied; despite a probable lack of marked geographic expansion among residents, they, too, experienced long-term growth. The second is that these initial growth periods are quite long, extending over much longer periods than we should expect populations to respond to opportunities for geographic expansion.

PSMC analyses excel at long time frames, but initial growth in residents and migrants suggests a strong role for ecological factors. Global climate reconstruction data do not suggest that there were similar long, steady trends during these time periods to match these thrush population size changes (Lisiecki & Raymo, 2005; Fig. S4). Perhaps geographic space does make a difference, on average, between the two life-history strategies. That these growth periods are long also suggest the importance of ecological factors, though. For example, the boreal and subarctic woodland breeding habitats of most C. ustulatus swainsoni and C. minimus expanded dramatically across North America as glaciers receded since the last glacial maximum (LGM), a period of < 20 Kyr (Pielou, 1991). A spatial response of such rapidity, less than 1% of the time of the initial historic growth phases of these two lineages (> 3 Myr; Table 1), leads us to infer that breeding range size alone is not a satisfactory explanation for our results. This reflects census size, and not effective population size (N_e), but it implies that a migratory Catharus lineage’s matching (or mapping) of N_e to average long-term available geographic space is not likely to take millions of years. Or, to put it another way, their comparatively recent rapid expansion into new space suggests a sort of preadaptation that is not apparent early in the lineages’ histories. Together, these aspects of our results suggest that while geographic factors are involved, ecological factors are also important. Teasing these factors apart might also explain the recent finding of Germain et al. (2023a) that migrants in general have had larger populations than residents during the last 0.67 my.

Identifying key ecological mechanisms affecting population sizes among migrants is difficult, and they will likely vary among lineages. Biotic interactions in general show greater effects on range limits at the warmer than the cooler edges (Paquette & Hargreaves, 2021). Long-distance migrants especially, in their extensive latitudinal movements, would receive a range expansion advantage in this respect regardless of continental landmass shape. Many ecological factors affect the population sizes of migrants versus residents, but the relationships between seasonal migration and predators, parasites, diseases, and competitors are extremely complex and key aspects remain poorly understood (Altizer et al., 2011; Newton, 2008).

Seasonal migration can be integrally related to predator avoidance, but among birds its effects are not understood in the full, circannual context relative to resident relatives and they are not generalizable (Brönmark et al., 2008; McKinnon et al., 2010; Newton, 2008; Rappole, 2013). Exposure to parasites and pathogens varies both spatially and seasonally. Although breeding range diminishment of these exposures can be important in some migrants (probably not among all lineages), full annual exposure to both disease and parasite diversity is likely higher, resulting in at least some lineages showing larger immune defense organs (Bennett & Fallis, 1960; Figuerola & Green, 2000; Jenkins et al., 2011; Møller & Erritzoe, 1998; Piersma, 1997; Ricklefs et al., 2017). Exposure to higher circannual parasite and pathogen diversity might not be correlated with susceptibility, however; migration creates a geographic break that can lower transmission, causing migrants overall to be less vulnerable to population declines (Hall et al., 2014).

Presently C. minimus shows an increased breeding and migratory prevalence and diversity of avian malaria, a blood parasite, but aspects of this relationship probably changed over the time scales involved (Pulgarín-R et al., 2019; Ricklefs et al., 2014). Avian migrants are more likely to be affected by blood parasites in ecological rather than evolutionary time, although the latter are also operating (Hellgren et al., 2007; Ricklefs et al., 2014). At the deeper times of avian families, seasonal migration is an important factor affecting blood parasite-host coevolution by breaking down the congruence of phylogenetic relationships between hosts and parasites, and there are also important differences among biting-fly vectors and their transmission of avian blood parasite hosts (Jenkins et al., 2011; Ricklefs et al., 2014). These issues suggest that reconstructing associations and their effects in deep time will be difficult.

Interspecific competition has long been studied between migrants and residents, but its potential role in affecting relative population sizes between the two groups is largely unknown (Rappole, 2013). As close relatives, these thrushes are likely to be each others’ closest competitors, so it is conceivable that the increased relative abundance of migrants might depress realized population sizes among residents. In this group, geographic ranges overlap broadly during the nonbreeding season, but the tropical residents tend to occupy higher elevations than nonbreeding migrants (Collar, 2005). Both interspecific competition and inland, montane, and mature-forest retreats among older taxa are major factors in taxon cycles (Ricklefs & Bermingham, 2002; Ricklefs & Cox, 1972; Wilson, 1961) and might also be operating in this system.

Conclusions

High-quality, genome-scale data offer important new insights into the effects of the major life-history trait of seasonal migration on long-term effective population sizes. Our results indicate larger and more variable migrant populations relative to residents, and these populations grow in a way suggesting an important role not only for seasonal migration and geography, but also for ancestry and evolutionary ecology. The differences we found are between group averages, and variation among lineages includes residents that are migrant-like (C. occidentalis; Table 1, Figs. 1, 2, S4) and vice-versa (e.g., C. u. ustulatus; Table 1, Fig. S4). Understanding why residents like C. occidentalis show migrant-like historic population characteristics will help us understand the mechanisms affecting these lineages. Such reasons could be as simple as realized range expansion (i.e., geographic and ecological opportunity), and this hypothesis appears to fit the deeply split eastern-western sister lineages within species in which smaller western ranges have smaller effective population sizes (C. guttatus and C. ustulatus; Fig. S4). Or the reasons might be as complex as a lineage dropping out of migration to become sedentary (C. occidentalis is sister to the migrant C. guttatus; Fig. S3).

What other patterns will emerge across geographic, temporal, and phylogenetic space? Further sampling among other closely related groups containing mixed lineages of migrants and residents are needed, as are complementary approaches to reconstructing demographic history, using more individuals per lineage and other population genomic methods applicable to that sampling (e.g., Beichman et al., 2017; Delmore et al., 2020; Marchi et al., 2021). Such studies will also help determine the relative contributions of seasonal migration and its range and ecological expansion components in the context of what appears to be a normal and phylogenetically influenced tendency for early lineage growth (Figs. 1, S4, Table 1).

Materials and methods

The thrush lineages in our study are from the genus Catharus and its sister, the single species in the genus Hylocichla. These comprise Nearctic-Neotropic seasonal migrants and Neotropical residents. Among the latter, some lineages (C. aurantiirostris, occidentalis, frantzii, and mexicanus) have some populations that exhibit some short-distance migration—likely partial migration or hard-weather movements, either in elevation or with extreme northern populations moving south (Collar, 2005). These lineages thus comprise two groups, which for simplicity we will refer to as ‘seasonal migrants’ and ‘residents’: seasonal migrants with extensive latitudinal movements, and residents that are either wholly sedentary or with some individuals moving rather limited distances within the Neotropics. The eight migrant lineages are: Hylocichla mustelina, Catharus fuscescens, C. guttatus E (eastern lineage), C. guttatus W (western lineage), C. minimus, C. ustulatus swainsoni, C. u. ustulatus, and C. bicknelli. (Note that four of these lineages represent two deeply divergent lineages each within what are currently considered two biological species, C. guttatus and C. ustulatus; Topp et al., 2013). The six resident lineages are: Catharus aurantiirostris, C. fuscater, C. frantzii, C. gracilirostris, C. mexicanus, and C. occidentalis. Over the full evolutionary history of this Hylocichla-Catharus clade, changes among lineages occurred in the trait of long-distance migration (Outlaw et al., 2003; Voelker et al., 2013; Winker and Pruett, 2006). Our study does not assume that this trait was constant within each lineage after the speciation events when these switches occurred. Instead, we consider that for a major trait like long-distance migration, losses and gains at the full lineage level are likely to be infrequent, and that current evidence of trait distribution on a phylogeny provides a powerful comparative framework in which to elucidate the effects of that trait on other lineage attributes.

DNA was extracted from high-quality voucher specimens or blood samples for 14 thrush lineages (Table S1). Whole-genome sequencing libraries were constructed using Nextera (Illumina) DNA Flex Library Prep kits, and libraries were sequenced on a NovaSeq S4 (Illumina, San Diego, CA) using paired-end 150 bp reads. Reads were aligned to the C. ustulatus genome (inland subspecies, Reference bCatUst1, NCBI PRJNA613294) using bwa (mem algorithm with default settings, Li and Durbin 2009). Resulting .sam files were converted to .bam format with samtools (Li et al., 2009); picardtools was used to clean, sort and add read groups to .bam files (http://broadinstitute.github.io/picard/). Between 97 and 98% of reads mapped to the reference for an average read depth of 23 bp (range 15.9-29.8).

We analyzed whole-genomic data from the .bam files using the pairwise sequential Markovian coalescent (PSMC) approach in the software package psmc (Li & Durbin, 2011). This analysis uses pairwise (allelic) estimates of coalescent times for each locus across the genome, developing a distribution of times to most recent common ancestor (TMRCAs), and the frequency of these coalescent events is inversely proportional to effective population size (Li & Durbin, 2011; Mather et al., 2020). We used a generation time of 2.5 yr, taken as an average among many of these species (Dellinger et al., 2020; Evans et al., 2020; Mack & Yong, 2020; Saether et al., 2005; Townsend et al., 2020; Whitaker et al., 2020), and a mutation rate of 2.3 × 10⁻⁹ mutations per site per year (Smeds et al., 2016). Single individuals represent lineage characteristics using this approach (Li & Durbin, 2011; Mazet et al., 2015; Nadachowska-Brzyska et al., 2015). Importantly for our focal questions, variations in these parameters do not affect the shape of the curve of population sizes (N_e) through time, but rather the values of the time and N_e estimates (Nadachowska-Brzyska et al. 2015). Being closely related and of similar size, we expect these taxa to have very similar generation times and mutation rates (e.g., Bird et al., 2020; Healy et al., 2014).

The PSMC method is not foolproof; for example, aspects of selection, inbreeding, and isolation can affect population size estimates (Mather et al., 2020). No method of genomic population demographic reconstruction is known to be consistently and reliably accurate (Marchi et al., 2021), but uniform application of a powerful method like PSMC to genome-scale data in a closely related group of lineages can provide new insights into the evolutionary effects of major life-history strategies.

For all variables, PSMC output for times more recent than 50 Kya was ignored, because recent estimates of effective population size are less reliable (Li & Durbin, 2011; Nadachowska-Brzyska et al., 2016). We consider this to be a conservative cutoff date. Our two focal variables, mean effective population size (N_e) and the coefficient of variation (SD/mean), were obtained from standard PSMC output. Our examination of early population growth involved deriving two additional variables from the PSMC output: degree of initial growth (1 - (N_trough/N_peak)), and rate of growth (degree/time period over which initial growth occurred; the latter we term deltaT; Fig. S2). We note that estimates of the dates of divergence events in Catharus evolution vary by study, genetic marker, methodology, and parameter estimates such as generation time and mutation rates (Outlaw et al., 2003; Topp et al., 2013; Voelker et al., 2013; Winker and Pruett, 2006). Date estimates in our study that differ from previous work reflect these methodological issues, and here it is the relative values among lineages that are most important.

We tested the assumption of phylogenetic independence of continuous characters obtained from the PSMC analyses using a phylogeny reconstructed using maximum likelihood from 1.5 Mb of autosomal sequence from chromosome 22 in MEGA X (Kumar et al. 2018), which produced the expected topology given prior work (e.g., Everson et al. 2019), and, for phylogenetic independence, the Abouheif-Moran approach (Abouheif, 1999; Münkenmüller et al., 2012) as implemented in the R package adephylo (Jombart et al., 2010). For those variables not showing phylogenetic signal, we tested our hypotheses using a one-tailed Mann-Whitney U-test. This nonparametric test has the advantage of not requiring that particular assumptions are met (e.g., homogeneity of variance), nor precise measurements (thus accommodating uncertainty in our N_e estimates; Sokal & Rohlf, 1995; Whitlock & Schluter, 2015).

Phylogenetic generalized least squares (PGLS) was used in a post-hoc analyses to test for differences in the single variable (deltaT) that showed phylogenetic signal and in followup analyses bringing in breeding latitude directly. PGLS can accommodate tests for relationships with discrete traits (here migratory-sedentary) and is basically a linear regression weighted in this case by detaT’s phylogenetic relationships (Freckleton, 2009; Symonds & Blomberg, 2014). Because this genus has experienced considerable interspecific gene flow (Everson et al. 2019), the accuracy of branch lengths is uncertain, so we ran analyses both using actual branch lengths and setting branch lengths to 1.0. We also extended our analysis of deltaT to examine its relationship to breeding latitude directly (using the same approach). Midpoints of breeding range latitudes were derived using ArcMap (ESRI 2021, ver.10.8.2) and range maps from BirdLife International and the Handbook of the Birds of the World (2021). PGLS analyses were done using the R packages phytools (ver. 0.7-80, Revell 2012) and nmle (ver. 3.1-152, Pinheiro et al., 2013) and using the corPagel model approach.

Supplementary Information

Table S1.

Table S2.

Figure S1.

Figure S2.

Figure S3.

Figure S4.