Biological evidence supports an early and complex emergence of the Isthmus of Panama.

The linking of North and South America by the Isthmus of Panama had major impacts on global climate, oceanic and atmospheric currents, and biodiversity, yet the timing of this critical event remains contentious. The Isthmus is traditionally understood to have fully closed by ca. 3.5 million years ago (Ma), and this date has been used as a benchmark for oceanographic, climatic, and evolutionary research, but recent evidence suggests a more complex geological formation. Here, we analyze both molecular and fossil data to evaluate the tempo of biotic exchange across the Americas in light of geological evidence. We demonstrate significant waves of dispersal of terrestrial organisms at approximately ca. 20 and 6 Ma and corresponding events separating marine organisms in the Atlantic and Pacific oceans at ca. 23 and 7 Ma. The direction of dispersal and their rates were symmetrical until the last ca. 6 Ma, when northern migration of South American lineages increased significantly. Variability among taxa in their timing of dispersal or vicariance across the Isthmus is not explained by the ecological factors tested in these analyses, including biome type, dispersal ability, and elevation preference. Migration was therefore not generally regulated by intrinsic traits but more likely reflects the presence of emergent terrain several millions of years earlier than commonly assumed. These results indicate that the dramatic biotic turnover associated with the Great American Biotic Interchange was a long and complex process that began as early as the Oligocene-Miocene transition.

The Isthmus of Panama is the narrow strip of land that connects North and South America and divides the Atlantic and Pacific oceans. The emergence of the Isthmus initiated one of the largest episodes of biological migration between previously disconnected landmasses, the Great American Biotic Interchange (GABI) (1–5), one of the best natural experiments on invasive species. The closure of the Central American Seaway (CAS, the oceanic pathway along the tectonic boundary between South America and the Panama Block) and rise of the Isthmus have been linked to the onset of both thermohaline oceanic circulation and northern hemisphere glaciation (6–8). Despite its broad importance, the formation of the Isthmus and its impact on the rich biodiversity of the Americas remains contentious (9). Therefore, a better understanding of when the formation of the Isthmus of Panama occurred has important implications in several scientific fields across multiple disciplines.

The timing of Isthmus formation has been assessed through different proxies. Previous studies have long suggested full closure by 3.5 Ma (7, 8, 10–16). More recent geological work has suggested a longer and more complex formation, where the initial collision between South America and the Panama Block occurred between 25 and 23 Ma (17). By 20 Ma the Panama Block is suggested to have been connected to North America (18–21) and the width of the CAS to be 200 km (19, 20). Full closure of the CAS occurred by 10 Ma, ending the exchange of deep and intermediate waters between the Caribbean and the Pacific (11, 19, 20, 22–24). However, the exchange of shallow waters between these oceans likely continued along pathways other than the CAS for many millions of years (7, 8, 10–16, 25).

Over the past two decades, hundreds of studies have assumed the Isthmus of Panama to have closed at 3.5 Ma (), causing the separation of widespread marine populations into distinct Pacific and Caribbean groups (vicariance) and the first possible dispersals between North and South America (with the exception of stochastic long-distance dispersals). Here, we address the following questions: Given the complexities and the recent evidence of a much older geological history of the Isthmus of Panama, is 3.5 Ma an adequate age for those events? Were the suggested water corridors across the Isthmus—even if shallow and narrow—indeed effective barriers against both dispersal of terrestrial organisms and conduits of marine ones until a full closure at 3.5 Ma? We address these questions using comprehensive biological data from living and fossil organisms, where the assumption is that any well-developed terrestrial corridor would lead to both more frequent biotic dispersal between North and South America as well as a division of widespread marine organisms into distinct Caribbean and Pacific lineages (26). Biological data provide a powerful tool for this purpose compared with geological evidence, which cannot inform on the subtle differences between exposed land and shallow waters that can be crucial to the dispersal of organisms.

Recent attempts to synthesize dispersal patterns across the Isthmus of Panama are based on dated phylogenies (27–29) that have included data that used a 3.5 Ma Isthmus closure as a calibration point (e.g., in calculating mutation rate). These assumptions lead to circularity if the goal is to examine timing of dispersal or vicariance across the Isthmus. Our molecular analyses differ from previous studies both qualitatively (by the elimination of the timing of the emergence of the Isthmus of Panama and the contrast between patterns in marine and terrestrial taxa) and quantitatively (735% more data points than a previous cross-taxonomic analysis) (27) and serve as a comparison with a comprehensive fossil dataset.

Another crucial aspect of the formation of the Isthmus of Panama and its impact on biological diversity entails the direction of terrestrial dispersal (e.g., north to south) and the ecology of dispersing organisms. In his seminal work on the GABI, Simpson (1) used the mammal fossil record to suggest that North American taxa had a competitive advantage over the South American fauna in that more North American taxa dispersed, survived, and diversified in southern ranges. Recent studies have suggested that ecological barriers to dispersal, such as dry savanna-like environments (30) or reduced rain forest cover (31), also prevented tropical South American taxa from migrating successfully to the north.

To explicitly address the timing of the formation of the Isthmus of Panama and its effect on geographical distributions of biota, we develop a Migration Rates Through Time (MRTT) model. Using this approach, in which migration is defined as dispersal and vicariance events collectively, we conduct an analysis of over 400 data points () of molecular divergence dates conferring dispersal events between North and South America or vicariance between the eastern Pacific and the Caribbean Sea. We compare results from extant data with over 23,000 records of American fossil mammals. Our results corroborate recent geological evidence that rather than a Pliocene, time-limited, single event, the formation of the Isthmus of Panama and the GABI were long and complex processes that began as early as the Oligocene–Miocene transition. We also show incongruent patterns in the direction of migration when comparing the molecular and fossil data, possibly reflecting differential diversification rates following dispersal. None of the ecological factors measured in the analyses can explain the statistical differences in the timing of migration across the Isthmus region among different taxa, and together with the MRTT results, this suggests that rather than intrinsic (biological) traits, extrinsic factors such as the presence of emergent land likely drove these patterns.

Results

Biological Data Reject a 3.5 Ma Assumption of Closure.

The most likely migration model for the molecular dataset as a whole demonstrates four statistically significant migration rate shifts. These are estimated at 41.1 Ma (95% confidence interval, 35.9–46.2), 23.7 (19.9–26.2), 8.7 (7.2–10.0), and 5.2 (5.1–6.0). The first shift at 41.1 Ma shows a sixfold rate increase from a negligible 0.001 (0–0.007) to 0.006 (0.003–0.008) migrations per family per million years, the shift at 23.7 Ma is a significant rate increase from 0.006 to 0.036 (0.029–0.044), whereas the shifts at 8.0 and 5.1 Ma represent drastic rate increases to 0.142 (0.106–0.207) and 0.371 (0.351–0.386) migrations per million years, respectively (Table 1 and ). All of the inferred migration rate shifts predate 3.5 Ma.

Table 1.

Migration models estimated from analysis of phylogenetic data

		Time of rate shift, Ma
Dataset, sample size	Model	First	Second	Third	Fourth	Log-likelihood	AICc	Delta AICc
All taxa, 426	Free model (five rates)	41.14	23.70	8.76	5.26	943.96	−1,869.48	0
	Standard model	3.5				670.72	−1,337.42	532.06
	Miocene model	20	10			876.85	−1,735.54	133.94
NA→SA, 144	Free model (three rates)	23.87	7.54			145.682	−280.928	0
	Standard model	3.5				78.94	−153.79	127.14
	Miocene model	20	10			128.66	−251.15	29.78
SA→NA, 154	Free model (three rates)	23.14	6.13			199.958	−389.511	0
	Standard model	3.5				102.85	−201.62	187.90
	Miocene model	20	10			128.66	−251.15	138.36
Marine, 86	Free model (four rates)	23.73	7.96	2.06		51.81	−88.17	0
	Standard model	3.5				5.41	−6.68	81.49
	Miocene model	20	10			41.07	−75.85	12.33

The results of migration rate analyses are shown for the full dataset, for the terrestrial taxa split by the direction of migration (North to South and South to North), and for the marine taxa. The free model—that is, the best-fitting model after optimizing the number and temporal position of the rate shifts—is compared against models with fixed times of shifts according to the “Standard” and the “Miocene” models. Rate estimates under the best model confidence intervals are provided in the MRTT plots (Fig. 2 and ).

Fig. 2.

Mammal fossil record in the Americas. The bar plots show mammal diversity through time based on fossil occurrences identified to the species level (reduced for clarity to the past 25 Ma, using time bins of 1 Ma) (see for more details). (A) Total number of fossil occurrences. (B) Diversity trajectories based on the first and last appearances of species (the thickness of curves reflects 95% confidence intervals of ages) (see also ) for North and South America. (C) Mean fraction of immigrant lineages—that is, North American lineages found in South America and vice versa (see for their definition)—plotted through time as a proportion of the total diversity sampled within each 1 Ma bin on each continent (confidence intervals are given in ). Approximate curves for global mean temperature (38) and eustatic sea level down-sampled to 0.1 Ma resolution (39) are plotted for reference.

To assess the effect of molecular data and uncertainty in divergence times on the estimated MRTT, the ages of the biogeographic events were randomized to simulate different proportions of error (). The dynamics of the MRTT obtained from the datasets with simulated errors were strikingly consistent with those estimated from the events’ mean ages, even when the relative error was more than twice as large as that observed from the data (50%; ). This indicates that our findings are robust even in the face of considerable amounts of potential error in the data, including ascertainment bias in taxonomic sampling, molecular clock calibration, and topological and branch length estimations (Fig. 1 and ).

Fig. 1.

MRTT estimated from molecular data. The plots show rescaled migration rate for dispersal (A–C) and vicariance (D) events across the Isthmus of Panama. The rates and their temporal dynamics were estimated using maximum likelihood and 1,000 randomizations to infer 95% confidence intervals of MRTT (shaded areas) (see also for a sensitivity analysis) and the time of each rate shift with their confidence intervals (circles and bars at the bottom of each panel). Dispersal between North and South America was calculated for (A) the entire dataset, (B) major taxonomic groups, and (C) direction of dispersal. (D) Comparison between vicariant events separating marine organisms in the Caribbean and the Pacific and the dispersal of terrestrial vertebrates (mammals, amphibians, nonavian reptiles, and birds) between North and South America. The number of rate shifts was selected by model testing via AICc (Table 1 and ). The dashed line indicates the generally accepted emergence of the Isthmus of Panama at 3.5 Ma.

We compared the Akaike Information Criterion (AIC) scores (32) of two migration models, one containing a single rate shift fixed arbitrarily at 3.5 Ma, consistent with the traditional interpretation of the formation of the Isthmus of Panama, and one encompassing two shifts, as recently suggested by Montes et al. (19, 20), at 25 Ma (initial collision of the Panama block and South America) and 10 Ma (closure of the CAS) ( and Table 1). In this explicit test of geological models for the emergence of the Isthmus, we found that although neither of the constrained migration models were as likely as those obtained under the best, unconstrained model, Miocene rate shifts were preferred over the previously accepted 3.5 Ma Pliocene shift. This is shown by an improvement in the second-order AIC by several orders of magnitude (Table 1). Results from these tests demonstrate that cross-Isthmian migrations occurred in incremental pulses over a long time period, presumably correlated with the stepwise formation and emergence of the Isthmus of Panama and in stark contrast to the assumption that a majority of events occurred at, or shortly after, a 3.5 Ma closure (dashed line in Fig. 1).

We split the molecular data to compare terrestrial dispersal patterns with those of marine vicariance events. We found that terrestrial migration exhibited significant increases at ∼20.23 (19.4–22.81) and 6.12 (5.1–7.89) Ma, whereas corresponding events separating marine organisms occurred at 23.73 (19.9–27.41) and 7.96 (7.75–8.96). Migration rates in marine organisms decreased by 50% after the third inferred shift at 2.06 (1.03–4.35) Ma.

We further compared our molecular results with inferences from fossil mammals, which provide the most abundant and best-studied fossil record. Although the number of fossil occurrences varies considerably through time due to differential sampling efforts and preservation rates (Fig. 2), there is no bias toward particular time periods at the species level, except for fossils 1 million years or younger (). The plots of the fossil data (Fig. 2 ) show that migration events occurred with increasing frequency since ca. 10 Ma and a drastic increase in North American taxa migrating south occurred over the last 3 Ma. The common explanation of this pattern is complete closure of the Isthmus and a competitive advantage of North American mammals over South American counterparts (1). However, we consider equally possible the scenario that the major drop in global temperatures in the last 3 Ma (Fig. 2, blue curve) acted as the triggering mechanism driving North American mammals southwards, as well as leading to the onset of northern hemisphere glaciations and rapid vegetation changes in North America. Dropping sea levels in the same period (Fig. 2, yellow curve) likely increased land exposure along coastal plains, which should have further facilitated migration.

Asymmetry in the Direction of Dispersal Events.

The directionality of migrations has long been an intriguing aspect of the GABI (2). In our analyses based on phylogenetic data, the asymmetric model (where migration rates and times of shift were considered as independent parameters) strongly outperformed the symmetric model [corrected ∆AIC (∆AICc) = 13.96; ]. This suggests that the migrations from South America to North America and in the opposite direction followed different dynamics. However, the estimated rates through time and the respective confidence intervals show that the rates in the two directions were largely overlapping across the long geological time and both underwent a similar rate shift around the transition between the Oligocene and Miocene (ca. 23 Ma; and Fig. 1). Migration rates only started to differ significantly after the most recent rate shift (around 6–7 Ma). During this time frame, lasting until the present, the overall migration rate from South America to North America is found to be around 30% higher than the rate from North America to South America.

Our analyses of the fossil mammal data recover a similar number of migrant lineages in South and North America between 11 and 3 Ma; however, these represent a higher fraction of the sampled diversity in South America (Fig. 2 and ). A large increase in the proportion of North American lineages in South America is observed over the last 3 Ma, resulting in about 45% of the sampled diversity of South American species (Fig. 2). This diversity is likely affected by our definition of migrant lineages and reflects the combined effects of migration and in situ diversification.

The asymmetry detected in the fossil record contrasts with the cross-taxonomic results from molecular phylogenies, where migrations in either direction were not found to be significantly different until ca. 6 Ma, when the migration rate from South to North America exceeds that in the opposite direction (Fig. 1 and ). The causes for this discrepancy are still uncertain. The fossil record of South American mammals is mostly derived from temperate latitudes and is heavily biased toward large body-sized grazer/browser animals, whereas the extant mammal data primarily include taxa inhabiting tropical regions. Perhaps fossil and molecular rates indicate the migration of two different biomes, temperate and tropical, respectively.

The Role of Ecology in Dispersal.

Organisms have different dispersal traits and life history strategies. Our results () support the view that plants, which have generally higher dispersal abilities compared with many animals, were the first organisms to migrate between North and South America (27, 33, 34). Freshwater fishes and amphibians, despite their requirements for survival, suggesting restricted dispersal ability, migrated relatively early compared with other organisms (). Birds began dispersing much later in the Pliocene, even though many are typically considered good dispersers (). The terrestrial vertebrate and marine data are highly congruent in their increase and decrease in the rate of exchange, respectively, between ca. 10 and 6 Ma (Fig. 1), suggesting that as the CAS closed and vicariance led to the split of species into distinct lineages in the Caribbean and the Pacific oceans, the Isthmus was formed and dispersal increased with connectivity between North and South America.

Despite the general trends detected among taxonomic groups, the Generalized Linear Mixed Model (GLMM) results showed that the variability in migration times is poorly explained by the metrics of dispersal and natural history we quantified (), and therefore the driving dynamics and mechanisms of migration remain elusive. Our results thus imply that migration events across the Isthmus of Panama were primarily determined not by intrinsic biological variables but potentially by extrinsic factors—including the availability of land, and sea and freshwater corridors, and the establishment of suitable climates and environments.

Discussion

The patterns inferred from the migration analyses of both molecular and fossil data (Figs. 1 and 2) are strikingly consistent and robust to uncertainties associated with age estimations (see confidence intervals in Figs. 1 and 2 and sensitivity analysis in ). Our results support an initial collision of the Panama Block and South America at 25–23 Ma with development of extensive terrestrial (although not necessarily fully connected) landscapes in Panama by 20 Ma, at least some 17 Ma earlier than generally assumed. Our results further imply that over the past 10 Ma there has been substantial dispersal and vicariance across the Isthmus and that pulses of dispersal of terrestrial organisms occurred over at least three periods during the last 30 Ma. The underlying causes for these discrete migration pulses remain elusive but may be associated with landscape formation, volcanism, climate change, and/or sea level fluctuations (25, 35).

North American taxa may have been preadapted to dispersal in that many of the successful lineages had migrated from Europe and Asia across northern land bridges, before their movement into South America. Despite this, our results on the North-to-South asymmetry detected in the fossil record contrast with the cross-taxonomic results from molecular phylogenies, where migrations in either direction were not found to be significantly different until ca. 6 Ma, when the migration rate from South to North America exceeds that in the opposite direction (Fig. 1 and ).

In summary, can we assume that “no vicariant date [3.5 Ma] is better dated than the Isthmus” (26)? Our results indicate that we cannot. We show that the GABI occurred over a much longer time period than previously proposed and comprised several distinct migrational pulses. Marine and terrestrial clades exhibit similar dispersal/vicariance pulses, occurring between ca. 23–20 Ma and 8–6 Ma. An earlier connection with North America, together with evidence of a southern connection with Antarctica until 30 Ma (36), challenges the long-standing idea of South America being an island continent that evolved in “splendid isolation” (1). These realizations impact our understanding of the temporal evolution and assembly of the American biota and urge a reevaluation of how the formation and emergence of the Isthmus of Panama impacted biological exchange, oceanic currents, atmospheric circulation, and global climate change.

Materials and Methods

Phylogenetic Data.

We performed an analysis of available molecular phylogenies with broad representation across taxonomic groups and habitats. We excluded studies that incorporated any assumption of the timing of the closure of the Isthmus (calibrations in dated phylogenies that were directly or indirectly derived from this event). We included phylogenetic studies with at least one unambiguous instance of dispersal across the Isthmus of Panama, for instance an ancestor with a North American distribution that had a descendant lineage in South America, as well as marine studies of sister species between the Caribbean and Pacific oceans. The dataset included 424 dated dispersal or vicariance events across the Isthmus of Panama (collectively referred to here as migration events). The data were compiled from 169 dated molecular phylogenies from across the tree of life and published in 29 peer-reviewed journals (, respectively). For each biogeographic event, we also recorded mean crown age, the lower and upper confidence intervals of ages, and direction of dispersal, as well as the taxon’s altitudinal range, dispersal capability, and biome of occurrence.

MRTT Analyses.

Biogeographic events with associated age uncertainties were analyzed in a maximum likelihood framework to estimate rates of migration () and their variation through time. Although, ideally, such migration rates should be estimated on a per-capita basis (as for other macroevolutionary rates), their estimation in this context would be affected by unobserved extinct lineages and further biases would likely arise from combining data from multiple independent phylogenies. We therefore modeled migrations as random events resulting from a stochastic Poisson process, with a rate parameter that describes the expected waiting time between successive events. As a mean to standardize the migration rates across datasets, we rescaled the estimated absolute Poisson rates by the number of families considered in each analysis. Thus, the migration rates shown in the MRTT plots indicate the expected number of migration events per Ma per family. We emphasize that this standardization is not intended as a correction for temporal biases within each dataset but rather as a tool to facilitate comparisons across different analyses. To account for deviations from a constant rate (homogeneous) process, we tested different nonhomogeneous Poisson processes with time-varying rates, including a model with exponentially increasing rates, which might capture potential biases owed to the exponential increase of lineages in extant-taxa phylogenies. To allow for temporal changes in the migration rates, we implemented a nonhomogenous Poisson process, in which rate shifts can occur through time. Our maximum likelihood algorithm involved (i) assessing the best-fitting number of rate shifts by a stepwise AICc procedure, (ii) the optimization of their temporal placement, and (iii) the estimation of the migration rates between shifts (). We used simulations to assess the most appropriate AICc thresholds, thus minimizing the risk of false positives in our analyses (). We generated and analyzed 1,000 datasets resampled from the uncertainty intervals for migration dates to generate the MRTT plots and calculate 95% confidence intervals.

In addition to jointly analyzing the full molecular dataset, we repeated the analyses on 14 subsets based on different criteria considering geographical, ecological, and taxonomic aspects. These tests allowed us to investigate the differential patterns of migrations linked to direction of migration and environment (terrestrial or marine) and to compare different taxonomic groups and their dispersal ability (Fig. 1 and ). Additionally, we compared the AIC scores of two migration models with a model containing a single rate shift fixed at 3.5 Ma and a model containing two shifts in the Miocene, at 25 Ma (initial collision of the Panama block and South America) and 10 Ma (closure of the CAS) ( and Table 1).

GLMM Analyses.

We used a GLMM approach on the molecular datasets of terrestrial/freshwater and marine taxa to test for intrinsic (biological, ecological) factors driving the differential timing of dispersal or vicariance. Ecological and biological variables (biome type, dispersal ability, dispersal direction, and elevation preference as derived from the literature and consultation with specialists) were treated as fixed effects, and we included a phylogenetic correction by treating taxonomic rank as a nested random effect (). Because the molecular dates in the dataset were calibrated using a range of approaches (e.g., fossils, molecular clocks), we included calibration type as an additional variable in the model to account for any potential bias associated with how the phylogenies were time-calibrated.

Fossil Data.

We compiled fossil occurrence data from all available publications of South American mammalian fossils and reviewed the Cenozoic American mammal fossil record in the Paleobiology Database (https://paleobiodb.org/#/), synthesizing data from across the entire American continent. The vetted dataset comprised 23,090 fossil records from 112 families and 3,589 species. Based on first and last appearances of each species in the fossil record, we plotted diversity trajectories through time, by counting the number of species within 1 Ma bins. To explicitly incorporate dating uncertainty of fossils, we randomized the age of each fossil 1,000 times within the temporal boundaries of the geological formation in which they were found (Fig. 2 and ). We classified each fossil occurrence as North or South American if the taxon or its ancestor was in either North or South America before 10 Ma, following Carillo et al. (37). We then plotted the amount of immigrant species as a proportion of the total diversity in each continent within 1 Ma time bins (Fig. 2).