dataesr/scientific-paragraphs-categorization huggingface dataset:
While future studies using phylogenomic approaches may address these issues, obtaining tissue samples for all these taxa will remain an obstacle for a long time. Although DAISIE does not directly use topological information (only divergence times are used), it is possible that the true topology for a clade may differ from that of the gene tree we have estimated and this could have an impact on our results by a) affecting colonisation and branching times (addressed in the paragraph below); or b) by altering the number of colonisation events. For instance, species that appear to have colonised an archipelago only once may have colonised multiple times. As with any phylogenetic study, we cannot rule out this possibility, but we assume that recolonisation of any of the archipelagos in our dataset by the same species is very rare, as these are all oceanic and isolated.
For archipelago lineages with cladogenesis (26 out of 502 lineages), alternative topologies could include non-monophyly of island radiations, with the corollary being that they would be the results of multiple colonisations. However, non-monophyly of these well-studied and highly isolated insular radiations is highly improbable. Regarding scaling of divergence times, we assessed how uncertainty in our estimated node ages could influence our results by running an analysis of 100 datasets. For each dataset we sampled the node ages (i. e. colonisation and branching times) at random from a uniform distribution centred on the posterior mean for that node in the Beast tree and extending twice the length of the highest posterior density (HPD) interval. For example, for a node with a 95% HPD interval of 2-3 Ma in our trees, the uniform distribution was set to between 1.
5 and 3. 5 Ma. The HPD interval will capture uncertainty under the selected phylogenetic and substitution models for the loci we used, but we conduct our sensitivity analysis over a broader range to accommodate the potential that the selected models and gene trees are inadequate. For cases where using this approach meant that the lower bound was lower than 0, we assigned a value of 0. 00001 Ma to the lower bound. We fitted the nine best models to the 100 datasets, using five initial starting parameters for each model (total 4,500 optimisations). We found that parameter estimates across the 100 datasets do not differ strongly from those in the main dataset (Table  XX ). Importantly, model selection was unaffected, with the M19 model being selected for all 100 datasets. This is because a lot of the information used for model selection is coming from the other sources of information DAISIE uses (island age, number of species, endemicity status) rather than colonisation/branching times.
To account for uncertainty in rates of molecular evolution, we repeated all BEAST dating analyses for markers that were not cyt-b using 1) the Weir and Schluter  50  cyt-b rate (dataset D1) and 2) marker specific rates estimated by Lerner et al. 41  , which are also widely used in the literature (dataset D2). Although the trees dated using the Lerner et al. rates provide younger ages, we found that the DAISIE results were very similar using either approach (same model preferred and similar parameters). Therefore, in the main text we only discuss the results of analyses of D1, i. e. applying the cyt-b rate of Weir and Schluter to all markers. For some taxa we did not use the stem age as the estimate of colonisation time, and instead used alternative nodes (see 'Colonisation and branching times' section).
To test whether our choice of nodes affects our main conclusions, we recoded all such taxa by extracting the stem ages, and used these ages as an upper bound for colonisation (DAISIE MaxAge option). We fitted all 28 models to this new dataset (D3) and found that the M19 model is preferred and that the parameters vary only slightly from those of the main analysis. We therefore conclude our results are robust to the node selection approach. If extinction has been high on the mainland, or if we failed to sample the closest relatives of the island taxa, this could lead to an overestimation of colonisation times when using the stem age as the precise time of colonisation. To investigate how this could have influenced our results, we ran analyses of datasets where we allow colonisation to have happened at any time since the stem age (i.
e. the time of divergence from the nearest relative of the taxon on the mainland). For this we use the DAISIE options Endemic_MaxAge or NonEndemic_MaxAge, which integrate over all possible ages between the given maximum age and the present (or the first branching event within the archipelago for cases where cladogenesis has occurred). We repeated this analysis coding all stem ages as maximum ages (D4), or coding only the 25% older stem ages as maximum ages (to account for the fact that on older stems there is potential for there to be more bias) (D5). We also ran analyses on 100 datasets (D6) for which we assigned younger ages by randomly selecting a value between the stem age and the present (or crown age for cladogenetic groups). For all these datasets (D4-D6) we found that the same model (M19) is preferred, but the initial values of the biogeographical rates (cladogenesis, extinction, colonisation and anagenesis) are estimated to be higher than in the main dataset.
Sensitivity to archipelago selection and isolation metrics The results of the following sensitivity analyses are presented in Supplementary Data 3 and the DAISIE R objects including these alternative datasets have been deposited in Mendeley. To test whether the inclusion of both true archipelagos and single islands in our dataset could affect the results, we repeated analyses excluding single island units and found that the same model is preferred. The estimated initial rate of cladogenesis (λ c 0) is higher if we exclude single islands, but this parameter is not different from a distribution of parameters estimated from datasets generated using a stratified-random sampling of both archipelagos and single islands. Alternative isolation metrics to Dm have been shown to explain varying and often higher amounts of variation in species richness on islands  60.
We tested two alternative metrics: distance to the nearest larger or equivalent sized-landmass (Db), and the mean between Dm and Db (metrics given in Supplementary Data 2). We found that the same DAISIE model with very similar parameters was preferred in both cases, and thus we used only the Dm metric, as this is more similar to the original model of MacArthur & Wilson. The Mascarenes (Mauritius Isl. , Reunion and Rodrigues) are often treated as a single biogeographical unit in analyses. We chose to analyse them as independent units because a) the distance between islands is much greater than our threshold for archipelago definition (more than 500 km between Mauritius Isl. and Rodrigues; more than 170 km between Reunion and Mauritius Isl. ); b) only two species of our target group are shared between the islands (Terpsiphone bourbonnensis found in Mauritius Isl.
S. E. by a  NOW   VICI grant ; M. M. by the  Portuguese Science and Technology Foundation  (Post-doctoral grant:  SFRH/BPD/100614/2014 ); S. M. C by the  National Geographic Society  ( CRE grant  #  9383-13 ); J. C. I. by a  GRUPIN research grant  from the  Regional Government of Asturias  (Ref. :  IDI/2018/000151 ); ); C. T. by the  'Laboratoire d'Excellence' TULIP  ( ANR-10-LABX-41 ).
100 datasets with variable ages Parameter Main dataset 2. 5 percentile 97. 5 percentile λ c 0 0. 040 0. 035 0. 043 d0 0. 026 0. 025 0. 027 μ0 1. 946 1. 840 2. 290 x 0. 150 0. 150 0. 172 γ0 67. 22 60. 99 69. 99 a 0. 294 0. 291 0. 297 λ a 0 0. 059 0. 028 0. 075 b 0. 383 0. 355 0. 495.
Games for entomological literacy and their relationships with sustainable development Using games for teaching and learning, while not the only solution, has several advantages, such as effective learning outcomes, personal and social development, and user-centred learning  [14, 15]. To get an overview of how research has approached the topic of games for entomological literacy, we conducted a search in the Web of Science on all the documents in all of the available databases of the platform for all the years until November 2019, using the following search terms: (game OR serious game OR game science) and (education OR learning OR teaching) and (insect OR entomology). After reading all the abstracts and removing erroneous entries and articles from the same study, we were left with only 12 records  [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27]  (Figure  1 ).
Advancing research into entomological literacy games to support sustainability While interest in research into games for learning has been increasing in recent years  [15] , we are interested to know how much of this interest is directed towards insect science. This may reveal the lack of research in this area, and therefore highlight the opportunity for further exploration. To this end, we conducted a search in the Web of Science platform to obtain the number of records on game research related to entomological literacy, and use that as our proxy indicator. Our search was on all the documents in all of the databases for all the years until November 2019 using the following search terms: (game OR serious game OR game science) and (education OR learning OR teaching). The query results yielded 51,448 records across 155 research areas, with 709 unclassified records (Figure  3 ).
Study Population and Characteristics Barcelona, the second largest city in Spain  (1, 650, 358  inhabitants in January 2019 by Municipal Population Register of Barcelona City Council), is located on the northeastern coast of the Iberian Peninsula. Barcelona has 10 districts, and a map of the 10 districts is provided in the Supplementary Materials. Barcelona's climate is Coastal Mediterranean. Winters are short and mild, while summers are long and hot. The average annual temperature is between 15 • C and 18 • C. The mortality data were collected from the mortality registry of Barcelona, which collects information from the mortality registry of the Department of Health of Catalonia and from the municipal population registry. The dataset includes individual all-cause (natural) mortality records for the residents in the city of Barcelona.
The individual mortality data were collected for the period of 1 January 1992 to 31 December 2015. Only the summer months' deaths were considered in the analysis (1 June-30 September). The detailed information of each individual death includes date of death (day, month and year), age (through date of birth) and sex (men and women). Then, as a measure of socio-economic position, we used education level and grouped them into two categories (without studies and with primary education or higher), the data for which were obtained from the population register of the city, after excluding the deaths from the youngest age groups (<13 years old). Also, residence of deaths was georeferenced, both address (latitude and longitude) and district of residence. Residents of Barcelona who died outside the city were excluded from the analysis.
Urban Climate Model (UrbClim) The UrbClim model was designed to simulate temperature and heat-stress at the city scale with a minimum amount of computation power  [21]. This model is composed of a simple urban physics, coupled with a 3-D atmospheric boundary layer module. In the land surface scheme, urban terrain is represented as an impermeable slab with appropriate parameter values for albedo, emissivity and aerodynamic and thermal roughness length, and also accounting for anthropogenic heat fluxes  [16]. UrbClim combines information about the urban structure (vegetation, soil sealing, typology, land use and land cover) through a simplified physics approach relevant at urban scales, and generates local climate data at very high spatial resolution of 100 m. The location and time were considered when climate data were extracted from the UrbClim model in order to merge with the georeferenced mortality dataset.
Exposure Definition and Lags Structures The exposure data were extracted using time and location information (geographical coordinates of residence address of the deceased and date of death) and calculated daily mean summer temperature percentiles (5th to 95th) for each death case and its controls. The exposure variable group was defined based on the 70th-100th percentile (preliminary analysis results showed adverse health effects at 70th percentile and above) of summer mean temperature and reference group was 0-30th percentile (minimum mortality percentile) of summer mean temperature. We categorized days in the 70th percentile and above as hot days and the days between the 0 to 30th percentile as a non-hot day. We have used the daily mean temperature because it is a better predictor of temperature-mortality relationships when compared with either maximum or minimum temperature  [23].
Results We collected mortality data by age, sex and education for the periods of 1 June-30 September. The dataset included 105,559 all-cause mortality (48. 84% among men and 51. 16% among women) given in Table  1. The characteristics of the temperature data of the UrbClim model presented in Table  2  and the 70th percentile of mean temperature (>27 • C) and 30th percentile (<24 • C) of mean temperature are compared.
There are a few limitations we faced in modelling the risk that merit discussion. Previous validation exercises have shown that UrbClim, coupled to reanalysis data, can sometimes show a little warm bias in the modeled air temperatures Therefore, it is important to recognize the specific limitations both in the analysis and the interpretation of results obtained by using model data  [21, 39]. The sociodemographic variables were derived from data available at a fixed point in time (the year 2012 until 2016). However, the underlying demographics data ordinances both changed over time, a process we were unable to capture using meta-regression approach. This introduces some uncertainty into the results, and future research should explore local-scale mortality patterns over both space and time using the Bayesian statistical approach.
Protected areas (PAs) are the cornerstones of global biodiversity conservation efforts, but to fulfil this role they must be effective at conserving the ecosystems and species that occur within their boundaries. Adequate monitoring datasets that allow comparing biodiversity between protected and unprotected sites are lacking in tropical regions. Here we use the largest citizen science biodiversity dataset -eBirdto quantify the extent to which protected areas in eight tropical forest biodiversity hotspots are effective at retaining bird diversity. We find generally positive effects of protection on the diversity of bird species that are forestdependent, endemic to the hotspots, or threatened or Near Threatened, but not on overall bird species richness. Furthermore, we show that in most of the hotspots examined this benefit is driven by protected areas preventing both forest loss and degradation.
Such location biases are inevitable because protected areas tend to be designated in regions of little economic interest (i. e. greater remoteness, higher altitudes, and lower agricultural suitability  6, 7  ), which are less likely to have suffered from human pressure both before and after protection. These differences can be statistically controlled for in counterfactual analyses of protected area effectiveness  [8] [9] [10]  , but this requires large datasets on the spatial distribution of the biodiversity features of interest across many protected and unprotected sites. Nowhere are effective protected areas more essential than in tropical regions, which host a disproportionately large share of the world's biodiversity  11  and face rapid habitat loss  11  and degradation  12, 13  , both major threats to biodiversity  [13] [14] [15].
Yet, evaluating protected area effectiveness in these regions is particularly challenging, given that the detailed biodiversity datasets required for counterfactual analyses are typically unavailable  16. Among the few analyses investigating biodiversity outcomes of tropical protected areas, most focused on protected area effects on land cover, finding that they mitigate both forest loss and forest degradation  10, [17] [18] [19]. While such analyses are possible using ex situ remote sensing data, investigating effectiveness in terms of species outcomes requires data collected in situ. Two global metaanalyses reviewed local-scale studies that had contrasted protected versus unprotected sites in terms of species diversity  9, 20. Both uncovered positive effects at the global scale, but-worryingly-weaker or mixed results in tropical regions, contrasting with reported positive effects of protected areas at reducing forest loss and degradation.
However, they contrast with what was known in tropical regions from two previous globalscale studies of protected area effectiveness, based on the metaanalysis of local-scale studies contrasting protected versus unprotected sites. One of these studies found higher species richness and abundances within protected areas in Africa and Asia but not in South America 20 ; the other found higher overall richness within protected areas, but no significant effects on the richness in rare and endemic species, including in the tropics  9. Nonetheless, the present study provides a stronger test of protected area effectiveness in tropical forests by focusing specifically on these biomes, using more comparable data (as emerging from a single, coherent dataset), and by exploring the underlying mechanisms of habitat loss and degradation.
Methods Study areas: biodiversity hotspots. We focused on eight biodiversity hotspots  21  : those with at least 25% of their extent within the "tropical and subtropical moist broadleaf forests" biome  44  and for which we obtained at least 1000 checklists from eBird (after applying the data selection procedure described below): Atlantic Forest, Tropical Andes, Tumbes-Chocó-Magdalena, and Mesoamerica (Americas); Eastern Afromontane (Africa); Western Ghats and Sri Lanka, Indo-Burma and Sundaland (Asia). Within each hotspot, we analysed only areas overlapping the "tropical and subtropical moist broadleaf forests" biome  44  (Fig. 1 , Supplementary Figs. 1, 4,  and 5 ), assumed to have been originally forested (see Supplementary Methods 4d). Data selection: eBird checklists. We obtained bird sightings from the eBird citizen science database  23.
The reporting system is based on checklists, whereby the observer provides: list of birds detected; GPS location; sampling effort (whether or not all detected species are reported; sampling duration; sampling protocol, e. g. , stationary point, travel, and banding; and distance travelled in case of travelling protocol); starting time of the sampling event; and number of observers. We used the eBird dataset released in December 2018  45  , focusing on records from 2005 to 2018, as data collected prior to 2005 were too scarce for analysis. We filtered this dataset to obtain high-quality checklists comparable in protocol and effort: we selected complete checklists only (i. e. in which observers explicitly declare having reported all bird species detected and identified); following either the "stationary points" or the "travelling counts" protocol; with durations of continuous observation of 0.
5-10 h; with observers travelling distances during the checklist < 5 km; only from experienced observers (≥10 checklists; ≥30 species per checklist on average; ≥100 different species in total); and removing potential duplicates (checklists made on the same day at the same place). We applied some taxonomical transformation to eBird data in order to fit with BirdLife International taxonomy (Supplementary Methods 1f). After data filtering (more details in Supplementary Methods 1), we obtained the final dataset used in the analyses, consisting of 66,777 checklists, covering 5467 species, from 6838 observers, in eight hotspots (Supplementary Figs. 4  and 5 ; Supplementary Table  2 ). Site characteristics. Our analyses included two types of sites: checklist sites, corresponding to the coordinates of each eBird checklist analysed (used in analyses I and III); and background sites, corresponding to the centre points of a regular grid of 2 × 2 km covering the whole area of each hotspot evenly (used in analysis II).
e. Vulnerable, Endangered, or Critically Endangered  25  ; more details in Supplementary Methods 2). Index of observer expertise. Heterogeneity in observers' birding skills, behaviours, and equipment increases data variability and potentially introduces biases to the analyses  55, 56. Heterogeneity is particularly high in citizen science datasets like eBird, where volunteers range from those only familiar with a few common local birds to experienced observers capable of detecting rare and cryptic species. As stated above, we only included checklists from relatively experienced observers. To account for the remaining variability in observer expertise, we calculated an observer expertise score (used as an explanatory variable in the statistical analyses), adapted from Kelling et al. 57  and from Johnston et al.
56  , and calculated separately for each continent. It estimates the variation in the number of species that observers are predicted to detect in similar conditions. To do so, we first ran a mixed general additive model (function gamm from "mgcv" R package  58  ) modelling species richness of checklists against potential confounding variables that are expected to affect either the number of species detected (sampling protocol; n. observers number of observers; duration of sampling; time of the day) or the true species richness (lat latitude; lon longitude; and Julian day), adding observer as a random effect gamm richness protocol þ n:observers þ s duration ð Þþs time ð Þ ð þ te lon; lat; day ð Þ þ random ¼ list observer 1 ð Þ Þ : ð1Þ The notation s() indicates that the variable was used as a smoothed term; and te() indicates that the variables have been used as interacting smooth terms, allowing here species richness to vary spatially during the year.
After fitting this model to each continental data subset, we used it to predict the logarithm of species richness that each observer would report for a fictive stationary point with all variables fixed to their median values. This resulted in an observer expertise score that we then assigned to all checklists; assigning the observer score of the observer with the highest expertise score in cases of multiple observers. This index ranged from 2. 2 to 4. 3 in Africa, from 2. 3 to 4. 4 in the Americas, and from 2. 8 to 4. 5 in Asia (more details in Supplementary Methods 3). Statistical analyses of protected area effectiveness. We investigated protected area effectiveness at retaining bird diversity through a set of three connected statistical analyses (Fig. 2 ), undertaken separately for each hotspot, using GAM models  58.