Temperature-related biodiversity change across temperate marine and terrestrial systems

Climate change is reshaping global biodiversity as species respond to changing temperatures. However, the net effects of climate-driven species redistribution on local assemblage diversity remain unknown. Here, we relate trends in species richness and abundance from 21,500 terrestrial and marine assemblage time series across temperate regions (23.5–60.0° latitude) to changes in air or sea surface temperature. We find a strong coupling between biodiversity and temperature changes in the marine realm, where species richness mostly increases with warming. However, biodiversity responses are conditional on the baseline climate, such that in initially warmer locations richness increase is more pronounced while abundance declines with warming. In contrast, we do not detect systematic temperature-related richness or abundance trends on land, despite a greater magnitude of warming. As the world is committed to further warming, substantial challenges remain in maintaining local biodiversity amongst the non-uniform inflow and outflow of ‘climate migrants’. Temperature-driven community restructuring is especially evident in the ocean, whereas climatic debt may be accumulating on land. Biodiversity time series from temperate regions reveal that marine communities in warmer places gain species but lose individuals with warming, but colder environments show weaker trends, whereas no systematic relationships between biodiversity and temperature change were detectable for terrestrial communities.

C limate change is driving a reorganization of ecological communities as species track changes in air and ocean temperatures [1][2][3] . However, global warming is not unfolding evenly across the planet, and this heterogeneity is layered over the uneven distribution of biodiversity. Populations of thermally restricted species may decline with warming as individuals die, fail to reproduce or move to more suitable locations [4][5][6] . Indeed, declines are typically expected for tropical species, since they have narrower thermal tolerances than temperate species and live closer to their upper thermal limits 5,[7][8][9][10] . In contrast, species benefiting from warming may exhibit abundance increases and expand their geographic ranges [1][2][3][4]11 . Thus, mid-to high-latitudes undergoing warming may provide suitable habitat for species expanding their ranges poleward 4,12 . Because the tropics hold the majority of the world's species 13 , lower-latitude warming temperate regions may experience larger increases in species richness and abundance compared to temperate locations at higher latitudes, due to a larger source pool of climate immigrants (Fig. 1). Biodiversity change may further depend on the baseline climate-that is, within latitudinal bands the effects of warming might differ between initially warmer versus colder regions 3,5,12,14 . For instance, warmer temperate regions may offer greater habitat suitability for climate immigrants from lower latitudes.
Warming-induced biodiversity change may also be stronger in the ocean than on land 3,15,16 . Marine species are highly responsive to temperature change and can track changing isotherms with fewer barriers to dispersal, compared to terrestrial species 3,[14][15][16][17][18][19] . Moreover, the availability of thermal microrefugia is limited in the ocean, while in terrestrial ecosystems organisms can seek shade or burrow in soil to buffer the effects of warming 17,20 (Fig. 1). Therefore, biodiversity responses to temperature change are expected to be more immediate and directly detectable for marine ecosystems; this expectation is supported by a growing literature which has quantitatively compared individual species' responses to temperature change between marine and terrestrial realms 3,[14][15][16] . However, the net effects of temperature-related species' movements and abundance changes on assemblage-level diversity have not yet been systematically investigated across realms.
Here, we quantify temperature-related species richness and total abundance change in marine and terrestrial assemblages across temperate regions of the planet (23.5-60.0° absolute latitude; Extended Data Fig. 1). Specifically, we test two predictions for the effects of temperature change on assemblage-level diversity: (1) species richness and total abundance will increase with warming, and such increases will be greatest across relatively warm regions that border the species-rich tropics; and (2) the coupling of assemblage and temperature change will be tighter in the ocean than on land. These predictions are informed by the interaction between the latitudinal gradients in species richness and warming tolerance, and by the differences between realms regarding environmental heterogeneity and species distribution shifts (Fig. 1).
We focus on local assemblage-level trends, rather than species-specific responses, and quantify changes in both total abundance and species richness. These two metrics can be decoupled from each other, and abundance is typically more responsive to environmental change than richness 21,22 . We further disentangle richness change into species gains and losses to better understand the underlying dynamics of temperature-related biodiversity change. To test our expectations, we used the largest database of assemblage time series, BioTIME 23 , which includes studies of plants, invertebrates, birds, mammals and fish. These assemblages consist of co-occurring species systematically sampled through time. Since Temperature-related biodiversity change across temperate marine and terrestrial systems spatial extent varies among studies in BioTIME, we harmonized the biodiversity observations to a common spatial resolution to minimize the influence of variation in spatial extent on our results 24 ; this allowed us to quantify the effect of temperature change at a standardized resolution across regions and realms. We first estimated trends in biodiversity and temperature separately, and then quantified the relationships between the two.
Specifically, for each study we allocated individual samples to 96-km 2 hexagonal grid cells based on their location (see Methods; 24 ); that is, each sample was assigned to a specific combination of study ID and grid cell based on its latitude and longitude, resulting in equal-extent assemblage time series with multiple samples across years. Because each time series contained samples from only one study, the integrity of sampling methods within each study was maintained. We used these spatially harmonized assemblage time series in our analysis, selecting data from temperate regions only (since these are the better sampled regions within BioTIME). We then selected time series with at least 5 years of sampling (mean = 9.2 years), yielding 21,500 assemblage time series across both realms (19,875 25 to equalize sampling effort among years and then quantified trends in richness, total abundance and number of species gained or lost. For the same locations and for the time spans matching the duration of the biodiversity monitoring periods, we extracted mean monthly temperature records from HadCRUT4 26,27 and estimated the corresponding annual rates of sea surface or air temperature change. We then quantified the relationships between changes in biodiversity and changes in temperature with meta-analytical Bayesian hierarchical models, allowing responses to vary among taxonomic groups. To test our expectations, we included an interaction term between temperature change and long-term average temperature (that is, baseline climate) or latitude in our models, fitted separately for the marine and terrestrial realms.

Results
Temperature trends were highly variable, with locations at similar latitudes exhibiting different directions and magnitudes of change (Fig. 2a). Yet, both sea surface and air temperatures increased on average across the locations and time periods of our study, even though the majority of our time series spanned less than 10 years (Extended Data Fig. 2). The warming signal was more pronounced on land than in the ocean ( Fig. 2b; the average mean temperature change rate was 0.022 °C per year on land, versus 0.012 °C per year in the ocean).
Biodiversity change was also highly variable among the assemblage time series (Fig. 2c). However, despite the variability in both temperature and biodiversity trends, coherent macroecological signals emerged in the marine realm ( Fig. 3 and Supplementary Table 2). We found an overall positive relationship between warming and species richness, but conditional on the baseline climate. Species richness increases were more pronounced in initially warmer locations (as indicated by the positive interaction effect) and were underpinned by higher rates of species gains, while there was no detectable signal for species losses (Figs. 3 and 4 and Extended Data Fig. 3). Additionally, warming coincided with losses of individuals in the warmest marine locations, whereas cooler locations tended to gain individuals with increasing temperature (Figs. 3 and 4 and Extended Data Fig. 3). In contrast, no systematic biodiversity responses emerged on land, where the 95% credible intervals overlapped zero for all biodiversity metrics included ( Fig. 3 and Supplementary Table 2).
Our analysis highlights the fundamental role of climate baselines in modulating biodiversity responses in the ocean, given that latitude showed no or very weak interacting effects with temperature change (Extended Data Fig. 4). Further tests revealed only a weak correlation between either baseline climate or latitude and temperature change (Spearman's rank correlation, ρ = 0.063 and 0.098, respectively). As such, the correlation structure of our covariates is unlikely to bias our results towards a positive effect of temperature change when combined with baseline climate rather than latitude.
Lower temperate latitudes are poised to gain more species than higher latitudes due to higher species richness in the tropics  These stem from the latitudinal gradients in species diversity (decrease with latitude) and species thermal tolerance breadths or thermal safety margins (increase with latitude) (a); and differences between marine and terrestrial realms in terms of environmental heterogeneity and geographical range shifts (b). Note that this is a schematic simplified representation of these large-scale mechanisms and patterns (for instance, thermal safety margins are not necessarily linear with latitude 15 and heat tolerance declines more steeply with latitude for marine ectotherms 60 ).
Our results were highly robust to a number of further sensitivity tests addressing biases in the BioTIME dataset, differences amongst covariates and model structure. We found similar responses in the marine realm when using subsets of the marine data that match the number of locations and latitudinal range of the terrestrial data, despite increased uncertainty in model estimates (Extended Data Fig. 5). Our results were consistent across different baseline temperature variables (long-term annual and maximum temperature, and annual mean temperature of the first year) and between different temperature databases for quantification of climate baselines (Extended Data Fig. 4; see Methods). Additionally, we did not detect any systematic effect of the number of years sampled, temporal duration or start year of the time series on the estimated rates of biodiversity change (Extended Data Figs. 6 and 7). Finally, there was no evidence for systematic differences in biodiversity responses among the different taxonomic groups (as estimated via the random slope model components; Extended Data Figs. 8 and 9; complete model outputs with posterior estimates are shown in Supplementary Tables 3-6).

Discussion
We reveal striking differences in warming-related biodiversity change between marine and terrestrial realms across temperate regions, with a much stronger signature of warming on marine assemblages. Our results are unlikely to be due to confounding factors, given that climate change is poorly spatially correlated with other anthropogenic drivers of change for both marine and terrestrial realms 28 . Furthermore, temperature is expected to be a strong driver of biodiversity change given its fundamental role in biological processes [29][30][31] .
The warming-related increase in local marine species richness is consistent with the expectation that as the Earth's climate warms, temperate regions undergoing warming will receive an influx of species tracking suitable temperatures, and increases in warm-affinity generalists 2,6,7,18,29 . This warming-related increase in species richness is probably, in part, underpinned by species from lower latitudes shifting their ranges poleward to avoid exceeding their upper thermal limits 4,12,14 . Indeed, projected rates of extirpation in response to recent and future warming are highest for tropical marine species 7,32 , which in turn are better able to track changing isotherms 14 . Here, we find that species gains outpaced losses in the warmest temperate locations in the ocean where temperature has also increased. This finding is consistent with asymmetrical responses between species range edges, with faster colonizations expected following climatic niche expansion, and with slower extirpations linked to delayed responses at species' trailing edges 3,4 . The prevailing influx of species with warming is likely to lead to novel biotic communities and interactions as species distributions are reshuffled, with potentially far-reaching consequences for ecosystems' functioning 3,7,32,33 .

Fig. 3 | Biodiversity responses to temperature change and its interaction with baseline climate.
Marine locations (blue) exhibited stronger responses compared to terrestrial locations (green), while baseline climate (that is, long-term average temperature) modulated these responses in divergent directions. The main effects of baseline climate were negligible, and are not shown. Bars represent the estimated parameters and whiskers indicate the 95% credible intervals from the Bayesian meta-analysis (note the different scales for the different metrics); estimated parameters were considered to represent signals in the responses when the credible intervals did not include zero (Supplementary Table 2; see Supplementary  Tables 3-6  Our results further highlight that loss of individuals can occur simultaneously with increases in richness, emphasizing that temperature-related changes in richness and abundance can be decoupled. Abundance declines may be occurring for more thermally restricted species, owing to reduced performance and population declines as critical thermal thresholds are crossed 4-6 , for instance related to the adverse effects of increasing temperatures on metabolic rates and primary production 29,34 . Influxes of climate immigrants can also drive local declines in populations due to greater grazing and predation rates 35,36 . The abundance declines across our warmest temperate locations in the ocean suggest warming-related destabilization of populations, possibly reflecting reductions in the carrying capacity of marine ecosystems.
We also find strong interacting effects of the baseline climate with temperature change on both abundance and richness change in the ocean, which responded in different directions in our warmest temperate locations. This highlights that rising temperatures in locations that are already warm can lead to greater assemblage restructuring. Our findings may reflect the patchiness in temperature regimes across similar latitudes-for example, due to altitudinal or depth gradients, proximity to the coast or ocean currents, and may explain why latitude is a poor predictor for assemblage-level changes in our analysis. Specifically, within the scope of our data, baseline climate better captured context-dependent spatial patterns and allowed us to detect biodiversity change across pockets of warmer and colder areas within latitudinal bands. Baseline climate therefore emerges as a major predictor of temperature-related biodiversity change in marine systems, as recently also found by others 14,18 . Specifically, marine species track isotherms more readily in initially warm waters 14 , which may be due to tropical species overall living closer to their upper thermal limits, and thus more promptly responding to warming trends.
Overall, and despite faster warming on land, we did not detect systematic responses in richness or abundance with temperature change for terrestrial assemblages. The stronger responses observed for marine assemblages are consistent with reported faster range shifts in the ocean and higher sensitivity of marine organisms to temperature change compared to terrestrial species 3,12,14-16 . Our findings are also consistent with warming-related local extirpations being twice as common in the ocean as on land 15 . The lack of systematic assemblage-level change associated with temperature change on land might be due to greater thermoregulation capacity and wider thermal safety margins of terrestrial taxa 8,15,16 . Additionally, compared to seascapes, higher landscape complexity enables terrestrial species to exploit thermal microhabitats, thus allowing for the persistence of local populations for longer periods 8,20 . Indeed, access to thermal refugia was reported to be a fundamental factor underlying the relatively low vulnerability of terrestrial ectotherms to warming, and emphasizes the potential deleterious combined effects of warming and land-use changes 15 . Finally, a weaker link between assemblage responses and temperature change on land may be due to other factors, such as land-use change and moisture availability, posing stronger constraints on local biodiversity. Nonetheless, the smaller magnitude or slower responses of terrestrial species to temperature change 11,12,15,24,37,38 , combined with the faster rates of warming on land, indicate a potentially higher accumulation of climatic debt (that is, response lags) among terrestrial versus marine taxa 3,12,14,19,37,[39][40][41] . Additional research with higher-resolution temperature data matching the scale of organisms' responses is needed to better quantify terrestrial assemblages' responses to temperature change, and these developments remain a major challenge for many different taxa.
Our results warrant discussing some caveats. First, as commonly noted for biodiversity data, both taxonomic and geographic biases exist 23,42 , with the majority of our data coming from the Northern Hemisphere and disproportionately representing marine fish and terrestrial birds. Nevertheless, we did not detect systematic differences in responses among the taxonomic groups in our data.
Additionally, our ability to detect systematic trends on land may have been affected by the smaller number of terrestrial studies analysed, as illustrated by credible intervals often overlapping zero when using smaller subsets of marine data to estimate biodiversity change. The apparent inconsistency between our results and previous range shifts studies on land can be further reconciled by the different ecological levels at which such changes are being evaluated. As marine species track changing isotherms more closely than terrestrial species 14 , range shifts are more immediately translated at the assemblage-level in the ocean. On land, systematic effects on richness and abundance can be harder to detect given the general higher thermal tolerance and slower responses of terrestrial species.
Our analyses represent a first step towards explaining divergent patterns of assemblage-level biodiversity change across the planet 24,43 . Overall, our results provide strong support for the expectation that divergent biodiversity trajectories may arise as a consequence of global warming, with polar and temperate regions probably acting as 'sinks' and tropical regions as 'sources' 5,7,32 . In turn, these responses could prompt a shift in the latitudinal diversity gradient towards higher latitudes, with faster rates of change in the ocean. While we focused here on temperate regions, tropical and polar biomes are predicted to undergo severe restructuring in response to temperature change, albeit along different trajectories 5,7,12 . However, lack of sufficient biodiversity monitoring data for tropical and polar systems 23,42 hampers a comprehensive assessment of assemblage-level responses to temperature change in these regions, which remains a key challenge for further research.
Future global warming impacts on biodiversity are likely to exceed and potentially diverge from the changes revealed here 7,32,34,[44][45][46] . Indeed, initial increases in richness and abundance in response to warming may be followed by losses if warming continues 6,47,48 . Declines in marine systems may affect food security and livelihoods of human populations that are dependent on the ocean 6,49 . Additionally, while a consistent signal of temperature change was not evident on land, future impacts on terrestrial assemblages are expected from continuing rising temperatures, extreme heat events, fires and lack of moisture 34,45,50 . Because the Earth is committed to further warming, a systematic reduction of greenhouse gas emissions alongside efforts to further prevent habitat loss and improve habitat connectivity will be fundamental to allow species to track suitable climates across increasingly impacted land-and seascapes, and to avoid severe biodiversity disruption and loss.

Methods
Biodiversity data and trends. BioTIME 23 is currently the largest global database of assemblage time series, including 386 individual studies (study ID; plus extended data sources) across different taxonomic groups, holding >12 million records of abundance for >45,000 species. For this analysis, we included only those studies reporting counts of individuals per species in terrestrial and marine systems. We excluded freshwater studies because these are too few to confidently analyse biodiversity trends across taxa and different regions.
Each study is comprised of distinct samples (that is, individual plots, transects, tows and so on, sampled at a given time), and the number of samples can vary among years within each study. Because the spatial extent varies among studies, we gridded those studies that had large extents and multiple sampling locations into hexagonal cells of ~96 km 2 ; many studies were not partitioned because they were contained within a single cell 24 . Specifically, each sample was assigned to different combinations of study ID and grid cell based on its latitude and longitude, resulting in new assemblage time series (each with multiple samples across years). These new time series were given a unique identifier that was the concatenation of the study ID and the grid cell reference number, and thus contained samples from only one study-that is, the integrity of each study and each sample was maintained. This process allowed us to relate biodiversity and temperature trends at a standardized resolution. To minimize the effect of unobserved species on estimates of biodiversity change, we calculated the abundance-based coverage 51 of each annual sample within each time series and removed all samples with coverage <0.85. To enable reliable estimation of biodiversity trends, we restricted our analysis to time series sampled in at least 5 years (not necessarily consecutive). Because the number of samples can vary among years, we used sample-based rarefaction 25 to standardize the number of samples among years for each time series before calculating the biodiversity metrics. Specifically, we identified the minimum number of samples taken in each year within each assemblage time series; this minimum was then used to randomly sample each year down to that number of samples. Finally, given the paucity of data representing polar and tropical regions, we excluded these regions (based on absolute latitudinal cut-offs at 60° and 23.5°, respectively). This process yielded 21,500 assemblage time series representing 156 original studies (Extended Data Fig. 1 and Supplementary  Table 1) between 1900 and 2016, across 19,875 marine and 1,625 terrestrial locations. The average number of years sampled across the time series was 9.2, with the longest time series spanning 97 years (Extended Data Fig. 2).
To quantify rates of biodiversity change, we calculated linear trends over time for species richness (logS), total abundance (logN) and number of species gains and species losses. Counts of gains and losses retained species identity information, and were quantified based on comparison with the first year sampled in each time series. For losses, a positive slope means that the number of species lost from a location is increasing through time; negative slopes represent time series where the magnitude in species losses decreased over time. We repeated the sample-based rarefaction process described above 199 times for each time series, recorded the values and took the median for each biodiversity metric in each year, to reduce the effect of any outlier samples on our estimates. We used ordinary least-squares regression because we were interested in the long-term direction and magnitude of the biodiversity trends, and to allow us to compare the rates of change among locations, realms and metrics. We retained the estimated slope and standard error for each time series for use in our second-stage meta-analytical models.
Temperature data and trends. We focus on temperature as a climate variable because of its influence on every level of biological organization, from individual metabolic rates to biological communities' dynamics and structure [29][30][31] . We extracted temperature records from HadCRUT4 26,27 , specifically the HadSST3 data for marine Sea Surface Temperature (SST) on a 1° resolution, and the CRUTEM4 data for air temperature on land on a 0.5° resolution. We did not harmonize the spatial resolution between the two data sources because we wanted to use the best available data in each realm. For the location of each biodiversity time series, we extracted monthly mean temperature records for the duration of the biodiversity monitoring period (Year start :Year end ), and estimated mean temperature trends using generalized additive models, including a temporally autocorrelated error structure (package mgcv 52 ). This also allowed us to assess whether accounting for seasonality within years would improve model performance. We used AIC to compare models with and without 'month' , selecting the best model for each time series. We extracted the linear slope from the model, which summarizes the trend for mean annual temperature change.
To test whether biodiversity responses to temperature change were modulated by the baseline climate at any given location, we extracted annual mean temperature data from the WorldClim 53 database for terrestrial time series, and from the Bio-ORACLE database 54,55 for marine time series (on a resolution of 0.01° and 0.1° for terrestrial and marine systems, respectively). For each realm, we standardized the long-term annual mean temperature across all the locations by subtracting the mean and dividing by s.d.

Meta-analysis.
Having estimated the trends for biodiversity and temperature independently for each individual time series, we assessed the effect of temperature change on the rates of change of each biodiversity metric in a second-stage analysis. We employed a meta-analytical Bayesian framework using the package brms 56,57 (v.2.6.0) and fitted generalized linear models to each realm separately, having initially evaluated that there was an effect of realm when fitting a full model. All models were created using the Stan computational framework (http://mc-stan. org/) accessed via brms. To determine whether the baseline climate modulated the biodiversity responses, models were fit with an interaction term between temperature change and the long-term average temperature at each location. Additionally, we fitted similar models using latitude. We used two random effect terms: one allowing for different slopes per taxonomic group (Taxon) and another allowing for different intercepts per study ID nested within Taxon. This allowed us to account for: (1) potentially different responses to temperature change among taxa; (2) differences in species richness among taxa, as well as different assemblage time series originating from the same study, and different studies monitoring the same taxonomic groups across the BioTIME database, respectively; and (3) spatial autocorrelation. The different taxonomic groups were informed by the original data sources metadata, and were: ' Amphibians' , 'Benthos' , 'Birds' , 'Fish' , 'Mammals' , 'Marine invertebrates' , 'Terrestrial invertebrates' , 'Terrestrial plants' and 'Multiple taxa' . The last of these refers to study IDs where the sampling methodology captured multiple taxa simultaneously. Additionally, taxonomic groups with only very few studies in BioTIME (for example, reptiles) were included in this group for practical reasons; we note that this category represented only 9% of the original studies and 5% of the standardized time series.
The overall model structure implemented for each realm was: ΔBiodiversity  0 þ ΔTemperature´Long�term average temperature þð0 þ ΔTemperaturejTaxonÞ þ ð1jTaxon=Study IDÞ: The meta-analytical framework implemented allows the response variables to be weighted by their uncertainty; specifically we included the estimated standard errors for the biodiversity trends using the function se() when specifying the formula for the models, using the brms syntax 56,57 : ΔBiodiversity | se(ΔBiodiversity, sigma = TRUE). Additionally, a global intercept equal to zero was specified to reflect a scenario where there is no biodiversity change when all the independent variables are zero; that is, the model assumes that no change in temperature would correspond to no change in the assemblage metrics for the average baseline climate (that is, at the mean value across our data).
Models were run using four chains, each with 8,000 iterations, with a warm-up of 4,000 and non-informative flat priors. Stan implements Hamiltonian Monte Carlo and its extension, the No-U-Turn Sampler (NUTS) algorithms, which converge quickly 57 . Convergence was assessed by visually examining trace plots and using Rhat values (the ratio of the effective sample size to the overall number of iterations, with values close to one indicating convergence) 57 . All analyses were run in R v.4.3.1 58 .

Sensitivity analysis.
To evaluate the robustness of potential interactions with the baseline climate, we additionally ran our models with two alternative baseline temperature variables. To that end, we extracted the variables 'Mean Temperature of Warmest Quarter' from WorldClim and 'Long-term maximum sea surface temperature' from Bio-ORACLE, as well as the average temperature in the first year sampled for each biodiversity time series from the same dataset that was used to quantify the trends (that is, the spatially less resolved HadCRUT4 dataset).
To evaluate whether uneven sampling could be driving the observed differences between the marine and terrestrial realms, we fitted models to subsets of the marine data that matched both the number of locations (1,625 time series) and the latitudinal range of the terrestrial data. We did not attempt to control for temperature change differences between realms because this is part of the signal to be modelled. We fitted the meta-analytical models to 100 random subsamples for each biodiversity metric, illustrating that the estimates for the marine realm are robust (Extended Data Fig. 5). This analysis also illustrated that a smaller sample size may have affected our ability to detect systematic biodiversity changes on land, given that confidence intervals overlapping zero occurred for several of the meta-analytical models fit to subsets of the marine data matching the number of terrestrial time series. We further evaluated that biodiversity responses did not show any clear pattern as a function of the number of years sampled, duration and start year of the time series. Although shorter time series tend to show higher rates of biodiversity change, these also tend to have higher standard errors, which are in turn accounted for within the meta-analytical framework. We further note there were no notable differences in the estimated rates of biodiversity change as a function of number of sampled years between marine and terrestrial time series. Thus, we believe the temporal heterogeneity of sampling in our data is very unlikely to drive our findings (Extended Data Figs. 6 and 7).
Reporting Summary. Further information on research design is available in the Nature Research Reporting Summary linked to this article.

Code availability
Code to reproduce the analysis is available on GitHub (https://github.com/ lauraantao/Temp_Biodiv_Change) and in an online archive at Zenodo 59 .

October 2018
Corresponding author(s): Laura H. Antão Last updated by author(s): Mar 16, 2020 Reporting Summary Nature Research wishes to improve the reproducibility of the work that we publish. This form provides structure for consistency and transparency in reporting. For further information on Nature Research policies, see Authors & Referees and the Editorial Policy Checklist.

Statistics
For all statistical analyses, confirm that the following items are present in the figure legend, table legend, main text, or Methods section.

n/a Confirmed
The exact sample size (n) for each experimental group/condition, given as a discrete number and unit of measurement A statement on whether measurements were taken from distinct samples or whether the same sample was measured repeatedly The statistical test(s) used AND whether they are one-or two-sided Only common tests should be described solely by name; describe more complex techniques in the Methods section.
A description of all covariates tested A description of any assumptions or corrections, such as tests of normality and adjustment for multiple comparisons A full description of the statistical parameters including central tendency (e.g. means) or other basic estimates (e.g. regression coefficient) AND variation (e.g. standard deviation) or associated estimates of uncertainty (e.g. confidence intervals) For null hypothesis testing, the test statistic (e.g. F, t, r) with confidence intervals, effect sizes, degrees of freedom and P value noted Give P values as exact values whenever suitable.

For Bayesian analysis, information on the choice of priors and Markov chain Monte Carlo settings
For hierarchical and complex designs, identification of the appropriate level for tests and full reporting of outcomes Estimates of effect sizes (e.g. Cohen's d, Pearson's r), indicating how they were calculated Our web collection on statistics for biologists contains articles on many of the points above.

Software and code
Policy information about availability of computer code Data collection All the data can be accessed through the BioTIME database on Zenodo (https://doi.org/10.5281/zenodo.1211105) or through the BioTIME website (http://biotime.st-andrews.ac.uk/). Studies that were not published with the database can be accessed via the primary references provided.

Data analysis
All analyses were performed in the R software v 4.3.1; specifically we used the brms package for the modelling (version 2.6.0). Code to reproduces the analysis is available on GitHub (https://github.com/lauraantao/Temp_Biodiv_Change) and deposited in Zenodo (http:// doi.org/10.5281/zenodo.3708095).
For manuscripts utilizing custom algorithms or software that are central to the research but not yet described in published literature, software must be made available to editors/reviewers. We strongly encourage code deposition in a community repository (e.g. GitHub). See the Nature Research guidelines for submitting code & software for further information.

Data
Policy information about availability of data All manuscripts must include a data availability statement. This statement should provide the following information, where applicable: -Accession codes, unique identifiers, or web links for publicly available datasets -A list of figures that have associated raw data -A description of any restrictions on data availability All the data used in the manuscript can be accessed through the BioTIME database on Zenodo (https://doi.org/10.5281/zenodo.1211105) or through the BioTIME website (http://biotime.st-andrews.ac.uk/) (plus extended data sources). Studies that were not published with the database can be accessed via the primary references provided.