Magdalena Falco
email: magdalena.falco@cima.fcen.uba.ar

Andrea F Carril

Claudio G Menéndez

Pablo G Zaninelli

Laurent Z X Li

Assessment of CORDEX simulations over South America: added value on seasonal climatology and resolution considerations

Keywords: Model assessment, Seasonal climatology, Added value, South America, CORDEX 1 Introduction

published or not. The documents may come   L'archive ouverte pluridisciplinaire

simulations. This has, in turn, prompted the need for understanding regional climate modeling strengths and weaknesses over the continent by establishing common experimental frameworks. Under these circumstances, [START_REF] Roads | International Research Institute/Applied Research Centers (IRI/ARCs) regional model intercomparison over South America[END_REF] were the first to analyze an ensemble of four RCMs driven by the NCEP/NCAR I global reanalysis [START_REF] Kalnay | The NCEP/NCAR 40-Year Reanalysis Project[END_REF], followed by CREAS regional program [START_REF] Marengo | Future change of climate in South America in the late twenty-first century: intercomparison of scenarios from three regional climate models[END_REF]. The CLARIS [START_REF] Menéndez | CLARIS Project: towards climate downscaling in South America[END_REF] and CLARIS-LPB Projects [START_REF] Boulanger | CLARIS-La Plata Basin: regional hydroclimate variability, uncertainties and climate change scenarios[END_REF] succeeded to create a large ensemble of regional simulations, involving different RCMs driven by ERA-Interim reanalysis [START_REF] Dee | The ERA-Interim reanalysis: configuration and performance of the data assimilation system[END_REF] and CMIP3 Global Climate Models [START_REF] Meehl | THE WCRP CMIP3 Multimodel Dataset: A New Era in Climate Change Research[END_REF]. The CLARIS-LPB Project comprehended two types of simulations for the recent past concerning the boundary forcing: the hindcast simulations and the historical simulations.

Simulating SAs climate is challenging due to the extension of the continent, and the local to regional climatic features associated to the South American Monsoon (SAM) System over tropical and subtropical latitudes, and also to frontal activities over mid-latitudes [START_REF] Garreaud | Present-day South American climate[END_REF]. While GCMs are generally capable of reproducing the main large-scale circulation characteristics, their performance deteriorates when evaluating regional aspects of climate [START_REF] Gulizia | Identification of the principal patterns of summer moisture transport in South America and their representation by WCRP/CMIP3 global climate models[END_REF][START_REF] Rusticucci | An intercomparison of modelsimulated in extreme rainfall and temperature events during the last half of the twentieth century. Part 1: mean values and variability[END_REF][START_REF] Silvestri | Evaluation of the WCRP-CMIP3 model simulations in the La Plata basin[END_REF][START_REF] Zazulie | Regional climate of the subtropical central Andes using high-resolution CMIP5 modelspart I: past performance (1980-2005)[END_REF]. For example, [START_REF] Gulizia | Identification of the principal patterns of summer moisture transport in South America and their representation by WCRP/CMIP3 global climate models[END_REF] and [START_REF] Zazulie | Regional climate of the subtropical central Andes using high-resolution CMIP5 modelspart I: past performance (1980-2005)[END_REF] found that CMIP5 GCM simulations [START_REF] Taylor | An Overview of CMIP5 and the Experiment Design[END_REF] do not reproduce adequately the moisture transport associated with the South American Low-Level Jet. Previous works have assessed regional climate simulations over SA within the context of CLARIS-LPB Project [START_REF] Carril | Performance of a multi-RCM ensemble for South Eastern South America[END_REF], 2016, López-Franca et al., 2016[START_REF] Sánchez | Regional climate modelling in CLARIS-LPB: a concerted approach towards twentyfirst century projections of regional temperature and precipitation over South America[END_REF][START_REF] Solman | Systematic temperature and precipitation biases in the CLARIS-LPB ensemble simulations over South America and possible implications for climate projections[END_REF][START_REF] Solman | Evaluation of an ensemble of regional climate model simulations over South America driven by the ERA-Interim reanalysis: model performance and uncertainties[END_REF][START_REF] Zaninelli | Explorando temperaturas máximas y mínimas en diferentes reanálisis. Parte 1: campos medios estacionales[END_REF], 2018[START_REF] Sanabria | Maps of wind hazard over South Eastern South America considering climate change[END_REF], Menéndez et al., 2016), indicating an overall good performance of RCMs in reproducing different aspects of the observed climate. However, none of these works focus entirely on the AV issue, and only a few of them [START_REF] Sánchez | Regional climate modelling in CLARIS-LPB: a concerted approach towards twentyfirst century projections of regional temperature and precipitation over South America[END_REF][START_REF] Solman | Systematic temperature and precipitation biases in the CLARIS-LPB ensemble simulations over South America and possible implications for climate projections[END_REF]) make a general comparison of the RCM with its driving data.

Apart from CLARIS-LPB Project, the CREMA experiment [START_REF] Coppola | Present and future climatologies in the phase I CREMA experiment[END_REF] succeeded to produce several simulations of RegCM4 RCM forced by different GCMs and reanalysis over SA. When comparing these simulations with their driving data, the results suggest an added value of RegCM4 in reproducing extreme events over SA [START_REF] Coppola | Present and future climatologies in the phase I CREMA experiment[END_REF], Giorgi et al., 2014). Also, the model partly reduces the dry bias of precipitation in the Amazon Basin when comparing to different CMIP5 GCMs [START_REF] Llopart | Climate change impact on precipitation for the Amazon and La Plata basins[END_REF], da Rocha et al., 2014). Some few more works have addressed the added value issue over SA, and there are no matching results: [START_REF] Roads | International Research Institute/Applied Research Centers (IRI/ARCs) regional model intercomparison over South America[END_REF] assessed the capabilities of RCMs to make skillful forecasts compared to its forcing data over SA, finding no noticeable improvement in precipitation when assessing the large-scale fields; [START_REF] Seth | RegCM3 regional climatologies for South America using reanalysis and ECHAM global model driving fields[END_REF] studied the performance of RCMs forced by reanalysis for different regions on tropical and extratropical latitudes in SA. In their work they concluded that regions where remote influences are strong (Northeast SA, Southeast SA, and Amazon regions) the AV of RCMs is low; by contrast, where local processes dominate and remote influences are weak (maximum monsoonal precipitation region) there is some potential for the regional model to add value. Finally, De Sales and Xue (2011) assessed the role of resolution on the South American climate simulation by performing a sensitivity test with the Andean topography represented by different heights. They suggested that height representation of the Andes Mountain Range topography impact on low-level moisture flux transport and therefore on precipitation in La Plata Basin.

Within the CORDEX framework [START_REF] Giorgi | Regional Dynamical Downscaling and the CORDEX Initiative[END_REF], several simulations with updated model versions have been recently incorporated to the CORDEX SA domain (Fig. 1), enabling the availability of state-of-the-art nested RCMs at a horizontal resolution of ∼0.44 • x∼0.44 • . Therefore, motivated by the availability of new simulations and the lack of studies about AV over the continent, this work aims to assess the AV of CORDEX simulations in South America. As a first step, we evaluate the capacity of RCMs to simulate seasonal mean precipitation and surface air temperature over SA through its comparison to their coarser-resolution driving data. We believe this is a fundamental step, not only for RCM evaluation but also due to the discrepancies in the literature whether RCMs can add value or not over large-scale mean values [START_REF] Llopart | Climate change impact on precipitation for the Amazon and La Plata basins[END_REF][START_REF] Roads | International Research Institute/Applied Research Centers (IRI/ARCs) regional model intercomparison over South America[END_REF]. The paper is organized as follows: in section 2 we first describe the observed and simulated data used in this study, followed by the introduction of model errors and observed uncertainties. In section 3 we assess the general and regional performance of individual and multi-model ensemble data in simulating the seasonal spatial patterns of the mean state, including a regional assessment of model errors and observational uncertainties. Finally, in section 4, we present a discussion of the most relevant results together with our main conclusions.

Am

Data and methods

CORDEX simulations and observed data

For our analysis, we selected the common time period of the integrations, covering from 1990 to 2004, to keep the same time-scale variability. The seasons considered were austral summer (December, January, and February) and winter (June, July, and August). All RCMs use roughly the same domain and horizontal resolution (∼0.44 • x∼0.44 • ), but different grid types. For a direct assessment and intercomparison of RCM simulations, we interpolated all original monthly precipitation and 2-meter-temperature data onto a regularly spaced latitude/longitude grid of 0.5 • x0.5 • . ERAi and CMIP5 GCM fields were also interpolated into a common regularly spaced grid of 1.5 • x1.5 • . Both variables were interpolated using a bilinear interpolation procedure, and a height correction was applied to the surface air temperature fields. The height correction consisted of adding the difference between the model native height (interpolated into the regular grid) and the observed height, multiplied by a standard lapse rate of -6,5 • /km. Different gridded data for precipitation and surface air temperature were used to evaluate the simulations. With the goal of measuring the uncertainty among observation data sets, we consider three station-based monthly data for precipitation (CRU, UDEL and CPC-UNI), over land on a 0.5 • x0.5 • resolution. We also considered two of the above data sets (CRU and UDEL) for surface air temperature. Complementary to the interpolation of model simulations, the observed datasets were also interpolated using the same procedure to the same 0.5 • x0.5 • regular grid, and also, to the 1.5 • x1.5 • grid to allow GCM evaluation. Throughout the work, the set of observed datasets is denoted as Obs, and the reference field for the assessment is calculated as the average across all the observed datasets.

Model error and observational uncertainties

The main advantage of RCMs is they are integrated with a higher spatial resolution than their boundary-condition providers and, therefore, add to the latter information on a regional scale. Thus, one can distinguish between two different spatial scales: the small scales that can be reproduced only by the RCM and the large scales reproduced by both simulations. However, when it comes to consciously comparing an RCM with its driving data, there is no optimal methodology to follow [START_REF] Prein | Precipitation in the EURO-CORDEX 0.11 • and 0.44 • simulations: high resolution, high benefits?[END_REF]. On the one hand, if we choose to interpolate the RCM and the observations to the resolution of the GCM, we would benefit the GCM, since we are removing the small-scale variability of the RCM. But if we interpolate the global model to the resolution of the regional model and the observations, then the GCM is disfavored since the reduction of the grid will generate non-consistent physical information. In this paper, we performed the evaluation of the experiments in their resolution, that is, in a 0.5 • grid for the regional simulations and a 1.5 • grid for the global simulations. However, for the evaluation to be adequate, we considered other sources of error associated with grid interpolation, resolution and observational uncertainty. [START_REF] Kanamitsu | The Added Value Index: A new metric to quantify the added value of regional models[END_REF] introduced a clear picture of model error or deviation from the observed reference by describing the whole process of model validation. Firstly, each grid point value in model represents an entire area of the field, and therefore the model will not be able to reproduce the sub-grid variability. This error or limitation can be called representativeness error and is a function of resolution of the model and variability of the field. Secondly, interpolation error can be identified as those changes of a field when interpolated from the model grid into a validation grid. The magnitude of this error is also influenced by the model's resolution and the variability of the field. Finally, the observations used for the evaluation are not perfect, and generally have errors associated with instrumental and data acquisition imperfections. Interpolation and representativeness errors should also be counted when using station-based gridded datasets. The concept of [START_REF] Kanamitsu | The Added Value Index: A new metric to quantify the added value of regional models[END_REF], expressed as an apparent deviation between the model and observation for a single location, can be qualitatively formulated as

M odel Error = [ε M ] + [ε R ] + ε I + ε obs , (1)
where [ ] is a spatial interpolation operator, ε M is the model deviation from observations, ε R is the representativeness error of model grid value, ε I is the interpolation error, and ε obs is the error associated with observation. Eq. ( 1) shows that model deviation seen in a validation operation is directly related to the performance of the model, i.e., the true model error, if and only if ε R , ε I , and ε obs are negligible. In real world, they are generally not negligible, and they need an appropriate evaluation, which allows us to have just an appreciation for the true performance (or error) of a model.

In our work, we revisit Eq. ( 1) and adapt it to our experimental design. Our reference observations for evaluation are station-based gridded datasets, so the grid spacing of the observational dataset defines the finest scale on which the comparisons are meaningful [START_REF] Prein | Precipitation in the EURO-CORDEX 0.11 • and 0.44 • simulations: high resolution, high benefits?[END_REF]. All observed datasets have a resolution of 0.5 • x0.5 • . This entails that it is impossible to quantify the sub-grid representativeness error of RCMs, as they have a finer resolution. Following Eq. ( 1), the deviation of both RCM experiments from the reference field will be affected by the following terms:

ε m RCM = f (ε m M , ε 0.5 • R , ε 0.5 • I , ε obs ), (2)
where m can be RCMs-fERAi or RCMs-fCMIP5, ε 0.5 • R is the representativeness error at a resolution of 0.5 • , that is neglected hereafter since its resolution is close to that of the observations, and ε 0.5 • I is the interpolation error at a resolution of 0.5 • . In the same way, the errors associated with GCMs will include the same terms but for a coarser resolution grid. Therefore, the deviation of both GCM experiments from the reference field will be affected by the following terms:

ε m GCM = f (ε m M , ε 1.5 • R , ε 1.5 • I , ε obs ), (3)
again, m can be ERAi or CMIP5, ε 1.5 • I is the interpolation error at a resolution of 1.5 • , and

ε 1.5 • R
is the representativeness error between the resolutions of 0.5 • and 1.5 • . In this case, the term associated with the representativeness error is not negligible, as degrading the observed field into the GCM resolution for evaluation implies losing details in the process of filtering.

In this work, we present a methodology that accounts for the different terms involved in the process of model evaluation. Despite strong constraints and difficulties in calculating precisely each term in Eqs. ( 2) and (3), we try to estimate their order of magnitude. We aim to find regions with significant contributions in order to determine if the model's deviations from observation are genuinely associated with model performance or not.

Starting with the model error, ε M , we cannot execute a temporal validation as historical simulations are not synchronized with the observed variability. Consequently, we perform a multimodel evaluation of the seasonal climatology by using the root-mean square (RMS) error metric. We define ε M as the RMS deviations of all individual members of an ensemble m from the reference data, that is

ε m M = 1 N m Nm i=1 ([F m,i ] res -[F ref ] res ) 2 (4)
with N m the number of associated simulations, and [ ] res the interpolation of the simulation and station-based observed datasets to a regular grid of resolution res (0.5 • for both RCMs and 1.5 • for both GCMs). It is worthy to note that the reference field is calculated as an average among available datasets in the resolution of each experiment.

The degradation of the observed field from 0.5 • to 1.5 • is done by aggregating the reference field at 0.5 • into a coarser-resolution grid of 1.5 • x1.5 • . The aggregation is performed by averaging all original sub-grid values inside each of the coarser grid boxes [START_REF] Di Luca | Potential for added value in precipitation simulated by high-resolution nested Regional Climate Models and observations[END_REF]. The aggregated values are then re-distributed into the original resolution, which will give an estimation of the details lost in the process of filtering. The error ε R is calculated as the RMS difference between the degraded and reconstructed fields

ε 1.5 • R = (F 0.5 • ref -F 1.5 • ref ) 2 =| F 0.5 • ref -F 1.5 • ref | . (5)
Eq. ( 6) permits to measure the information contained in the range of resolutions from 0.5 • to 1.5 • . Unlike [START_REF] Kanamitsu | The Added Value Index: A new metric to quantify the added value of regional models[END_REF], we won't neglect ε I and ε obs as both terms can be significant over some regions. While there is no easy way to calculate the interpolation error, we estimate it by interpolating the mean observed climatological field into a regularly-spaced but staggered grid, and then interpolating it back to the original grid. For the 0.5 • x0.5 • resolution we consider a displacement of 0.25 • in both latitude and longitude directions and for the 1.5 • x1.5 • a displacement of 0.75 • . The interpolation error is then estimated as half the difference between the double interpolated field and the original field

ε res I = ([[F res ref ]] -F res ref ) 2 2 = | [[F res ref ]] -F res ref | 2 . ( 6
)
Finally, we need to point out that a main part of observation errors comes from the representativeness error and interpolation error when station-based information is converted into a regular grid. But its estimation is out of the scope of our study. Nevertheless, we can estimate the observational uncertainty as the spreading among individual datasets. For that, we substitute m by Obs and res by 0.5 • in Eq. ( 4), and obtain

ε obs = 1 N Obs N Obs i=1 ([F Obs,i ] 0.5 -[F ref ] 0.5 ) 2 .
(7)

3 Results

Seasonal spatial patterns of the mean state

The challenge for regional climate modelers in representing the South American climate comes from the extent of the computational domain and the complexity of the physical phenomena involved. During the austral summer, precipitation maximizes over tropical and subtropical latitudes and is largely associated to the SAM System (Fig. 2a). All multi-model ensemble climatologies can reproduce large-scale features of the summer precipitation, but the presence of deviations from the reference indicates some systematic biases (Figs. 2b-2e). There are some common biases in all experiments, such as an overestimation of precipitation intensity along the Andes Cordillera, a dry bias in northern and central Brazil and a wet bias in Northeast Brazil (NeB). This last behavior is accentuated in historical experiments, especially in RCMs-fCMIP5.

All the three ensembles, RCMs-fERAi, RCMs-fCMIP5 and CMIP5, show a wet bias over the South Atlantic Convergence Zone (SACZ) and a dry bias over La Plata Basin (LPB). On the contrary, ERAi has opposite biases in these two regions. The spatial distribution of surface air temperature biases (Figs. 2g-2j) reveals a cold bias over most of the continent in all experiments. This is not true over the Amazon basin, where RCMs-fERAi and CMIP5 ensembles show small positive values. Also, over LPB, most ensembles except ERAi product have a warm bias. Given that results are presented as ensemble mean, it is important to assess the multi-model dispersion over each ensemble. Following [START_REF] Rinke | Evaluation of an ensemble of Arctic regional climate models: spatiotemporal fields during the SHEBA year[END_REF] and [START_REF] Solman | Evaluation of an ensemble of regional climate model simulations over South America driven by the ERA-Interim reanalysis: model performance and uncertainties[END_REF], we can highlight those areas where the dispersion among datasets is larger than the interannual variability of the field (hatched areas in Fig. 2). In the case of summer precipitation, all three ensemble mean climatologies have large inter-model spread over most of tropical and subtropical SA (Fig. 2b, 2d, and 2e), suggesting that the ensemble climatology is not very representative of individual climatologies. The absence of hatched areas for observed climatologies (Fig. 2a) indicates that the inter-dataset spread is small over most of SA. Regarding surface air temperature, only CMIP5 ensemble has large hatched areas that denote a large dispersion among simulations over most of the continent. Winter precipitation climatology (Fig. 3a) exhibits two maxima: the first over the northern part of the continent, north of the equator, associated to the North American Monsoon (NAM) System; and the other one over the meridional Andes, south of 40 • S, associated to baroclinic activities. While all experiments show the northern precipitation maximum associated to the NAM System (not shown), they fail to capture the observed intensity by largely underestimating it (Figs. 3b-3e). Even ERAi fails to simulate the southward extension of the NAM influence, which is noticed by the wet bias at the southern edge of the region. All ensembles present a significant multi-model spread over the tropics, indicating a wide variation in the representation of the NAM. In the extratropics, all multi-model ensembles overestimate precipitation over the western slope of the Andes. The bias is maximal in RCMs, as they also present an unrealistic extension over northern Chile. Regarding 2-meter temperature (Fig. 3, bottom panels), we note again that most experiments have a cold bias over most of the continent, except RCMs-fERAi multi-model ensemble over the southern Amazon basin and CMIP5 multi-model ensemble over northern Brazil and LPB. During this cold season, unlike in summer, both RCM ensembles have large inter-model spread for surface-air temperature climatology over the tropical and subtropical latitudes, and also over Patagonia in the case of RCMs-fCMIP5. CMIP5 simulations continue to show large inter-model spread over the same latitudes as in summer.

Taylor diagrams [START_REF] Taylor | Summarizing multiple aspects of model performance in a single diagram[END_REF] are multi-parameter-based diagrams that allow intercomparison of models performance in reproducing a reference pattern. They are a useful and visual tool to summarize how well climate models reproduce the observed climate. We first assess models performance over the South American low-level terrain (lower than 1500m) in Fig. 4. Results suggest that there is a clear distinction between models ability to reproduce 2-meter temperature climatology and precipitation climatology (Fig. 4). The spatial correlation between observed and modeled climatology is always higher for 2-meter temperature (ranging from 0.9 to 0.99) than for precipitation (from 0.5 to 0.9). It is interesting to note that in all cases the ensemble mean improves the field when comparing to each individual simulation, as suggested by other authors [START_REF] Carril | Performance of a multi-RCM ensemble for South Eastern South America[END_REF], López-Franca et al., 2016). Surface air temperature tends to be closer from each other in the Taylor diagrams compared to precipitation points, which indicates a larger inter-model spread for precipitation. During summer, ERAi reanalysis and RCMs-fERAi multi-model ensemble have both the best performance in reproducing the observed summer precipitation climatology (Fig. 4a). Their overall performance is similar in terms of correlation but not in terms of normalized standard deviation (nSTD). When compared with the reference, ERAi has lower spatial variability (nSTD<1) and RCMs has higher spatial variability (nSTD>1). These characteristics are also valid when evaluating the ensemble of RCMs-fCMIP5 and their driving data ensemble. During this season, all RCMs-fCMIP5 individual simulations show worse performance, compared to RCMs-fERAi individual simulations, probably due to the quality of the different boundary conditions. Nevertheless, all RCM individual simulations have larger nSTD when compared to their driving data, suggesting there is an excessive spatial variability. Regarding surface air temperature, the ensemble of RCMs-fERAi is the one showing the best match with the observed climatology pattern by presenting the highest spatial correlation coefficient. In both experiments, there is an added value from RCM ensembles in reproducing the referenced surface air temperature, each individual simulation being also the case. Again, as noted in Fig. 2, CMIP5 simulations tend to have a significant spread during this season.

RCMs-fCMIP5 ensemble is the most skillful multi-model ensemble in simulating the winter precipitation climatology (Fig. 4b). RCMs-fERAi has similar nSTD than its driving data but lower spatial correlation. In the historical experiments, RCMs-fCMIP5 has higher correlation and a spatial variability nearer the observed when comparing to its driving data. Alike summer, RCMs-fERAi continue having the best performance in reproducing the winter 2-meter temperature climatology, but all ensemble and individual simulations have high correlation values and near to 1 nSTD, suggesting an overall good performance.

3.2 Regional performance of individual simulations and multi-model ensembles Regional-scale climate assessment is fundamental for understanding the potential added value of RCMs since AV is expected only over those regions where stationary climatic features are associated to small-spatial scale processes [START_REF] Di Luca | Potential for added value in precipitation simulated by high-resolution nested Regional Climate Models and observations[END_REF]. For this reason, we emphasize on the evaluation of regional performance by selecting several key areas over the whole continent (Fig. 1). It is clear that RCMs ensembles can add value to the simulated climate when compared to their driving data (Figs. 2 and3), but it is worthy to note that there is a large spread among simulations (Figs. 2, 3 and 4). We need to evaluate the performance of individual RCM simulations against their driving data, which can reveal the robustness of results. To do so, we use an AV index defined as the ratio between the absolute error (AE) of the RCM and the AE of the GCM. The AE is calculated as the absolute difference between the simulation climatology with the reference, followed by an average over each individual region. Hence, if the index is smaller than one, the RCM has a lower AE than its driving data, suggesting a better representation of the mean climatology over a region. The index is also calculated for the multi-model ensemble climatology to see how representative it is.

Table 2 Added Value Index (AVI) of individual regional climate simulations against their driving data for austral summer (DJF) and for each study region. The index is calculated as the ratio of the area-averaged absolute error (AE) of RCM over that of GCM. Bold figures highlight cases with AVI smaller than 1, indicating added value in RCM. When analyzing the AV index over all summer regions (table 2), we find that even though the ensemble of RCMs-fERAi represents the mean value of precipitation better than its driver does, individual simulations show different statistics. Over the five regions, only three of them (SAM, SACZ and LPB) show added value within at least one individual simulation: three simulations (over five) in the case of SAM region, two in the case of SACZ and only one in LPB. For RCMs driven by CMIP5, the RCMs add value over the CMIP5 driver only in Am region, coincident with results in Fig. 9, since 3 of 5 individual simulations have lower bias than its driver. Other regions where individual simulations add value in precipitation are SAM with two simulations, NeB with one simulation, SACZ with one simulation and LPB with three simulations.

Experiment

Summer surface air temperature shows a more positive outlook for the RCMs-fERAi ensemble, as in Am (5 simulations), SAM (5 simulations), NeB (3 simulations) and SACZ (4 almost all simulations, together with the ensemble mean, show a lower AE than ERAi. RCMs-fCMIP5 adds value over two regions: Am (4 simulations) and LPB (3 simulations). For Am, the ensemble mean did not represent individual performance. This counter-intuitive result is probably caused by the large spreading among CMIP5 models over this region (Fig. 2j). For winter precipitation, Table 3 shows that AV can only be found in some individual simulations over NAM region (2 simulations of CMIP5-driving RCMs and one of ERAi-driving RCMs), four simulations in LPB (two in each experiment) and only one simulation in SSA. Regarding surface air temperature, most models show an improvement over NAM region for both experiments (5 in ERAi-driving RCMs and 3 in CMIP5-driving RCMs), but in the case of CMIP5-driving RCMs the higher skill is not represented by the ensemble mean. For RCMs-fERAi, there are two individual simulations showing added value in LPB and also two simulations in SSA with an AV index of 1.

Box-and-whisker plots are an excellent non-parametric tool for evaluating quartile distribution, allowing us to visually analyze and compare various L-estimators beyond the mean value. Assessing the distribution of rainfall and surface air temperature allows us not only to assess the biases in the mean climate but also to evaluate the interannual variability of the simulations. Generally speaking, models are far from correctly reproducing the observed variability (Figs. 5 and6). However, in most regions, data sets broadly reproduce the observed statistical distribution by intersecting at least part of the inter-quartile range (IQR; that is 25 th -75 th quartile distance) with the corresponding observed one.

Most of the experiments tend to have a systematic cold bias over most regions during summer (Fig. 5). This is especially true for ERAi reanalysis and RCMs-fCMIP5 ensemble. Precipitation is Every box represents the inter-quartile range (25 th -75 th percentile), the median is indicated with a white line and the mean with a green triangle. The whiskers indicate the 95 th and 5 th percentiles. The regional time series of each simulation is obtained by spatially averaging all values within the region in every time step. Also, all percentiles of the multi-model ensembles are calculated as the average of individual percentiles. The left panels correspond to precipitation (mm/day; blue) and the right panels to 2-meter temperature ( • C; red).

well simulated in all experiments and over most regions, an exception is RCMs-fCMIP5 ensemble which has a systematic wet bias over SAM, NeB and SACZ. Concerning the AV of RCMs in precipitation, there is a little improvement compared to its driving data. RCMs-fERAi improves ERAi by better reproducing the observed median value and the interannual variability over Am region. There are also some improvements over LPB and SACZ for the median value but not so much for the interannual variability. RCMs-fCMIP5 only slightly improve its driving data for precipitation distribution over LPB. Regarding surface air temperature, the added value is more evident, especially for the case of hindcast simulations. RCMs-fERAi improves the reproduction of the observed distribution over Am, SAM, NeB and SACZ region. On the other hand, RCMs-fCMIP5 ensemble does not seem to reduce the bias of CMIP5 simulations, but they rather improve the representation of the interannual variability; this is the case for Am, SAM, SACZ and LPB regions.

During winter, there is a general disability of climate models in reproducing the observed statistics over NAM, by systematically underestimating both temperature and precipitation (Fig. 6). AV hardly manifests, but RCMs-fERAi reduces the large cold bias that ERAi has over the NAM region and the same region sees an improvement in RCMs-fCMIP5 which reduces the dry bias of CMIP5 simulations. Winter precipitation statistics over SSA indicate a general systematic overestimation of rainfall by all models. This is especially true for both RCM ensembles, which lowers its performance when comparing to their driving data.

Assessment on model errors and observational uncertainties

The assessment performed in the previous sections combines different types of error, which makes difficult to understand how the RCM simulations add value (or not), compared to its coarserresolution driver. We follow the resolution-aware methodology described in the methods section to identify the type of error that should not be neglected during the assessment. A regional comparison of all the terms presented in Eqs. ( 2) and ( 3) is performed over each region for both temperature and precipitation. We present the average value of each term over each region for DJF and JJA respectively (Fig. 7). The observational uncertainty, ε obs , is considered as a reference value, as no assessment can be done when errors are smaller than this term. Results show that model error (ε M ) is usually the largest one. We recall that it generally represents the model performance that we actually want to evaluate. However, ε R , and ε I are not always negligible. The spatial distribution of these errors (see supplementary material) shows they are particularly large over complex terrain regions, especially for surface air temperature. This result is also evident in the regional diagram (Fig. 7), where the coarse-resolution errors, i.e., ε 1.5 I and ε 1.5 R , are relatively large over SACZ in summer and NAM and SSA in winter. ε R is usually larger than ε I , being the latter the one that presents the smallest values when comparing to the other terms. ε I is also dependent of resolution, being the interpolation error associated to a grid of 1.5 • resolution (ε 1.5 I ) larger than the one from a grid of 0.5 • resolution (ε 0.5 I ). In summer, all ε 1.5 R and ε 1.5 I are of a similar order of magnitude, close to ε obs for both temperature and precipitation. ε 1.5 R often exceeds ε obs value, as is the case in NeB and SACZ for precipitation and in NeB, SACZ and LPB for surface air temperature (Fig. 7, left panels). When analyzing the regional performance of models, ERAi has the best performance in reproducing regional precipitation climatology, since ε ERAi 2) and (3). Top panels are for precipitation (mm/day), and bottom panels for 2-meter temperature ( • C). Red, yellow, light-blue and green colors represent our estimations in RCMs forced by reanalysis, reanalysis, RCMs forced by CMIP5 and CMIP5, respectively. Black markers represent our estimation for observational uncertainty. Grey empty markers correspond to interpolation (square markers) and representativeness (circle markers) errors.

During winter, ε obs is large for precipitation over the NAM region (Fig. 7, right panels), as was already seen with the spread among datasets in Fig. 3a. The uncertainty is also relatively large (larger than 0.5 • C) over NAM and SSA for surface air temperature, probably because of the selection of domains placed over complex terrain regions. ε 1.5 R is larger than the observational uncertainty for surface air temperature, by exceeding it over all three winter regions (i.e., NAM, LPB and SSA). ε 1.5 I is also larger than ε obs in NAM for surface air temperature. Regarding ε M , RCMs usually have similar or larger RMSE in precipitation than its driving data; except for the case of NAM region where ε RCM s-f CM IP 5 M has a lower value than ε CM IP 5 M . The aforementioned also represents the case for 2-meter temperature, except for RCMs-fERAi over LPB.

Conclusions and discussion

Although many CORDEX simulations over SA became available recently, there was a general lack of works investigating the AV of RCMs over the continent. It is thus very relevant for us to assess these regional simulations, which can provide a useful guide for researchers or end users of these simulations. The particularity of our approach was in the joint evaluation of RCM simulations and their coarser-resolution driving data. We first evaluated their capability in reproducing mean climate conditions over the continent. We used two types of RCM experiments with different nature in lateral boundary conditions: hindcast simulations with RCMs driven by ERA-Interim re-analysis and historical simulations with RCMs driven by CMIP5 global models. Our database included five hindcast simulations and five other historical simulations, as well as their driving data. The multi-model ensembles and individual simulations were evaluated against two or three observed gridded datasets for 2-meter temperature and total precipitation, respectively. We performed the analysis on a seasonal basis, for austral summer and winter, over a common period ranging from 1990 to 2004.

Generally speaking, all multi-model ensembles can reproduce main features of summer and winter climatologies. However, individual simulations exhibit a large inter-model spread, especially for precipitation, suggesting the need for assessing models individually. The ensembles also have some common model biases, as it is the overestimation of precipitation intensity along the Andes Cordillera. An important bias is found during the warm season since the ensembles differ substantially in their representation of location and intensity of the maximum monsoon rainfall. Regarding surface air temperature, models tend to simulate lower temperatures than the observed, as shown by in the bias fields and the box-and-whisker plots. This cold bias is maximized by ERAi and RCMs-fCMIP5 over the Amazon Basin region. Over the tropical and subtropical latitudes, CMIP5 models show a large spread among simulations in reproducing the surface air temperature climatology for both winter and summer seasons. [START_REF] Solman | Evaluation of an ensemble of regional climate model simulations over South America driven by the ERA-Interim reanalysis: model performance and uncertainties[END_REF] assessed the capability of a set of 7 hindcast RCMs, performed in the framework of the CLARIS-LPB Project, in reproducing the mean climate conditions over the South American continent. If we compare the CLARIS-LPB simulations with the ones shown here for CORDEX, we find that the distribution and intensity of biases are very similar (see Figs. 2 and 3 from Solman et al.). However, some improvements are found for CORDEX such as a reduction of the precipitation bias over LPB during winter, and a reduction of the warm bias over central SA in the same season. Also, our results are consistent with [START_REF] Solman | Evaluation of an ensemble of regional climate model simulations over South America driven by the ERA-Interim reanalysis: model performance and uncertainties[END_REF] in terms of spread among simulations, both showing that the inter-model uncertainty is larger than the observational uncertainty.

We performed a regional assessment in several key areas of the whole continent in summer and winter. In summer, RCMs-fERAi ensemble shows a better reproduction of observed precipitation and surface air temperature when compared to ERAi. However, some opposite results can be found for precipitation when conducting the comparison individually. There is a clear warning message on the practice of only analyzing the multi-model ensemble, especially over regions with large inter-model spread. Surface air temperature exhibits some encouraging results, with a general improvement when compared to ERAi. Regarding historical simulations with RCMs driven by CMIP5, results are inconclusive, since added value is only found in certain regions and in limited simulations. It is important to highlight that the complexity to be considered is increased when assessing the AV of RCMs in their historical simulations, as it should include GCM/RCM physical compatibility. We need to consider different GCM/RCM configurations to understand if biases are inherited from GCM through boundary conditions. During winter, RCMs do not noticeably improve the representation of precipitation or surface air temperature. As a matter of fact, they even tend to degrade the winter climatology.

We are aware that the direct comparison between the modeled field and the observed field include some additional sources of error independent of model performance. Thus, model error is associated to model performance only in those regions where ε R , ε I , and ε obs are negligible. In general, model error, ε M , is the largest one over full SA continent, allowing a direct evaluation of climate models. However, we have identified some areas where other sources of error are not at all negligible. ε R has its highest values over regions with strong surface forcing, such as the Andes Cordillera. ε I also showed the same behavior. ε obs is also large over complex-terrain regions, especially for surface air temperature. This last result is probably due to the scarcity of surface meteorological stations over isolated mountainous regions and to different interpolation techniques used in each dataset.

There is a general disagreement over literature regarding added value of RCMs for large-scale mean values, as some works find an improvement on representing the mean climatic features [START_REF] Kerkhoff | Assessment of Bias Assumptions for Climate Models[END_REF], Llopart et al., 2014[START_REF] Nikulin | Precipitation Climatology in an Ensemble of CORDEX-Africa Regional Climate Simulations[END_REF][START_REF] Torma | Added value of regional climate modeling over areas characterized by complex terrain-Precipitation over the Alps[END_REF][START_REF] Fotso-Nguemo | On the added value of the regional climate model REMO in the assessment of climate change signal over Central Africa[END_REF][START_REF] Veljovic | Regional climate modeling: Should one attempt improving on the large scales? Lateral boundary condition scheme: Any impact?[END_REF][START_REF] Prein | Precipitation in the EURO-CORDEX 0.11 • and 0.44 • simulations: high resolution, high benefits?[END_REF], while others find a deterioration [START_REF] Coppola | Present and future climatologies in the phase I CREMA experiment[END_REF][START_REF] Feser | Enhanced Detectability of Added Value in Limited-Area Model Results Separated into Different Spatial Scales[END_REF], Giorgi et al., 2014[START_REF] Roads | International Research Institute/Applied Research Centers (IRI/ARCs) regional model intercomparison over South America[END_REF][START_REF] Mishra | CMIP5 vs. CORDEX over the Indian region: how much do we benefit from dynamical downscaling[END_REF]. Our results are consistent with this general uncertainty since they suggest that added value of RCMs depends on the different driving fields (i.e., hindcast or historical), the surface properties of the area, the season and the variable analyzed. The dependency of AV to diverse factors was already noted by [START_REF] Torma | Added value of regional climate modeling over areas characterized by complex terrain-Precipitation over the Alps[END_REF]. The methodology used throughout this paper leaves an important message about comparing RCMs with its coarser-resolution driving data, as the model assessment should not only consider model performance but also the resolution of the model and the uncertainty in observations. As stated in the Introduction, the present work is a first step to assessing the AV of state-of-the-art RCMs over SA. However, to understand the potential added value of RCMs, further investigation should include the analysis of high-frequency temporal scales and small-scale spatial scales.

Fig. 1

 1 Fig. 1 CORDEX region 1 domain, with ETOPO2 (National Geophysical Data Center, 2006) terrain height (m) indicated in shading. Boxes in solid-black lines denote areas selected for area-mean assessment.

Fig. 2

 2 Fig. 2 Summer mean precipitation (mm/day; top panels) and 2-meter surface air temperature ( • C; bottom panels) of the observed fields and the multi-model ensembles. The first column is the average of observed datasets and the second, third, forth and fifth columns are the climatologies of RCMs-fERAi, ERAi, RCMs-fCMIP5, and CMIP5, respectively. Hatched areas indicate regions where the inter-model standard deviation is larger than the natural variability, represented by the mean interannual standard deviation of all simulations involved.

Fig. 3

 3 Fig. 3 As in Fig. 2 but for the winter season.

Fig. 4

 4 Fig.4Taylor diagrams showing the correlation, normalized standard deviation and centered RMSE for individual simulations (empty markers) and for multi-model ensembles (full markers). The left panel corresponds to summer and the right one to winter. The observation ensemble is considered as the reference data. Red, yellow, light-blue and green colors represent the RCMs forced by reanalysis, the reanalysis itself, RCMs forced by CMIP5 and CMIP5 models themselves, respectively. Black markers represent individual observed datasets. Circles correspond to seasonal precipitation climatology and triangles to seasonal surface air temperature climatology. The calculation excludes grid points with topography higher than 1500m.

Fig. 5

 5 Fig.5Box-and-whisker plots showing the mean temporal summer distribution statistics in each study region. Every box represents the inter-quartile range (25 th -75 th percentile), the median is indicated with a white line and the mean with a green triangle. The whiskers indicate the 95 th and 5 th percentiles. The regional time series of each simulation is obtained by spatially averaging all values within the region in every time step. Also, all percentiles of the multi-model ensembles are calculated as the average of individual percentiles. The left panels correspond to precipitation (mm/day; blue) and the right panels to 2-meter temperature ( • C; red).

MMFig. 6

 6 Fig. 6 As in Fig. 5 but for the winter season.

Fig. 7

 7 Fig. 7 Diagrams summarizing different errors or deviations obtained with our best estimation for each of our regions of interest. Left panels are for austral summer DJF, and right ones for austral winter JJA. Values plotted in the form of different symbols represent the spatial average of different terms in Eqs. (2) and (3). Top panels are for precipitation (mm/day), and bottom panels for 2-meter temperature ( • C). Red, yellow, light-blue and green colors represent our estimations in RCMs forced by reanalysis, reanalysis, RCMs forced by CMIP5 and CMIP5, respectively. Black markers represent our estimation for observational uncertainty. Grey empty markers correspond to interpolation (square markers) and representativeness (circle markers) errors.

Table 1

 1 Overview of Regional Climate Models used in the present study and their boundary conditions.

	Experiment Acronym	RCM id / Reference	Resolution /	Boundary condition id /	Boundary
			Number of	Reference	resolution
			grid points
	HadRM3P-fERAi	MOHC-HadRM3P Gordon et al. (2000)	0.44 • x0.44 • 146x167	ECMWF-ERAINT Dee et al. (2011)	1.5 • x1.5 •
	REMO2009-fERAi	MPI-CSC-REMO2009 Jacob et al. (2012)	0.44 • x0.44 • 143x167	ECMWF-ERAINT	1.5 • x1.5 •
	RCA4-fERAi	SMHI-RCA4 Kupiainen et al. (2014)	0.44 • x0.44 • 146x167	ECMWF-ERAINT	1.5 • x1.5 •
	RegCM4-fERAi	ICTP-RegCM4-3 Giorgi et al. (2012)	0.44 • x0.44 • 189x199	ECMWF-ERAINT	1.5 • x1.5 •
	WRF341I-fERAi	UCAN-WRF341I Skamarock et al. (2005)	0.44 • x0.44 • 146x167	ECMWF-ERAINT	0.75 • x0.75 •
	RCA4-fEC-EARTH	SMHI-RCA4 Kupiainen et al. (2014)	0.44 • x0.44 • 146x167	EC-EARTH Hazeleger et al. (2010)	1.25 • x1.25 •
	REMO2009-fMPI	MPI-CSC-REMO2009 Jacob et al. (2012)	0.44 • x0.44 • 143x167	MPI-ESM-LR Stevens et al. (2013)	1.875 • x1.875 •
	WRF341I-fCanESM2	UCAN-WRF341I Skamarock et al. (2005)	0.44 • x0.44 • 146x167	CCCma-CanESM2 Chylek et al. (2011)	1.875 • x1.875 •
	RegCM4-fHadGEM2	ICTP-RegCM4-3 Giorgi et al. (2012)	0.44 • x0.44 • 189x199	MOHC-HadGEM2-ES Collins et al. (2011)	1.875 • x1.25 •
	LMDZ4-fIPSL-CM5A	IPSL-LMDZ4 Hourdin et al. (2006)	∼0.48 • x∼0.48 • 184x180	IPSL-CM5A Hourdin et al. (2012)	3.75 • x1.9 •

Table 3

 3 Same as in Table2but for austral winter (JJA).

	Experiment	NAM		LPB		SSA
	Acronym	pr	tas	pr	tas	pr	tas
	RCMs-fERAi	1.05	0.41	1.25	1.06	1.69	0.94
	HadRM3P-fERAi	1.26	0.37	0.70	0.70	2.25	1.00
	REMO2009-fERAi	0.77	0.27	0.93	0.65	2.34	1.03
	RCA4-fERAi	1.61	0.50	2.73	1.42	1.64	1.00
	RegCM4-fERAi	1.24	0.41	1.08	1.20	1.26	1.27
	WRF341I-fERAi	1.64	1.08	1.58	2.68	1.58	1.17
	RCMs-fCMIP5	0.71	1.73	1.08	2.08	2.23	1.08
	RCA4-fEC-EARTH	1.08	0.84	1.94	1.98	1.62	1.04
	REMO2009-fMPI	0.65	0.71	1.11	1.37	2.31	0.84
	WRF341I-fCanESM2	0.69	1.20	0.74	3.89	1.55	1.35
	RegCM4-fHadGEM2	1.73	1.04	1.53	1.43	0.82	0.91
	LMDZ4-fIPSL-CM5A	1.10	0.94	0.70	2.09	2.43	1.18

Acknowledgements This work was supported by the Argentinean agencies CONICET (PIP 11220150100402CO) and ANPCyT (PICT 2014-0887, PICT-2015-3097), and by the French National Program LEFE/INSU (AO2015-876370). The first author was supported by CONICET and by the Saint Exupry Program (MED-MAEDI). The authors would like to thank the two anonymous reviewers for comments that helped improve the manuscript.


Catherine Damerval

Wafa Ben Othman

Domenica Manicacci

Florian Jabbour

Distribution area of the two floral morphs of Nigella damascena L. (Ranunculaceae): A diachronic study using herbarium specimens collected in France

Keywords: Collection bias, collecting effort, floral dimorphism, love-in-a-mist, museum time-series, natural history

Perianth dimorphism in Nigella damascena (Ranunculaceae) has been documented since four centuries, and its genetic basis has started to be studied almost one century ago. Based on herbarium collections, we reconstructed the evolution of the distribution area of the wild-type (perianth composed of a calyx and nectariferous petals) and mutant (apetalous) floral morphs in mainland France during the last two centuries. We showed that the collection localities of the wild-type morph are mostly concentrated in the southern coastal areas of France, whereas the collection localities of the mutant morph were less numerous but more scattered in France.

The collection effort of both morphs peaked in the second half of the 19 th century and continuously decreased afterwards. Results are interpreted taken into account the biases in the constitution of a herbarium collection. Knowledge about the spatial distribution and the frequency of occurrence of wild-type and mutant populations through time will help botanists and evolutionary biologists investigate the origin and ecological grounds for the maintenance and propagation of the apetalous morph in the wild.

Introduction

Floral dimorphism has intrigued botanists and evolutionary biologists since centuries. Darwin devoted a book to this topic, focusing on heterostylous Primula flowers [START_REF] Darwin | The different forms of flowers on plants of the same species[END_REF].

Nigella damascena L. (Ranunculaceae), or love-in-a-mist, is a herbaceous species mainly distributed in the Mediterranean Basin, whose ancient use by human populations is well documented (Heiss and Oeggl 2005). In this species, two types of flowers are known at least since the end of the 16 th century [START_REF] Clusius | Melanthium pleno flore[END_REF]. The wild-type flower has a perianth composed of a pentamerous petaloid calyx and a series of 5-10 pouch-shaped nectariferous petals. By contrast, the unipartite perianth of the mutant type consists of numerous petaloid tepals [START_REF] Deroin | Floral vascular patterns of the double-flowered and wild-type morphs of Nigella damascena L. (Ranunculaceae)[END_REF][START_REF] Deroin | Floral vascular patterns of the double-flowered and wild-type morphs of Nigella damascena L. (Ranunculaceae)[END_REF]. The monogenic control of the dimorphism was first established by [START_REF] Toxopeus | Erblichkeitsuntersuchungen an Nigella damascena L[END_REF] based on segregation analyses, and its molecular basis was unraveled almost one hundred year later [START_REF] Gonçalves | An APETALA3 homolog controls both petal identity and floral meristem patterning in Nigella damascena L. (Ranunculaceae)[END_REF][START_REF] Zhang | Disruption of the petal identity gene APETALA3-3 is highly correlated with loss of petals within the buttercup family (Ranunculaceae)[END_REF].

For simplicity, we will adopt the terminology of [START_REF] Gonçalves | An APETALA3 homolog controls both petal identity and floral meristem patterning in Nigella damascena L. (Ranunculaceae)[END_REF] to designate the wildtype (homozygous dominant or heterozygous for the combination of alleles controlling the dimorphism) and mutant (homozygous recessive) floral morphs, called [P] and [T], respectively. Human selection for plants with conspicuous flowers harbouring supernumerary petaloid organs likely played an important role in the long-term maintenance of the [T] phenotype, but [T] populations are also expected to occur in the wild because individuals are self-fertile and plants can easily escape from gardens and become naturalized in favorable environments.

Herbarium collections are repositories of past and current plant biodiversity. Compiling the data held by all the herbarium specimens of a species allows reconstructing the evolution of its distribution area through time. In the case of N. damascena, knowledge about the spatial distribution and the frequency of occurrence of [P] and [T] populations through time would help botanists and evolutionary biologists investigate the origin and ecological grounds for the maintenance and propagation of the [T] morph in the wild. It should be possible to examine whether [T] populations can persist per se in the wild or if they are consistently found in the neighborhood of towns and villages, suggesting constant gene flow from gardens; alternatively, [T] individuals could consistently be found in the vicinity of [P] individuals, suggesting that the presence of nectariferous flowers in [P] plants may favor pollination and reproductive success of [T] plants with flowers lacking pollinator reward, and/or that the mutation is slightly deleterious and preferentially maintained in the heterozygous state. A last use of distribution maps would be to identify localities where N. damascena individuals from each morph could be currently collected.

To get insight into these topics, we built a quasi-exhaustive database of the world collection of N. damascena herbarium specimens. For this special issue of Botany Letters, we present the results concerning mainland France, explaining the database structure and the protocol to select the herbarium specimens for georeferencing and building the temporal sequence of distribution maps. We observe different distributions of the two morphs in mainland France, and discuss possible causes.

Materials and Methods

Data collection and database construction

We built a database gathering the information found on herbarium specimens of N. damascena (Fig. 1) from 39 herbaria, from the Mediterranean Basin or from other regions but housing large amounts of Mediterranean species or species growing in areas with a Mediterranean climate. The herbaria sent specimens on loan to P (Paris Herbarium), where they were studied. They were photographed for the record. Alternatively, the herbaria sent images of their collection of N. damascena.

The data gathered from the biological specimen (Fig. 1A-B and1D) and from the label (Fig. 1C) are organized into 16 fields: (A) herbarium code; (B) barcode or reference in the herbarium; (C) collector's name; (D) collection number; (E-G) day, month, and year of collection; (H) country of origin, transcribed from the label or interpreted; current country names are provided; (I) collection locality; (J) other notes on the label (in the original language); (K) label mentioning that the specimen was in cultivation at the time of collection (yes/no); (L) assumption that the specimen was in cultivation at the time of collection (yes/no); (M) presence of flower(s) (yes/no); (N) flower morph (P/T, Fig. 1B and1D); (O) presence of capsule(s) (yes/no); (P) GPS coordinates in decimal degrees, in the WGS 84 system. In some cases, GPS coordinates were available directly from the label. In most cases, localities were georeferenced by us using Google Maps (Google Inc., Mountain View, CA, USA). The collection locality could be precisely specified on the label (e.g. "Gorges d'Estrabols, commune de Salvagnac, canton de Villeneuve"), allowing precise GPS coordinates to be determined (in red in Supplementary Table 2). Less precisely described collection sites, for example just mentioning the closest city, were also attributed GPS coordinates (in blue in Supplementary Table 2). In other cases, e.g. when only the country of origin is specified or when the plant was in cultivation when collected with no mention of the seed origin, localization data was considered unavailable and (P) cells were marked with a double-slash (//). Except for (H) and (L), all data are transcribed from the label without interpretation or come from the observation of the specimen. Label data was provided as a csv file by the following herbaria: BC, BEI, BM, E, FR, HUJ, and MA; it was double-checked before being included in our database.

Diachronic study of the distribution of Nigella damascena in France

For the present study, we focused on data relative to all specimens collected in France (654 specimens; Supplementary Table 1, Fig. 1), recorded from 20 herbaria (B (Berlin), BC (Barcelona), BEI (Beirut), BM (London), CAI (Cairo), E (Edinburgh), FI (Florence), FR (Frankfurt), G (Geneva), HUJ (Jerusalem), K (Kew), MA (Madrid), MEL (Melbourne), MPU (Montpellier), P (Paris), PAL (Palermo), S (Stockholm), TO (Torino), UPS (Uppsala), W (Vienna)); herbarium codes follow [START_REF] Thiers | Index Herbariorum: A global directory of public herbaria and associated staff, New York Botanical Garden's Virtual Herbarium[END_REF]. For the analysis of the distribution of the two morphs over time, we did not consider specimens without flowers (or with flowers compressed in such a way that the floral morph was not distinguishable), duplicate specimens (belonging to the same collection, i.e. same collector and collection number) and specimens with no collection year. We considered two datasets extracted from the database:

(i) the georeferenced dataset excludes specimens with a double slash in (P), and gathered 409 specimens from 18 herbaria (all specimens from BEI and CAI herbaria were discarded). Among these specimens, 12 were composed of two different collections (different combination [collector + collection number]). After splitting these accessions into two, the final dataset consisted of 421 entries. Supplementary Table 2 shows the collections that were used to build the distribution maps. We classified the entries into four time intervals ([1800-1849], , , ), and into two categories ([P] and [T] morph) within each interval. When the collection date is an interval, the latest year was considered.

After the sorting process described above, no specimen collected in France after 2000 was included in the dataset. Maps were built using SimpleMappr [START_REF] Shorthouse | SimpleMappr, an online tool to produce publication-quality point maps[END_REF].

(ii) the [P]/[T] dataset was built to compare the number of collected specimens from each morph among each 50-year period. It includes the cultivated specimens and the specimens that could not be georeferenced, and comprises 469 entries.

Results

More than two centuries of biogeographical data captured in herbarium specimens

The oldest and most recent French specimens of N. damascena included in our database were collected in 1784 by an unknown collector and in 2015 by F. Jabbour, respectively. These specimens were both cultivated (in "Jardin de Sally" and in Gif-sur-Yvette (see [START_REF] Deroin | Floral vascular patterns of the double-flowered and wild-type morphs of Nigella damascena L. (Ranunculaceae)[END_REF], respectively), and hence not georeferenced. The whole dataset consists of 654 collections. The three herbaria that contributed most in terms of number of collections were P (Paris, 313 collections), MPU (Montpellier,105) and G (Geneva,43). The floral morph could be identified on 469 specimens (371 [P] and 98 [T]) and 421 specimens (352 [P] and 69 [T]) could be georeferenced. The three most active collectors of N. damascena specimens in France were Elisée Reverchon (1835Reverchon ( -1914;;11 specimens), Florian Jabbour (1983-;9 specimens), andBrother Sennen (Etienne Marcellin Granié-Blanc, 1861-1937;8 specimens).

The collection localities of the [P] morph are mostly concentrated in the southern coastal areas of France

A peak of collection in the second half of the 19 th century

The collection of both morphs was the highest between 1850 and 1899 (233, about half of the specimens of the dataset), before decreasing continuously afterwards (Fig. 3A). The peak also coincides with more individuals from both morphs collected from the northern half of France.

The collection of cultivated specimens, even though less numerous, follows the same temporal trend as the wild ones (Fig. 3B). In the wild collections, [P] morph individuals are more frequent than [T] morph ones in all time intervals (2.6 to 10 times more). The proportion of [T] morph individuals was the highest in the first half of the 20 th century, reaching 28%. In contrast to the wild collection, more [T] than [P] individuals were found in the cultivated collection at each time period (Fig. 3B). During the last two centuries, 28 of the c. 300 collectors of N. damascena collected both floral morphs, but only 15 of them collected both morphs in the same locality.

Discussion

Caution required when analyzing herbarium collection records

A collection of natural history is by nature biased. Not all the biodiversity can be captured in time and space, and the representativeness of a collection depends on the spatio-temporal uniformity of the collection effort [START_REF] Auderset | Occurrence of Characeae in Switzerland over the last two centuries (1800-2000)[END_REF]. Many collectors collect objects of natural history in the places where they live or travel, and they can be biased towards or against unusual objects (e.g., with an unusual morphology or spatial occurrence). They can collect in priority individuals with unusual phenotypes in order to have a record in their own collection, or they can ignore them because they deviate from the standard morphology. Indeed, botanists can purposely ignore plants considered as weeds and plants that are locally common. Hence a herbarium collection can paradoxically lack plants that are common or present at a given locality for long. Knowing these limitations, herbarium collection data should not be over-interpreted. They give a picture of biodiversity at certain times and places, filtered by collectors' interest that may vary according education and societal context. For a given species, a way to partly circumvent bias is to consider the most exhaustive set of herbarium collections from the highest number of herbaria. This is what we intended to do for N. damascena.

The evolution of Nigella damascena collection effort in France

In France, the collection effort concerning N. damascena has not been uniform during the last two centuries. After a peak in the second half of the 19 th century, the number of specimens added to collections of natural history decreased continuously. Does it reveal that N. damascena distribution area has expanded and then shrunk down, or that the collection effort has been lower during the last decades? We tend to favor the second hypothesis. Indeed, we cannot figure out any natural cause that could explain a decrease in abundance in the wild, since this species has been common in France (see SIFlore and INPN databases;Just et al. 2015;Muséum national d'Histoire Naturelle 2003[START_REF] Thiers | Index Herbariorum: A global directory of public herbaria and associated staff, New York Botanical Garden's Virtual Herbarium[END_REF]Fig. 4). However, the decrease in specimen number collected in the first half of the 20 th century could be due to the impact of the two world wars on the activity of the scientists (Le Bras et al. 2015). Furthermore, we witness nowadays a general decline in collecting due to less funding for research in the field of natural history and the societal perception that collecting local flora is unnecessary [START_REF] Renner | Is plant collecting in Germany coming to an end?[END_REF].

The fact that N. damascena has often been collected next to areas with high human population density could be due to the local presence of many collectors and to the naturalization of populations after escaping from gardens. Whilst this explanation could stand for the Parisian region, we do not think it is valid for the whole French territory, as N. damascena has not been heavily collected from many densely inhabited areas in northern France. Thus, we believe that the distribution of the species at the Atlantic coast, at the Rhône Valley and at the Mediterranean coast represent its true distribution area. The data made available online by SIFlore (Just et al. 2015), gathering distributional data for vascular plants in France from herbarium records but also from the literature and data from botanical conservatories and covering the last five centuries (from 1545 to 2014) and by the INPN (Muséum national d'Histoire Naturelle 2003[START_REF] Thiers | Index Herbariorum: A global directory of public herbaria and associated staff, New York Botanical Garden's Virtual Herbarium[END_REF] are in accordance with our results, considering that these datasets distinguishes neither the [P] from the [T] morph nor wild from cultivated specimens (Figure 4).

The morph [T] has been less collected than the [P] morph over the last couple of centuries suggesting that such populations are less common in the wild than the [P]-morph populations (72/430 = 17%). One-tenth of the collectors of N. damascena collected both floral morphs, and more than half of them found both morphs in the same locality. This suggests that at least 21% (15/72) of the [T] phenotype appears in mixed [P]-[T] populations.

Additionally, the two morphs are not equally distributed in France. Relatively more [T] specimens have been collected in the northern part of France and there is no reason to think that the collectors mostly living or botanizing in the northern half of France were more sensitive to unusual phenotypes. The presence of [T] populations further away from the Mediterranean and the Atlantic coasts could be due to the naturalization of populations in the vicinity of towns and villages after escaping from the gardens. The fact that most cultivated specimens, which are enriched in [T] morphs, were also collected from the northern half of France is consistent with this hypothesis.

A dataset valuable for diachronic evo-devo studies

Museum time-series (set of specimens from the same species collected at different times) are the basic material for floras and conservation studies. Thanks to the information gathered on the collection label, the distribution area of a single species can be traced over time (e.g., [START_REF] Auderset | Occurrence of Characeae in Switzerland over the last two centuries (1800-2000)[END_REF]. Such diachronic studies can be conducted at the scale of some centuries and do not result from any inference-based analysis. They differ from studies of historical biogeography in which older (thousand-to million-year old) dispersal events are inferred from the distribution areas of the extant taxa and the topology of a phylogenetic tree.

The distribution maps we built are based on the georeferencing of the collection localities of 421 individuals across two hundred years. Among this set of individuals, only three could be georeferenced with high precision, meaning the maps presented in Fig. 2 mostly show broad localities. However, comparing maps obtained from [P] and [T] specimens can allow us targeting localities where populations of both morphs could currently coexist. Such natural populations would be required for in situ long-term ecological studies.

Ancient DNA (i.e. extracted from relatively old herbarium specimens) has been sequenced to track invasion routes (Matsuhashi et al. 2016) or analyse the evolution of the genetic diversity of a plant pathogen (Li et al. 2007) through time. In the case of N. damascena, mining our database for old and modern specimens from different origins will make it possible to test whether the molecular basis for the [T] phenotype is identical at the scale of the Mediterranean Basin as a whole, during the last two centuries. This research is currently under way in our team.

Disclosure statement

No potential conflict of interest was reported by the authors.

Notes on contributors

Catherine Damerval is senior researcher at Centre national de la Recherche Scientifique. Her research deals with evo-devo of the flower in Ranunculales.

Wafa Ben Othman is an agronomist, recently graduated from AgroParisTech. Domenica Manicacci is an assistant professor at Université Paris-Sud. She is a population geneticist working on different angiosperm models such as maize and Nigella damascena. 1. Herbarium specimens of Nigella damascena collected in France, label data transcribed and interpreted. For each specimen, the data is organized into 16 fields: 2. This table corresponds to Supp. Table 1 from which all the entries not taken into consideration (see our sorting procedure in the Materials and Methods section) for building the distribution maps were excluded.

Supplementary Table

Supplementary Table

Figure 2

 2 Figure2shows the evolution of the distribution area of both floral morphs of N. damascena

Florian

  Jabbour is an associate professor at the Muséum national d'Histoire Naturelle and a curator at Paris Herbarium. His research focuses on the evolution of plant taxa and plant shapes. He is editor for Botany Letters. For this publication, CD, WBO, DM and FJ gathered data, CD and FJ analyzed the results and wrote a first version of the manuscript, and all coauthors contributed to the last version of the manuscript. Heiss, A. G., and K. Oeggl. 2005. The oldest evidence of Nigella damascena L. (Ranunculaceae) and its possible introduction to central Europe. Vegetation History and Archaeobotany 14 (4): 562-570. doi:10.1007/s00334-005-0060-4. Jabbour, F., M. Udron, M. Le Guilloux, B. Gonçalves, D. Manicacci, S. Nadot, and C. Damerval. 2015. "Flower development schedule and AGAMOUS-like gene expression patterns in two morphs of Nigella damascena (Ranunculaceae) differing in floral architecture." Botanical Journal of the Linnean Society 178 (4): 608-619. doi:10.1111/boj.12297. Just, A., J. Gourvil, J. Millet, V. Boullet, T. Milon, I. Mandon, and B. Dutrève. 2015. "SIFlore, a dataset of geographical distribution of vascular plants covering five centuries of knowledge in France: Results of a collaborative project coordinated by the Federation of the National Botanical Conservatories." PhytoKeys 56: 47-60. doi: 10.3897/phytokeys.56.5723. Le Bras, G., M. Pignal, M. L. Jeanson, S. Muller, C. Aupic, B. Carré, G. Flament, M. Gaudeul, C. Gonçalves, VR. Invernón, F. Jabbour, E. Lerat, P. P. Lowry, B. Offroy, E. P. Pimparé, O. Poncy, G. Rouhan, and T. Haevermans. 2017. "The French Muséum national d'histoire naturelle vascular plant herbarium collection dataset." Scientific Data 4: 170016. doi:10.1038/sdata.2017.16. Li, W., Q. Song, R. H. Brlansky, and J. S. Hartung. 2007. "Genetic diversity of citrus bacterial canker pathogens preserved in herbarium specimens." Proceedings of the National Academy of Sciences 104 (47): 18427-18432. doi:10.1073/pnas.0705590104. Matsuhashi, S., H. Kudoh, M. Maki, M. Cartolano, M. Tsiantis, T. Itagaki, and S. Sakai. 2016. "Invasion history of Cardamine hirsuta in Japan inferred from genetic analyses of herbarium specimens and current populations." Biological Invasions 18 (7): 1939-1951. doi:10.1007/s10530-016-1139-9.

Figure legends Figure 1 .

 legends1 Figure legends

Figure 2 .

 2 Figure 2. Distribution area of Nigella damascena in France during the last two centuries,

Figure 3 .

 3 Figure 3. Frequency of [P] and [T] morphs of Nigella damascena in a Museum time-series

Figure 4 .

 4 Figure 4. Distribution maps of Nigella damascena (both morphs, wild and cultivated) based

(

  A) Herbarium code. (B) Barcode or reference in the herbarium. (C) Collector's name. (D) Collection number. (E-G) Day, month and year of collection. (H) Country of origin, transcribed from the label or interpreted. Current country names are provided. (I) Collection locality. (J) Other notes on the label. (K) Is there anything on the label specifying that the specimen was in cultivation at the time of collection? (L) Assumption that the specimen was in cultivation at the time of collection. (M) Presence of flower(s). (N) Flower morph. (O)Presence of capsule(s). (P) GPS coordinates, transcribed from the label or inferred. Except for (H) and (L), all data are transcribed from the label without interpretation or come from the observation of the specimen.

Acknowledgments

We thank the curators and staff members of B (Robert Vogt and Ulrike Starck), BC (Neus Nualart Dexeus), BEI (Nada Sinno Saoud), BM (Mark Carine and Jacek Wajer), CAI (Hasnaa Hosni), E (Lesley Scott), FI (Lorenzo Cecchi), FR (Georg Zizka and Rainer Doering), G (Cyrille Chatelain), HUJ (Jotham Ziffer-Berger), K (Marie-Hélène Weech), MA (Charo Noya Santos), MEL (Helen Barnes), MPU (Caroline Loup), PAL (Gianniantonio Domina), S (Arne Anderberg), TO (Laura Guglielmone), UPS (Mats Hjertson), W (Ernst Witek), for sending herbarium specimens on loan and/or sending digital images. We thank the staff members of P Herbarium for their technical help and the Institute for Systematics, Evolution and Biodiversity for funding the internship of W. Ben Othman. We also thank two anonymous reviewers for their critical and helpful comments and suggestions.