dataesr/scientific-paragraphs-categorization huggingface dataset:
Once there, biologists on board registered the location (GPS) and time of the sighting, as well as species, group size and composition, behaviour and other observations. As research was not the main purpose of the activity, collection of standardized data, always registered manually, was restricted so as to interfere as little as possible with the main goal, tourism. Therefore, routes were not recorded, but starting and finishing times of each trip were recorded. Only those sightings for which the species identification was confirmed (genus accepted for beaked whales) and a reliable location recorded (longitude, latitude) were considered for the analyses (Fig  2 ). Following this process, around 93% of the sightings data were retained for analysis. The study period ran from May 2008 to December 2014.
Effort As animals were spotted firstly from land and sightings at sea are thus not associated with the amount of boat-based search effort during a trip, it is not possible to calculate real search effort or blue whale sightings per unit effort for this dataset. However, we recorded the number of trips, trip duration and time spent at sea each season and year as a rough indicator of the seasonal and year-to-year distribution of effort. As described below, we used sightings of other cetacean species as "absences" for blue whales.
Environmental variables The main environmental variables used for the analyses are summarized in Table  1. We obtained depth data from the General Bathymetric Chart of the Oceans (GEBCO, 2012) with a resolution of 30 arc-seconds. Slope was derived in degrees from the depth raster using the Surface Analysis tool from IDRISI Kilimanjaro. Each sighting was linked to its corresponding depth and slope values using Matlab R2012a. Distance to the coast from each sighting was calculated with ArcGIS 9. 3, using as a reference high resolution coastline data provided by the Instituto Hidrográfico de Portugal in 2011. The distribution of a given species usually depends on the distribution of its prey. Oceanographic variables such as temperature or chlorophyll can act as proxies for primary production and prey availability or distribution.
The oceanographic dynamic variables considered were prepared at three different temporal scales: daily, weekly and monthly. The products were chosen according to the availability of homogeneous data over the studied period. All the dynamic variables were processed with Matlab R2012a. Data on Sea Surface Temperature (SST) were gathered from two different products, one with high spatial resolution to get a detailed view of the study area and better account for the local oceanographic features (Fig  3 ), and one with coarser resolution that simplifies data processing and reduces data loss. The latter was derived from Operational SST and Sea Ice Analysis (OSTIA), which is run at the UK Meteorological Office on a daily basis. It combines high resolution products from infrared and microwave satellite instruments and in situ SST data, resulting in a foundation SST which avoids diurnal variability.
A climatological mean (daily, weekly and monthly mean) was also calculated over the area for the 7 years of study to account for seasonality. A SST anomaly was calculated as the difference between the SST at the sighting location and the corresponding climatological mean for each observation. This allows us to account for local spatial or temporal variations. Furthermore, location of thermal fronts were calculated from SST images (data gaps interpolated) using a Canny edge detection function  [57]  with the upper threshold set at 1e-5˚C/m (~1˚C/ 100 km). Distance, gradient and SST of the nearest front point were calculated for every sighting. Chlorophyll-a concentration is considered as a proxy for primary production. As blue whales feed mainly on krill, primary production should directly affect their prey availability with a relatively short time-lag.
Krill development varies with environmental conditions, but juvenile stages can be reached in less than three (or even two) months after spawning, which usually happens after the phytoplankton bloom  [58]. Chlorophyll data were obtained from GlobColour (http://globcolour. info) L3 products that combine, when available, data from Sea-WIFS, MODIS, MERIS and VIIRS under a weighted average merging method, with weightings based on the sensor/product characterisation. As chlorophyll concentration is a dynamic variable, and its distribution can change in a short period of time, we decided to use daily, 8-day and monthly (processed from the 8-day ones) composites with 1km resolution. For each temporal resolution, we extracted a value for the sighting location, and a mean value and standard deviation over the Azores archipelago (35-42˚N, 33-23˚W), São Miguel (37-38.
Wind data were obtained from the European Centre for Medium-Range Weather Forecasts (ECMWF) ERA-Interim 0'5 6-hourly surface analysis, which combines an assimilation system, constrained by the available in situ observations, and atmospheric modelling re-analysis. Its gridded resolution is 0. 5˚(~54 km) but the effective resolution is around 79 km  [60]. Wind data were obtained only for the day and location of the sighting, in order to avoid misleading temporal and spatial averages.
Variable selection After preparing all the variables for the different spatial and temporal scales we had six main datasets: low spatial resolution monthly (47 variables), weekly  (57)  and daily  (59) , and high spatial resolution monthly (41 variables), weekly  (60)  and daily  (60). A detailed data exploration was carried out to assess collinearity among explanatory variables. We computed Pearson correlation matrices for each dataset in order to investigate the linear dependence between the predictors. These were visualized as dendrograms to better illustrate the similarity among the variables. The bigger the vertical distance between two variables, the more independent they are. We retained all the variables with a distance bigger than 0. 8 from other variables  [48, 61]. The variables with distances under this value were selected firstly according to their horizontal distribution, retaining at least one of each cluster of correlated variables.
When two highly collinear variables had similar ecological meanings, the one with most missing values was removed. For daily datasets, chlorophyll mean and standard deviation values from south of São Miguel were kept instead of those from around São Miguel to better account for local events. After this first selection, collinearity was assessed again based on the Variance Inflation Factor (VIF), removing all the variables whose VIF was higher than 5, one by one, starting with the highest value and recalculating VIFs at each step  [62]. Finally, some additional variables were excluded from the high resolution selections (according to the dendrograms) due to greatly reducing (by hundreds and even thousands) the sample size available. After this process our datasets comprised of 19 (low spatial resolution daily), 18 (low spatial resolution weekly), 17 (low spatial resolution monthly), 17 (high spatial resolution daily), 19 (high spatial resolution weekly) and 14 (high spatial resolution monthly) explanatory variables.
Temporal distribution To assess the temporal distribution of whale sightings, taking into account the survey effort, we calculated an Encounter Rate (ER) as the number of sightings of blue whales divided by the total number of trips over each month or season. There were normally two trips per day (morning and afternoon). We consider the season as the period that included 99% of the blue whale sightings, hence March to June every year. We did not use data from 2008 in seasonal averages since data onboard were recorded only in surveys conducted during six months in 2008: from May to September and December.
Habitat preferences analysis Descriptive analyses were carried out in order to visualize available habitat conditions within the study area and understand distribution of the blue whale sightings in relation to habitat variables. The available area is considered as the Minimum Convex Polygon containing 100% of the sightings calculated with the "sp" and "adehabitatHR" packages for R  [63, 64, 65]. It accounts for a total extension of 5780 km 2. We determined the range of available values for the static variables in the estimated available area, although it should be noted that the true surveyed area is probably slightly larger. To model species distribution, available environmental conditions need to be represented in the analysis  [66]. Usually these "background" data are chosen at random to represent the range of environmental conditions in the study area, but this process can fail to account for sample bias, as sightings are generally recorded in more accessible areas  [67].
Furthermore, the changing marine environment can present different conditions in the exactly same location at different times, therefore creating a different habitat. We decided to use a pseudo-absence approach in order to have effort-related 'non-occurrence points' and to reduce possible spatial bias  [66, 68]. To do so, we use the presence of the other cetaceans recorded during the study period as pseudo-absence points for blue whales  [69] , making sure to exclude occasions when blue whales were recorded at the same time as the other species. We thus ensured that those locations were sampled and that blue whales were not seen there at the given time, while accounting for, or at least reducing, sampling bias. The approach used here could be considered conceptually closer to the background sampling approach, as the number of "pseudoabsences" is, in this study, much bigger than the number of presences and the pseudo-absence dataset well characterizes the environmental conditions present across the study region  [66].
Nevertheless, environmental conditions that were unsuitable for launching the boat (more prevalent in winter) would not have been included in the dataset or, therefore, our models. Statistical analyses were carried out with the open-source software R 2. 15. 0 (R Development Core Team, 2012) and the mgcv library  [70]. We applied Generalized Additive Models with a binomial distribution and a logit link function to each of the six datasets previously prepared in order to investigate habitat preferences of blue whales. The number of splines (k -knots which measure the complexity of the fitted curve-) was set to 4 for all explanatory variables to prevent overfitting and avoid unrealistically complex relationships. We used a backwards stepwise procedure to select the best model for each dataset, with removal of variables at each step being based on individual significance of explanatory variables and on the overall goodness of fit as measured by the Akaike Information Criterion (AIC) value.
Starting from the first model with all possible variables, the least significant one (normally with highest p-value) was dropped out each time. If the resulting model had a lower AIC or a bigger sample size (and in the latter case was not a markedly poorer model according to the AIC value), it was retained. In some cases confidence bands were rather wide, even in the model with the lowest AIC, so extra model selections steps were needed  [71]. In order to avoid a number of predictors higher than m/10, where m is the number of observations of the least represented category (here m = 89)  [72] , we retained in the final model only the variables significant at p<0. 01. The best model for each dataset was achieved and the deviance explained was noted. We applied a temporal k-fold cross-validation for each of the resulting models  [73] , using sightings from six of the seven available years as the training dataset, and the remaining one for evaluation (all but 2008 because no blue whales were recorded that year during the six months of data available).
Then we compared the performance of each final model based on the Area Under the Curve (AUC) of the Receiving Operating Characteristic plot computed with pROC package for R  [74]. AUC indicates how well the model adjusts to our presence/ absence distribution without accounting for overfitting, so that more complex models usually have better AUC  [75]. Hence, as described above, we attempted to avoid very complex models to carry out validation. We obtained an AUC for each year (i. e. for each option for the testing dataset) and an average AUC value and its corresponding standard deviation for each model, allowing us to compare the different models. We consider that a random model has a 50% chance to correctly distinguish between presence and non-presence locations, so we assume that any model with AUC higher than 0.
88) and the AUC standard deviation the highest (0. 100) (Table  2 ), which means that for some years this monthly model fits quite well, while for others it fits poorly, leading to substantial differences in AUC among the years considered within the model. Distance to the coast was selected in all models, both high and low resolution, always with the same shape of the fitted curve, which indicates a preference for waters further than 10km from the coast (Fig  9 ). To put this in context, observers on land could detect whale blows up to around 30 km offshore but all sightings used in the analysis were from the whale watching boats, which normally travelled up to 30 km from the shore, recording most of the sightings (77. 7%) within 15 km from the coast (  Fig 10). Therefore, although data might present bias towards shorter distances from the coast due to the initial searching for animals being done from land, blue whales sightings are among those recorded furthest away from land.
Discussion Modelling results can be evaluated from two main perspectives: statistical and ecological performance. The former seeks a robust empirical description of the existence, strength and form of relationships between the response variable (whale presence) and the explanatory variables (environmental conditions); the latter looks for the ecological explanations underlying the observed empirical relationships. Ideally, and as is the case here, statistical and biological performance is linked, because putative explanatory variables are selected on the basis of hypothesised (or at least plausible) relationships between whales and their environment, especially in relation to the distribution of feeding opportunities. Arguably the biggest challenge relates to the opportunistic nature of the dataset.
Initial land search for the animals could bias sightings towards shorter distances from the shore. However, this bias would affect both blue whale sightings and pseudo-absences (sightings of other species), which will tend to counteract the original bias. What is likely is that the number of data points will decline with increasing distance from the shore, resulting in wider confidence limits and a reduced ability to detect preferences. Finally, as our models have a descriptive purpose (and not predictive), we analyze habitat preferences of blue whale in a certain study area, limited in this case by the distance observed from the coast. The lack of well-quantified effort information directly affects abundance estimations and hampers calculation of encounter rates per unit of area or time. We calculated sightings per number of trips to provide a rough overall abundance index, but did not use this in the models.
Acknowledgments We thank all the biologists and crew onboard the whale watching boats of Futurismo Azores Adventures, particularly  Miranda van der Linde  for her efforts to ensure the dataset effective and useful;  Roberto Soares , who spent very long hours looking for the animals. We also warmly thank  Marc Fernandez Morro  ´n for his advice on the model evaluation and  Cristina Gonzalez Haro , for her support on the final manuscript.
This work was supported by the  Applied Physics Department of the Universidad de Vigo (Spain) ; a predoctoral grant from  Consellerı ´a de Cultura, Educacion e Ordenacion Universitaria, Xunta de Galicia [Modalidade B]  ( PRE/2013/409 ) for the  Marine Science, Technology and Management (DOMAR)  PhD Program (www. edu. xunta. gal/);and a grant from  LabexMER action  for young researchers in 2016 (https://www. labexmer. eu/). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
In this study, after investigating whether V. v. nigrithorax is at equilibrium in its invasion range, we used more than 10,000 European invasion occurrences recorded between 2004 and 2015 to test whether occurrences from the later stage of invasion could have been predicted by models calibrated using occurrences from the earlier stage of invasion. We also took advantage of having independent validation data to test whether SDMs would have a better predictive accuracy of the ongoing invasive range if native data were accounted for, thereby responding to another strong question regarding the use of SDMs in invasion biology.
Presence data of V. v. nigrithorax in its native and invaded ranges Presence data from the native Asian range was obtained by gathering information on museum specimens, published records and hornet sampling performed in China  [50]  (see S1 Fig in Supporting Information). As for the invaded range in Europe, the species was mostly observed in France, where it was first seen in 2004. Data from the French part of the invaded range came from the INPN database that aggregates all validated French records (https://inpn. mnhn. fr/), including nests but also presences of workers in regions where nests have not been found yet. To this French database, we added the recent locations reported in other European countries (Spain, Portugal, Italy, Belgium and Germany)  [44] [45] [46] [47] [48] 51]. Overall, we obtained 10,395 records in the European invaded range observed between 2004 and 2015.
Climate data We used the same eight climatic variables as in previous studies for the niche modelling of V. v. nigrithorax  [50, 52]. We considered: (1) annual mean temperature, (2) mean temperature of the warmest month, (3) mean temperature of the coldest month, (4) temperature seasonality, (5) annual precipitation, (6) precipitation of the wettest month, (7) precipitation of the driest month and (8) precipitation seasonality. The seasonality is the coefficient of variation of the monthly means. Current data were downloaded from the worldclim database  [53]  (http:// www. worldclim. org/) as 2. 5 arc-min grids (subset of the 19 bioclim variables). These data are interpolations from observed data representative of current climatic conditions.
Climate suitability modeling Climate suitability of V. v. nigrithorax was modeled by running eight different modeling techniques implemented within the biomod2 package (3. 3-7 version)  [56]  in R  [57] : three regression methods (GLM, GAM and MARS), two classification methods (CTA and FDA) and three machine learning methods (ANN, BRT and RF). We built two sets of models: one that accounted for presence data from the invaded range only and one that accounted for presence data from both the native and the invaded range. In order to evaluate whether the ongoing invasion in Europe can be predicted by climate suitability modeling, we used only presence data from the earlier stage of the invasion (2004-2010), so that records from the later stage of the invasion could be used as evaluation data.
Model predictive accuracy was evaluated by assessing how well data recorded during the later stage of the invasion (2011-2015, hereafter called evaluation data) was predicted by the models that were calibrated with data from the early stage of the invasion only (2004-2010, hereafter called calibration data). To this end, we extracted the predicted climate suitability values for all invasive records from the evaluation data. Although the species displays one of the fastest invasive rate, with founder queens in flight mill experiments able to travel over 40 km a day  [42]  or fly an average of 18 km per day covering up to 200 km over 10 days  [43] , dispersal (in all cardinal directions) remains likely to be a limiting factor to the natural spreading. V. v. nigrithorax is therefore more likely to colonize areas close to where it first invaded.
Acknowledgments This work was supported by the  INVACOST  project ( ANR & Fondation BNP Paribas ),  Biodiversa Eranet (FFII) ,  France AgriMer  ( Programme communautaire pour l'Apiculture ,  2008-2011  and  2013-2016 ), the  MEDDE (Ministère de l'Ecologie, du Developpement Durable et de l'Energie)  and the  MAAF (Ministère de l'Agriculture, de l'Agroalimentaire et de la Forêt ). We are most thankful all persons and organizations that provided records of hornet nests in France and in Europe (their list is available on the INPN website).
All the data used in this study are available from the INPN database, an aggregated and heterogeneous database about biodiversity. Aggregated data (number of nests by 10x10 km grid cells) are freely available, but the precise GPS records cannot be made publicly available as they come from citizen Data requests can be made by filling in the request QR), the authors had special access to the data. In addition to the restrictions explained above, the authors cannot provide these data freely due to French regulations. However, the authors confirm that interested, qualified researchers will be able to access the data without restriction after filling in the request form and specifying that the data will not be used for commercial purposes, agreeing to provide attribution in any published papers, and that they will not be distribute these data to third parties.
Spatial heterogeneity In the last decade, there was an almost exponential increase in studies addressing ES value mapping  (Schägner et al. 2013). ES value mapping combines results from economic valuation studies and multiple spatial datasets on, for example, ES indicators. Accounting for spatial heterogeneity, i. e. , spatial variation of ES monetary values, as well as opportunity and intervention costs, is important for the design and the cost- effectiveness of public regulations, e. g. , PES schemes. To facilitate the use of economic studies in mapping exercises, it is important that the spatial scale and context be explicitly described and accounted for when valuating ES  (Ninan and Inoue 2013). This also includes being explicit about the availability of alternative local sites, making it possible to account for substitution in the provision of services.
The need for long historical data series to identify long-term decisions in forestry and long, uncertain, difficultly observable ecological outcomes explains the lack of empirical results about competition or synergies in the joint production of FES. Moreover, the complex biological processes in forest ecosystems where dimensions of space and time play a crucial role in their development should encourage economists to develop new models of FES provision. These models established in the traditional input-output framework could also be reconsidered. Accounting for spatial and temporal aspects also challenges the existing approaches of valuation and gives rise to a new research field: ES mapping, where multiple spatial datasets are linked to ecological and valuation data. As highlighted by  (Dobbs et al.
Morphological analyses Genus or species of each mat sample were identified using direct microscopy. The phenotypic identification was performed using a Zeiss Axioskop 40 microscope on which a Zeiss camera (AxicCam ICc1) was mounted. Using the computer program Axiovision (Carl Zeiss, Jena, Germany), we determined morphological features such as sheath, calyptras, end cell morphology, cell sizes (n = 50 per sample) and proportions and degree of constriction at cross-walls. Subsequent taxonomic identification of the cyanobacterial species was based on the identification keys of Komarek and Anagnostidis (2005) and Komarek and Hauer (2010). The relative abundance of each genus or species within each sample was scored on a semiquantitative scale (0: 0%, 1: 1-19%, 2: 20-39%, 3: 40-69%, 4: 70-89%, 5: 90-100%) after screening 3-5 subsamples.