Global and Planetary Change 227 (2023) 104153 Contents lists available at ScienceDirect Global and Planetary Change journal homepage: www.elsevier.com/locate/gloplacha Macroevolutionary and macroecological response of Iberian rodents to late Neogene climatic oscillations and events Jan A. van Dam a, b, *, Pierre Mein c, Miguel Garcés d, Ronald T. van Balen e, f, Marc Furió b, g, Luis Alcalá h a Department of Earth Sciences, Utrecht University, Vening Meinesz Building A, Princetonlaan 8a, 3584 CB Utrecht, the Netherlands b Institut Català de Paleontologia Miquel Crusafont, c/ Columnes s/n Campus de la UAB 08193, Cerdanyola del Vallès, Barcelona, Spain c A.R.P.A Université Lyon 1, Bâtiment Géode 69622, Villeurbanne Cedex, France d Departamento de Dinàmica de la Terra i de l’Oceà, Facultat de Ciències de la Terra, Universitat de Barcelona, Barcelona, Spain e Department of Earth Sciences, Vrije Universiteit Amsterdam, De Boelelaan 1085, 1081 HV Amsterdam, the Netherlands f TNO-Geological Survey of the Netherlands, Princetonlaan 6, 3584 CB Utrecht, the Netherlands. g Serra Húnter Fellow, Universitat Autònoma de Barcelona, Geology Department, Cerdanyola del Vallès (08193), Barcelona, Spain. h Parque de las Ciencias de Andalucía, Av. de la Ciencia, s/n, Granada (18006), Spain. A R T I C L E I N F O A B S T R A C T Editor: Dr. Howard Falcon-Lang Biozones are routinely used for chronological purposes, but their nature is seldomly questioned. Here, we attempt to find the evolutionary-ecological basis for Iberian rodent biozones for the interval 8.5–2 Ma based on Keywords: currently available paleontological, stratigraphic and paleoclimatic information. Our comparison of biozone Neogene boundary age uncertainty intervals to records of marine SST, δ18O and δ13C, terrestrial hydrological proxies and Spain astronomical parameters suggests an orbitally forced climatic origin for the majority of biozones. Zone boundary Rodentia ages during the late Miocene are mostly associated with humidity changes during 1.2-Myr obliquity nodes and Climate Precipitation 405-kyr eccentricity minima. Tectonics and strong regional cooling may additionally explain faunal change Milankovitch cycles during the Messinian (7–5 Ma). A partly reversed pattern characterizes most Pliocene zone boundaries, which are Macroevolution mostly associated with wetter conditions during obliquity nodes and 405-kyr eccentricity maxima. A third Macroecology configuration culminates in the early Pleistocene and consists of the combination of strong obliquity amplitude Migration maxima and 2.4-Myr eccentricity minima. It is suggested that larger-scale, European Neogene mammal (MN) Clade replacement units have an astronomical basis as well, with an important role of 2.4-Myr eccentricity cycle and the 1.2-Myr obliquity cycle, with the latter becoming especially important after the mid-Miocene cooling (~14 Ma). Whereas rodent events on the Iberian Peninsula during the late Miocene are shaped by replacements within resident communities dominated by dry-adapted clades of true mice (Murinae) and hamsters (Cricetinae), transitions during the cooler Pliocene involve invasions of newly emerging Eurasian clades of ‘microtoid ham­ sters’ and voles (Arvicolinae) preferring wetter and cooler environments. A new model of stepwise clade displacement is proposed in which the combination of long-period Milankovitch cycles and gradual long-term climatic change allows for the periodic functioning of migration corridors, along which increasingly pervasive dispersal events by members of new clades cause the gradual extinction of old clades. 1. Introduction processes that control ecological communities (macroecology) and into the biotic and abiotic factors that trigger species extinction, origination Since long, stratigraphers have been using biozones for correlating and migration (macroevolution; Vrba, 1995; Barnosky, 2001; Fraser sediments across a wide range of temporal and spatial scales. The et al., 2021). Apart from tectonic or eustatic sea-level events that oc­ evolutionary-ecological basis of such zones is often not well understood, casionally lead to openings or closures of migration routes between however. To comprehend the observed persistence of relatively fixed formerly separated (sub-)continents, climatic events would be among sets of taxa over very long periods, insight is required into the long-term the first to be held responsible for controlling million year-scale faunal * Corresponding author at: Department of Earth Sciences, Utrecht University, Vening Meinesz Building A, Princetonlaan 8a, 3584 CB Utrecht, the Netherlands. E-mail address: j.a.vandam@uu.nl (J.A. van Dam). https://doi.org/10.1016/j.gloplacha.2023.104153 Received 13 February 2023; Received in revised form 16 May 2023; Accepted 18 May 2023 Available online 24 May 2023 0921-8181/© 2023 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/). J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 dynamics. formulation of an overarching, Iberian-scale rodent stratigraphic Rodents are the most important organisms for subdividing Cenozoic framework (Table 1), on the basis of which the age accuracy of many continental stratigraphy. Their rapid evolution has helped geologists to sites could be increased. define local biozones with durations in the order of several hundred In this study we use species relative abundances (molar frequencies) thousands of years at the local basin scale (van de Weerd, 1976; Kälin for reconstructing both temporal and geographical gradients in com­ and Kempf, 2009; van der Meulen et al., 1999) and bioprovincial scale munity structure. Focusing on the dominant family of Muridae, our (van Dam et al., 2023) to 1–2 Myr at the continental scale (Mein, 1990; analysis of these gradients leads to a better understanding of patterns of Hilgen et al., 2012; Speijer et al., 2020). With abiotic factors, such as migration and evolution and their combined effect during ‘clade climate change, and biotic factors, such as competition, continuously replacement’ (the substitution of one higher taxonomic unit by another, affecting communities at ecological time scales, biozone-scale commu­ e.g., Benton, 1996; Fraser et al., 2021). nity stasis should be qualified as a remarkable phenomenon that de­ serves an explanation (Ivany, 1996). 2. Material and methods In a previous analysis covering a large part of the Neogene, we found support for a climatic mechanism explaining the million-year scale 2.1. Chronology clustering of rodent turnover, with faunal change being concentrated during low-amplitude intervals (nodes) of the Earth’s Milankovitch cy­ Age uncertainty ranges for Iberian assemblage biozone boundaries cles of eccentricity and/or obliquity (2.4-Myr and 1.2-Myr cycles; van and individual mammal sites are taken from van Dam et al. (2023), who Dam et al., 2006). In that study, which was based on records from combined Iberian rodent biostratigraphy with local magneto­ central Spain, increased turnover during such nodes was hypothesized to stratigraphy and interpolations based on sedimentation and/or evolu­ work in two ways: 1) remotely, via high-latitude global cooling events tionary rates (Table 1). Iberian rodent biozones are coded by ‘Ib-N’ and their global effects, and 2) locally/regionally, via the breaking down (Iberian Neogene) and followed by a number (11, 12, etc.) that roughly of the ‘regular’ evolutionary-ecological response to short-period cycles corresponds to the time-equivalent MN unit, and a letter to indicate the (i.e., short-distance migration and re-assembly of incipient species with actual Zone (Ib-N12a, Ib-N12b, etc.). All paleontological and climatic mother species). records included in this paper are calibrated to the newest Geomagnetic To study the response to medium-period cycles such as the ~100-kyr Polarity Time Scale (GPTS; Raffi et al., 2020), which, in turn, is cali­ and 405-kyr eccentricity cycles, a high temporal resolution of fossil lo­ brated to the Astronomical Time Scale (ATS; Laskar et al., 2004; Hilgen calities is indispensable. Furthermore, to better understand the spatial et al., 2020). dimension of macroevolutionary and macroecological processes, the inclusion of a geographical component is highly desirable. Here, we 2.2. Relative abundances focus on the climatically dynamic late Miocene to early Pleistocene (8.5–2 Ma), studying rodent communities from a series of Iberian basins Relative abundances for 27 murid lineages consisting of true mice (Fig. 1), for several of which new rodent and chronological data have (Murinae), hamsters (Cricetinae, Cricetodontinae), vole-like hamsters recently become available. Our recent revision of all paleontological and (Baranomyinae, Trilophomyinae), and voles (Arvicolinae) were stratigraphic evidence (van Dam et al., 2023) has resulted in the compiled per locality and per zone. Several rare lineages and Fig. 1. Map of the Iberian Peninsula with Cenozoic basins used in this study. 2 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 Table 1 Iberian rodent zones with boundaries and ages (after van Dam et al., 2023). Iberian Assemblage Zone Zone boundary Age (Ma) Duration (Myr) Age (Ma) Duration (Myr) Name Code Maximum Minimum Midpoint Uncertainty Midpoint Average top Ib-N17b 1.90 1.80 1.85 0.10 Mimomys realensis Ib-N17b 1.98 0.25 Ib-N17a-b 2.15 2.05 2.10 0.10 Mimomys medasensis Ib-N17a 2.33 0.45 Ib-N16b-17a 2.56 2.54 2.55 0.02 Kislangia ischus Ib-N16b 2.80 0.50 Ib-N16a-b 3.06 3.03 3.05 0.03 Mimomys hassiacus Ib-N16a 3.12 0.14 Ib-N15b-16a 3.19 3.18 3.19 0.01 Dolomys adroveri Ib-N15b 3.58 0.78 Ib-N15a-b 4.01 3.92 3.97 0.09 Mimomys davakosi Ib-N15a 4.10 0.27 Ib-N14c-15a 4.27 4.19 4.23 0.08 Trilophomys castroi Ib-N14c 4.42 0.38 Ib-N14b-c 4.71 4.51 4.61 0.20 Paraethomys baeticus Ib-N14b 4.75 0.28 Ib-N14a-b 4.91 4.87 4.89 0.04 Apocricetus barrierei Ib-N14a 5.15 0.51 Ib-N13c-14a 5.47 5.33 5.40 0.14 Paraethomys meini Ib-N13c 5.83 0.86 Ib-N13b-c 6.26 6.25 6.26 0.01 Apocricetus alberti Ib-N13b 6.34 0.16 Ib-N13a-b 6.56 6.27 6.42 0.29 Castromys inflatus Ib-N13a 6.69 0.55 Ib-N12b-13a 7.02 6.90 6.96 0.12 Parapodemus meini Ib-N12b 7.05 0.18 Ib-N12a-b 7.18 7.09 7.14 0.09 Huerzelerimys turoliensis Ib-N12a 7.40 0.53 Ib-N11-12a 7.68 7.65 7.67 0.03 Huerzelerimys vireti Ib-N11 8.18 1.03 base Ib-N11 8.75 8.63 8.69 0.12 unidentified forms not assignable to any of the main lineages were consuming invertebrate prey) in the same micromammal community. combined with jumping mice, gerbils and myocricetodontines into a Eq. (2) differs slightly from the one published by van Dam and single rest group (Muridae plus Dipodidae rest group). The remaining Utescher (2016), which contains a minor error (resulting in a difference taxa (Sciuridae, Gliridae, Castoridae, Hystricidae) were aggregated into of maximally 1.5 mm). A second modification regarding the latter work a second rest group. Abundances were based the total number of molars concerns the ecology of the squirrel genus Csakvaromys (which includes (M/m123). In case only counts of measurable specimens were available ‘Spermophilinus’; Sinitsa et al., 2022) for which we assume an arboreal from the literature, those counts were used as an alternative for total way of life based on a study of its microwear and cranial shape (laterally counts. Abundance data were restricted to zones or localities with n ≥ 25 directed orbits; Gusovsky and Sinitsa, 2021; Sinitsa et al., 2022). specimens and to localities with age uncertainties ≤200 kyr. Prediction intervals for individual sites typically are ±350, ±20, and ± 60 mm wide for MAP, LMP and HMP, respectively (van Dam and Utescher, 2016). 100 specimens was taken as a minimum for Eqs. (1)–(3) 2.3. Precipitation and lake level proxies to be applied. In the absence of quantitative information on eulipoty­ phlans and lagomorphs, a minimum of 70 rodent specimens was used for We used the Climate-Diversity Approach of van Dam (2006) and van a locality to be included (85 in case one group was missing). In the total Dam and Utescher (2016) to predict precipitation levels. The following absence of any frequency data, we considered 10 species as the mini­ equations for mean annual precipitation (MAP) and driest- and wettest- mum for a locality to be included. Because the lagomorph material of month precipitation (LMP, HMP) were used: some localities has sometimes remained unstudied, we used the lago­ morph species richness and micromammal specimen counts of time- MAP = 9.466 + 45.119 ln2 (PA + 1) + 12.743 PI (1) equivalent localities to approximate the number of lagomorph species ( R2 = 0.83 ) as follows (with nR+Li = number of rodent plus lipotyphlan specimens): one ochotonid species if nR+Li ≥ 10, one leporid species in Ib-N14a-15b LMP = 1.689 + 10.801 ln(PA + 1) (2) if nR+Li ≥ 100 and in Ib-N16a-17b if nR+Li ≥ 10, and zero leporid species in Ib-N13c or lower (unless reported). Furthermore, unidentified (indet.) ( R2 = 0.79 ) Eulipotyphla, Soricidae, Talpidae, Sciuridae were assumed to be terrestrial, Arvicolinae and Mimomys indet. were assumed to be terres­ HMP = 8.108 + 1.674 ln3 (PA + 1) + 1.896 PI (3) trial up and including Zone Ib-N15a, and (semi-)aquatic for Zones Ib- N15b-17b. ( R2 = 0.74 ) The proportion f(aq) of (semi-) aquatic micromammal species in the fauna (beavers, watermoles, aquatic shrews and water voles) was used with as a proxy for the presence of water bodies. We note that the extrapo­ PA = percentage of arboreal species in the non-aerial, non-(semi) lation of a semi-aquatic lifestyle of several modern voles (Arvicola, aquatic micromammal community (rodents, eulipotyphlans, Ondatra) to their supposed Pliocene relatives (Mimomys hassiacus - M. lagomorphs). pliocaenicus, M. medasensis and Dolomys) is hypothetical, although PI = percentage of invertivorous species (species dominantly generally agreed upon by specialists (see van Dam, 2006). A second, 3 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 more qualitative hydrological indicator consists of the stratigraphic et al., 2023). Similar peak-like occurrences of Apocricetus occur in Ib- distribution of lacustrine/palustrine facies complexes (presence of car­ N13a and 13c at 6.6 and 6.0 Ma (Crevillente 14 and 6; Freudenthal bonate/marl lithologies) in the Teruel Basin, where the sedimentary et al., 1998). record is the most complete. 3.2. Precipitation and lake level 3. Results Precipitation predictions were inferred for 78 localities resulting in 3.1. Relative abundance trends and gradients mean values of 380, 9 and 59 mm for MAP, LMP and HMP, respectively (Fig. 3D; Table S3). Highest average and wettest-month values (MAP = 3.1.1. Zone-level comparisons 831, HMP = 114 mm) were estimated for Concud 3 (Ib-N12a, Teruel Lineage relative abundances per zone and region are shown in Basin, Central Spain). A driest-month maximum was recorded for Fig. 2A (data in Table S1). Three major patterns can be discerned: Crevillente 8 (Ib-N12a, Bajo Segura Basin, LMP = 34 mm). Very dry annual estimates of only ~10 mm/yr were estimated for several 1) A dominance shift from Murinae and Cricetinae in the late Miocene southeastern Ib-N13a-14c localities in the Fortuna Basin (Sifón de Lib­ to vole-like hamsters (Baranomyinae, Trilophomyinae) and Arvico­ rilla 1, 61, 413 and P), Cabriel Basin (La Bullana 2B, Fuente del Viso), linae in the Pliocene and Early Pleistocene. Granada Basin (La Mina 2) and Guadix-Baza Basin (Gorafe 4/A). These 2) An exclusive presence or a presence with larger abundances of new extreme values, which would point to truly desertic conditions, are a lineages in the North (i.e., Central Spain) compared to the South. direct reflection of the absence of (insectivorous) eulipotyphlans (PI = 0; This pattern applies to Stephanomys in Ib-N13a, Blancomys sanzi in Ib- Eq. (1)), a group with standard presence in Iberian micromammal N13a, Celadensia in Ib-N14a, Promimomys in Ib-N14b, Trilophomys faunas. Our method is effectively breaking down for these sites with from Ib-N14c onwards, Mimomys davakosi in Ib-N15a, M. gracilis in MAP/12 values exceeding HMP estimates (Table S3). This discrepancy Ib-N15b, M. hassiacus in Ib-N16a, and Kislangia ischus in Ib-N16a. An occurs, because even the driest Sahara desert sites in the modern cali­ exception is the arvicoline Dolomys adroveri in Ib-N15b, which has a bration data set for Eqs. (1)–(3) have at least one eulipotyphlan species slightly larger proportion in the South (Guadix-Baza Basin) than in present (van Dam, 2006). Regardless of the possible existence of some the North (Teruel Basin). On the other hand, the scarcity of arvico­ sort of taphonomic bias against eulipotyphlans, it is safe to say that lines in the same zone in the Jerez Basin, the most south-western eastern/southern Spain was occasionally extremely arid during the lat­ region in this study for which there are data, confirms the overall est Miocene to earliest Pliocene. Even today, this area is the driest part of pattern. Superposed on the latitudinal pattern is a temporal pattern the Iberian Peninsula (e.g., Tabernas Desert; province of Almeria). with younger immigrants such as Trilophomys, but especially Dol­ The overall trend in precipitation consists of a gradual decrease from omys, Mimomys and Kislangia penetrating further south into the the late Miocene to early Pleistocene, with a distinct maximum occur­ Peninsula than older immigrants such as Blancomys sanzi, Celadensia ring during Ib-N12a (7.65–7.15 Ma) and smaller maxima concentrated and Promimomys, which fail to colonize the South. As a net result, the around the Ib-N14a-b boundary and during Ib-N15a and 16a (Fig. 3C). quantitative replacement in the South of the older true mouse and In addition, a trend towards seasonally more homogeneous levels hamster clades by the vole-like and vole clades lags the equivalent characterizes Ib-N15b. Maxima in the occurrence of lacustrine/palus­ replacement in the north by several millions of years (Fig. 2B, trine sediments in the Teruel Basin (indicating lake level, Fig. 3F) are Table S2). approximately coeval with precipitation maxima. 3) A late Miocene east-west gradient in true mice versus hamsters. Whereas the former are more abundant in the north (Teruel and 4. Discussion: The origin of rodent zones Madrid Basins) and extreme south (Granada Basin), the latter are especially abundant in the east (Bajo Segura, Cabriel and Alcoy Ba­ Our search for mechanisms behind rodent zone transitions primarily sins Basin). This eastern part also functions as the major colonization consists of a systematic investigation of time-correlative and potentially area for the murine Paraethomys, which is regarded as an African perturbing extrinsic events. Three types of events are chronologically immigrant. explored: 1) high-latitude climatic events as reflected in North Atlantic δ18O and sea-surface temperature (SST) records (ODP drilling sites), 2) 3.1.2. Locality-level comparisons regional, mid-latitude/Mediterranean climate- and tectonics related At the locality level, the following patterns can be noted (Fig. 2C): events as observed in marine δ18O, δ13C and SST records and terrestrial hydrological proxy records, 3) critical climatic episodes as predicted 1) The most important abundance changes occur at boundaries. As a from astronomical parameters (eccentricity/precession and obliquity). result, lineages may show high and relatively constant abundance These chronological explorations will be followed by a summary of levels across a specific zone, and much lower levels below and above those events that most successfully predict faunal change across zone it (Fig. 2C: left part). For example, Castillomys shows a more or less transitions. The section finishes with the formulation of a new model of constant level of temporary raised abundances in the Teruel Basin clade displacement and a discussion on the origin of European-scale during Zone Ib-N15b. Similarly, although it is relatively rare, mammal chronological (MN) units. Blancomys has its maximum abundance in Ib-N16a. Comparable examples from the Bajo Segura Basin include Cricetulodon lucentensis 4.1. High-latitude climatic events across Ib-N12a, Neocricetodon across Ib-N13a, and the Parapodemus barbarae – P. meini – Castromys lineage across both zones (Fig. 2C). The Ib-N11/12a transition between 7.68 and 7.65 Ma is not associ­ 2) Peak events involving abundances as high as 70–80% at the base of a ated with any distinct positive oceanic benthic δ18O excursion (cooling/ zone, followed by a decline (Fig. 2C: upper right part). Examples are glaciation) event in the NE Atlantic (Fig. 3B) or more distant regions (e. Mimomys davakosi, Dolomys adroveri and M. hassiacus in the Teruel g., equatorial Atlantic and Pacific sites ODP-926 and IODP-U1337; Basin. Shackleton and Hall, 1997; Drury et al., 2017). Alkenone-based SST 3) Occasionally, even more extreme peak events occur during which a records do show some cooling, however, with values starting to drop at sudden abundance change is immediately neutralized. Examples are ~7.9 Ma at site ODP-982 and culminating in two small cooling peaks at short-lasting maxima of Parapodemus gaudryi (Fig. 2C: lower right 7.57 and 7.52 Ma (Fig. 3A). A prominent positive benthic δ18O excursion part) and Apocricetus aff. plinii at the Ib-N11-12a boundary in the of ~0.5‰ at 7.56 Ma ± 0.04 Myr was nevertheless reported in the Bajo Segura Basin (Crevillente 4; de Bruijn et al., 1975; van Dam southeastern Atlantic (ODP-1085; Westerhold et al., 2005). In turn, this 4 J.A. van Dam et al. 5 Global and Planetary Change 227 (2023) 104153 Fig. 2. Rodent relative abundance patterns. (A) Lineage patterns per basin and per Ib Zone (Ib-N11-N17b, specimen n ≥ 25); long red bars: zone boundary uncertainty intervals; blue, green and purple colors: true mice (Mur., M. = Murinae; lineages 1–11); red, yellow and brown colors: hamsters (Cric., C. = Cricetinae/Cricetodontinae; lineages 12–17); khaki and black colors: vole-like hamsters (mc = ‘microtoid cricetids’; lineages 18–20); grey and white colors: voles (Arv., A. = Arvicolinae; lineages 21–27); pink colors: rest groups (28–29). Individual lineages: 1) Huerzelerimys vireti - H. turoliensis, 2) Progonomys hispanicus – Occitanomys sondaari – O. adroveri - O. brailloni, 3) Parapodemus gaudryi – Apodemus gudrunae – A. gorafensis, 4) Parapodemus lugdunensis - Apodemus atavus, 5) Parapodemus barbarae - P. meini - Castromys juniensis - C. littoralis - C. inflatus, 6) Stephanomys ramblensis – S. dubari – S. cordii – S. margaritae – S. vandeweerdi – S. minor – S. prietaensis – S. balcellsi – S. progressus, 7) Stephanomys donnezani - S. thaleri, 8) Occitanomys alcalai, 9) Paraethomys meini, 10) Paraethomys baeticus - P. abaigari - P. aff. jaegeri - P. jaegeri, 11) Castillomys gracilis – C. crusafonti – C. rivas, 12) Neocricetodon occidentalis - N. seseae, 13) Cricetulodon meini - C. lucentensis, 14) Apocricetus aff. plinii - A. plinii - A. alberti - A. barrierei - A. angustidens, 15) Hispanomys aragonensis - H. peralensis – H. freudenthali – Ruscinomys schaubi – R. lasallei – R. europaeus, 16) Hispanomys adroveri - H. aff. adroveri - H. romeroi, 17) Blancomys sanzi - B. neglectus - B. meini, 18) Celadensia nicolae, 19) Trilophomys castroi - T. vandeweerdi, 20) Trilophomys depereti - T. sp. A., 21) Mimomys davakosi, 22) Dolomys adroveri - Kislangia cappettai - K. aff. cappettai, 23) Kislangia ischus - K. gusii, 24) Mimomys gracilis - M. stehlini, 25) Mimomys hassiacus - M. polonicus - M. praepliocaenicus - M. pliocaenicus, 26) Mimomys medasensis, 27) Mimomys realensis. Rest groups: 28) remaining lineages in Muridae and Dipodidae, E = Eozapus, 29) remaining lineages in Castoridae, Gliridae, Hystricidae and Sciuridae. For data see Table S1. (B) Trends in the proportion of ‘old clades’ (Murinae, Cricetinae, Cricetodontinae) versus ‘new clades’ (‘microtoid cricetids’, Arvicolinae) in central Spain (Teruel Basin) and southern Spain (Granada and Guadix-Baza Basins); black squares represent zone averages, with values for southern Spain representing average proportions for the two basins; localities with ambiguous zone assignment are excluded. For data see Table S2. (C) Trends for selected taxa (see text); included are sites with n ≥ 25 and age uncertainty ≤200 kyr; red bars represent zone boundary uncertainty intervals. For data see van Dam et al. (2023). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 6 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 Fig. 3. Correlation diagram of Ib-N Zone boundaries with Mediterranean and Atlantic climate records, estimated precipitation (Spain), proportions (semi-)aquatic species and major sedimentary units (Teruel Basin), event was correlated to a (low-resolution) isotopic event as observed in Atlantic (1085; Westerhold et al., 2005), N. Adriatic Sea (Monte del Casino; Kouwenhoven and van der Zwaan, 2006), Mediterranean (Pliocene stack; Wang et al., 2010), and global oceans (Pliocene stack; Lisiecki and Raymo, 2005). (C) precipitation levels estimated from micromammals (this paper) for localities with age uncertainty ≤200 kyr; squares: average monthly level, left triangle: driest month level, right triangle: wettest 2004); astrochronozone numbering after Hilgen et al. (2020); black line: eccentricity, grey line: envelope (100-kyr maxima). (E) proportion (semi-)aquatic species in non-aerial micromammal community (f(aq)), and chronology by Lawrence et al. (2009) and Drury et al. (2017) (ODP-982) and by van der Schee et al. (2016) (IODP-1387, 3.30–5.36 Ma); only magnetostratigraphically calibrated part younger than 7.3 Ma included for 907 (Channell et al., 1999). (B) Benthic (BE) and planktic (PL) δ18O and δ13C records from the N. Atlantic (ODP-982; Hodell et al., 2001; re-calibrations: Lawrence et al. (2009), Drury et al. (2017)), S. and orbital parameters (eccentricity, obliquity). Horizontal red bars: zone boundary uncertainty ranges, horizontal black arrows: important events that correlate to zone boundaries (see text). (A) Alkenone (Uk’37)-based month level, blue arrows: trends, black: Teruel Basin, green: Bajo Segura Basin, red: Granada Basin, magenta: Cabriel and Júcar Basins, pink: Alcoy Basin, yellow: Guadix-Baza Basin. (D) Eccentricity record (Laskar et al., major sedimentary units in Teruel Basin (Mein et al., 1990; Opdyke et al., 1997; Alcalá et al., 2000; van Dam et al., 2001; this paper); Units I-V according to Alcalá et al. (2000), blue arrows: trends f(aq), blue bars: lacustrine/palustrine facies. (F) Obliquity record (Laskar et al., 2004); astrochronozone numbering after Hilgen et al. (2020). (For interpretation of the references to colour in this figure legend, the reader is referred to SST records from N. Atlantic (ODP-907, 982, 1387) and Mediterranean (Herbert et al., 2016; Tzanova et al., 2015; Tzanova and Herbert, 2015) re-calibrated according to Conte et al. (2006) and with re-calibrated northern Atlantic site ODP-608 (MTi4; Abreu and Anderson, 1998). With site ODP-1085 potentially affected by Antarctic Circumpolar (in addition to North Atlantic) Deep Water, this positive excursion may have been triggered by one of a series of major glaciation pulses on Antarctica. Cross-equatorial atmospheric forcing with northward ITCZ shifting would then have to be invoked as a mechanism for transporting the effects to the Northern Hemisphere (NH) (as similarly hypothesized for the Middle Miocene; e.g., John et al. (2003), Holbourn et al. (2010)). A peak in Southern Hemisphere (SH) glaciation may thus have led to a climatic perturbation in the northern part of the NH subtropical-dry region (for instance as a drying event at the Ib-N11-12a boundary, Section 4.2). Besides the hypothetical nature of this teleconnection, there is a small discrepancy in timing (7.67–7.66 Ma for the fauna, maximally 7.60 Ma for MTi4) that still needs to be explained. No major well-known open ocean δ18O events are associated with the Ib-N12a-b and 12b-13a transitions at ~7.15 and ~ 6.95 (Fig. 3B). As can be seen in Fig. 3B, it is not until 6.4 Ma that strongly positive values are recorded (‘C3An.1n events’, Fig. 3B). These events (which are preceded by a negative shift in δ18O culminating at 6.43 Ma) herald the onset of the so-called ‘Messinian glaciation’. Cooling is also observed in the Atlantic SST record at this time (Fig. 3A). Because of the large age un­ certainty associated with the Ib-N13a-b boundary, it is not possible to establish whether faunal change is associated with this cooling or with the preceding warming. The Ib-N13b-c boundary at 6.26–6.25 Ma correlates to the end of the ‘C3An.1n-events’ period and to low SST values between 6.37 and 6.24 Ma (e.g., C3An.18O.14/16 at 6.28 and 6.24 Ma; Hodell et al., 1994). The entry of the African rodent immigrant Paraethomys as recorded in sites Venta del Moro and Sifón de Librilla 52 (~6.25 Ma) may have occurred during the second of these isotopic events. If this was indeed the case, faunal change at the Ib-N13b-c boundary would qualify as an example of high-latitude forcing (via sea level). At the same time, this case is unique, as the Messinian is the only interval during which this type of direct immigration from Africa would have been possible. (We note that besides Paraethomys also other lineages such as Myocricetodon and Debruijnimys almenarensis crossed the Gibraltar barrier. However, because of their rarity and poorly constrained ages, their first appear­ ance cannot be dated with accuracy and could even be coeval with that of Paraethomys.) Perhaps surprisingly, both the onset of the Messinian Salinity Crisis (5.96 Ma) and the two well-known glacial events TG20/22 and TG12 (5.75–5.80 and 5.55 Ma; Shackleton et al., 1995; Fig. 3B) are not asso­ ciated with any zone boundary. Similarly, no major open ocean climatic event (δ18O and/or SST) is recorded at the Ib-N13c-14a boundary (range 5.47–5.33 Ma). If any climatic link to the Ib-N13c-14a boundary had to be made, it would be the occurrence of the ‘interglacial’ interval at 5.5–5.3 Ma (Fig. 3A,B). The Ib-N14a-b boundary (4.91–4.87 Ma), which involves the disap­ pearance of Celadensia and the simultaneous expansion of Paraethomys (two lineages) and Stephanomys, correlates to well-known glacial event Si6 at 4.89 Ma (Shackleton et al., 1995). The possibility cannot be excluded that also Promimomys, the first arvicoline immigrant in Europe, entered the Iberian Peninsula during Si4–6. However, in contrast with south-eastern Europe (where it is already well-represented from 5.3 Ma onwards; Hordijk and de Bruijn, 2009), Promimomys had a very marginal presence on the Iberian Peninsula with occurrences in only three not the web version of this article.) very accurately dated Ib-N14b faunas of the Teruel Basin (La Gloria 4, Celadas 9/12; Adrover et al., 1993; van Dam et al., 2023). The lack of any sign of SST cooling until 4.0 Ma west of Gibraltar (U1387, record starting 5.4 Ma; Tzanova and Herbert, 2015) may be due to low sam­ pling resolution (Fig. 3A). In Greece and on Sicily, on the other hand, a drying effect of event Si6 has been inferred based on an expansion of herbs including Artemisia (Kloosterboer - van Hoeve Kloosterboer - van Hoeve, 2000). Whereas the Ib-N14b-c and Ib-N14c-15a transitions (4.71–4.51, 7 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 4.27–4.19 Ma) are not associated with any well-known high-latitude attained by van Dam (2006)). The increase in precipitation around 7.6 climatic event (Fig. 3B), the Ib-N15a-b boundary (4.01–3.92 Ma) Ma is consistently associated with an increase of the average proportion temporally correlates to glacial event Gi20 (4.00 Ma; Lisiecki and of (semi-)aquatic micromammals from 0.01 in Ib-N11 to 0.08 in Ib- Raymo, 2005). On the other hand, the long Zone Ib-N15b (~4–3.2 Ma) N12a. bridges well-known glacial events such as Gi2/4 (3.67, 3.63 Ma) and M2 Precipitation increase at the transition from Ib-N11 to N12a is (3.30 Ma). consistent with lake expansion in the Teruel Basin as indicated by an The lack of a clear expression of the global δ18O event M2 at 3.3 Ma increase in carbonate/marl lithologies (Fig. 3E). In fact, the age of the not only characterizes the Spanish rodent record, but also the Mediter­ boundary (7.68–7.65 Ma) exactly coincides with a major lithological ranean SST record (Fig. 3A). Likewise, the astronomically calibrated switch from alluvial to lacustrine sediments as recently dated for the record of Punta Piccola (Sicily, representing the 3.7–2.6 Ma part of the base of a 30 m thick cyclic, limestone/marl unit in the Concud area (van Mediterranean composite SST record of Herbert et al., 2016; Fig. 3A) is Dam et al., 2023). This unit correlates to various other local carbonate- not showing any glacial peak at this time. Similarly, the ensuing Ib- dominated bodies in the basin that have collectively been informally N15b-16a transition at 3.19–3.18 Ma and the Ib-N16a-b transition at named ‘Páramo 1’ (e.g., Godoy et al., 1983; Unit II of Alcalá et al., 2000; 3.06–3.03 Ma do not correspond to any well-known glacial event. The genetic unit TN3 of Ezquerro et al., 2020; Fig. 3E). For instance, 3–4 km apparent ‘immunity’ of Mediterranean terrestrial environments to these north and northeast of Concud this unit is 15-20 m thick and likewise events seems to contrast with biotic responses in more northern areas; associated with N12a/b sites (Las Casiones and Lomas de Casares sec­ for instance, pollen records from Norway show an expansion of (similar tions: localities la Casiones 2 and Barranco de Cuevas 1–4; Mein et al., to present-day) alpine ‘cool’ vegetation around 4.2–4.0 and 3.4–3.3 Ma 1990; van Dam et al., 2001). Lacustrine sedimentation started somewhat (Panitz et al., 2016), intervals that include well-known cooling events earlier in the more central and central-eastern parts of the basin Gi20 and M2, respectively. (Fig. 3E), where carbonate/marl series are recorded just above the Ib- The Ib-N16b-17a boundary at 2.56–2.54 Ma is temporally very close, N11 sites Masada Ruea 2 (8.5 Ma), Los Aguanaces 1/3, Alfambra (8.3 although slightly preceding isotope stage MIS 100 (2.52 Ma). Together, Ma), Tortajada A (8.2 Ma) and El Búnker 6/7 (= Valdecebro 4, 7.9 Ma) stage 100 and stages 98 and 96 at 2.49 and 2.43 Ma represent the first (van Dam et al., 2001). Precession cycle counting in the Concud area has important glacial phase of the Early Pleistocene (e.g., Lisiecki and yielded an age uncertainty range of 7.05–7.00 Ma for the top of the Raymo, 2005; Fig. 3B). At 2.52–2.51 and 2.48 Ma, Mediterranean SST Páramo 1 unit. Whether this top is isochronous across the basin is un­ drops ~5 ◦ C and Atlantic SST decreases as well (ODP-982, Fig. 3A). The sure; for instance, in the El Búnker section its age uncertainty interval is next rodent boundary (Ib-N17a-b at 2.15–2.05 Ma) also coincides with estimated at 7.00–6.71 Ma (van Dam et al., 2023). the onset of a series of strongly 18O-enriched peaks (stage MIS 78 at 2.12 The observed correspondence between precipitation increase and Ma, followed by stage 82 at 2.07 Ma), and with another Mediterranean lake expansion suggests an important role for climate in addition to SST drop (2.08–2.06 Ma). The top of Ib-N17b (1.9–1.8 Ma), on the other tectonic subsidence as an agent for lake formation in the Teruel basin. In hand, is not associated with a single specific, well-known δ18O event. the Bajo Segura Basin, a MAP increase is observed as well, with mean In summary, the evidence for a triggering of Iberian faunal transi­ values increasing from 388 to 574 mm between Ib-N11 and N12a (n = 4 tions by high-latitude climatic/glaciation events is relatively thin for the and 2 sites, respectively). The return of the hamster Cricetulodon in this latest Miocene and Pliocene (8.7–3 Ma, exceptions being events Si4–6 at basin in Ib-N12a (Freudenthal et al., 1998) is consistent with a humidity 4.9 Ma and Gi-20 at 4.0 Ma), but relatively strong for the earliest increase, as can be inferred from its earlier disappearance from the Pleistocene (2.5–2 Ma). Iberian record at 9.5 Ma during a drying episode associated with the ‘Vallesian Crisis’ (van Dam, 2006; Casanovas-Vilar and Agustí, 2007). 4.2. Regional events Finally, the large mammal record is consistent with the small mammal record in showing wetter conditions in Ib-N12a-correlative sites as In this section we compare the timing of Iberian rodent transitions to revealed by the analysis of microwear and diet in ruminants (Teruel regional, Mediterranean-wide marine proxy records (SST, planktic δ13C Basin; De Miguel et al., 2019). and δ18O) as well as to our own Iberian paleoprecipitation and lake The recent study of time-equivalent sediments in the Jumilla-La Celia expansion/contraction records. Several important tectonics-related ob­ basin in the northeastern Betics (~7.6–6.9 Ma; van Dam et al., 2014; van servations will be included into this discussion as well. Balen et al., 2015; van Dam et al., 2023; Fig. 1) has shown that large- scale lacustrine conditions did not develop, but that small lakes were 4.2.1. Late Miocene climatic and tectonic events present at 7.6–7.5 Ma and 7.3–7.2 Ma as evidenced by the presence of Because both North Atlantic and Mediterranean SST records fail to carbonate and gypsum deposits (Los Barracones section). A parallel show any clear sign of cooling around 7.6 Ma (Fig. 3A,B) other expla­ section (La Celia section) consists almost entirely of siliciclastics, but the nations for rodent turnover should be explored. Here, we favor a link of excess production of paleosol carbonate reveals at these same times Iberian records to the sudden 1.5–2‰ negative shift in surface water reveals the occurrence of wetter conditions also there. A specific dry δ13C starting at 7.63 Ma and culminating at 7.60–7.58 Ma, as reported event as recorded in the in the Bajo Segura Basin around 7.65 Ma has for northern Italy (Monte del Casino section; Romagnan Appenines; thus far not been encountered, but future study may reveal its existence. Kouwenhoven et al., 1999). According to Kouwenhoven and van der At ~7.65 Ma, the composition of the Ib-N11-12a transitional fauna Zwaan (2006) this shift could well be related to increased runoff from Crevillente 4 in the Bajo Segura Basin differs markedly from typical Ib- the continental (and hence 13C- depleted) hinterland (sensu Gudjonsson N11 as well as Ib-N12a compositions by peak occurrences of Apocricetus and van der Zwaan, 1985), in turn attributed to the emergence of Italian and the Parapodemus gaudryi - Apodemus gudrunae lineage well before and Alpine landmasses. The planktic δ18O record from Monte Casino these forms become regular inhabitants (in Ib-N13b). A low diversity for shows a shift towards more depleted values at 7.6 Ma as well, which, CR4 (only six terrestrial rodents and one eulipotyphlan) points to arid given the lack of any major global SST change at that time, can best be conditions, which is confirmed by a predicted mean MAP of only 220 interpreted as an increase of the contribution of fresh water (e.g., Roh­ mm/yr (Table S3). Consistently, A. mystacinus, the supposed living ling, 1999). An overall increase of precipitation and runoff in northern descendant of the latter lineage (Martín Suárez and Mein, 1998) is Italy coincides with an increase of precipitation in the Teruel Basin, as currently restricted to the drier region of southwestern Europe and the evidenced by the transition from a drier Ib-N11 Zone (average MAP of Middle East. Finally, the unique late presence of crocodiles at this site 475 mm; driest month with zero precipitation) to an overall wetter Ib- (Böhme and Ilg, 2003) further indicates special conditions with tem­ N12a Zone (average MAP of 633 mm, driest month with non-zero pre­ peratures being higher than before or after this time. cipitation in ~60% of the sites; Table S3, Fig. 3C; similar results were Interestingly, the Ib-N11-12a boundary overlaps with the ‘Tortonian 8 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 Salinity’ crisis, a short interval of evaporite deposition in the eastern With the final disappearance of Huerzelerimys and Cricetulodon at the Betics in SE Spain, which has been dated at 7.8–7.6 Ma (Krijgsman et al., Ib-N12a-b transition (7.18–7.09 Ma), Iberian faunal structure starts to 1999). Although no indications for a specific drying event can be echo that of Ib-N11, as evidenced by a low faunal diversity, the re- inferred from the sea surface isotopic records of the Adriatic Sea (which expansion of Occitanomys and the re-occurrence of Parapodemus lugdu­ only indicate the shift to wetter conditions, Fig. 3B), there is recent age nensis (in the Teruel Basin, van Dam et al., 2023). The resemblance of information for a significant drying event (‘Khersonian drying’) between Zone Ib-N12b to N11 is consistent with the very rapid return to higher 8.02 and 7.84 Ma in the Paratethys domain that ends at 7.65 Ma Mediterranean SST temperatures of 24–25 ◦ C at the base of Ib-N12b (‘Maeotian transgression’; Palcu et al., 2019, 2021; Butiseaca et al., (7.13 Ma, and to a lesser degree also at 7.18 and 7.09 Ma, Fig. 3A). 2021) time-correlative with the Ib-N11-12a boundary. Possibly this On the other hand, it should be noted that drying is observed to start drying event was not even restricted to Europe, as suggested by the already high in Ib-N12a (~7.25 Ma) as indicated by the disappearance documentation of a late Miocene low in (wet season) precipitation in in the Teruel Basin (Tortajada, Concud Cerro de la Garita 5, Aljezar B) of China at 7.6 Ma as well (Hui et al., 2023). Whether Eurasian drying was aquatically adapted and arboreal species (with the latter loss translating leading in bringing aridity to the peripheral Iberian Peninsula, or all into zero precipitation in the driest-month, Fig. 3C). These seasonally areas responded to a single external, even more remote event (e.g., dry conditions, which normally apply to the summer in western parts of Antarctic glaciation, Section 4.1) cannot be established with certainty. continents (including the Mediterranean region; Henderson-Sellers and In any case, a tectonic explanation for evaporite deposition in southern Robinson, 1994; Bosmans et al., 2015a), basically persist throughout the Spain (uplift leading to the formation of a submarine sill and partial remaining part of the Miocene (van Dam, 2006). On the other hand, the isolation of marginal basins in the Eastern Betics; Krijgsman et al., 1999) continuing presence of invertivorous groups indicates reasonably high could be supplemented by a climatic explanation. precipitation levels during the rest of the year, as reflected in mean MAP After a long period of stable and high SST levels of ~27–28 ◦ C or values exceeding 500 mm in these higher placed Ib-N12a sites higher between 9 and 7.5 Ma, a surprisingly low SST in the Mediterra­ (Table S3), and also in Zone Ib-N12b (La Roma 5A, Valdecebro 5), where nean Sea (Adriatic region) is recorded between 7.5 and 6.6 Ma (Tzanova values of 600 mm are occasionally reached (not in Table S3, because age et al., 2015). During two successive peaks, SST dramatically drops to uncertainty >200 kyr). ~20 ◦ C, i.e., even lower than at North Atlantic site ODP-982. (A com­ A low level of faunal similarity between Iberian Basins during Ib- parable SST drop has recently been recorded in the eastern Mediterra­ N12b and across the Ib-N12a-b transition (van Dam et al., 2023) nean; Kontakiotis et al., 2022). The first cooling peak starts at 7.28 Ma points to the existence of important barriers, for which mountain ranges and then culminates at 7.19 and 7.14 Ma. After a short reversal, SST constitute the best candidates. The idea of enhanced tectonic activity drops again with strong minima at 6.92 and 6.85–6.84 Ma. Both cooling leading to topographic barriers at this time is supported by several lines phases would qualify as potential triggers for faunal change, with the of evidence, including the intensification of volcanic activity in south­ first event corresponding to the Iberian Ib-N12a-b boundary (between eastern part of the Iberian Peninsula at ~7.7 Ma (Fortuna Basin; Kuiper 7.18 and 7.09 Ma), and the second event falling in the uncertainty range et al., 2006) with a marine to continental transition occurring at 7.6 Ma of the Ib-N12b-13a boundary (7.02–6.90 Ma). Alternatively, the latter (Garcés et al., 2001). Volcanic rocks in association with diapirism of transition could correspond to the preceding SST rebound to 27.0 ◦ C at Triassic gypsum in the Jumilla-La Celia Basin (Fig. 1) have been radio­ 6.97 Ma, or to a combination of the two events. metrically dated at ~7.6–7.2 Ma (Nobel et al., 1981; van Balen et al., The nature of the two Mediterranean cooling peaks centered at 7.2 2015) and magnetostratigraphically dated at 7.6–7.5 (tuffaceous marl and 6.9 Ma is rather enigmatic, as a clear expression in open ocean re­ and basalt) and 7.05 Ma (ash and clasts) (van Dam et al., 2023). It may cords is lacking and both events precede the SST drops in the open ocean not be coincidental that a major phase of diapirism of Triassic evaporites associated with the major phase of the ‘Messinian glaciation’ starting at is recorded in the Teruel Basin as well around this time (Tortajada 6.3 Ma (Fig. 3A). A small SST drop is recorded at ODP-982 at 7.2 Ma, but Formation, van de Weerd, 1976; Alonso-Zarza and Calvo, 2000). Spe­ not at 6.9 Ma, whereas the ODP-907 record shows no sign for strong cifically, we report a transition from limestone to gypsum occurring a cooling events at all. A moderately strong, positive step in the δ18O re­ few meters above the Ib-N12a site Tortajada Pueblo (van Dam et al., cord in ODP-982 occurs at 7.14–7.135 Ma (Fig. 3B), but at this time 2023). Also, other local gypsum bodies in the basin are associated with Mediterranean SST starts to rise again. The relatively low resolution of Ib-N12a sites, for example at Los Mansuetos (Calvo et al., 1999) and Los the samples used for SST reconstruction (especially for site 982) may Aljezares (see: López-Martínez, 1989). Somewhat higher in the sequence play a role, however. (locality Modorras, Zone Ib-N13b or N13c) gypsum and alluvial beds The occurrence of a strong regional cooling during the early Messi­ interfinger. The age of the top of the Tortajada Formation is poorly nian is consistent with the ‘siphon event’ hypothesis (Hodell et al., 1989; constrained but appears to project at least into Zone Ib-N13b (above Benson et al., 1991) involving a current reversal event through the sites Masada del Valle 7 and Regajo 5; van Dam et al., 2001), for which a Rifian Corridor by which cool waters from below the Atlantic thermo­ maximum age of 6.56 Ma was inferred (Table 1). From the viewpoint of cline are injected into the Mediterranean. The planktic record of the general basin formation, however, there are no major indications for a Adriatic Sea may reflect this event (Kouwenhoven et al., 1999) with major tectonic pulse in the Teruel Basin around this time, with Zone various strong positive δ18O spikes at 7.14 Ma, and at 6.97, 6.94 and Ib-N12a positioned in the middle of a major tectonic unit (TN3, Ezquerro 6.85 Ma correlating to the above-mentioned cooling events. The et al., 2020) and relief rejuvenation and alluvial fan development simultaneous occurrence of positive δ13C excursions in the Adriatic re­ occurring ~1 Myr later (TN3–4 transition). In a broader context, how­ cord (Fig. 3B) seems more difficult to explain, however, because the ever, the ~7.6–6.8Ma interval corresponds to a period of increased influx of below-thermocline Atlantic source waters is expected to lead to tectonic activity in the Mediterranean (leading up to the Messinian more negative δ13C values. Here we suggest that stronger, SST- Salinity Crisis; Krijgsman et al., 1999, 2018) and Alps (Valla et al., dependent carbon isotopic equilibrium fractionation effects during air- 2012). sea CO2 exchange played an important role (~0.1‰/◦ C, Zhang et al., Contemporary with a short period of fluctuating Mediterranean SST 1995) with potentially strongly counteracting non-equilibrium effects (Fig. 3A), the Iberian fauna was truly reshuffled at the Ib-N12b-N13a (e.g., see Schmittner et al., 2013) being subdued in the restricted Med­ transition (between 7.02 and 6.90 Ma). New lineages (Stephanomys iterranean setting. With significant productivity increase being ruled out and a second, small-sized Occitanomys) enter, existing lineages re- (Kouwenhoven et al., 1999) the only other explanation for these δ13C expand (Parapodemus gaudryi - Apodemus, Neocricetodon) or rapidly maxima would be a temporary reduction of continental runoff. How­ transform (P. barbarae - Castromys). At this point it is worthwhile to note ever, such a reduction is not expected in an otherwise active tectonic that the improved stratigraphic resolution for the 7.4–6.8 interval (van regime (van der Meulen et al., 1999). Dam et al., 2023) now shows that faunal changes previously considered 9 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 to belong to one event (SM13 event of van Dam (1997), MME 1 event of waterbodies (Fig. 3E). Zone Ib-N14a faunas in the Teruel Basin also Agustí et al. (2006)) are distributed across two phases that correspond to witness the re-appearance of the jumping mouse Eozapus (Dipodidae), a zone boundaries Ib-N12a-b and N12b-13a. form well represented in Central Europe (with occurrences at latitudes The only well-documented fauna in the Teruel Basin just after the Ib- as least as high as the level of Poland), and the temporary presence of a 12b-13a transition (Gloria 6, Ib-N13a) has an associated estimated mean lineage of the murine Micromys, a form especially well represented in the MAP of ~500 mm, i.e., comparable to the preceding upper Ib-N12a and eastern Mediterranean (Greece, Turkey; Hordijk and de Bruijn, 2009). N12b faunas. The fauna nevertheless witnesses the return of the Although these latter two genera are rare in Spain (they belong to one of arboreally-adapted dormouse Muscardinus, pointing to a temporary in­ our rest groups), their entry together with that of Celadensia in the crease in summer precipitation. By contrast, the Ib-N13a fauna Crev­ central region can be regarded as the start of an important series of Plio- illente 22 (age between 6.78 and 6.60 Ma) in the Bajo Segura Basin has Pleistocene Central-Eastern European rodent invasions into the Iberian an estimated mean MAP of only 350 mm, which is lower than that of Ib- Peninsula. (The entry of the rare hypsodont hamster Blancomys sanzi at N12b faunas in that area. As already mentioned in Section 3.2, estimates the Late Miocene Ib-N13a-b transition in Central Spain may be regarded for the more eastern Fortuna Basin (Sifón de Librilla) are extremely low as an early precursor of these events.) reflecting the almost complete absence of insectivorous species. Possibly The Ib-N14a-b, N14b-c and N14c-15a transitions around ~4.9, 4.6 these low levels could be linked to two strong, concerted positive δ13C and 4.2 Ma are not associated with distinct changes in SST at Mediter­ and δ18O excursions at 6.7–6.6 Ma (still within Ib-N13a, Fig. 3A) as ranean latitudes (including the Gulf of Cadiz, Fig. 3A). However, espe­ recorded in Adriatic surface waters (Kouwenhoven and van der Zwaan, cially the latter two transitions are associated with low values in both 2006; Fig. 3A) indicating a reduction of freshwater runoff (aridity). Mediterranean planktic δ13C and δ18O (Wang et al., 2010; Fig. 3B). With an average MAP of ~400 mm/yr, Ib-N13b faunas from the Freshwater input is the most widely accepted process to explain these Teruel Basin are associated with somewhat drier conditions than Ib- intervals with low δ13C (Gudjonsson and van der Zwaan, 1985) and δ18O N12b-13a faunas (sites not shown in Fig. 3C, because of large age un­ (Lourens et al., 1992). These intervals are also associated with an certainties, but see van Dam (2006)). Consistently, Ib-N13b witnesses increased deposition of organic-rich beds (equivalent to younger sap­ the expansion of the Parapodemus gaudryi - Apodemus gudrunae lineage ropels) of which the genesis is explained by an increase in circum- that also peaked at the dry Ib-N11-12a ‘boundary event‘ (Section 3.2). Mediterranean humidity (Rossignol-Strick, 1985). Consistently, the The entry in Ib-N13b of Blancomys, a hamster with robust, hypsodont successive colonizations of the Iberian Peninsula by the vole-like ham­ molars, is also consistent with drying. Although the Mediterranean SST ster Trilophomys and the true vole Mimomys davakosi at the Ib-N14b-c record shows high values until 6.4 Ma (start of a gap in this record, and N14c-15a boundaries fit a humidity-related scenario, as both Fig. 3A), the association of Ib-N13b to warm sea water conditions cannot forms are common at more northern latitudes in both Western and be demonstrated unequivocally, because of the large age uncertainty Eastern Europe where precipitation levels are higher (van Dam, 2006; surrounding the Ib-N13a-b boundary (6.56–6.27 Ma, Table 1). van Dam and Utescher, 2016). Consistent with the absence of a clear MAP values are very low for the Granada Basin around the Ib-N13b-c correlation to Mediterranean δ13C and δ18O (Fig. 3B) and a pollen-based transition (Fig. 3C), which would be consistent with the expansion of the indication for aridity during glacial event Si6 (4.89 Ma, Section 4.1), a terrestrial, supposedly dry-adapted African rodent immigrant Para­ second, larger-sized Iberian lineage of the African-derived potentially ethomys at this time (6.24–6.23 Ma, MME 2 event of Agustí et al., 2006). dry-adapted genus Paraethomys appears and expands at the Ib-N14a-b Values for the more eastern Fortuna Basin are even drier (Sifón de transition (4.91–4.87 Ma). Very low (Fuente del Viso) or moderately Librilla). low (Alcoi 0) precipitation levels as predicted for eastern Spain during Ib-N14b (Fig. 3C) are consistent with these observations. 4.2.2. Pliocene and early Pleistocene climatic events An average MAP of ~330 mm (range: 200–600 mm) characterizes The Ib-N13c-14a boundary (between 5.47 and 5.33 Ma) correlates to the well-dated and well-covered upper Ib-N14c to lower N15b interval the final phase of the Messinian Salinity Crisis between 5.55 and 5.33 Ma in the Teruel Basin. Two estimates for the southern Guadix-Baza Basin (Upper Gypsum/Lago Mare; Krijgsman et al., 2018). Differences be­ (Gorafe 1 in Ib-N14c and Gorafe 2 in Ib-N15b) reveal very low values of tween ‘Late Miocene’ Ib-N13c and ‘Early Pliocene’ Ib-N14a faunas in the 100-200 mm. The entry of the flying squirrel Pliopetaurista in Ib-N15a southeastern basins (Alcoy, Fortuna, Granada Basins) are surprisingly (Villalba Alta Rio 2, 4.16 Ma) indicates the start of period of small, however (Fig. 2A), which seems at odds with the extreme increasing tree density and non-zero precipitation during the driest paleogeographic changes that would accompany a (partly) desiccated month (Fig. 3C). On the other hand, winter rainfall decreases implying a Mediterranean. More to the north (Teruel Basin), on the other hand, reduction in rainfall seasonality. Whereas most of the Ib-N15a localities faunal change between Ib-N13c and 14a is more marked because of the still lack Pliopetaurista, its presence becomes more permanent after 4.1 immigration of Celadensia, the first ‘vole-like hamster’ to enter the Ibe­ Ma (Ib-N15b). A more homogeneous rainfall regime appears to char­ rian Peninsula. Celadensia has more hypsodont and prismatic molars acterize the Iberian Peninsula as a whole; for instance, the Ib-N16a fauna compared to typical hamsters, which can be seen as an adaptation to a Tollo de Chiclana 1B in the Guadix-Baza Basin (which also contains more herbivorous diet and to the contamination by grit and soil asso­ Pliopetaurista) has a comparable precipitation level as the more northern ciated with more open landscapes (Jardine et al., 2012). (This type of sites in the Teruel Basin. dental morphological change foreshadows a trend that is further per­ A new phase of lake expansion in the Teruel Basin that is especially fected in arvicolines, see below.) Celadensia is morphologically similar to well recorded at Orrios (Fig. 3E) starts at ~4.5 Ma and lasts until ~3 Ma. the eastern genus Microtodon as initially defined on material from Inner Based on the more peripherally located sedimentary sequences at Vil­ Mongolia (Fejfar et al., 2011). Its subsequent entry in Greece is dated at lalba Alta Río and Escorihuela, we infer that lake area was maximal at 5.40–5.23 Ma (Hordijk and de Bruijn, 2009), an interval that overlaps 4.27–4.01 Ma (during Ib-N15a), between 3.57 and 3.46 Ma (an interval with the uncertainty range surrounding the entry of Celadensia in Spain. within Ib-N15b, but without localities) and ~ 3.3–3.0 Ma (during Ib- As this entry occurs during the period of Mediterranean desiccation N16a), which is in line with the results of recent studies near Orrios (5.59–5.33 Ma, Krijgsman et al., 1999) dispersal to the Iberian Peninsula and Villalba Alta (Ezquerro et al., 2022; Gao et al., 2022). The fact that may have been eased by the creation of new migration routes. the entry of the vole Dolomys coincides with lake contraction (at 4.0 Ma) Given its Asian affinities, we hypothesize that a preadaptation to seems at odds with its assumed semi-aquatic lifestyle. This latter continental (cooler, more seasonal) conditions allowed Celadensia to assumption, which is based on Repenning et al. (1990) hypothesis of a preferentially colonize the northern, interior and higher-altitude parts of relationship with the living muskrat Ondatra (tribe Ondatrini) has, in our study area (Teruel Basin). Precipitation values in this area during Ib- fact, been doubted by other workers (Martin, 2008). Moreover, loco­ N14a show a small rise, which agrees with the temporary expansion of motor styles are flexible in voles, as evidenced by the existence of a 10 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 Table 2 Iberian rodent zone boundaries with associated orbital configurations, hydrological changes, main murid dispersal events, glacial events and tectonically relevant events. Astrochronozone numbering after Hilgen et al. (2020). Zone Age (Ma) Obliquity Eccentricity Hydrology Dispersal Correlations boundary amplitude events Maximum Minimum (1.2-Myr (2.4-Myr (405-kyr (Precipitation) (Teruel Basin (New murid (Isotopic events, tectonics) cycle) cycle) cycle) lakes) clades) top Ib- min. Med. surface δ13C- N17b 1.90 1.80 min (3 → 2) max (2) max (5) Allophaiomys δ18O similar/ Ib-N17a-b 2.15 2.05 min (3 → 2)? min (2 → 1) min (6 → 5) uncertain MIS 82,78 Ib-N16b- similar/ 17a 2.56 2.54 max (3) min (2 → 1) max (7) uncertain MIS 100–96 min. Med. surface δ13C- Ib-N16a-b 3.06 3.03 min (4 → 3) min (2 → 1) max (8) wetter→drier contraction δ18O Ib-N15b- seasonality Mimomys min. Med. surface δ13C- 16a 3.19 3.18 min (4 → 3) min (2 → 1) max (8) increase expansion hassiacus δ18O min (11 → 10) / max seasonality Ib-N15a-b 4.01 3.92 min (5 → 4) min (3 → 2) (10)? decrease contraction Dolomys Gi20 Ib-N14c- Mimomys min. Med. surface δ13C- 15a 4.27 4.19 min (5 → 4) min (3 → 2) max (11) drier→wetter expansion davakosi δ18O min. Med. surface δ13C- Ib-N14b-c 4.71 4.51 intermediate min (3 → 2) max (12) uncertain expansion Trilophomys δ18O contraction min (13 → after Ib-N14a-b 4.91 4.87 max (5) min (3 → 2) 12) wetter→drier expansion Si4,6 max (14) / Ib-N13c- min (14 → tectonics?→ 14a 5.47 5.33 min (6 → 5) max (3) 13)? drier→wetter Celadensia Messinian salinity crisis tectonics / C3An18O.14,16 → Af.- Ib-N13b-c 6.26 6.25 intermediate max (3) max (16) drier? Eur. landbridge min (17 → 16) / max Ib-N13a-b 6.56 6.27 min (7 → 6) max (3) (16)? wetter→drier contraction Ib-N12b- min (18 → similar/ tectonics / regional 13a 7.02 6.90 max (7) min (4 → 3) 17) uncertain contraction SST anomalies tectonics / regional Ib-N12a-b 7.18 7.09 max (7) min (4 → 3) max (18) wetter→drier similar SST anomalies Ib-N11- min (20 → dry event, MTi4, Paratethys drying 12a 7.68 7.65 min (8 → 7) min (4 → 3) 19) drier→wetter expansion and expansion terrestrial to semi-aquatic spectrum within living species of Arvicola and N17a-17b transitions are associated with low values of Mediterra­ (Nowak, 1999; Chevret et al., 2020). Alternatively, Dolomys may have nean δ13C and δ18O (Fig. 3B), suggesting temporary higher humidity preferred small rivers rather than lake margins. Given the relatively low levels. A small raise of MAP around the Ib-N16b-17a boundary in the levels of wettest month precipitation during Zone Ib-N15b (Fig. 3C), Teruel Basin may also reflect this increase (Fig. 3D). lake reduction between 3.99 and 4.01 Ma could be related to the reduction of winter precipitation, with a certain base level of humidity 4.3. Orbital cycles still secured by (non-zero) summer precipitation. Although more well-dated localities are needed to be certain, pre­ Although many of the climatic events described above find their cipitation levels for Ib-N16b are relatively low. We note that a high ultimate origin in astronomical variations, orbital extremes do not al­ percentage (44%) of (semi-)aquatic micromammals (two voles and two ways lead to proportionally strong climatic extremes (Zachos et al., eulipotyphlans) in Ib-N16b fauna Concud Estación 1 (Fig. 3F) stems 2001a). Because more subtle but persistent orbitally forced changes in from a low richness of the total fauna, with estimated MAP amounting to local insolation patterns, precipitation (including seasonality) and 260 mm only. Precipitation for site Concud Estación 2 is somewhat vegetation are expected to affect rodent communities via ecological- higher (390 mm). The development of an aeolian dune field near scale processes (e.g., reproduction, survival and competition), a direct Escorihuela between 2.9 and 2.6 Ma (Fig. 3E) suggests dry conditions comparison between rodent transitions and astronomical curves may (Rodríguez-López et al., 2012). reveal additional insights into relevant forcing mechanisms. As noted in Section 4.1, the timing of the Ib-N16b-17a and N17a-17b Given the 105–6 year-scale of biozones, our focus is not directed at boundaries (2.56–2.54, 2.15–2.05 Ma) corresponds to that of glacial individual ~21-kyr precession and 41-kyr obliquity cycles, but on events MIS 100 and 78, respectively. It is around 2.5 Ma but especially at longer-scale cycles of eccentricity and obliquity. It should be noted, 1.85 Ma, that the Mediterranean SST record finally shows a permanent however, that short-period cycles may affect proxy results based on drop after having stayed at a stable high level during the Pliocene, fossil content, such as our precipitation proxies. This bias may take the Fig. 3A). This cooling trend around the Plio-Pleistocene boundary form of aliasing: because mammal sites are bound to certain lithologies compares well to trends as inferred from pollen (Garraf core, north- within small-scale cycles (usually organic-rich marl beds deposited near eastern Spain; Fauquette et al., 1998), although these appear to occur aquatic margins), reconstructions for (sets of) sites may be biased to­ slightly earlier (2.9–2.7 Ma). Finally, our data for Ib-N17a and N17b wards either drier (shallow facies) or wetter conditions (deeper facies). faunas in eastern and southern Spain (Júcar and Guadix-Baza Basins) On the other hand, the relatively fixed structure of assemblages as show a continuation of dry conditions without any significant summer observed across entire biozones suggests that this type of bias may not be precipitation. Like most Pliocene zone transitions, also the Ib-N16b-17a extreme. 11 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 4.3.1. The overall pattern subtle obliquity forcing must be assumed, with the sustained absence of The average duration of rodent zones Ib-N11-17b is 0.46 Myr. 41-kyr obliquity maxima during 1.2-Myr nodes playing a more impor­ Looking more closely (Fig. 3), a clustering of zones into three 1.4–1.5 tant role. The approximate time equivalence of the ‘Messinian glacia­ Myr intervals (7.6–6.2, 5.4–4.0 and 3.2–1.9 Ma) can be discerned tion’ (6.3 to 5.5 Ma), our long zone Ib-N13c (6.25 to 5.5–5.3 Ma) and the (Fig. 3D, F). These compound intervals are separated from each other by lower-amplitude obliquity cycle between ~6.5–5.5 Ma (cycle 6, Fig. 3F) long intervals of ~0.8–0.9 Myr corresponding to a single zone only (Ib- suggests obliquity forcing presumably via NH high-latitude climate. We N13c and N15b, at 6.3–5.4, and 4.0–3.2 Ma, respectively). The oldest will now proceed by discussing orbital relationships in more detail per zone (Ib-N11) has a long duration as well (1.1 Myr). Excluding the three epoch. long zones, average zone duration reduces to 0.35 Myr. Given the pre­ viously established correlation of the long 2.4 Myr eccentricity cycle to 4.3.2. Late Miocene orbital forcing Neogene rodent turnover (van Dam et al., 2006), an effect of the next According to classical theory (Milankovitch, 1930), the strongest shorter eccentricity cycle (405 kyr) is to be expected. The 460 and 350 expected correlations between Mediterranean climate and lower- kyr durations as calculated for average zone duration (including and frequency orbital parameters involve eccentricity (as a modulator of excluding long zones) strongly suggest an imprint of this cycle on the precession), because precession rather than obliquity is the dominant timing of the rodent zones. component in insolation at low- to middle latitudes. The temporal cor­ Whereas the three long zones Ib-N11, N13c and N15b are associated respondence of the Ib-N11-12a boundary at 7.68–7.65 Ma with an with strongly fluctuating eccentricity (middle parts of 2.4-Myr cycles extreme in the 405-kyr eccentricity cycle (minimum between cycles 2–4, Fig. 3D; cycle numbers after Hilgen et al., 2020), the compound 20–19, Fig. 3D) thus fits the theory. On the other hand, an additional intervals are associated with intervals with low average eccentricity (2.4 control of obliquity on Iberian hydrology could help to explain the Myr minima, Fig. 3D, Table 2). Several late Miocene zone boundaries occurrence of a dry event at this boundary (as recorded in the Bajo associated with strong turnover, such as the Ib-N11-N12a boundary at Segura Basin, Section 4.2.1) during the rare moment of overlap of a node 7.66–7.65 Ma and the Ib-12b-13a boundary between 7.02 and 6.90 Ma in the 1.2-Myr obliquity cycle (~8.2–7.6 Ma) and minima in the 2.4-Myr are additionally associated with low 405-kyr eccentricity values. eccentricity cycle (~7.7–7.2 Ma; Figs. 3D,F) and 405-kyr eccentricity (Although the uncertainty interval of the latter boundary includes a 100- cycle (centered at ~7.65 Ma). Although a remotely working high- kyr maximum at 6.96 Ma, it is associated with a minimum at the 400-kyr latitude mechanism could have played a role as well at this moment scale). An age of the Ib-N13a-b boundary close to that of the next 405- (SH glaciation, Section 4.1.1), the persistence of weak 41-kyr obliquity kyr eccentricity minimum at 6.5 Ma would fit this pattern, but the maxima until ~7.6 Ma would be consistent with reduced regionally considerable age uncertainty for this boundary (6.56–6.27 Ma, Table 1) (Atlantic)-derived (winter) precipitation (see Bosmans et al., 2020). In makes this correlation uncertain. During the Plio-Pleistocene, the 405- any case, the response to the 405-kyr minimum at ~7.65 Ma (dry con­ kyr eccentricity-biozone relationship appears to be largely reversed ditions as indicated by small mammals) differs from the one during the with most zone boundaries associated with maxima (at 5.4, 4.6, 4.2, 3.2, subsequent 405-kyr maximum (cycle 19), as well as that of the next 3.0, 2.5, 1.8 Ma) and only two (at 4.9 and 2.1 Ma) to minima. (The minimum at ~7.25 Ma (both with relatively wet conditions; Fig. 3D, E). correlation at 4.0 Ma is uncertain.) It is thus conceivable that stronger obliquity variation between 7.5 and Van Dam et al. (2006) noted a deviating Pliocene pattern of orbital 7.0 Ma (cycle 6, Fig. 3F) is preventing average conditions from becoming forcing compared to the rest of the Neogene, with subdued rodent too dry in this interval. response and lake reduction instead of expansion during the 2.4-Myr The observed increase in humidity between 7.5 and 7.0 Ma cannot be eccentricity minima between ~4.9–4.4 Ma. Here, we suggest that this linked to the eccentricity record in a simple and straightforward way. Pliocene ‘reversal’ in rodent turnover is due to a shift in the response to The suppression of the ~100-kyr cycle by the 2.4-Myr eccentricity cycle the 405-kyr cycle, with maxima not associated with central parts of 2.4- in this interval leaves the 405-kyr cycle as the only eccentricity cycle Myr minima turning into hydrologically critical intervals as well. As with a clear expression (cycle 19, Fig. 3D). From the viewpoint of further discussed below, this phase reversal could best be explained by a insolation, the ‘damping’ of precession-related climatic extremes under slowly proceeding underlying change in the average global climate (late low eccentricity leads to the prolonged absence of both extreme warm Cenozoic cooling) resulting in an equatorward migration of the NH summers and cool NH winters (during high-seasonality precession temperate belt. minima) alternating with the absence of cool NH summers and warm NH In addition to eccentricity, the distribution of biozone ages bears the winters (during low-seasonality precession maxima). The prediction of imprint of obliquity. The long ~0.8 Myr zones are associated with a such ‘intermediate’ conditions seems inconsistent with the observed specific, low-amplitude expression of the 1.2-Myr obliquity cycle (cycles increase in precipitation and lake area during in the Teruel Basin this 6 and 4; Fig. 3F). The oldest zone (Ib-N11, 1.1 Myr duration), corre­ interval (Section 4.2.1, Fig. 3D). sponds to such a low-amplitude obliquity cycle as well (cycle 8). (The The emergence of a distinct climatic signal during an interval that is relationship of these long zones with both long-period eccentricity and otherwise orbitally ‘balanced’ implies the occurrence of non-linearities obliquity is not surprising, as these long cycles are astronomically in the interaction between climatic and/or biotic/sedimentary compo­ linked; Laskar et al., 2004). As mentioned above (Section 4.1), zone nents. A well-known high-latitude example of non-linear climate-related boundaries and high-latitude glacial events are coeval twice during the response is the increase of polar ice volume at long-period orbital nodes, Pliocene (Si6 at 4.9 and Gi20 at 4.0 Ma) and twice during the Early which depends on the reduction of summer melting rather than winter Pleistocene (MIS 100 at 2.5 and MIS 78, 82 at 2.1 Ma). During this latter growth during strong ~21-kyr precession minima (NH) or maxima (SH) interval NH ice sheets became permanent with obliquity becoming the or ~ 41-kyr obliquity maxima (Zachos et al., 2001b; Lourens and Hilgen, dominant orbital parameter controlling ice volume and, in turn, global 1997). Other examples include lake expansion associated with 2.4-Myr climate (Ruddiman et al., 1989; Lisiecki and Raymo, 2005). Pliocene eccentricity minima or maxima (van Dam et al., 2006), minima (Abels events Si4–6 and Gi20 are associated with obliquity nodes, whereas MIS et al., 2009, 2010), maxima (Valero et al., 2014) and peatland expansion 100 is associated with a strong obliquity antinode. This latter configu­ during 2.4-Myr minima (Valero et al., 2016). ration with obliquity varying maximally around 2.5 Ma (1.2-Myr cycle Model results with near-modern boundary conditions (for tempera­ 3, Fig. 3F) figures in the hypothesis that especially strong 41-kyr ture, CO2, topography, etc.) predict net more humid conditions during obliquity minima associated with reduced summer insolation accel­ eccentricity maxima associated with strongly expressed precession erate NH polar ice growth at this time (non-linear response; e.g., Lourens minima (NH), with positive effects on both Mediterranean summer and and Hilgen, 1997; see also Section 4.3.2). For earlier periods (late winter precipitation (Kutzbach et al., 2014; Bosmans et al., 2015a). The Miocene and Pliocene), on the other hand, a different type of more observed precipitation increase and lake expansion during low- 12 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 eccentricity interval Ib-N12a (Fig. 3D) seems at odds with these results. As discussed in Section 4.2.1, the age of the top of the ‘Páramo 1’ at Here, we propose that it is the average global climate at the time of study Concud in the Teruel Basin (between 7.05 and 7.00 Ma; van Dam et al., (during coolhouse - icehouse climatic evolution; Westerhold et al., 2023) post-dates the first strong 100-kyr eccentricity maximum at 7.16 2020) that controls the amplitude and phase of sedimentary and hy­ Ma after the low-eccentricity interval 7.7–7.2 Ma by at least 110 kyr drological responses to long-term orbital cycles. Before ~7 Ma global (Fig. 3F). A similar delay in timing of large-scale (formation-level) li­ temperatures were still between 5 and 10 ◦ C higher on the Iberian thologies with respect to long-period eccentricity has been described for Peninsula (van Dam and Reichart, 2009; Herbert et al., 2016), which older parts of the basin (Cascante area), with a 400-kyr/2.4-Myr ec­ was positioned well into the drier subtropical belt (average MAP in the centricity minimum at 9.7 Ma leading the establishment of full lacus­ studied Iberian Basins rises from ~450 mm in Ib-N11 to ~600 mm in Ib- trine conditions by ~130 kyr (Abels et al., 2008). Comparable lags have N12a; driest-month precipitation increase from zero to non-zero values; been detected with respect to the two next older (middle Miocene) long- Fig. 3D) south of a much wetter (warm temperate) central Europe (van term eccentricity minima (Madrid Basin; Abels et al., 2010). In this latter Dam, 2006). With temperatures in Spain and western Europe remaining case, lags in the timing of large-scale lithological changes with respect to high until ~7 Ma, a decrease in summer evaporation during weakly orbital parameters have been explained by a ‘memory effect’ in such a expressed precession minima would help in explaining the existence of a way that groundwater reservoirs must first be filled before a lake can more positive water budget. In addition to this potential decrease in form (Abels et al., 2010). Likewise, it will take a considerable time for evaporation, driest-month precipitation (LMP) estimates increase the lake to finally disappear (e.g., at 7.0 Ma) because reservoirs will still (reflecting the expansion of the arboreally adapted Csakvaromys and to a feed the lake long after the surface water budget has turned negative. On lesser degree Muscardinus). The mismatch between our proxy results and the other hand, if obliquity is believed to be the dominant parameter model results based on near-modern boundary conditions results could affecting lake area, there would be less need to invoke a ‘memory effect’- imply that the latter need some modification in order to properly model type of explanation given the approximate match of the maxima in lake much warmer intervals such as the Miocene and earlier periods (see also area and obliquity amplitude (Fig. 3F). Herold et al., 2010). The broad age uncertainty interval of the next transition (Ib-N13a-b, An alternative orbital explanation for the observed wetter conditions 6.56–6.27 Ma) includes both a 405-kyr eccentricity minimum and between 7.6 and 7.2 Ma involves obliquity. With this interval corre­ maximum, making the correlation uncertain. The fifth transition (Ib- sponding to the middle part of the well-expressed 1.2 Myr obliquity N13b-c) coincides with a 405-kyr eccentricity maximum (cycle 16) but cycle (cycle 7, ~7.9–6.7 Ma, Fig. 3F), strong obliquity maxima should differs from all other transitions by the existence of a unique migration then be held responsible for an increase of Mediterranean winter pre­ route between Africa and Europe. The terminal Miocene Ib-N13c-14a cipitation, in line with some recent studies (Bosmans et al., 2015a, boundary uncertainty interval (5.47–5.33 Ma) includes a 100-kyr ec­ 2015b). Because it can be assumed that winter rainfall is quantitatively centricity maximum at 5.43 Ma, but also extends into the first of three more important than summer rainfall in filling groundwater reservoirs, a 405-kyr minima (between cycles 14–13, 13–12 and 12–11 in Fig. 3D) link between obliquity and long-term lake expansion/contraction may that accompany/fringe the 2.4-Myr minimum centered at 4.7 Ma. In this be hypothesized. (On shorter time scales, however, precession would respect, the orbital configuration is reminiscent of the one at 7.65 Ma still be the orbital parameter controlling shallow lake level variations; e. (Ib-N11-12a transition), which likewise is associated with a minimum in g., van Dam et al., 2023). An obliquity-dependent scenario would be 1.2-Myr obliquity amplitude variation (transition from cycles 8 to 7 in consistent with the clear expression of the 1.2-Myr obliquity cycle in Fig. 3F). records from the Chinese loess plateau between 7.5 and 4.5 Ma (Qin et al., 2022: fig. 9), where a higher magnetic susceptibility is associated 4.3.3. Pliocene orbital forcing with wetter conditions during obliquity antinodes (see Gao et al., 2022, Compared to the late Miocene, the Pliocene shows a dominantly and reference therein). As further discussed below (Section 4.6), pacing reversed phase relationship between Iberian rodent transitions and ec­ by the 1.2-Myr cycle would then be temporally decoupled from forcing centricity, with zone boundaries at 4.6, 4.2, 3.2 and 3.0 Ma associated by regular 41-kyr cycles, which, before the Messinian glaciation (van der with high eccentricity. Three of these boundaries (at 4.6, 4.2 and 3.0 Laan et al., 2012), are mainly restricted to intervals of low eccentricity Ma) correlate to centers of 405-kyr eccentricity maxima (where the 100- (2.4-Myr minima) in marine (e.g., Hilgen, 1991; Lourens et al., 1996) kyr cycle has a maximum expression; cycles 12, 11 and 8; Fig. 3D), while and terrestrial settings, including the Teruel Basin (Abdul Aziz et al., one boundary (3.2 Ma) correlates to the second strongest 100-kyr 2004; Abels et al., 2008). maximum within 405-kyr maximum 8. (The other strong 100-kyr The next zone boundary (Ib-N12a-b, 7.18–7.09 Ma) post-dates the maximum correlates to the 3.0 Ma boundary.) A similar configuration preceding eccentricity minimum at 7.23 Ma by 0.05–0.14 Myr. Putting could characterize the Ib-N15a-b boundary at 4.0–3.9 Ma (association aside alternative, tectonics-related explanations (Section 4.2.1), the with the 100-kyr maximum at 3.93 Ma), although in this case the timing of the Ib-N12a-b boundary could thus theoretically be affected by boundary age uncertainty is too large to exclude an alternative corre­ the resumption of strong 100-kyr eccentricity cycles associated with lation to the 100/405-kyr minimum at 3.99 Ma (Table 1). In contrast to decreasing humidity after a long interval of low eccentricity corre­ the first four boundaries mentioned above, the Ib-N14a-b boundary at sponding to the transition between 2.4-Myr cycles 4 and 3. (Obliquity 4.9 Ma shows a late Miocene-type of association to a 405-kyr/2.4-Myr variation is still large during Ib-N12b.) eccentricity minimum. As mentioned in Section 4.3.1, the Ib-N12b-N13a boundary The Pliocene zone boundaries at 4.6, 4.2 and 3.0 Ma that correspond (7.02–6.90 Ma) is again coeval with a 405-kyr minimum (between cycle to 405-kyr eccentricity maxima 12, 11 and 8 also correlate to minima in 18 and 17, Fig. 3D). Except for a strong reduction of lake area around both Mediterranean planktic δ13C and δ18O, for which, in turn, signifi­ this time, no drying event based on fauna like the one at ~7.65 Ma could cant coherences with the 405-kyr eccentricity cycle were demonstrated be found. A low chronological resolution of micromammal sites could (Wang et al., 2010; Fig. 3B). The inference of wetter conditions during play a role here. On the other hand, the weakly seasonal climate asso­ eccentricity maxima as inferred from the marine record (Section 4.2.2) ciated with the negative eccentricity ‘spike’ at 6.90 Ma could temporary fits the observation that rodent immigrants entering Spain during these have led to a relaxation of one or more ecological barriers between re­ transitions originate from more northern, temperate latitudes in both gions, allowing new lineages (Stephanomys, Castromys, Occitanomys Western and Eastern Europe (Trilophomys, Mimomys davakosi, alcalai) to take over dominance. The second (potentially unrelated) M. hassiacus and Kislangia at 4.6, 4.2, 3.2 and 3.0 Ma, respectively, and Mediterranean SST cooling as recorded by Tzanova et al. (2015) at 6.92 possibly Dolomys at 4.0 Ma, although the phase relationship with the Ma (and later, at 6.85–6.84 Ma, Section 4.2.1) could have boosted the 405-kyr cycle is uncertain; Fig. 3, Table 2). The association of these vole- faunal instability generated by orbitally forced changes. like hamsters and voles with relatively humid conditions is supported by 13 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 a comparison of time-equivalent lacustrine and karstic localities in Spain and a 405-kyr eccentricity maximum (cycle 5, Fig. 3D) as well as with a (van de Weerd and Daams, 1978), showing that these groups abound in distinct Mediterranean δ13C-δ18O depletion event (Table 2). The the former type of environments, while true mice have larger abun­ occurrence of more humid conditions associated with these older Plio­ dances in upland sites. The correspondence between 100-kyr and 405- cene transitions also characterizes the 1.9–1.8 Ma interval, which at kyr eccentricity maxima and higher Pliocene humidity is supported by least partly corresponds to the ‘Tiglian’ interglacial (Zagwijn, 1992); top recent work in the Teruel Basin on magnetic susceptibility variation dated at ~1.8 Ma; Pillans and Gibbard, 2012). Immediately afterwards, (4–3.6 Ma; Gao et al., 2022) and carbonate oxygen isotope ratios SST starts to fluctuate strongly around 1.85–1.8 Ma in the North Atlantic (Ezquerro et al., 2022), respectively, although in the latter study tem­ and particularly in the Mediterranean, reaching peaks significantly perature is expected to play a role as well. Finally, our interpretation of cooler than before (Fig. 3A). Because precise age control is lacking, it is the climatic effects of the 405-kyr cycle during the Pliocene is consistent uncertain whether the immigration of rootless voles (including Allo­ with marine and terrestrial cyclostratigraphic studies in the Gulf of phaiomys) that defines the top of Ib-N17b occurred in association with Cadiz (Sierro et al., 2000) and Guadix-Baza Basin (Agustí et al., 2001) relatively warm or cold conditions (for a discussion: see Repenning that indicate overall drier conditions in southern Spain during Pliocene et al., 1990). Regardless of the climatic mode, the evolutionary novelty 405-kyr minima. of ever-growing molars has allowed these voles to increase the per The correlation of slightly enhanced precipitation rates and espe­ lifetime consumption of calory-poor bulk food and to be well prepared cially lake expansions in the Teruel Basin at 4.3–4.0 and 3.3–3.0 Ma (as for a survival in the harsher environments that accompany Pleistocene indicated by distribution maxima of the ‘Páramo 2’ lithostratigraphic glacials. unit of Godoy et al., 1983; Unit IV of Alcalá et al., 2000; upper part TN4- lower part TN5 of Ezquerro et al., 2020; Fig. 3F) to periods of minimum 4.4. Synthesis: Biozone boundaries and dominant forcing mechanisms obliquity amplitude variation (1.2-Myr cycle) once more suggests a control of obliquity on the regional water cycle (Fig. 3D,F). As in the case The previous discussion shows that the timing of Iberian rodent zone of eccentricity, the phase relationship with lake level is opposite to the boundaries can best be explained as the result of a series of changes in ones before that time, when lacustrine facies complexes reached their regional/Mediterranean humidity. Most of these changes are astro­ maximum size during periods of maximum obliquity amplitude (be­ nomically controlled and can be linked to extreme configurations of the tween ~7.5–7.0 Ma (Fig. 3F) and between ~9.7–9.0 (van Dam et al., 405-kyr and 2.4-Myr eccentricity cycle and the 1.2-Myr obliquity cycle 2001; Abels et al., 2008)). The reversal in the expression of the 1.2-Myr (Table 2). Long biozones (i.e., ‘no change’, Ib-N11, N13c and N15b) obliquity cycle between the Miocene and Pliocene fits recent results occur during 2.4-Myr eccentricity maxima that concur with weak from the Chinese loess plateau (Qin et al., 2022), where a typical late maxima in obliquity variation (8, 6, 4; Fig. 3F). During the late Miocene Miocene correlation between 1.2-Myr obliquity antinodes and wet and earliest Pliocene (up to ~4.8 Ma) eccentricity appears to be the conditions (as indicated by magnetic susceptibility) disappears around dominant parameter, with zone boundaries restricted to 405-kyr 4.5 Ma. The recovery of individual 41-kyr obliquity cycles in magnetic minima at 7.7, 6.9 and 4.9 Ma, and possibly also 6.5 and 5.4 Ma. susceptibility between 3.3 and 3.18 Ma in the Teruel basin (Gao et al., These five boundaries are all positioned within broad 2.4-Myr eccen­ 2022) fits observations of a globally increasing influence of obliquity tricity minima (characterized by weak 100-kyr eccentricity cycles) be­ during the Pliocene (Lisiecki and Raymo, 2005) and of stronger obliq­ tween cycles 4 and 3 (8–6.5 Ma) and between 3 and 2 (after 5.3 Ma; uity expression during low eccentricity intervals (Hilgen, 1991; 405-kyr Fig. 3D). The role of obliquity seems subsidiary at this time; although the minimum between cycles 9 and 8, Fig. 3D), but may not have been 7.7, ~6.5 and 5.4 Ma zone boundaries coincide with 1.2-Myr nodes, no understood fully yet, as 41-kyr obliquity variation itself is at a minimum boundaries occur at the onset or middle of the ~8.2–7.7 Ma node (cycle at this time (Fig. 3F). 8–7) or onset or middle of strong variation in the next antinode at 7.5–7.3 Ma (cycle 7). The typical climatic condition associated with late 4.3.4. Early Pleistocene orbital forcing Miocene 405-kyr and 1.2-Myr minima is one of reduced humidity. Comparable to the early Pliocene N14b-c boundary at ~4.6 Ma, the After ~4.8 Ma, faunal change is still associated with low eccentricity Ib-N16b-17a boundary at ~2.55 Ma is sitting in the middle of a 2.4-Myr as dictated by 2.4-Myr eccentricity minima (weak 100-kyr cyclicity eccentricity minimum (between cycles 1–2) close to a 405-kyr maximum between 5.4 and 4.0 and 3.2–1.8 Ma) but the phase relationship with the (cycle 7, Fig. 3D). The boundary is also associated with maximum 405-kyr cycle is reversed, with zone boundaries at 4.6, 4.2, possibly 3.9, obliquity variation (1.2-Myr obliquity cycle 3, Fig. 3F). Strong 41-kyr 3.2, 3.0 Ma and 1.8 Ma associated with (wetter) maxima instead of obliquity minima associated with reduced summer insolation are (drier) minima. The 405-kyr eccentricity cycle thus affects the Iberian considered critical moments for NH polar ice growth (e.g., Lourens and hydrological cycle in a consistent way at least from the late middle Hilgen, 1997). A strong 41-kyr obliquity minimum of 22.06◦ at 2.49 Ma, Miocene onwards with dry conditions occurring during minima (Sierro which slightly postdates the Ib-N16b-17a transition (2.54–2.56 Ma), et al., 2000; Abdul Aziz et al., 2003; Hui et al., 2023). However, the correlates to MIS 100, which, together with the ensuing MIS 98 and 96 correlation between Plio-Pleistocene zone boundaries and obliquity is represent the first set of strong δ18O events associated with NH glacia­ much more pronounced than for late Miocene zone boundaries, with tion. The occurrence of strong glacial events during intervals of each of the three 1.2-Myr nodes around 4, 3 and 2 Ma corresponding to maximum obliquity variation seems to be foreshadowed by the occur­ two separate zone boundaries, and the strong antinodes at ~5.0 and ~ rence of events Si4–6 at 4.9 Ma during the previous strong obliquity 2.5 Ma associated with a zone boundary as well. Unlike the 405-kyr antinode (cycle 5, Fig. 3F), when NH ice sheets were still in an eccentricity cycle, but like the 2.4-Myr eccentricity cycle, the 1.2-Myr ephemeral state (Zachos et al., 2001a). obliquity cycle changes its phase relationship with hydrology around The orbital configuration at the subsequent Ib-17a-17b boundary the Miocene-Pliocene transition. (between 2.15 and 2.05 Ma) reiterates orbital configurations at 6.95 and Although we regard climate change as the dominant agent in the 4.0 Ma, with the boundary age uncertainty interval including a 405-kyr pacing of terrestrial fauna, tectonic changes will additionally have eccentricity minimum at the younger fringe of a 2.4-Myr minimum. affected the communities at times (Table 2). Especially during the Here, we consider a correlation to the strong 100/405-kyr minimum Tortonian-Messinian transition (7.5–6.7 Ma), major tectonic activity between cycles 6 and 5 at 2.07 Ma more probable than the correlation to (uplift and volcanism related to the alpine orogenesis; Krijgsman et al., the preceding, sub-maximal 100-kyr maximum at 2.12 Ma (Fig. 3D). 2018; Valla et al., 2012; Section 4.2.1) will have created environmental The orbital configuration at the top of Ib-N17b (1.9–1.8 Ma), on the gradients and barriers preventing faunal exchange, subsequently lead­ other hand, resembles that of older Pliocene boundaries, such as Ib- ing to increased endemism and extinction (e.g., Zone Ib-N12b, van Dam N14c-15a and N16a-b in its association with a 1.2-Myr obliquity node et al., 2023). In addition, new physical migration routes were created as 14 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 a result of tectonic activity in the western peri-Mediterranean region, i.e., before Antarctica was more or less fully glaciated) or NH ice sheets facilitating faunal exchange between Africa and Iberia (e.g., Para­ (latest Pliocene-Pleistocene) were very dynamic (Zachos et al., 2001a; ethomys entry at base Ib-N13c) or possibly across desiccated areas Rohling et al., 2021). During the early-middle Miocene, the effects of (Celadensia at base Ib-N14a during the Messinian Salinity Crisis). large fluctuations of the Antarctic ice sheet could have affected Northern Here, we suggest that the shift in key orbital configurations affecting Hemisphere climate via northward climate belt shifting (e.g., John et al., the fauna is linked to the late Neogene cooling trend controlling the 2003; Holbourn et al., 2010), with the most extreme responses occurring position of the subtropical high-pressure zone in between the low- during the Mi glacial/cooling events. Low-amplitude eccentricity and pressure tropical and temperate zones. With high-latitude climate forc­ obliquity variation around 7.6 Ma may thus have been one of the last ing becoming more and more important during the Neogene (West­ expressions of this typical Oligo-Miocene ‘teleconnection mechanism’ erhold et al., 2020; De Vleeschouwer et al., 2020) the Mediterranean (MTi4 glacial event; Westerhold et al., 2005; Fig. 3B). region is expected to undergo a net transition from a subtropical warm and dry Miocene climate via a more temperate and less arid Pliocene 4.5. Stepwise clade displacement climate to a cooler Pleistocene climate characterized by glacial (dry)- and interglacial (less dry) alternations (e.g., Tzedakis, 2007). Our From a macroevolutionary point of view, the fossil record consists of observed shift from 405-kyr eccentricity minima to maxima as key phylogenetic units (clades) that are continuously replacing each other moments for faunal change correlates in time to the reversal in phase (e.g., see: Fraser et al., 2021). Competing views on the underlying pro­ around 4.5 Ma between global deep-sea δ13C and δ18O at the 100 and cesses range between passive mechanisms rooted in the opportunistic 405-kyr eccentricity scale as observed for the Atlantic Ocean (transition diversification of surviving groups after extinction (‘clade replacement’; between climate – carbon cycle intervals II-III of De Vleeschouwer et al. Gould and Calloway, 1980; Benton, 1996) to active mechanisms rooted (2020)). This reversal is attributed to a shift from net terrestrial carbon in ecological-time, inter-species competition (‘clade displacement’; storage to net loss during eccentricity minima as the result of the global Sepkoski, 1996; Kemp, 2006). To judge whether clades competed or not, increase of less productive high-latitude biomes. In line with this sce­ specific test criteria, such as stratigraphic overlap, inverse correlation in nario of global terrestrial productivity decrease and net aridification diversity and abundance (‘double wedge’), and utilization of common (with most change occurring at high latitudes and mid-continent areas), resources can be used (Benton, 1996; Krause, 1986). subtropical belts (which include the Iberian Peninsula) may have In our case, true mice and hamsters (indicated as ‘clade group A’) are become temporarily more humid (Pliocene), because of the equatorward replaced by vole-like hamsters and voles (indicated as ‘clade group B’, movement of the temperate belt. Fig. 2A,B). A and B share two major ecological properties (e.g., Nowak, At first sight, the correspondence of long periods of strongly fluctu­ 1999) that allow us as to regard these groups as potential competitors: 1) ating 100-kyr eccentricity (during 2.4-Myr eccentricity maxima; cycles a similar locomotor style (terrestrial) and body size, although also 4, 3, 2; Fig. 3D) to long biozones (Ib-N11, N13c, 15b) is counterintuitive. aquatic adaptations occur in B, and 2) a dominantly vegetarian diet, On the one hand, obliquity could play a role, because extreme 41-kyr although A is more specialized towards granivory and B towards her­ minima or maxima are not reached during these intervals (weak 1.2- bivory. In terms of preferred climate, A is optimally adapted to the warm Myr amplitude cycles 8, 6, 4; Fig. 3F). However, the question remains and dry subtropics, and B to the cooler and wetter temperate zone why strong 100-kyr eccentricity cycles and associated precessional (Pliocene) or colder and drier boreal zone (Pleistocene). The observed changes had no disruptive effects on the fauna, despite a larger contri­ stratigraphic overlap in combination with negatively correlated di­ bution of precession to insolation compared to obliquity at mid-latitudes versity/abundance patterns of the four major murid groups (Murinae, (Milankovitch, 1930; Berger, 1978). The lack of such response either Cricetinae, Arvicolinae and ‘microtoid cricetids’, Fig. 2A,B) points to means that environments did not change on shorter time scales of 100 active competitive displacement. kyr or less (which is hard to believe and disproven by the literature on But how is competition scaling down to the level of individual line­ astronomically forced climate) or that species were somehow resilient to ages? On the one hand, intra-clade competition is expected to be such shorter-lasting climate swings. stronger than inter-clade competition because members of the same In the case of rodents, we previously suggested that buffering by clade are closer phylogenetically. On the other hand, members of a short-distance (within-bioprovince) migration and recurrent amalgam­ newly arriving clade may possess key innovative features that work to ation of incipient species with mother species could prevent long-term their evolutionary-ecological advantage. In our case we suggest that a change (van Dam et al., 2006). Here, we extend this theoretical frame­ slow environmental change as ultimately guided by late Cenozoic work by assuming that immigration of new species from other bio­ cooling constituted the crucial background process driving clade provinces constitutes the key disruptive factor in explaining faunal displacement. During this gradual change, clade group B replaces group change at longer time scales. In this view, successful dispersal is only A in a stepwise fashion, with the steps being dictated by extremes of possible as a result of environmental changes that have sufficient long-period Milankovitch cycles. The new group B takes over, because magnitude and duration across areas sufficiently large and homoge­ its clade members are better adapted to the changing conditions that neous (in terms of vegetation) in order to be crossed successfully. We more and more resemble those in their area of origin. postulate that only longer cycles such as the 405-kyr and 2.4-Myr ec­ Due to the persistence of a warm and dry Miocene climate and a lack centricity cycle and the 1.2-Myr obliquity cycle lead to the creation of of suitable migration corridors, well-adapted resident groups such as such semi-permanent migration corridors that can be crossed by new true mice and hamsters could retain their dominant role in the late immigrant species in large enough numbers in order to successfully Miocene Iberian fauna. Because of their dominance, most lineage events outcompete the sitting species. could only take place within these groups during any upcoming The relative insignificance of high-latitude forcing of the Iberian perturbation. When paleoclimatic conditions changed in the Pliocene fauna during the latest Miocene and parts of the Pliocene should not be (with mid-latitude Atlantic sea surface temperatures dropping from 27 seen as evidence of for an absence of such forcing in older periods. For to 28 to ~20 ◦ C; Fig. 3A), southern Europe climatically ‘migrated’ out of instance, most middle to late Miocene oxygen isotope (‘Mi’) events the warm subtropics into a mixed subtropical dry - temperate wet zone. primary have been correlated to extreme points of long-period cycles in We suggest that during humid 405-kyr eccentricity maxima and 1.2-Myr obliquity and eccentricity via SH glaciation (e.g., Lourens and Hilgen, obliquity minima (Table 2) new, periodically functioning migration 1997; Turco et al., 2001; Zachos et al., 2001b) and to Iberian rodent corridors were created across central and southern Europe. These routes turnover (van Dam et al., 2006). Our observations confirm previous acted as ‘filters’ across areas that previously were too arid to be suc­ ideas that the latest Miocene - Pliocene constituted a special interval in cessfully traversed by northern and eastern forms (see Brown and between a period when either SH ice sheets (Oligocene to late Miocene, Lomolino, 1998, for a discussion on dispersal barriers). The immigration 15 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 Fig. 4. Correlation between European Mammal Neogene (MN) and Quaternary (MQ) units and orbital parameters (obliquity, eccentricity; Laskar et al., 2004; astrochronozone numbering after Hilgen et al., (2020)). MMCT: middle Miocene climate transition. (A) Obliquity record. (B) MN units (Mein, 1975, 1990) and MQ units with boundaries corresponding to selected First Historical Appearance events (Hilgen et al., 2012). (C) Eccentricity record. waves of eastern European microtoid hamsters and voles during the precedes immigration and radiation (the ‘incumbent replacement’ Iberian Pliocene (Celadensia, Trilophomys, Mimomys davakosi, Dolomys, model of Rosenzweig and McCord, 1991), does not apply very well to Mimomys hassiacus) could be explained this way. (As noted above, our case. This model, which functions particularly well in the context of Mediterranean desiccation may have additionally helped Celadensia in mass extinctions, may not be very suitable to explain ‘mainstream’ (sub) reaching Spain). Several of these waves have the actual form of an in­ family-level replacements that shape the fossil record at the 106–7- vasion, i.e., as an abundance peak event followed by a decline (for million-year scale. In our model, ‘incumbency’ (i.e., the power of ‘simply instance, Mimomys davakosi, Dolomys adroveri, M. hassiacus; Fig. 2C, being there first‘) is still an important concept, however. At the lineage Section 3.1.2). As one would expect, central parts of the Iberian level, this concept applies to the long intervals of rapidly oscillating Peninsula were colonized earlier and more extensively than the southern ~21-kyr precession and 41-kyr obliquity, during which successful im­ parts by most of these immigrants. migrations of new forms largely fail. At the clade-level, incumbency With global temperatures dropping another time in the early Pleis­ applies to the middle part of our group A clade history, when most tocene, Iberian voles were well prepared for these new conditions, given lineages involved in replacements simply belong to the same clades (true their area of origin (Asia) with its continental climate and low (winter) mice and hamsters). Here, “privileges of incumbency (including sub­ temperatures. Having become the new main residents in the Iberian stantial geographic range, large numbers of organisms and multiple Peninsula, evolutionary/radiational events started to outweigh immi­ species, heavy use of resources preventing new taxa from gaining a grations also in the arvicoline clade, as evidenced by the development of foothold) buffer existing taxa against the slowly changing conditions.” endemic lineages (Mimomys medasensis, M. realensis, and Kislangia gusii). (Gilinsky and Bambach, 1987: p444). However, in the post-bloom part The short interval of endemic vole evolution corresponding to Zones Ib- of our clade A (Plio-Pleistocene), geographical range boundaries are N16b to 17b, was terminated by another immigration event that was breached during critical phases (400-kyr eccentricity maxima; 1.2-Myr part of a major expansion of a new type of voles (rootless voles, Allo­ obliquity minima) and faunal composition changes. Metaphorically, phaiomys event) across the Holarctic at 1.9–1.8 Ma. one could say that it needs multiple ‘sieges’ at the right moments to The above model can easily be extended to other terrestrial mam­ finally ‘conquer’ the ‘enemy’s’ territory. As Gilinski and Bambach note mals, as it only requires two groups adapted to different climate zones, a from the perspective of the conquered clade: “Thus, as time passes, the climatic trend, and ever-present Milankovitch cycles. For instance, a group’s genetic and developmental system becomes increasingly comparable case of clade ‘displacement’ (i.e., active replacement) has removed from the context in which it arose, and increasingly fewer been described involving the substitution of pelycosaurs by therapsids branching events lead to successful establishment of new species and during the Permian (Kemp, 2006). Although Milankovitch cycles (which families”. also existed hundreds of million years ago with short as well as long periods; Olsen and Kent, 1999) were not included in that study, resource 4.6. Astronomical pacing of European Mammal Neogene (MN) units utilization and shifting climate zones are thought to have played a crucial role. Fig. 3 suggests an astronomical origin for larger-scale bio­ An alternative clade replacement model, in which extinction chronological units as well, with MN12, MN14 and MN16–17 16 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 corresponding to periods of low average eccentricity (2.4-Myr cycle) and Mediterranean humidity and vegetation. Pacing of larger-scale, Neogene strong obliquity amplitude maxima (1.2-Myr cycle), and MN11, MN13 European mammal chronological (MN) units by million-year-scale ec­ and MN15 associated with high average eccentricity and weak obliquity centricity and obliquity cycles is strongly suggested as well. Phase re­ amplitude maxima. lationships between orbital parameters and environment/fauna are not When the complete set of Neogene MN units is inspected (Fig. 4), the constant, however, but depend on the state of the background climate. association with long-period orbital cycles becomes even more evident. With climate belts moving equatorward during gradual late Neogene Whereas the sequence of early to early middle (pre-14 Ma) Miocene cooling, more and more rodent forms with a more temperate / higher- units MN1–5 temporally correlates to 2.4-Myr eccentricity cycles 10–7 latitude origin (‘microtoid hamsters’ and voles) are successively intro­ on an almost one-to-one basis (the placement of the short MN4 between duced into the existing Mediterranean fauna dominated by true mice the longer MN3 and MN5 is an exception), the sequence of post-14 Ma and hamsters. A new model of stepwise clade displacement is thus units from MN6 onwards correlates on a one-to-one basis to the 1.2-Myr proposed in which a series of increasingly pervasive, orbitally forced obliquity cycles 11–1. The only exception here is the temporal corre­ dispersal events of the new clade gradually disrupts the supremacy of spondence of two units (MN16 and 17) to one cycle (no. 3), of which the the old clade and leads to its gradual extinction. maximum amplitude point corresponds to their boundary as well as to Although an analysis of actual micromammal lineage turnover the Neogene-Quaternary boundary. (sensu van Dam et al., 2006) is expected to give more detailed results, The tight correlation between astronomical cycles and MN units is all our correlations at the biozone level are specific enough to formulate the more remarkable, as the latter were defined about fifty years ago several important implications for macroevolutionary/macroevolu­ based on faunal content only (Mein, 1975, 1990). At that time, cali­ tionary theory: brations to the geological time scale were still poor and there was no knowledge of orbital forcing of faunas on a million-year time scale. 1) Oscillations of long-term climatic parameters (notably the 405-kyr, Although the definition of MN units is not without problems (Hilgen 1.2-Myr and 2.4-Myr envelopes of shorter-term ~21-kyr, 41-kyr et al., 2012), discussions have mainly focused on boundaries rather than and ~ 100-kyr Milankovitch cycles) constitute an attractive mech­ on contents or number of units. In this respect it should be noted that the anism for triggering background species-level and (sub)family-level boundaries between MN units as drawn in Fig. 4 should not be treated as evolution. Besides exerting a more direct influence via global cooling moments of sudden fauna-wide turnover. Rather, most of the boundary events (e.g., Pleistocene), these long cycles affect the connectivity definitions depend on choices for specific appearance events, which may between bioprovinces (by the expansion and shrinking of biome or may not be isochronous with other (dis)appearance events. (A con­ vegetation), affecting overall rates of animal migration and inter­ crete example is the choice of the Mimomys entry as base for MN15, specific competition. resulting in a very short MN14 unit, Fig. 4; Hilgen et al., 2012.) In fact, 2) Evolutionary rates on the 105–6-year scale are higher during some the gradual change from one continental-scale fauna to another is intervals than others, because of the concentration in time of faunal consistent with the slow change in long-period orbital parameters, to events including the addition of new, competing species from remote which species respond differently (although it is expected that unit areas (dispersal) and from nearby areas (local speciation following transitions last much shorter than the units themselves). local extinction). Previous hypotheses (e.g., ‘Turnover-pulse’ and It thus appears that the two major Neogene global climate transitions ‘Court-Jester’ hypotheses; Vrba, 1995; Barnosky, 2001) have (mid-Miocene cooling transition (MMCT) with an expanding SH ice stressed the importance of climatic triggers (including Milankovitch sheet around 14 Ma and the onset of NH-derived glacial-interglacial cycles) for evolution, focusing particularly on extreme global cool­ cycles around 2.5 Ma, Fig. 4; Zachos et al., 2001a; Westerhold et al., ing/warming events. Our model is less demanding, stressing the 2020; Rohling et al., 2021) also represent dividing points in the mac­ importance of changes in regional insolation patterns as controlled roecological and macroevolutionary history of European mammals. The by long-period Milankovitch cycles. In addition, the model includes absence of a clear imprint of the 1.2-Myr obliquity cycle during the long-term resilience to short-period climate change as a central early-middle Miocene (before ~14 Ma, but especially before 17 Ma) is component (e.g., Ivany, 1996), as a result of which permanent intriguing, as this cycle has repeatedly been encountered in the ecosystem changes only occur after relatively long (> ~ 100-kyr) geological archive, including older Oligocene to middle Miocene benthic intervals with deviating conditions. δ18O records across the World’s oceans, where obliquity nodes have 3) The stepwise clade displacement model is relevant to theories of been found to correspond to glacial episodes on Antarctica (Zachos species radiation, because it predicts a concentration of speciation et al., 2001b; Pälike et al., 2006). Apparently, European-scale faunal events in time (e.g., Eldredge and Gould, 1972). change during the early Miocene was mostly dependent on regional, 4) The faunal dynamics model proposed here is relevant to macro­ eccentricity-precession-controlled climatic change, with high-latitude ecological theory, because it provides a time scale for species-level forcing (obliquity) playing a minor role. After the early-middle stasis (in aspects such as diversity, abundance, distribution; Brown Miocene transition (17–16 Ma), and even more so after the and Lomolino, 1998; Rosenzweig and McCord, 1991), while mid-Miocene cooling event (14 Ma), a stronger obliquity influence on acknowledging the importance of bioprovincial exchange and the Atlantic-Eurasian climate made itself felt. Given the amplification of dispersal. insolation variations on land (Short et al., 1991; Seltzer et al., 2023) and the increasing importance of aridity and wet-dry variations Admittedly, our study is limited to one group (rodents), abundance continent-inward, the faunal replacement model as proposed here for class (mostly common species), area (Iberian Peninsula), time frame late Neogene rodent clades could easily be extended to mammal faunas (late Neogene) and data resolution type (no data at 103–4-year time as a whole. scale). On the other hand, the approximately similar, 105–6-year-level duration of biozones in other marine and terrestrial realms across 5. Conclusions geological history (Gradstein et al., 2020) suggests that a similar type of forcing was at work. Obviously, much more testing will be necessary A comparison of late Miocene to early Pleistocene Iberian rodent before astronomically forced climate change can be accepted as a prime biozone boundary ages with records of marine SST, δ18O and δ13C, cause for species and community turnover in the fossil record. terrestrial hydrological proxies, and astronomical parameters suggests an orbital background for most zones, with specific configurations Declaration of Competing Interest (minima, maxima) of the 1.2-Myr obliquity cycle and the 405-kyr/2.4- Myr eccentricity cycles controlling faunal changes via oscillations in The authors declare that they have no known competing financial 17 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 interests or personal relationships that could have appeared to influence Bosmans, J., Drijfhout, S., Tuenter, E., Hilgen, F.J., Lourens, L.J., Rohling, E., 2015a. Precession and obliquity forcing of the freshwater budget over the Mediterranean. the work reported in this paper. Quat. Sci. Rev. 123, 16–30. Bosmans, J., Hilgen, F., Tuenter, E., Lourens, L.J., 2015b. Obliquity forcing of low- Data availability latitude climate. Clim. Past 11, 1335–1346. Bosmans, J., van der Ent, R.J., Haarsma, R., Drijfhout, S.S., Hilgen, F.J., 2020. Precession- and obliquity-induced changes in moisture sources for enhanced Research data are attached in the Supplementary precipitation over the Mediterranean Sea. Paleoceanogr. Paleoclimatol. 36, Tables accompanying this paper e2019PA003655 https://doi.org/10.1029/2019PA003655. Brown, J.H., Lomolino, M.V., 1998. Biogeography. Sinauer Associates, Inc., Sunderland, Massachusetts. Acknowledgments Butiseaca, G.A., Vasiliev, I., van der Meer, M.T.J., Krijgsman, W., Palcu, D.V., Feurdean, A., Niedermeyer, E.M., Mulch, A., 2021. Severe late Miocene droughts We thank Raef Minwer Barakat for allowing us to use his quantitative affected western Eurasia. Glob. Planet. Chang. 206, 103644 https://doi.org/ 10.1016/j.gloplacha.2021.103644. data for Guadix-Baza Basin rodents and Tanja Kouwenhoven and Alex­ Calvo, J.P., Alcalá, L., Alonzo Zarza, A., van Dam, J.A., Gutiérrez Santolalla, F., 1999. andrina Tzanova for their data from Mediterranean Sea records. Frits Estratigrafía y structura del área de Los Mansuetos (Cuenca de Teruel). Precisiones Hilgen and Gert-Jan Reichart are thanked for discussion and the re­ para la definición del estratotipo del Turoliense. Geogac. 25, 55–58. Casanovas-Vilar, I., Agustí, J., 2007. Ecogeographical stability and climate forcing in the viewers for their comments. This publication is part of project I+D+i late Miocene (Vallesian) rodent record of Spain. Palaeogeogr. Palaeoclimatol. PID2020-117289GBI00 funded by MCIN/AEI/10.13039/5011000110 Palaeoecol. 248, 169–189. 33. Our research has also been supported by CERCA Programme/Gen­ Channell, J.E.T., Amigo, A.E., Fronval, T., Rack, F., Lehman, B., 1999. Magnetic stratigraphy at Sites 907 and 985 in the Norwegian-Greenland Sea and a revision of eralitat de Catalunya, AEI–FEDER EU (to M.F.), AGAUR (Consolidated the Site 907 composite section. In: Raymo, M.E., Jansen, E., Blum, P., Herbert, T.D. Research Group, 2017 SGR 960 to M.F. and 2017 SGR 596 to M.G.), and (Eds.), Proc. Ocean Drill. Program Sci. Res., 162 Ocean Drilling Program, College the Department of Science, University and Society of Knowledge of the Station, Texas, pp. 131–148. Chevret, P., Renaud, S., Helvaci, Z., Ulrich, R.G., Quéré, J.P., Michaux, J.R., 2020. Government of Aragon (Reference Research Group E04_20R FOCON­ Genetic structure, ecological versatility, and skull shape differentiation in Arvicola TUR funding to Fundación Dinópolis -L.A.). water voles (Rodentia, Cricetidae). J. Zool. Syst. Evol. Res. 58, 1323–1334. Conte, M., Sicre, M.-A., Ruhlemann, C., 2006. Global temperature calibration of the alkenone unsaturation index (UK’37) in surface waters and comparison with surface Appendix A. Supplementary data sediments. Geochem. Geophys. Geosyst. 7, Q02005. https://doi.org/10.1029/ 2005GC001054. Supplementary data to this article can be found online at https://doi. de Bruijn, H., Mein, P., Montenat, C., van de Weerd, A., 1975. Correlations entre les gisements de rongeurs et les formations marines du Miocène terminal d’Espagne org/10.1016/j.gloplacha.2023.104153. meridionale; provinces d’Alicante et de Murcia. Proc. K. Ned. Akad. Wet. Ser. B 78, 282–313. References De Miguel, D., Azanza, B., Morales, J., 2019. Regional impacts of global climate change: a local humid phase in Central Iberia in a late Miocene drying world. Paleontol. 62, 77–92. Abdul Aziz, H., Sanz-Rubio, E., Calvo, J.P., Hilgen, F.J., Krijgsman, W., 2003. De Vleeschouwer, D., Drury, A.J., Vahlenkamp, M., Rochholz, F., Liebrand, D., Pälike, H., Paleoenvironmental reconstruction of a middle Miocene proximal alluvial fan to 2020. High-latitude biomes and rock weathering mediate climate-carbon cycle feed- cyclic shallow lacustrine depositional system in the Calatayud Basin (NE Spain). backs on eccentricity timescales. Nat. Commun. 11, 1–10. Sedimentol. 50, 211–236. Drury, A.J., Westerhold, T., Frederichs, T., Tian, J., Wilkens, R., Channell, J.E.T., Abdul Aziz, H., van Dam, J.A., Hilgen, F.J., Krijgsman, W., 2004. Astronomical forcing in Evans, H., Johns, C.M., Lyle, M., Rohl, U., 2017. Late Miocene climate and time scale Upper Miocene continental sequences: implications for the Geomagnetic Polarity reconciliation: Accurate orbital calibration from a deep-sea perspective. Earth Time Scale. Earth Planet. Sci. Lett. 222, 243–258. Planet. Sci. Lett. 475, 254–266. Abels, H.A., Abdul Aziz, H., Calvo, J.P., Tuenter, E., 2008. Shallow lacustrine microfacies Eldredge, N., Gould, S.J., 1972. Punctuated equilibria: analternative to phyletic document orbitally paced lake-level history in the Miocene Teruel Basin (NE Spain). gradualism. In: Schopf, T.J.M. (Ed.), Models in paleobiology. W. H. Freeman, San Sedimentol. 56, 399–419. Francisco, pp. 82–115. Abels, H.A., Aziz, H.A., Ventra, D., Hilgen, F.J., 2009. Orbital climate forcing in mudflat Ezquerro, L., Luzón, A., Simón, J.L., Liesa, C.L., 2020. Segmentation and increasing to marginal lacustrine deposits in the Miocene Teruel Basin (Northeast Spain). activity in the Neogene-Quaternary Teruel Basin rift (Spain) revealed by J. Sediment. Res. 79, 831–847. morphotectonic approach. J. Struct. Geol. 135, 104043 https://doi.org/10.1016/j. Abels, H.A., Aziz, H.A., Krijgsman, W., Smeets, S.J.B., Hilgen, F.J., 2010. Long-period jsg.2020.104043. eccentricity control on sedimentary sequences in the continental Madrid Basin Ezquerro, L., Muñoz, A., Liesa, C.L., Simón, J.L., Luzón, A., 2022. Late Neogene to early (middle Miocene, Spain). Earth Planet. Sci. Lett. 289, 220–231. Quaternary climate evolution in southwestern Europe from a continental Abreu, V.S., Anderson, J.B., 1998. Glacial eustasy during the Cenozoic: sequence perspective. Glob. Planet. Chang. 211, 103788 https://doi.org/10.1016/j. stratigraphic implications. AAPG Bull. 82, 1385–1400. gloplacha.2022.103788. Adrover, R., Mein, P., Moissenet, E., 1993. Roedores de la transición Mio-Plioceno de la Fauquette, S., Guiot, J., Suc, J.-P., 1998. A method for climatic reconstruction of the región de Teruel. Paleontol. Evol. 26-27, 47–84. Mediterranean Pliocene using pollen data. Palaeogeogr. Palaeoclimatol. Palaeoecol. Agustí, J., Oms, O., Remacha, E., 2001. Long Plio-Pleistocene terrestrial record of climate 144, 183–201. change and mammal turnover in southern Spain. Quat. Res. 56, 411–418. Fejfar, O., Heinrich, W.D., Kordos, L., Maul, L.C., 2011. Microtoid cricetids and the early Agustí, J., Garcés, M., Krijgsman, W., 2006. Evidence for African-Iberian exchanges history of arvicolids (Mammalia, Rodentia). Palaeontol. Electron. 14, 27A. https:// during the Messinian in the Spanish mammalian record. Palaeogeogr. palaeo-electronica.org/2011_3/6_fejfar/index.html. Palaeoclimatol. Palaeoecol. 238, 5–14. Fraser, D., Soul, L.C., Tóth, A.B., Balk, M.A., Eronen, J.T., Pineda-Munoz, S., Alcalá, L., Alonso-Zarza, A., Álvarez Sierra, M.A., Azanza, B., Calvo, J.P., Cañaveras, J.C., Shupinski, A.B., Villaseñor, A., Barr, W.A., Behrensmeyer, A.K., Du, A., Faith, J.T., van Dam, J., Garcés, M., Krijgsman, W., van der Meulen, A.J., Morales, J., Peláez- Gotelli, N.J., Graves, G.R., Jukar, A.M., Looy, C.V., Miller, J.H., Potts, R., Lyons, S.K., Campomanes, P., Pérez González, A., Sánchez Moral, S., Sancho, R., Sanz Rubio, E., 2021. Investigating biotic interactions in deep time. Trends Ecol. Evol. 36, 61–75. 2000. El registro sedimentario y faunístico de las cuencas de Calatayud-Daroca y Freudenthal, M., Mein, P., Martín Suárez, E., 1998. Revision of late Miocene and Pliocene Teruel. Evolución paleoambiental y paleoclimática durante el Neógeno. Rev. Soc. Cricetinae (Rodentia, Mammalia) from Spain and France. Treb. Mus. Geol. Barc. 7, Geol. Esp. 13, 323–343. 11–93. Alonso-Zarza, A.M., Calvo, J.P., 2000. Palustrine sedimentation in an episodically Gao, P., Nie, J., Breecker, D.O., Gallagher, T., Serach, L., Alonso- Zarza, A.M., 2022. subsiding basin: the Miocene of the northern Teruel Graben (Spain). Palaeogeogr. Similar magnetic enhancement mechanisms between Chinese loess and alluvial Palaeoclimatol. Palaeoecol. 160, 1–21. sediments from the Teruel Basin, NE Spain, and paleoclimate implications. Geophys. Barnosky, A.D., 2001. Distinguishing the effects of the Red Queen and Court Jester on Res. Lett. 49, e2021GL096977 https://doi.org/10.1029/2021GL096977. Miocene mammal evolution in the northern Rocky Mountains. J. Vertebr. Paleontol. Garcés, M., Krijgsman, W., Agustí, J., 2001. Chronostratigraphic framework and 21, 172–185. evolution of the Fortuna basin (Eastern Betics) since the late Miocene. Basin Res. 13, Benson, R.H., Rakic-El Bied, K., Bonaduce, G., 1991. An important current reversal 199–216. (influx) in the Rifian Corridor (Morocco) at the Tortonian-Messinian boundary: the Gilinsky, N.L., Bambach, R.K., 1987. Asymmetrical patterns of origination and extinction end of Tethys Ocean. Paleoceanogr. 6, 165–192. in higher taxa. Paleobiol. 13, 427–445. Benton, M.J., 1996. On the nonprevalence of competitive replacement in the evolution of Godoy, A., Olivé, A., Moissenet, E., 1983. Mapa y memoria explicativa de la Hoja 567 tetrapods. In: Jablonski, D., Erwin, D.H., Lipps, J.H. (Eds.), Evolutionary (Teruel) del Mapa Geológico de España, a escala 1:50.000 (2a Serie). IGME Madrid, Paleobiology. Chicago University Press, Chicago, pp. 185–210. 1-70. Berger, A.L., 1978. Long-term variations of caloric insolation resulting from the Earth’s Gould, S.J., Calloway, C.B., 1980. Clams and brachiopods: ships that pass in the night. orbital elements. Quat. Res. 9, 139–167. Paleobiol. 6, 383–396. Böhme, M., Ilg, A., 2003. Database on fossil fishes, amphibians and reptiles: fosFARbase. http://www.wahre-staerke.com. 18 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 Gradstein, F.M., Ogg, J.G., Schmitz, M.D., Ogg, G.M., 2020. The Geologic Time Scale Laskar, J., Robutel, P., Joutel, F., Gastineau, M., Correia, A.C.M., Levrard, B., 2004. 2020, 2. Elsevier, Boston, MA. A long-term numerical solution for the insolation quantities of the Earth. Astron. Gudjonsson, L., van der Zwaan, G.J., 1985. Anoxic events in the Pliocene Mediterranean: Astrophys. 428, 261–285. stable isotope evidence of run-off. Proc. K. Ned. Akad. Wet. Ser. B 88, 69–82. Lawrence, K.T., Herbert, T.D., Brown, C.M., Raymo, M.E., Haywood, A.M., 2009. High- Gusovsky, V.V., Sinitsa, M.V., 2021. Dietary preferences and autecology of a basal amplitude variations in North Atlantic Sea surface temperature during the early Holarctic ground squirrel Spermophilinus (Rodentia: Sciuridae: Xerinae), inferred Pliocene warm period. Paleoceanogr. 24, PA2218. from dental microwear analysis. Hist. Biol. 33, 54–61. Lisiecki, L.E., Raymo, M.E., 2005. A Pliocene-Pleistocene stack of 57 globally distributed Henderson-Sellers, A., Robinson, P.J., 1994. Contemporary Climatology. John Wiley & benthic δ18O records. Paleoceanogr. 20, PA001071. https://doi.org/10.1029/ Sons, New York. 2004PA001071. Herbert, T.D., Lawrence, K.T., Tzanova, A., Peterson, L.C., Caballero-Gill, R., Kelly, C.S., López-Martínez, N., 1989. Revisión sistemática y biostratigráfica de los lagomorfos 2016. Late Miocene global cooling and the rise of modern ecosystems. Nat. Geosci. 9, (Mammalia) el Terciario y Cuaternario inferior de España. Diputación General de 843–849. Aragón, Zaragoza. Herold, N., Muller, R.D., Seton, M., 2010. Comparing early to middle Miocene terrestrial Lourens, L.J., Hilgen, F.J., 1997. Long-periodic variations in the earth’s obliquity and climate simulations with geological data. Geosph. 6, 952–961. their relation to third-order eustatic cycles and late Neogene glaciations. Quat. Int. Hilgen, F.J., 1991. Extension of the astronomically calibrated (polarity) time scale to the 40, 43–52. Miocene/Pliocene boundary. Earth Planet. Sci. Lett. 107, 349–368. Lourens, L.J., Hilgen, F.J., Gudjonsson, L., Zachariasse, W.J., 1992. Late Pliocene to early Hilgen, F.J., Lourens, L.J., van Dam, J.A., 2012. The Neogene Period. In: Gradstein, F.M., Pleistocene astronomically forced sea surface productivity and temperature Ogg, J.G., Schmitz, M.D., Ogg, G. (Eds.), A geological time scale 2012. Elsevier, variations in the Mediterranean. Mar. Micropaleontol. 19, 49–78. Amsterdam, pp. 923–978. Lourens, L.J., Hilgen, F.J., Antonarakou, A., Van Hoof, A.A.M., Vergnaud-Grazzini, C., Hilgen, F.J., Lourens, L.J., Pälike, H., Research Support Team, 2020. Should unit- Zachariasse, W.J., 1996. Evaluation of the Pliocene to early Pleistocene astronomical stratotypes and astrochronozones be formally defined? A dual proposal (including time scale. Paleoceanogr. 11, 391–413. postscriptum). Newsl. Stratigr. 53, 19–39. Martin, R.A., 2008. Arvicolidae. In: Janis, C.M., Gunnell, G.F., Uhen, M.D. (Eds.), Hodell, D.A., Benson, R.H., Kennett, J.P., Rakic-El Bied, K., 1989. Stable isotope Evolution of Tertiary Mammals of North America, vol. 2. Cambridge University stratigraphy of latest Miocene sequences in Northwest Morocco: the Bou Regreg Press, Cambridge, pp. 480–497. Section. Paleoceanogr. 4, 467–482. Martín Suárez, E., Mein, P., 1998. Revision of the genera Parapodemus, Apodemus, Hodell, D.A., Benson, R.H., Kent, D.V., Boersma, A., Bied, K.R.E., 1994. Rhagamys and Rhagapodemus (Rodentia, Mammalia). Geobios 31, 87–97. Magnetostratigraphic, biostratigraphic, and stable-isotope stratigraphy of an upper Mein, P., 1975. Résultats du Groupe de Travail de Vertébrés. In: Senes, J. (Ed.), Report on Miocene drill core from the Salé-Briqueterie (northwestern Morocco): a high- Activity of the RCMNS Working Groups (1971–1975). RCMNS, Bratislava, resolution chronology for the Messinian stage. Paleoceanogr. 9, 835–855. pp. 78–81. Hodell, D.A., Curtis, J.H., Sierro, F.J., Raymo, M.E., 2001. Correlation of late Miocene to Mein, P., 1990. Updating of MN zones. In: Lindsay, E., Fahlbusch, V., Mein, P. (Eds.), early Pliocene sequences between the Mediterranean and North Atlantic. European Neogene Mammal Chronology. NATO ASI Series A, pp. 73–90. Paleoceanogr. 16, 164–178. Mein, P., Moissenet, E., Adrover, R., 1990. Biostratigraphie du Neogène supérieur du Holbourn, A., Kuhnt, W., Regenberg, M., Schulz, M., Mix, A., Andersen, N., 2010. Does bassin de Teruel. Paleontol. Evol. 23, 121–139. Antarctic glaciation force migration of the tropical rain belt? Geology 38, 783–786. Milankovitch, M., 1930. Mathematische Klimalehre und astronomische Theorie der Hordijk, K., de Bruijn, H., 2009. The succession of rodent faunas from the Mio/Pliocene Klimaschwankungen. In: Köppen, W., Gieger, R. (Eds.), Handbuch der Klimatologie I lacustrine deposits of the Florina-Ptolemais-Servia Basin (Greece). Hell. J. Geosci. (A), pp. 1–176. 44, 21–103. Nobel, F.A., Andriessen, P.A.M., Hebeda, E.H., Priem, H.N.A., Rondeel, H.E., 1981. Hui, Z., Liu, J., Chevalier, M., Wei, X., Chen, P., Zhang, J., Peng, T., Zhou, X., 2023. Isotopic dating of the post-alpine Neogene volcanism in the Betic Cordilleras, Multiple forcing on late Miocene East Asian Summer Monsoon Precipitation Southern Spain. Geol. Mijnb. 60, 209–214. Variability in NE Tibetan Plateau. Catena 221 A, 106752. https://doi.org/10.1016/j. Nowak, R.M., 1999. Walkers’ Mammals of the World. Johns Hopkins University Press, catena.2022.106752. Baltimore. Ivany, L.C., 1996. Coordinated stasis or coordinated turnover? Exploring intrinsic vs. Olsen, P.E., Kent, D.V., 1999. Long-period Milankovitch cycles from the late Triassic and extrinsic controls on pattern. Palaeogeogr. Palaeoclimatol. Palaeoecol. 127, early Jurassic of eastern North America and their implications for the calibration of 239–256. the early Mesozoic time-scale and the long-term behaviour of the planets. Philos. Jardine, P.E., Janis, C.M., Sahney, S., Benton, M.J., 2012. Grit not grass: concordant Trans. R. Soc. Lond. A 357, 1761–1786. patterns of early origin of hypsodonty in Great Plains ungulates and Glires. Opdyke, N., Mein, P., Lindsay, E., Pérez-Gonzáles, A., Moissenet, E., Norton, V.L., 1997. Palaeogeogr. Palaeoclimatol. Palaeoecol. 365, 1–10. Continental deposits, magnetostratigraphy and vertebrate paleontology, late John, C.M., Mutti, M., Adatte, T., 2003. Mixed carbonate-siliciclastic record on the North Neogene of Eastern Spain. Palaeogeogr. Palaeoclimatol. Palaeoecol. 133, 129–148. African margin (Malta) - coupling of weathering processes and mid Miocene climate. Palcu, D.V., Vasiliev, I., Stoica, M., Krijgsman, W., 2019. The end of the Great Khersonian Geol. Soc. Am. Bull. 115, 217–229. Drying of Eurasia: magnetostratigraphic dating of the Maeotian transgression in the Kälin, D., Kempf, O., 2009. High-resolution stratigraphy from the continental record of Eastern Paratethys. Basin Res. 31, 33–58. the Middle Miocene Northern Alpine Foreland Basin of Switzerland. N. Jahrb. Geol. Palcu, D.V., Patina, I.S., Sandric, I., Lazarev, S., Vasiliev, I., Stoica, M., Krijgsman, W., Paläontol. Abh. 254, 177–235. 2021. Late Miocene megalake regressions in Eurasia. Sci. Rep. 11, 11471. Kemp, T.S., 2006. The origin and early radiation of the therapsid mammal-like reptiles: a Pälike, H., Norris, R.D., Herrle, J.O., Wilson, P.A., Coxall, H.K., Lear, C.H., Shackleton, N. palaeobiological hypothesis. J. Evol. Biol. 19, 1231–1247. J., Tripati, A.K., Wade, B.S., 2006. The heartbeat of the Oligocene climate system. Kloosterboer - van Hoeve, M., 2000. Cyclic changes in the late Neogene vegetation of Science 314, 1894–1898. northern Greece. LPP Contrib. Ser. 12, 1–131. Panitz, S., Salzmann, U., Risebrobakken, B., De Schepper, S., Pound, M.J., 2016. Climate Kontakiotis, G., Besiou, E., Antonarakou, A., Zarkogiannis, S.D., Kostis, A., Mortyn, P.G., variability and long-term expansion of peat lands in Arctic Norway during the late Moissette, P., Cornée, J.-J., Schulbert, C., Drinia, H., Anastasakis, G., Karakitsios, V., Pliocene (ODP Site 642, Norwegian Sea). Clim. Past 12, 1043–1060. 2022. Decoding Sea surface and paleoclimate conditions in the eastern Pillans, B., Gibbard, P., 2012. The Quaternary Period. In: Gradstein, F.M., Ogg, J.G., Mediterranean over the Tortonian-Messinian transition. Palaeogeogr. Schmitz, M.D., Ogg, G. (Eds.), A Geological Time Scale 2012. Elsevier, Amsterdam, Palaeoclimatol. Palaeoecol. 534, 109312. pp. 979–1010. Kouwenhoven, T.J., van der Zwaan, G.J., 2006. A reconstruction of late Miocene Qin, J., Zhang, R., Kravchinsky, V.A., Valet, J.-P., Sagnotti, L., Li, J., Xu, Y., Anwar, T., Mediterranean circulation patterns using benthic foraminifera. Palaeogeogr. Yue, L., 2022. 1.2 Myr band of Earth-Mars obliquity modulation on the evolution of Palaeoclimatol. Palaeoecol. 238, 373–385. cold late Miocene to warm early Pliocene climate. J. Geophys. Res. Solid Earth 127, Kouwenhoven, T.J., Seidenkrantz, M.-S., van der Zwaan, G.J., 1999. Deep-water e2022JB024131. https://doi.org/10.1029/2022JB024131. changes: the near-synchronous disappearance of a group of benthic foraminifera Raffi, I., Wade, B.S., Pälike, H., 2020. Chapter 29 - the Neogene Period. In: Gradstein, F. from the late Miocene Mediterranean. Palaeogeogr. Palaeoclimatol. Palaeoecol. 152, M., Ogg, J.G., Schmitz, M.D., Ogg, G.M. (Eds.), The Geologic Time Scale 2020, vol. 2. 259–281. Elsevier, Boston, MA. Krause, D.W., 1986. Competitive exclusion and taxonomic displacement in the fossil Repenning, C.A., Fejfar, O., Heinrich, W.D., 1990. Arvicolid rodent biochronology of the record: The case of rodents and multituberculates in North America. In: Flanagan, K. Northern Hemisphere. In: Fejfar, O., Heinrich Wolf, D. (Eds.), International M., Lillegraven, J.A. (Eds.), Vertebrates, Phylogeny, and Philosophy. University of Symposium on Evolution, Phylogeny and Biostratigraphy of Arvicolids (Rodentia, Wyoming, Laramie. Mammalia). Geol. Survey, Prague, pp. 385–418. Krijgsman, W., Hilgen, F.J., Raffi, I., Sierro, F.J., Wilson, D.S., 1999. Chronology, causes Rodríguez-López, J.P., Liesa, C.L., van Dam, J., Lafuente, P., Arlegui, L., Ezquerro, L., de and progression of the Messinian salinity crisis. Nature 400, 652–655. Boer, P.L., 2012. Aeolian construction and alluvial dismantling of a fault-bounded Krijgsman, W., Capella, W., Simon, D., Hilgen, F.J., Kouwenhoven, T.J., Sierro, F.J., intracontinental aeolian dune field (Teruel Basin, Spain); a continental perspective Tulbure, M.A., van den Berg, B.C.J., van der Schee, M., Rachel Flecker, R., 2018. The on late Pliocene climate change and variability. Sedimentol. 59, 1536–1567. Gibraltar Corridor: watergate of the Messinian Salinity Crisis. Mar. Geol. 403, Rohling, E.J., 1999. Environmental control on Mediterranean salinity and δ18O. 236–246. Paleoceanogr. 14, 706–715. Kuiper, K.F., Krijgsman, W., Garcés, M., Wijbrans, J., 2006. Revised isotopic (40Ar/39Ar) Rohling, E.J., Yu, J., Heslop, D., Foster, G.L., Opdyke, B., Roberts, A.P., 2021. Sea level age for the lamproite volcano of Cabezos Negros, Fortuna Basin (Eastern Betics, SE and deep- sea temperature reconstructions suggest quasi-stable states and critical Spain). Palaeogeogr. Palaeoclimatol. Palaeoecol. 238, 53–63. transitions over the past 40 million years. Sci. Adv. 7, eabf5326. https://doi.org/ Kutzbach, J.E., Chen, G., Cheng, H., Edwards, R.L., Liu, Z., 2014. Potential role of winter 10.1126/sciadv.abf5326. rainfall in explaining increased moisture in the Mediterranean and Middle East Rosenzweig, M.L., McCord, R.D., 1991. Incumbent replacement: evidence for long-term during periods of maximum orbitally-forced insolation seasonality. Clim. Dyn. 42, evolutionary progress. Paleobiol. 17, 202–213. 1079–1095. 19 J.A. van Dam et al. Global and Planetary Change 227 (2023) 104153 Rossignol-Strick, M., 1985. Mediterranean Quaternary sapropels, an immediate response van Dam, J.A., Reichart, G.J., 2009. Oxygen and carbon isotope signatures in late to African monsoons to variations of insolation. Palaeogeogr. Palaeoclimatol. Neogene horse teeth from Spain and application as temperature and seasonality Palaeoecol. 49, 237–263. proxies. Palaeogeogr. Palaeoclimatol. Palaeoecol. 274, 64–81. Ruddiman, W.F., Raymo, M.E., Martinson, D.G., Clement, B.M., Backman, J., 1989. van Dam, J.A., Utescher, T., 2016. Plant- and micromammal-based paleoprecipitation Pleistocene evolution of Northern Hemisphere climate. Paleoceanogr. 4, 353–412. proxies: comparing results of the Coexistence and Climate-Diversity Approach. Schmittner, A., Gruber, N., Mix, A.C., Key, R.M., Tagliabue, A., Westberry, T.K., 2013. Palaeogeogr. Palaeoclimatol. Palaeoecol. 443, 18–33. Biology and air-sea gas exchange controls on the distribution of carbon isotope ratios van Dam, J.A., Alcalá, L., Alonzo Zarza, A., Calvo, J.P., Garcés, M., Krijgsman, W., 2001. (δ13C) in the ocean. Biogeosci. 10, 5793–5816. The upper Miocene mammal record from the Teruel-Alfambra region (Spain). The Seltzer, A.M., Blard, P.-H., Sherwood, S.C., Kageyama, M., 2023. Terrestrial MN system and continental stage/age concepts discussed. J. Vertebr. Paleontol. 21, amplification of past, present, and future climate change. Sci. Adv. 9, eadf8119. 367–385. https://doi.org/10.1126/sciadv.adf8119. van Dam, J.A., Abdul Aziz, H., Álvarez Sierra, M.Á., Hilgen, F.J., van den Hoek Sepkoski, J.J., 1996. Competition in macroevolution: The double wedge revisited. In: Ostende, L.W., Lourens, L.J., Mein, P., van der Meulen, A.J., Peláez-Campomanes, P., Jablonski, D., Erwin, D.H., Lipps, J.H. (Eds.), Evolutionary Paleobiology. Chicago 2006. Long-period astronomical forcing of mammal turnover. Nature 443, 687–691. University Press, Chicago, pp. 211–255. van Dam, J.A., Furió, M., van Balen, R.T., 2014. Re-interpreting the biochronology of the Shackleton, N.J., Hall, M.A., 1997. The late Miocene stable isotope record, Site 926. In: La Celia and Los Gargantones mammal sites (late Miocene, Murcia, Spain). Geobios Shackleton, N.J., Curry, W.B., Richter, C., Bralower, T.J. (Eds.), Proc. Ocean Drill. 47, 155–164. Program Sci. Res 154. Ocean Drill. Program, College Station, Texas, pp. 367–373. van Dam, J.A., Mein, P., Garcés, M., van Balen, R.T., Furió, M., Alcalá, L., 2023. A new Shackleton, N.J., Hall, M.A., Pate, D., 1995. Pliocene stable isotope stratigraphy of site rodent chronology for the late Neogene of Spain. Geobios 76, 53–74. 846. In: Proc. Ocean Drill. Program Sci. Res. 138. Ocean Drill. Program, College van de Weerd, A., 1976. Rodent faunas of the Mio-Pliocene continental sediments of the Station, Texas, pp. 337–355. Teruel-Alfambra region, Spain. Utrecht Micropaleontol. Bull. Spec. Publ. 2, 1–185. Short, D.A., Mengel, J.G., Crowley, T.J., Hyde, W.T., North, G.R., 1991. Filtering of van de Weerd, A., Daams, R., 1978. Quantitative composition of rodent faunas in the Milankovitch cycles by Earth’s geography. Quat. Res. 35, 157–173. Spanish Neogene and paleoecological implications. Proc. K. Ned. Akad. Wet. Ser. B Sierro, F.J., Ledesma, S., Flores, J.-A., Torrescusa, S., del Olmo, W.M., 2000. Sonic and 81, 448–473. gamma-ray astrochronology: cycle to cycle calibration of Atlantic climatic records to van der Laan, E., Hilgen, F.J., Lourens, L.J., de Kaenel, E., Gaboardi, S., Iaccarino, S., Mediterranean sapropels and astronomical oscillations. Geology 28, 695–698. 2012. Astronomical forcing of Northwest African climate and glacial history during Sinitsa, M.V., Čermák, S., Kryuchkova, L.Y., 2022. Cranial anatomy of Csakvaromys the late Messinian (6.5-5.5 Ma). Palaeogeogr. Palaeoclimatol. Palaeoecol. 313, bredai (Rodentia, Sciuridae, Xerinae) and implications for ground squirrel evolution 107–126. and systematics. J. Mamm. Evol. 29, 149–189. van der Meulen, M.J., Kouwenhoven, T.J., van der Zwaan, G.J., Meulenkamp, J.E., Speijer, R.P., Pälike, H., Hollis, C.J., Hooker, J.J., Ogg, J.G., 2020. The Paleogene Period. Wortel, M.J.R., 1999. Late Miocene uplift in the Romagnan Apennines and the In: Gradstein, F.M., Ogg, J.G., Schmitz, M.D., Ogg, G. (Eds.), Geological time scale detachment of subducted lithosphere. Tectonophys. 315, 319–335. 2020. Elsevier, Amsterdam, pp. 1087–1140. van der Schee, M., Sierro, F.J., Jiménez-Espejo, F.J., Hernández-Molina, F.J., Flecker, R., Turco, E., Hilgen, F.J., Lourens, L.J., Shackleton, N.J., Zachariasse, W.J., 2001. Flores, J.A., Acton, G., Gutjahr, M., Grunert, P., García-Gallardo, Á., Andersen, N., Punctuated evolution of global climate cooling during the late Middle to early late 2016. Evidence of early bottom water current flow after the Messinian Salinity Crisis Miocene: High-resolution planktonic foraminiferal and oxygen isotope records from in the Gulf of Cadiz. Mar. Geol. 380, 315–329. the Mediterranean. Paleoceanogr. 16, 405–423. Vrba, E.S., 1995. On the connections between paleoclimate and evolution. In: Vrba, E.S., Tzanova, A., Herbert, T.D., 2015. Regional and global significance of Pliocene Sea Denton, G.H., Partridge, T.C., Burckle, L.H. (Eds.), Paleoclimate and Evolution, with surface temperatures from the Gulf of Cadiz (Site U1387) and the Mediterranean. Emphasis on Human Origins. Yale University, New Haven, pp. 24–45. Glob. Planet. Chang. 133, 371–377. Wang, P.X., Tian, J., Lourens, L.J., 2010. Obscuring of long eccentricity cyclicity in Tzanova, A., Herbert, T.D., Peterson, L., 2015. Cooling Mediterranean Sea surface Pleistocene oceanic carbon isotope records. Earth Planet. Sci. Lett. 290, 319–330. temperatures during the late Miocene provide a climate context for evolutionary Westerhold, T., Bickert, U., Röhl, U., 2005. Middle to late Miocene oxygen isotope transitions in Africa and Eurasia. Earth Planet. Sci. Lett. 419, 71–80. stratigraphy of ODP site 1085 (SE Atlantic): new constrains on Miocene climate Tzedakis, P.C., 2007. Seven ambiguities in the Mediterranean palaeoenvironmental variability and sea-level fluctuations. Palaeogeogr. Palaeoclimatol. Palaeoecol. 217, narrative. Quat. Sci. Rev. 26, 2042–2066. 205–222. Valero, L., Garcés, M., Cabrera, L., Costa, E., Sáez, A., 2014. 20 Myr of eccentricity paced Westerhold, T., Marwan, N., Drury, A.J., Liebrand, D., Agnini, C., Anagnostou, E., lacustrine cycles in the Cenozoic Ebro Basin. Earth Planet. Sci. Lett. 408, 183–193. Barnet, J.S.K., Bohaty, S.M., De Vleeschouwer, D., Florindo, F., Frederichs, T., Valero, L., Cabrera, L., Sáez, A., Garcés, M., 2016. Long-period astronomically-forced Hodell, D.A., Holbourn, A.E., Kroon, D., Lauretano, V., Littler, K., Lourens, L.J., terrestrial carbon sinks. Earth Planet. Sci. Lett. 444, 131–138. Lyle, M., Pälike, H., Röhl, U., Tian, J., Wilkens, R.H., Wilson, P.A., Zachos, J.C., Valla, P.G., van der Beek, P.A., Shuster, D.L., Braun, J., Herman, F., Tassan-Got, L., 2020. An astronomically dated record of Earth’s climate and its predictability over Gautheron, C., 2012. Late Neogene exhumation and relief development of the Aar the last 66 million years. Science 369, 1383–1387. and Aiguilles Rouges massifs (Swiss Alps) from low-temperature thermochronology Zachos, J., Pagani, M., Sloan, L., Thomas, E., Billups, K., 2001a. Trends, rhythms, and modeling and 4He/3He thermochronometry. J. Geophys. Res. 117, F01004. https:// aberrations in global climate 65 Ma to present. Science 292, 686–693. doi.org/10.1029/2011JF002043. Zachos, J., Shackleton, N.J., Revenaugh, J.S., Pälike, H., Flower, B.P., 2001b. Climate van Balen, R.T., Forzoni, A., van Dam, J.A., 2015. Active faulting and folding along response to orbital forcing across the Oligocene-Miocene boundary. Science 292, Jumilla Fault Zone, northeastern Betics. Spain. Geomorphol 237, 88–97. 274–278. van Dam, J.A., 1997. The small mammals from the upper Miocene of the Teruel- Zagwijn, W.H., 1992. The beginning of the ice age in Europe and its major subdivisions. Alfambra region (Spain): paleobiology and paleoclimatic reconstructions. Geol. Quat. Sci. Rev. 11, 583–591. Ultraiect. 156, 1–204. Zhang, J., Quay, P.D., Wilbur, D.O., 1995. Carbon isotope fractionation during gas-water van Dam, J.A., 2006. Geographic and temporal patterns in the late Neogene (12-3 Ma) exchange and dissolution of CO2. Geochim. Cosmochim. Acta 59, 107–114. aridification of Europe: the use of small mammals paleoprecipitation proxies. Palaeogeogr. Palaeoclimatol. Palaeoecol. 238, 190–218. 20