François Duchenne
email: françoisduchenne@mnhn.fr

Gabrielle Martin

Emmanuelle Porcher

European plants lagging behind climate change pay a climatic debt in the North, but are favored in the South

Keywords: range shifts, occupancy trends, climate change, plant distribution, global change drivers

For many species, climate change leads to range shifts that are detectable, but often insufficient to track historical climatic conditions. These lags of species range shifts behind climatic conditions are often coined "climatic debts", but the demographic costs entailed by the word "debt" have not been demonstrated. Here we used opportunistic distribution data for ~4,000 European plant species to estimate the temporal shifts in climatic conditions experienced by these species and their occupancy trends, over the last 65 years. The resulting negative relationship observed between these two variables provides the first piece of evidence that European plants are already paying a climatic debt in Alpine, Atlantic and Boreal regions. In contrast, plants appear to benefit from a surprising "climatic bonus" in the Mediterranean. We also find that among multiple pressures faced by plants, climate change is now on par with other known drivers of occupancy trends, including eutrophication and urbanization.

Introduction

Climate change is recognized as a major threat for biodiversity [START_REF] Sala | Global Biodiversity Scenarios for the Year 2100[END_REF]. Species have two main ways to persist under such change: they can track their climatic optimum in space (geographic range shift) or they can respond adaptively to survive and reproduce in a altered environment [START_REF] Chevin | Adaptation, Plasticity, and Extinction in a Changing Environment: Towards a Predictive Theory[END_REF], i.e. modify their climatic optimum. Concerns about the capacity of species to achieve one or the other quickly enough have increased [START_REF] Orme | caper: Comparative Analyses of Phylogenetics and Evolution in R[END_REF][START_REF] Chevin | Adaptation, Plasticity, and Extinction in a Changing Environment: Towards a Predictive Theory[END_REF][START_REF] Hoffmann | Climate change and evolutionary adaptation[END_REF], with many studies showing that living organisms are currently moving poleward and upward in response to climate warming [START_REF] Orme | caper: Comparative Analyses of Phylogenetics and Evolution in R[END_REF][START_REF] Kelly | Rapid shifts in plant distribution with recent climate change[END_REF][START_REF] Devictor | Differences in the climatic debts of birds and butterflies at a continental scale[END_REF][START_REF] Lenoir | Species better track climate warming in the oceans than on land[END_REF]). Yet, most studies so far have shown that range shifts are rarely as fast as climate change [START_REF] Menéndez | Species richness changes lag behind climate change[END_REF][START_REF] Devictor | Differences in the climatic debts of birds and butterflies at a continental scale[END_REF][START_REF] Vanderwal | Focus on poleward shifts in species' distribution underestimates the fingerprint of climate change[END_REF][START_REF] Lenoir | Species better track climate warming in the oceans than on land[END_REF], i.e. climate variables move faster through space than most species do. This lag of species movements behind climate change is often coined a "climatic debt" [START_REF] Devictor | Differences in the climatic debts of birds and butterflies at a continental scale[END_REF][START_REF] Monsinjon | The climatic debt of loggerhead sea turtle populations in a warming world[END_REF][START_REF] Lenoir | Species better track climate warming in the oceans than on land[END_REF]. It can be evidenced by a temporal change in the so-called "climatic niche" of a species, as measured by the average of one or several climatic variables throughout its range [START_REF] Vanderwal | Focus on poleward shifts in species' distribution underestimates the fingerprint of climate change[END_REF], hereafter "species climatic indices" (SCIs).

The consequences of species shifting their ranges slower than the movement of suitable climate conditions, however, remain to be quantified. First, the resulting temporal change in the climatic conditions experienced by a species need not necessarily translate into a "debt": most species can and do also respond to climate change via adaptive plastic or evolutionary trait changes, which could be sufficient to sustain populations despite changing climatic conditions [START_REF] Hoffmann | Climate change and evolutionary adaptation[END_REF], albeit the most recent meta-analysis to date suggests otherwise [START_REF] Radchuk | Adaptive responses of animals to climate change are most likely insufficient[END_REF]. Thus, to verify the existence of a true climatic debt, one needs to demonstrate that the spatial lag of species behind changing climate conditions results in decreased individual fitness or population growth; this has not yet been done. Second, most species face additional threats beyond climate change, e.g. habitat loss or pollution, which may be the dominant drivers of current species loss [START_REF] Maxwell | Biodiversity: The ravages of guns, nets and bulldozers[END_REF]. Thus, a climatic debt, should it exist, might be of little consequence for population trends compared to other drivers. However, few studies so far have compared quantitatively the impact of different drivers across a large number of species, and none has included the climatic debt.

Costs associated with the putative climatic debt remain poorly investigated likely because the climatic debt concept has been developed at the community level mostly [START_REF] Bertrand | Changes in plant community composition lag behind climate warming in lowland forests[END_REF][START_REF] Devictor | Differences in the climatic debts of birds and butterflies at a continental scale[END_REF], while cost estimation is easier at species level, e.g. via an assessment of species persistence. By shifting the concept of climatic debt from community to species level, one can correlate delayed spatial responses (i.e. temporal trends in SCIs) with species persistence. Under climate warming, if both species movements and adaptive responses are insufficient, we expect (1) an increase in the temperature SCI of a given species over time (limited spatial response) and (2) a negative relationship between the temporal trend in SCI for temperature and the species persistence (limited adaptive response, [START_REF] Radchuk | Adaptive responses of animals to climate change are most likely insufficient[END_REF]. The latter observation only is suggestive of costs at individual and population levels, i.e. a climatic debt. Alternatively, uncorrelated SCI trends and species persistence would indicate an absence of climatic debt, which could be explained by adaptive responses buffering limited spatial responses or by species insensitivity to temperature [START_REF] Rodríguez-Sánchez | Uncertainty in thermal tolerances and climatic debt[END_REF]. This logic applies similarly to other climatic variables beyond temperature, such as precipitation.

Here we examine the temporal shifts in the climatic conditions experienced by a species throughout its range, as measured by SCI trends, and occupancy trends, for more than 4,000 European plant species over the last 65 years, using a large dataset of opportunistic distribution records. These two elements allowed us to verify the existence of a climatic debt, which we estimate via the relationship between occupancy and SCI trends, and to compare its strength to other potential drivers of occupancy trends besides climate change via a trait-based approach.

Methods

Data collection; trends in SCIs and species occupancy

Plant database

We focused on the most common vascular plant species in geographical Europe, i.e. with at least 500 records in the Global Biodiversity Information Facility (GBIF) database (https://www.gbif.org) between 1950 and 2014 within a rectangle bounded by longitudes [-13°, 34°] and latitudes [34°,75°]. We downloaded all species sightings during the time period from 1951 through 2014, excluding 1950 as this year contains data not precisely dated, but corresponding to the mid-twentieth century. The DOIs associated with the extraction are presented in the Supplementary methods. We considered only records from the European mainland, stopping at 34° of longitude because there were too few data to the East of this meridian (Fig. S1). Note that the area also includes the western part of Turkey. This yielded a dataset containing 111,549,494 occurrence records, characterized by a species name, a location and a date. Of these, we analyzed temporal trends in species climatic indices (SCIs) and occupancy only for the species observed in at least 20 years between 1951 and 2014 and with at least one record between 1951 and 1980. We removed crop species and considered invasive species separately (see Supplementary methods), because the drivers of their occupancy trends are likely different from those for non-invasive species. This selection resulted in 4,120 native and naturalized plant species (listed in Table S1), plus 58 invasive species.

Bioclimatic variables

Climate change is not limited to increases in annual mean temperature; hence we characterized the climatic conditions with three bioclimatic variables related to temperature and three bioclimatic variables related to precipitation, because temperature and precipitation are strong predictors of plant distribution [START_REF] Franklin | Modeling plant species distributions under future climates: how fine scale do climate projections need to be?[END_REF]. We used previously published European time series [START_REF] Fréjaville | The EuMedClim Database: Yearly Climate Data (1901-2014) of 1 km Resolution Grids for Europe and the Mediterranean Basin[END_REF] to extract annual mean temperature, maximum temperature of the warmest month, minimum temperature of the coldest month, annual precipitation, precipitation of the wettest month, and precipitation of the driest month.

These bioclimatic variables are the same as in Worldclim (bio1,5,6,12,13,14) but with annual data. For computational reasons we aggregated 1kmx1km raster cells, which decreased spatial resolution from 1 to about 100km² (~10km×10km).

Calculation of annual species climatic indices and their temporal trends

Species climatic indices (SCIs) are often calculated as the average of a climatic variable, e.g. temperature, across a species range [START_REF] Devictor | Differences in the climatic debts of birds and butterflies at a continental scale[END_REF]. However, the heterogeneity of opportunistic datasets can bias this index [START_REF] Loiselle | Predicting species distributions from herbarium collections: does climate bias in collection sampling influence model outcomes?[END_REF][START_REF] Beck | Spatial bias in the GBIF database and its effect on modeling species' geographic distributions[END_REF]). Because we are interested in inter-annual comparisons, we corrected a temporal bias in the spatial distribution of sampling pressure: the average latitude and longitude of the GBIF records both decrease significantly with time. To reduce such bias we defined annual SCI as the mean climatic variable of all 100km² cells occupied by a given species, while weighting the contribution of each cell by the ratio of the number of records of this species on the number of records of all plant species for the given year and grid cell. Such method enables estimations of the SCIs that are more independent from the sampling pressure than a weighting by the number of records of the given species only (Fig. S2).

Following the method above, we calculated one SCI for each bioclimatic variable, species and year. We then assessed the temporal trend in each SCI and species separately, using the following linear model:

𝑆𝐶𝐼 𝑘 = 𝜇 + 𝛽 × 𝑦𝑒𝑎𝑟 𝑘 + 𝜀 𝑘 (1)
where 𝑆𝐶𝐼 𝑘 is the species climatic index of year k, 𝜇 is the grand mean (intercept), 𝛽 is the year effect and 𝜀 𝑘 is an error term (independent and identically distributed, following N(0,σ²)).

Observations (𝑆𝐶𝐼 𝑘 ) are weighted by the square root of the number of grid cells included in their calculation for each year, which gives greater weight to years with more data. We finally estimated the phylogenetic signal in SCI trends (Supplementary Methods).

Occupancy trends

To model the temporal occupancy trends, we first discretized the dataset spatially and temporally, to define areas occupied or not by a species for a given time period. Such discretization then allowed us to estimate variation in occupancy probability among time periods by taking the sampling pressure into account. As for SCI calculation, we used a grid cell of about 10km×10km to discretize the dataset spatially, and we aggregated records temporally by years. For each year, a grid cell is considered as visited by at least an observer if it contains at least one plant species record. Non-visited cells are discarded. A visited cell is considered as occupied by a species if it contains at least one record of the given species in the given year, and unoccupied otherwise. We obtained a dataset composed of annual presences and pseudo-absences in each 10km×10km grid cell.

Before analyzing the data, we discarded all grid cells visited only one year, to improve occupancy estimations by decreasing the confusion between grid cell and year effects [START_REF] Isaac | Statistics for citizen science: extracting signals of change from noisy ecological data[END_REF]. To save computing time, for each species we removed non-informative grid cells, i.e.

cells with no record of the species over the whole study period. Finally, to estimate a yearly occupancy probability (p), we explained remaining presences and pseudo-absences for each species separately using the following binomial generalized linear model with a logit link:

log ( 𝑝 𝑖𝑘 1-𝑝 𝑖𝑘 ) = 𝜇 + 𝛽 1 × 𝑦𝑒𝑎𝑟 𝑘 + 𝛽 2 × log(𝑆𝐿 𝑖𝑘 ) + 𝜑 𝑖 (2)
Where 𝑝 𝑖𝑘 is the occupancy probability of grid cell i for year k, 𝜇 is the intercept of the model, 𝛽 1 is the year effect (i.e. the occupancy trend), and 𝜑 𝑖 a random grid cell effect. Finally, 𝛽 2 is the effect of the logarithm of the species list length (𝑆𝐿 𝑖𝑘 , i.e. the number of species observed in a given grid cell and year), used as a proxy for the sampling pressure [START_REF] Isaac | Statistics for citizen science: extracting signals of change from noisy ecological data[END_REF].

Those models were implemented using the R package glmmTMB [START_REF] Brooks | glmmTMB Balances Speed and Flexibility Among Packages for Zeroinflated Generalized Linear Mixed Modeling[END_REF]. We estimated occupancy trends for all 4,120 native and naturalized species, as well as for the 58 invasive species. This latter step allowed us to verify that occupancy trend estimates were large and positive for invasive species (Fig. S3), confirming that the data and statistical methods to estimate occupancy trends yield results consistent with known species trends. As for SCI trends, we estimated the phylogenetic signal in occupancy trends (Supplementary Methods).

Additional potential drivers of occupancy trends

To compare the strength of the climatic debt with that of other potential drivers of occupancy trends, we tried to include all additional drivers in the analysis of the relationship between occupancy trends and SCI trends (see below). We used the plant and global change literature to identify potential drivers of plant trends, which we took into account via species traits: historical climatic niche [START_REF] Martin | Short-term climateinduced change in French plant communities[END_REF], lifespan [START_REF] Martin | Short-term climateinduced change in French plant communities[END_REF], habitat affinity [START_REF] Aronson | A global analysis of the impacts of urbanization on bird and plant diversity reveals key anthropogenic drivers[END_REF][START_REF] Buse | Conservation of threatened habitat types under future climate change -Lessons from plant-distribution models and current extinction trends in southern Germany[END_REF], nitrophily [START_REF] Sala | Global Biodiversity Scenarios for the Year 2100[END_REF][START_REF] Bobbink | Global assessment of nitrogen deposition effects on terrestrial plant diversity: a synthesis[END_REF], moisture preferences [START_REF] Moeslund | Topographically controlled soil moisture drives plant diversity patterns within grasslands[END_REF]) and pollinator dependency [START_REF] Biesmeijer | Parallel Declines in Pollinators and Insect-Pollinated Plants in Britain and the Netherlands[END_REF]. We detail all calculations for these traits in Supplementary Methods.

Species traits, SCI and occupancy trends are available in Table S1.

Evidencing the climatic debt: relationship between SCI trends and occupancy trends

Analysis at continental scale

As explained above, a climatic debt can be revealed by a negative relationship between species occupancy trends and SCI trends. To assess this relationship between SCI trends and occupancy trends, while controlling for the species traits cited above, we used linear models.

We first checked the correlations among the six SCI trends and six historical climatic niche indices. We noticed high correlations (r>0.7) among historical temperature indices, among historical precipitation indices, among SCI trends related to temperature and among SCI trends related to precipitation (Fig. S4). In order to avoid multicollinearity issues, we retained only the two most used and integrative bioclimatic variables, both for the historical climatic niche and for SCI trends: annual mean temperature and annual precipitation.

In summary, we considered the following correlates of occupancy trends: SCI trends related to annual mean temperature and annual precipitation, historical climatic index related to annual mean temperature and annual precipitation, nitrophily and moisture Ellenberg indicator values, pollinator dependency, habitat affinity and lifespan. To be able to compare the strength of relations across explanatory variables, we scaled them before using the Phylogenetic Generalized Least Squares regression, except for habitat affinity and lifespan, as we wanted to extract the effect on the intercept for the first and as the second is a qualitative variable. We used the caper R package [START_REF] Orme | caper: Comparative Analyses of Phylogenetics and Evolution in R[END_REF] to implement models controlling for phylogenetic signal in the residuals:

∆𝑂 𝑠𝑗 = 𝜇 + 𝛽 1 × ∆𝑆𝐶𝐼 𝑏𝑖𝑜1 𝑠 + 𝛽 2 × ∆𝑆𝐶𝐼 𝑏𝑖𝑜12 𝑠 + 𝛽 3 × 𝐵𝐶𝐼 𝑏𝑖𝑜1 𝑠 + 𝛽 4 × 𝐵𝐶𝐼 𝑏𝑖𝑜12 𝑠 + 𝛽 5 × 𝑀𝐸 𝑠 + 𝛽 6 × 𝑁𝐸 𝑠 + 𝛽 7 × 𝑃𝑜𝑙𝑙 𝑠 + 𝜃 𝑗 + ∑ 𝛽 ℎ𝑎𝑏𝑖𝑡𝑎𝑡 𝑖 × 𝐻𝑎𝑏𝑖𝑡𝑎𝑡 𝑖 ℎ 𝑖=1 + 𝜀 𝑠𝑗 (3)
Where ∆𝑂 𝑠𝑗 is the occupancy trend of species s with lifespan class j, 𝜇 is the grand mean (intercept), 𝛽 1 and 𝛽 2 are the slopes of SCI trends, 𝛽 3 and 𝛽 4 are the effects of historical climatic niche indices, 𝐵𝐶𝐼 𝑏𝑖𝑜1 and 𝐵𝐶𝐼 𝑏𝑖𝑜12 , related to annual mean temperature and annual precipitation respectively. 𝛽 5 and 𝛽 6 are the effects of Ellenberg indicator values for moisture and nitrogen respectively. 𝛽 7 is the effect of pollinator dependency, 𝜃 𝑗 is the qualitative lifespan effect and 𝛽 ℎ𝑎𝑏𝑖𝑡𝑎𝑡 𝑖 is the effect of affinity to habitat i, with six habitat classes (h = 6, Table S2), woodland being the reference habitat. Finally, 𝜀 𝑠𝑗 is an error term, after correction by Pagel's λ value at the likelihood maximum.

This model included all species matching the phylogeny and with full trait data (n = 2,013).

To add the 67 species that were not in the phylogeny but for which all traits were available, we also performed a linear mixed-effect model similar to the phylogenetic regression but including random taxonomic effects of class (𝜑 𝑐 ) and of family nested in class (𝜑 𝑓 ) instead of a phylogeny: + 𝜑 𝑐 + 𝜑 𝑓 + 𝜀 𝑠𝑗𝑐𝑓 (4)

∆𝑂 𝑠𝑗𝑐𝑓 = 𝜇 +
This linear mixed-effect model was weighted by the inverse of the standard errors associated with the occupancy trends.

Analysis by biogeographic region and time period

To examine the spatial variation in climatic debt, we also conducted the same analysis within biogeographic regions. We focused on the five biogeographic regions with at least 1,000 species: Alpine, Atlantic, Boreal, Continental and Mediterranean regions. For each region, as for the main analysis, we retained only species with at least 20 years of data.

We re-calculated SCI and occupancy trends within each region independently, using the same method as above, but including only plant records from the focal biogeographic region.

With these new estimates, we assessed regional climatic debts as the relationship between SCI and occupancy trends, using the same models as in the European analysis (equations ( 3) and ( 4)) but removing the random effect of taxonomic class 𝜑 𝑐 , to avoid singularities. While SCI and occupancy trends were calculated within biogeographic regions, we considered other species traits as constant throughout Europe.

Finally, to confirm that the costs of the shifts in experienced climatic conditions occur only after the acceleration of climate change (1980-2000, Fig. S5) we performed the same set of analyses on the earliest data, as detailed in Supplementary Methods. Results are shown in Fig. S6.

Combining the effects of precipitation and temperature SCI trends to estimate an overall climatic debt

To combine the effects of annual precipitation and temperature SCI trends on occupancy trends for each biogeographic region, we multiplied each unscaled effect of SCI trends by the observed change in the corresponding climatic indices. This yielded a measure of the effective cost/bonus due to a given SCI trend, i.e. to the lag of species behind climate change. For each region, we then summed the values of effective cost/bonus for temperature and precipitation to measure the overall climatic debt/bonus over the study period. We also calculated the relative contribution of precipitation and temperature SCI trends to the overall debt/bonus by dividing each by the sum of the absolute values of both.

Results

Temporal trends in Species Climatic Indices

During the study period , all SCIs change significantly, in direct relation with climate change. All bioclimatic variables related to temperature and precipitation increase on average over the study area, but with substantial spatial heterogeneity for precipitation (Fig. S5). Consistently with the trends in bioclimatic variables related to temperature, we show that the temperature SCIs increase over time for a large majority of species (Fig. 1). Precipitation

SCIs also increase over time on average (Fig. 1, Table S3), but the distribution of precipitation SCI trends is closer to zero, with more numerous negative trends than for temperature SCIs.

Temperature and precipitation SCI trends are not significantly correlated (Fig. S4), probably due to the fact that annual precipitation and temperature trends exhibit contrasting spatial distributions (Fig. S5c-d). We also find a significant phylogenetic signal in all SCI trends (Fig. S7).

Occupancy trends and their drivers

The average trend in occupancy over all native and naturalized species is slightly positive: 0.0048 ± 3.29 10 -4 year -1 (mean ± SE), but the number of species with significant increase in the occupancy estimates (1,721 species, 42%) is comparable to the number of species with significant decline in occupancy (1,519 species, 37%). The average positive trend over all species (Fig. 2a) is mainly explained by the skewed distribution of trends towards positive values, with a couple of species exhibiting strong increases even after exclusion of invasive species (Fig. S3). Furthermore, we find that plant occupancy trends exhibit a strong The analysis of the correlates of species occupancy trends reveals that plant species pay a climatic debt, but only in some parts of Europe. While at the continental scale, the negative relationship between SCI trends and occupancy trends is not significant (Fig. 2b), analysis by biogeographic regions reveals significant correlations, with a strong heterogeneity among regions and across bioclimatic variables. Regarding temperature, in the two coldest biogeographic regions (i.e. Boreal and Alpine regions), a temporal increase in temperature throughout a species range, a consequence of an insufficient range shift to keep pace with historical climatic conditions, is associated with more negative occupancy trends over time (Fig. 2c). Surprisingly, in the warmest Mediterranean region, the opposite pattern is observed:

species that have experienced a temperature increase throughout their range tend to increase (Fig. 2c), suggesting that climate change elicits a bonus instead of a debt in this area. Regarding precipitation, plant species occupancy trends are negatively related with precipitation SCI trends in the Alpine and Atlantic regions with the highest rainfall, while this relationship tends to be positive in the drier Boreal, Mediterranean and Continental regions (Fig. 2c).

These contrasting consequences of lagging behind climate change come in addition to expected effects of the historical climatic niche (Fig. 2b), that however vary in space. At the continental scale, species from rainy and warm historical niches exhibit higher occupancy trends. However, within biogeographic regions, the advantage of warm historical niches was observed only in cooler parts of Europe (Boreal, Atlantic and Alpine regions, Fig. S8).

Similarly, the benefit of rainy niches is seen in Boreal and Continental regions only (Fig. S8),

where rainfall has increased the most (Fig. S5), while species with dry niches seem to be favored in the Alpine region with decreasing precipitation (Fig. S8).

The overall consequences of climate change, combining both temperature and precipitation variables, are a climatic debt in Alpine, Atlantic and Boreal regions (Fig. 3a), but a climatic bonus in the Mediterranean (Fig. 3a), driven mostly by temperature changes but with a significant contribution of precipitation locally (Fig. 3b,c). In the Continental region the observed lags in range shifts did not have any overall significant effect on plant persistence (Fig. 3a), although when combining all the non-significant effects, lags in range shifts tended to benefit plants there, similarly to the Mediterranean region (Fig. S9). The examination of the relative contributions of temperature and precipitation SCI trends to the climatic debt/bonus shows that temperature is generally the major driver (Fig 3c), except in the Atlantic region, where precipitation shifts are the only significant driver of estimated climatic costs.

Importantly, SCI and occupancy trends are not significantly related during 1951-1990 (Fig. S6), when climate was relatively stable.

Finally, plant occupancy trends are also expectedly related to other drivers beyond climate change. At a continental scale, nitrophily and urban affinity are significant correlates of plant occupancy trends (Fig. 2b). We also find a negative but non-significant effect of pollinator dependency on occupancy trends (Fig. 2b). A majority of the remaining variables, such as most habitat affinities except urban affinity and moisture preferences, have contrasting effects on occupancy trends across biogeographic regions (Fig. S8).

Discussion

Our study of the consequences of climate change for plant species first confirms a spatial lag in species responses to climate change, evidenced by an increase in both temperature and precipitation SCIs (Fig. 1) suggesting that species are not moving fast enough to track their historical climatic conditions (i.e. to keep constant SCIs). Those SCI trends are phylogenetically structured, which could be explained by the already known climatic niche conservatism in plants (Prinzing [START_REF] Andreas | The niche of higher plants: evidence for phylogenetic conservatism[END_REF][START_REF] Preston | Adaptation to seasonality and the winter freeze[END_REF][START_REF] Hawkins | Community phylogenetics at the biogeographical scale: cold tolerance, niche conservatism and the structure of North American forests[END_REF][START_REF] Liu | Strong phylogenetic signals and phylogenetic niche conservatism in ecophysiological traits across divergent lineages of Magnoliaceae[END_REF] but also by the phylogenetic structure in the ability of plant species to track their optimum spatially via colonization [START_REF] Baeten | Disentangling dispersal from phylogeny in the colonization capacity of forest understorey plants[END_REF]. Furthermore, we also find that temperature and precipitation exhibit contrasting temporal trends in Europe, inducing uncorrelated SCI changes and possibly leading to a spatial trade-off for European plants between tracking precipitation and temperature historical conditions, as has been shown along the elevation gradient [START_REF] Crimmins | Changes in Climatic Water Balance Drive Downhill Shifts in Plant Species' Optimum Elevations[END_REF].

Analyses within biogeographic regions reveal that the lag in species response to temperature change translated into a climatic debt in the North but a surprising climatic "bonus" in the South, resulting in no overall significant signal for a climatic debt at the European scale.

However, differences between northern and southern Europe are consistent with a climatic debt at a continental scale: species that track climate change spatially, i.e. with SCI trends close to zero, have positive occupancy trends on their leading edge (northern margins) but negative trends on their trailing edge (southern margins). In addition, among-region differences also match regional patterns of correlations between occupancy trends and historical climatic niche, suggesting that climate could have distinctive effects among regions.

In Boreal and Alpine regions, the effect of temperature SCI trends comes in addition to an effect of the historical climatic niche, with larger occupancy trends for species from historically warmer area. This pattern is consistent with previous results on French plants [START_REF] Martin | Short-term climateinduced change in French plant communities[END_REF], and its significance only in the cooler parts of Europe (Fig. S8) confirms the well-known stronger effect of climate warming at higher latitudes [START_REF] Parmesan | Influences of species, latitudes and methodologies on estimates of phenological response to global warming[END_REF]. Occupancy trends correlate significantly both with temperature SCI trends and temperature of the historical niche in cooler parts of Europe, suggesting that climate is an important driver of species persistence in these areas. In contrast, the surprising climatic bonus in the Mediterranean region is consistent with the absence of correlation between the historical climatic niche and occupancy trends there, suggesting that climate is currently not a strong driver of plant occupancy trends in this area, as previously shown for colonization patterns [START_REF] Normand | Postglacial migration supplements climate in determining plant species ranges in Europe[END_REF]. This unexpected climatic bonus, which is generally overlooked, could be caused by changes in competitive interactions, an important driver of species responses to climate change [START_REF] Alexander | Novel competitors shape species' responses to climate change[END_REF], although we cannot exclude a role of other types of interspecific interactions, such as facilitation or herbivory [START_REF] Descombes | Novel trophic interactions under climate change promote alpine plant coexistence[END_REF]. Plants with limited or no northward shift (i.e. plants with an increasing temperature SCI) may benefit from competitive release associated with the range shift of more mobile species, without being in competition with novel competitors from southern regions, because of the numerous geographic barriers limiting plant colonization in the Mediterranean region [START_REF] Normand | Postglacial migration supplements climate in determining plant species ranges in Europe[END_REF]). These apparent benefits may however disappear when focusing on the whole Mediterranean region: here we ignored the southern margin of many Mediterranean plants, located in Northern Africa, a region with few plant records. Moreover, this climatic bonus is likely to be reversed by sustained climate change on the longer term, when climatic conditions exceed the climatic tolerance of species.

In addition to these effects of temperature, the climatic debt can also be driven by changes in precipitation, albeit to a lesser extent. When we combine the effects of temperature and precipitation SCI trends, we show that the inability of plant species to track their historical climatic conditions has been costly in the Alpine, Atlantic and Boreal regions, but beneficial in the Mediterranean region. These patterns substantiate further the notion of climatic debt in the former areas, and confirm the climatic bonus in the Mediterranean, although lags behind climate are most often interpreted as a climatic debt there [START_REF] Bertrand | Ecological constraints increase the climatic debt in forests[END_REF]. The effects of lagging behind changing precipitation are variable however. In relatively dry biogeographic regions, a decrease in the annual precipitation SCI of plant species over the past decades is associated with negative, or less positive, occupancy trends, which suggests that climate change causes water-deficit stress with detrimental consequences for plant population dynamics, a well-known phenomenon [START_REF] Breshears | Regional vegetation die-off in response to global-change-type drought[END_REF][START_REF] Allen | A global overview of drought and heat-induced tree mortality reveals emerging climate change risks for forests[END_REF][START_REF] Zhao | Drought-Induced Reduction in Global Terrestrial Net Primary Production from 2000 Through 2009[END_REF]. This applies to the Mediterranean region, in which precipitation shifts can be as important a driver as temperature changes, although it is widely overlooked in climatic debt assessments.

In contrast, in relatively wet areas, plant occupancy trends seem to be hindered by an increase in annual precipitation SCI, which suggests water-excess stress, via e.g. waterlogging.

Such consequence of climate change, via an increase in precipitation, is less documented but has been shown to drive downhill shifts in plant species elevation against temperature changes in mountain areas [START_REF] Crimmins | Changes in Climatic Water Balance Drive Downhill Shifts in Plant Species' Optimum Elevations[END_REF]. Their general contribution to the climatic debt relative to temperature is however moderate, except in the Atlantic region.

Beyond the effects of climate change, our results also strongly suggest that nitrogen deposition and urbanization are important disturbances for plants [START_REF] Aronson | A global analysis of the impacts of urbanization on bird and plant diversity reveals key anthropogenic drivers[END_REF][START_REF] Bobbink | Global assessment of nitrogen deposition effects on terrestrial plant diversity: a synthesis[END_REF]. However, while nitrogen deposition is sometimes cited as the first driver of changes in plant species composition [START_REF] Bobbink | Global assessment of nitrogen deposition effects on terrestrial plant diversity: a synthesis[END_REF], our results challenge this statement: by assessing response traits simultaneously, we find stronger links of occupancy trends with historical climatic niche or urban affinity than with nitrophily. Hence, our results are consistent with the recent acceleration of climate change in Europe and suggest that climate warming has caught up with urbanization and nitrogen deposition to become an important driver of plant persistence. We thus provide further evidence that biodiversity is often affected by multiple global change drivers rather than by single threats [START_REF] Brooks | glmmTMB Balances Speed and Flexibility Among Packages for Zeroinflated Generalized Linear Mixed Modeling[END_REF]. Consistent with previous results and with the pollinator decline [START_REF] Biesmeijer | Parallel Declines in Pollinators and Insect-Pollinated Plants in Britain and the Netherlands[END_REF], we find a negative effect of pollinator dependency on occupancy trends, but the latter is non-significant. This lack of signal for an effect of pollinator loss on plant occupancy trends may be attributable to contrasting plant trends depending on the group of pollinators [START_REF] Biesmeijer | Parallel Declines in Pollinators and Insect-Pollinated Plants in Britain and the Netherlands[END_REF]).

Here we show that plants are under multiple pressures from global change, and that plant occupancy trends exhibit a strong phylogenetic signal, which entails a risk of important evolutionary history losses associated with the forecasted extinctions. In particular, in some regions plant persistence is already affected by climate change and the resulting climatic debt/bonus, while these climate-related costs/benefits are often considered long-term. The climatic debt/bonus evident here is an integrative measure of all ecological and evolutionary costs/benefits associated with climate change, which we are not able to partition. For example, the costs we observe in Northern Europe could be due to insufficient adaptive response to buffer a spatial lag, to the arrival of novel competitors [START_REF] Alexander | Novel competitors shape species' responses to climate change[END_REF], and/or to the demographic cost of an ongoing adaptive response [START_REF] Lynch | Evolution and extinction in response to environmental change[END_REF] buffering the spatial lag. As plant adaptation to climate change opens the door to a possible evolutionary rescue for species that track their climatic optimum poorly in space [START_REF] Gonzalez | Evolutionary rescue: an emerging focus at the intersection between ecology and evolution[END_REF], assessing the contribution of ecological and evolutionary mechanisms of the climatic debt or bonus is a

remaining key challenge to predict future effects of climate change on plants.

Although this study tackled two dimensions of the species climatic niches simultaneously, our estimation of the climatic debt faces some limitations. First, we examined average climatic conditions only, thereby probably underestimating the cost of climate change, which also includes changes in variability, such as in the frequency of extreme events. Moreover, other environmental conditions beyond climate have changed during the last decades, such that their effects are difficult to unravel. Our climatic debt/bonus measure is thus integrative, likely including costs related to interactions between climate change and other drivers. For example landscape alteration can hinder species tracking of their historical climatic conditions [START_REF] Bertrand | Ecological constraints increase the climatic debt in forests[END_REF][START_REF] Gaüzère | Where do they go? The effects of topography and habitat diversity on reducing climatic debt in birds[END_REF]. We used a correlative approach on the basis of species traits, while taking phylogeny into account: the remaining phylogenetic signal in plant occupancy trends points to a likely omission either of some drivers or of synergistic effects among drivers.

These could be investigated through local studies.

Finally, we present the first overview of plant occupancy trends at continental scale. This was made possible by the use of opportunistic data, which are often the only data source to obtain long time-series at large spatial extent [START_REF] Biesmeijer | Parallel Declines in Pollinators and Insect-Pollinated Plants in Britain and the Netherlands[END_REF][START_REF] Bartomeus | Historical collections as a tool for assessing the global pollination crisis[END_REF], together with statistical methods aiming to correct the potential biases associated with those data [START_REF] Isaac | Statistics for citizen science: extracting signals of change from noisy ecological data[END_REF]. The fact that we find strong positive trends for invasive species suggests that trends estimated from GBIF data provide an accurate picture of actual changes in species occupancy. However, finding independent datasets and methods that allow turning the clock back and studying past effects of global change on biodiversity is a major challenge to confirm our results and anticipate future threats for biodiversity.

Figure 1 :

 1 Figure 1: Linear trends in species climatic indices (SCIs) for all species (n=4,120). The first row shows the trends in temperature SCIs and the second row the trends in precipitation SCIs over time. The WorldClim abbreviation for bioclimatic variables is indicated in parentheses atop each panel. The vertical red line indicates zero (no change across years) while the vertical dashed blue lines show the average values across the 4,120 species. Filled bars represent the count of species with significant trends (p-value<0.05) whereas open bars represent the count of species with non-significant trend (p-value>0.05).Occupancy trends and their drivers

  phylogenetic signal both before (Pagel's λ = 0.62, p-value<0.01; n=2,785) and after (Pagel's λ = 0.45, CI95% = [0.32,0.57]; n=2,013) considering species traits.

Figure 2 :

 2 Figure 2: Occupancy trends of native and naturalized European plant species and their correlates. (a) Histogram of occupancy trends, on a logit scale y -1 . The red vertical line indicates zero while the vertical dashed yellow line shows the average value (n=4,120). Filled bars represent the count of species with significant trends (p-value<0.05) whereas open bars represent the count of species with non-significant trends (p-value>0.05). (b) The three left panels represent the estimates (±CI95%) from phylogenetic regressions (PGLS, n=2,013) and linear mixed-effect model (LME, n= 2,080, see Methods) explaining occupancy trends with temporal trends in species climatic indices (SCIs) and other species traits. The two right panels show predicted averaged occupancy trends (±CI95%) for annual species with complete affinity for each habitat (habitat affinity score = 1 and lifespan = annual), and for lifespan categories, predicted at the average of all other variables. Red symbols represent significant effects while black symbols represent non-significant correlations. (c) Estimates (±CI95%) from PGLS and LME of the effect of standardized temperature and precipitation SCI trends on occupancy trends, for each biogeographic region.The analysis of the correlates of species occupancy trends reveals that plant species pay a

Figure 3 :

 3 Figure 3: Climatic debt/bonus in Europe and its climatic drivers. (a) Climatic debt/bonus averaged over all species over the last 65 years. The gradient from white to red indicates a climatic debt (cost of climate change in terms of species occupancy), while the gradient from white to blue indicates a climatic bonus (benefits of climate change in terms of species occupancy); white represents no cost of range shift lags on average for plants. Relative contribution of trends in species climatic indices (SCIs) related to (b) temperature and (c) precipitation to the climatic debt/bonus, in percentage. Black regions are biogeographic regions with too few data. The maps were generated using only predictions for effects of SCI trends that are significant in both the phylogenetic regression and the linear mixed-effect model averaged over these two models.

  𝛽 1 × ∆𝑆𝐶𝐼 𝑏𝑖𝑜1 𝑠 + 𝛽 2 × ∆𝑆𝐶𝐼 𝑏𝑖𝑜12 𝑠 + 𝛽 3 × 𝐵𝐶𝐼 𝑏𝑖𝑜1 𝑠 + 𝛽 4 × 𝐵𝐶𝐼 𝑏𝑖𝑜12 𝑠 + 𝛽 5 × 𝑀𝐸 𝑠 + 𝛽 6 × 𝑁𝐸 𝑠 + 𝛽 7 × 𝑃𝑜𝑙𝑙 𝑠 + 𝜃 𝑗 + ∑ 𝛽 ℎ𝑎𝑏𝑖𝑡𝑎𝑡 𝑖 × 𝐻𝑎𝑏𝑖𝑡𝑎𝑡 𝑖

	ℎ 𝑖=1

Acknowledgments

We thank the GBIF institution and GBIF data providers for building the database, making it publicly available and populating it. Statistical computation was partly performed at the HPCaVe center of Sorbonne Université. This project was funded by the Ministère de la Transition Écologique et Solidaire as part of the project "What is the sensitivity of pollinators to global warming in France" (convention n° SJ 3-17) led by Christophe Daugeron and Colin Fontaine, by the Institut de la Transition Écologique, Sorbonne Université, as part of the project Yapludsaison and by a grant from the Program Hubert Curien Alliance 2020.

Data availability statement

R scrips used to perform the analysis are available here: https://github.com/fduchenne/European-plants-lagging-behind-climate-change-. All data supporting the analysis, excepting row data from the GBIF, can be downloaded here: http://doi.org/10.5281/zenodo.4550500.

Author Contributions

FD designed the study, extracted the data and performed all the statistical analyses. GM and FD extracted and compiled species traits. FD and EP wrote the paper with contributions from GM.

Competing interests

The authors declare no competing interests.


Nathalie M Grob

Roger Schibli

Martin Behe

Ibai E Valverde

Thomas L Mindt

,5-disubstituted 1,2,3-triazoles as amide bond isosteres yield novel tumor-targeting minigastrin analogs

Keywords: 1, 2, 3-Triazoles, Peptidomimetics, Structure-Activity Relationships, Radiopharmaceuticals, Tumor Targeting, Cancer

1,5-disubstituted 1,2,3-triazoles (1,5-Tz) are considered bioisosteres of cis-amide bonds. However, their use for enhancing the pharmacological properties of peptides or proteins is not yet well established. Aiming to illustrate their utility, we chose the peptide conjugate [Nle 15 ]MG11 (DOTA-DGlu-Ala-Tyr-Gly-Trp-Nle-Asp-Phe-NH2) as a model compound since it is known that the cholecystokinin-2 receptor (CCK2R) is able to accommodate turn conformations. Analogs of [Nle 15 ]MG11 incorporating 1,5-Tz in the backbone were synthesized, radiolabeled with lutetium-177, and their pharmacological properties (cell internalization, receptor binding affinity and specificity, plasma stability, and biodistribution)

were evaluated and compared with [Nle 15 ]MG11 as well as their previously reported analogs bearing 1,4-disubstituted 1,2,3-triazoles. Our investigations led to the discovery of novel triazole-modified analogs of [Nle [START_REF] Rostovtsev | A Stepwise Huisgen Cycloaddition Process: Copper(I)-Catalyzed Regioselective Ligation of Azides and Terminal Alkynes[END_REF] ]MG11 with nanomolar CCK2R-binding affinity and twofold increased tumor uptake. This study illustrates that substitution of amides by 1,5-disubstituted 1,2,3-triazoles is an effective strategy to enhance the pharmacological properties of biologically active peptides.

LETTER

Five-membered N-heterocycles have attracted attention as amide bond mimics for the last twenty years. Due to the similar size, polarity, planarity, and their ability to create hydrogen bonds, imidazoles, oxa(dia)zoles, tetrazoles, and triazoles, among others, have been used to replace amide bonds in order to improve the pharmacological characteristics of biologically active compounds. [START_REF] Kumari | Amide Bond Bioisosteres: Strategies, Synthesis, and Successes[END_REF] These replacements are performed to overcome the inherently low proteolytic stability of amides or to probe binding conformations since heterocycles can mimic conformationally locked amide bonds. [START_REF] Meanwell | Synopsis of Some Recent Tactical Application of Bioisosteres in Drug Design[END_REF][START_REF] Sun | Applications of amide isosteres in medicinal chemistry[END_REF][START_REF] Bonandi | The 1,2,3-triazole ring as a bioisostere in medicinal chemistry[END_REF][START_REF] Davis | Synthesis of sansalvamide A peptidomimetics: triazole, oxazole, thiazole, and pseudoproline containing compounds[END_REF][START_REF] Rani | CuAAC-ensembled 1,2,3-triazole-linked isosteres as pharmacophores in drug discovery: review[END_REF] Despite the large number of examples illustrating the potential of replacing amide bonds with small heterocycles in drug-like compounds with low molecular weight, [START_REF] Davis | Synthesis of sansalvamide A peptidomimetics: triazole, oxazole, thiazole, and pseudoproline containing compounds[END_REF][START_REF] Rani | CuAAC-ensembled 1,2,3-triazole-linked isosteres as pharmacophores in drug discovery: review[END_REF] the use of this approach in peptides is not very common yet.

Disubstituted tetrazoles as well as 1,2,4-and 1,2,3-triazoles are among the heterocycles that allow a precise replacement of an amide bond (Figure 1). [START_REF] Yu | Synthesis and chemical properties of tetrazole peptide analogs[END_REF][START_REF] Hitotsuyanagi | A cis amide bond surrogate incorporating 1,2,4-triazole[END_REF][START_REF] Valverde | 1,2,3-Triazoles as Amide-bond Surrogates in Peptidomimetics[END_REF] Triazole synthesis can be performed in solution or on solid support, the latter making their use particularly appealing for peptide chemists. [START_REF] Boeglin | Solution and solidsupported synthesis of 3,4,5-trisubstituted 1,2,4-triazole-based peptidomimetics[END_REF][START_REF] Tischler | Braces for the Peptide Backbone: Insights into Structure-Activity Relationships of Protease Inhibitor Mimics with Locked Amide Conformations[END_REF][START_REF] Aucagne | Towards the Simplification of Protein Synthesis: Iterative Solid-Supported Ligations with Concomitant Purifications[END_REF] In this context, 1,2,3-triazoles (Tz) are the only heterocycles that, while having a similar dipolar moment, are capable of mimicking both trans-and cis-amide bonds, depending on their 1,4-or 1,5-substitution pattern, respectively, while being synthetically accessible (Figure 2). [START_REF] Pedersen | 2,3-Triazoles in Peptidomimetic Chemistry[END_REF] 1,2,3-triazoles have been seldomly used as peptidomimetics until the discovery of the copper(I)-catalyzed azide-alkyne cycloaddition (CuAAC), which selectively gives access to 1,4disubstituted 1,2,3-Triazoles (1,4-Tz). [START_REF] Tornøe | Peptidotriazoles on Solid Phase: [1,2,3]-Triazoles by Regiospecific Copper(I)-Catalyzed 1,3-Dipolar Cycloadditions of Terminal Alkynes to Azides[END_REF][START_REF] Rostovtsev | A Stepwise Huisgen Cycloaddition Process: Copper(I)-Catalyzed Regioselective Ligation of Azides and Terminal Alkynes[END_REF] Since then, 1,4-Tz have been quickly adopted as bioisosteres of amide bonds in peptidomimetics to modify the pharmacological properties of peptides or proteins. [START_REF] Galibert | Substrate-derived triazolo-and azapeptides as inhibitors of cathepsins K and S[END_REF][START_REF] Horne | Probing the bioactive conformation of an archetypal natural product HDAC inhibitor with conformationally homogeneous triazole-modified cyclic tetrapeptides[END_REF][START_REF] Decourt | A synthetic kisspeptin analog that triggers ovulation and advances puberty[END_REF][START_REF] Valverde | Synthesis of a Biologically Active Triazole-Containing Analogue of Cystatin A Through Successive Peptidomimetic Alkyne-Azide Ligations[END_REF][START_REF] Haj Salah | How are 1,2,3-triazoles accommodated in helical secondary structures?[END_REF] The systematic replacement of amide bonds by 1,4-Tz (termed a triazole scan) has been used as a tool to probe the contribution of the backbone to the stability and the activity of biologically relevant peptides similar to an N-methyl amide scan. [START_REF] Rečnik | 1,4-Disubstituted 1,2,3-Triazoles as Amide Bond Surrogates for the Stabilisation of Linear Peptides with Biological Activity[END_REF] So far, this approach has been applied to enkephalin, [START_REF] Proteau-Gagné | Systematic replacement of amides by 1,4-disubstituted[1,2,3]triazoles in Leu-enkephalin and the impact on the delta opioid receptor activity[END_REF] angiotensin, [START_REF] Vrettos | Single peptide backbone surrogate mutations to regulate angiotensin GPCR subtype selectivity[END_REF] bombesin, [START_REF] Valverde | 1,2,3-Triazoles as Amide Bond Mimics: Triazole Scan Yields Protease-Resistant Peptidomimetics for Tumor Targeting[END_REF][START_REF] Valverde | Probing the Backbone Function of Tumor Targeting Peptides by an Amide-to-Triazole Substitution Strategy[END_REF][START_REF] Valverde | Radiolabeled antagonistic bombesin peptidomimetics for tumor targeting[END_REF] neurotensin, [START_REF] Mascarin | 1,2,3-Triazole Stabilized Neurotensin-Based Radiopeptidomimetics for Improved Tumor Targeting[END_REF][START_REF] Mascarin | Radiolabeled analogs of neurotensin (8-13) containing multiple 1,2,3-triazoles as stable amide bond mimics in the backbone[END_REF] and minigastrin. [START_REF] Grob | Triazolo-Peptidomimetics: Novel Radiolabeled Minigastrin Analogs for Improved Tumor Targeting[END_REF][START_REF] Grob | Design of Radiolabeled Analogs of Minigastrin by Multiple Amide-to-Triazole Substitutions[END_REF] Examples of peptidomimetics based on 1,5-disubstituted 1,2,3triazoles (1,5-Tz) are scarce. 1,5-Tz, obtained by the ruthenium-catalyzed azide-alkyne cycloaddition (RuAAC), are useful to probe the bioactive conformation of peptides, since they can mimic constrained structures such as b-turns or hairpins, [START_REF] Horne | Probing the bioactive conformation of an archetypal natural product HDAC inhibitor with conformationally homogeneous triazole-modified cyclic tetrapeptides[END_REF][START_REF] Tietze | Development of a Functional cis-Prolyl Bond Biomimetic and Mechanistic Implications for Nickel Superoxide Dismutase[END_REF][START_REF] Pokorski | Introduction of a triazole amino acid into a peptoid oligomer induces turn formation in aqueous solution[END_REF] and their use as valuable amide bioisosteres in peptide backbones has been demonstrated both theoretically [START_REF] Marion | Amber-Compatible Parametrization Procedure for Peptide-like Compounds: Application to 1,4-and 1,5-Substituted Triazole-Based Peptidomimetics[END_REF][START_REF] Kracker | 1,5-Disubstituted 1,2,3-Triazole-Containing Peptidotriazolamers: Design Principles for a Class of Versatile Peptidomimetics[END_REF] and experimentally. [START_REF] Tischler | Braces for the Peptide Backbone: Insights into Structure-Activity Relationships of Protease Inhibitor Mimics with Locked Amide Conformations[END_REF][START_REF] Horne | Probing the bioactive conformation of an archetypal natural product HDAC inhibitor with conformationally homogeneous triazole-modified cyclic tetrapeptides[END_REF][START_REF] Lai | Cyclic tetrapeptide HDAC inhibitors as potential therapeutics for spinal muscular atrophy: Screening with iPSC-derived neuronal cells[END_REF] We have recently reported the triazole scan of [Nle 15 ]MG11 (compound 1, Figure 3) labeled with lutetium-177 ( 177 Lu). [START_REF] Grob | Triazolo-Peptidomimetics: Novel Radiolabeled Minigastrin Analogs for Improved Tumor Targeting[END_REF][START_REF] Grob | Design of Radiolabeled Analogs of Minigastrin by Multiple Amide-to-Triazole Substitutions[END_REF] This conjugate consists of a peptide with high affinity for the CCK2R coupled to the macrocyclic chelator DOTA to enable labeling with radiometals and shows potential in nuclear medicine for targeting a variety of tumors such as medullary thyroid carcinoma and small cell lung cancer. We have previously shown that single or multiple 1,4-Tz can be incorporated between residues DGlu 10 to Trp 14 of compound 1 while, in the majority of the examples, the nanomolar affinity of the resulting peptidomimetics toward the target receptor is preserved and an increased proteolytic stability was observed. In addition, several linear and cyclic peptides and peptidomimetics mimicking b-turns have shown high affinities toward the CCK2R. [START_REF] De Tullio | Recent advances in the chemistry of cholecystokinin receptor ligands (agonists and antagonists)[END_REF][START_REF] Blommaert | Cholecystokinin peptidomimetics as selective CCK-B antagonists: design, synthesis, and in vitro and in vivo biochemical properties[END_REF][START_REF] Blommaert | Structure-based design of new constrained cyclic agonists of the cholecystokinin CCK-B receptor[END_REF] Based on these reports, we set out to study the influence of amide-to-1,5-Tz switches on the CCK2R-binding properties of compound 1 in positions likely to be involved in b-turn conformation. Triazole-containing analogs of compound 1 were synthesized, radiolabeled with [ 177 Lu]Lu 3+ , and evaluated in vitro (cell internalization, receptor affinity and specificity, metabolic stability) and in vivo (biodistribution in mice bearing CCK2Rpositive xenografts).

With the aim of identifying 1,5-Tz-modified derivatives of compound 1 with retained tumortargeting properties, bonds between Gly 13 -Trp [START_REF] Tornøe | Peptidotriazoles on Solid Phase: [1,2,3]-Triazoles by Regiospecific Copper(I)-Catalyzed 1,3-Dipolar Cycloadditions of Terminal Alkynes to Azides[END_REF] (5), Tyr 12 -Gly 13 (6), Ala 11 -Tyr 12 (7) were substituted with a 1,5-Tz moiety. Modifications between Trp 14 -Nle [START_REF] Rostovtsev | A Stepwise Huisgen Cycloaddition Process: Copper(I)-Catalyzed Regioselective Ligation of Azides and Terminal Alkynes[END_REF] , and DGlu 10 -Ala 11 were discarded since these bonds are distant from the suspected location of the turn. [START_REF] Grob | Triazolo-Peptidomimetics: Novel Radiolabeled Minigastrin Analogs for Improved Tumor Targeting[END_REF][START_REF] Grob | Design of Radiolabeled Analogs of Minigastrin by Multiple Amide-to-Triazole Substitutions[END_REF] The heterocycles were inserted in the peptide sequence following the strategy exemplified for conjugate 7 in Scheme 1. Azido benzyl esters of amino acids were synthesized via diazo transfer using 1H-imidazole-1-sulfonyl azide. [START_REF] Goddard-Borger | An Efficient, Inexpensive, and Shelf-Stable Diazotransfer Reagent: Imidazole-1-sulfonyl Azide Hydrochloride[END_REF] Amino alkynes were obtained following published procedures (see Supporting Information (SI)). [START_REF] Valverde | 1,2,3-Triazoles as Amide Bond Mimics: Triazole Scan Yields Protease-Resistant Peptidomimetics for Tumor Targeting[END_REF][START_REF] Mascarin | 1,2,3-Triazole Stabilized Neurotensin-Based Radiopeptidomimetics for Improved Tumor Targeting[END_REF] Amino alkyne mimics of Fmoc-Ala-OH and Fmoc-Tyr(tBu)-OH were isolated in a L:D ratio of 66:34 and 90:10, respectively (see SI). Cp*RuCl(PPh3)2 and reaction temperatures of 60 °C40, 41 afforded the cycloaddition products in solution within 2-8 h (see SI). RuAAC was followed by deprotection of the carboxylic acid by hydrogenation using Pd/C as catalyst. The obtained pseudodipeptide building blocks were coupled to the growing peptide sequence on solid support using HATU as a coupling reagent. The remaining amino acids were coupled using standard Fmoc/tBu solid-phase peptide synthesis (SPPS) until completion of the DOTA-peptide conjugates (Scheme 1). After cleavage from solid support, deprotection by TFA, and separation of the diastereomers by RP-HPLC, the three DOTA-substituted peptidomimetics 5-7 were obtained in satisfactory yields (30-50%) and characterized by analytical RP-HPLC and HR-MS (see SI).

Compounds 5 to 7 were radiolabeled with [ 177 Lu]LuCl3 to obtain radioconjugates [ 177 Lu]Lu-5-7 (Table 1) in high radiochemical yields and purities of ≥95%, with specific molar activities between 25-50 MBq/nmol (not optimized). [START_REF] Grob | Triazolo-Peptidomimetics: Novel Radiolabeled Minigastrin Analogs for Improved Tumor Targeting[END_REF][START_REF] Grob | Design of Radiolabeled Analogs of Minigastrin by Multiple Amide-to-Triazole Substitutions[END_REF] Representative chromatograms of quality controls by g-HPLC can be found in the SI.

Receptor-mediated internalization of the radioconjugates was evaluated using A431 cells stably transfected with the CCK2R (A431-CCK2R cells) [START_REF] Aloj | In vitro and in vivo evaluation of 111In-DTPAGlu-G-CCK8 for cholecystokinin-B receptor imaging[END_REF] and is expressed as the percent of added radioactivity per 0.85 million cells (Table 1 and Figure 4). Blocking experiments (performed by addition of a 5000-fold excess of minigastrin) decreased the internalization of radiolabeled compounds to less than 1% in all cases, hence demonstrating specific receptormediated interaction (see SI). With the exception of compound [ 177 Lu]Lu-5, the other 1,5-Tz analogs of compound 1 were able to bind and internalize into receptor-positive A431-CCK2R cells. Conjugates [ 177 Lu]Lu-6 and [ 177 Lu]Lu-7 reached similar or higher cell internalizations than reference compound 1 (51 to 30% versus 32%, respectively), and the internalization rates and kinetics resembled closely the ones of their 1,4-substituted counterparts [ 177 Lu]Lu-3 and -4, respectively (54 and 30% respectively, Table 1, Figure 4).

The affinity (IC50, Table 1) of the novel triazolopeptides toward the CCK2R was determined by a competition-binding assay. Briefly, the established CCK2R-ligand [ 177 Lu]Lu-PP-F11N (DOTA-(DGlu)6-Ala-Tyr-Gly-Trp-Nle-Asp-Phe-NH2) [START_REF] Ritler | Elucidating the Structure-Activity Relationship of the Pentaglutamic Acid Sequence of Minigastrin with Cholecystokinin Receptor Subtype 2[END_REF] was added to A431-CCK2R cells and displaced by increasing concentrations of non-radioactive 175 Lu-labeled peptidomimetics 5-7 for 1 h. The cells were washed and lysed to determine the bound fraction of [ 177 Lu]Lu-PP-F11N as a function of the concentration of added competitors (see Table 1 conformation from the C-terminus until Gly [START_REF] Pedersen | 2,3-Triazoles in Peptidomimetic Chemistry[END_REF] . This observation confirms our previous hypothesis that Gly [START_REF] Pedersen | 2,3-Triazoles in Peptidomimetic Chemistry[END_REF] separates the C-terminal part of [Nle 15 ]MG11, which is buried in the receptor cavity, from the N-terminal part that is pointing toward the extracellular domain of the receptor. [START_REF] Grob | Triazolo-Peptidomimetics: Novel Radiolabeled Minigastrin Analogs for Improved Tumor Targeting[END_REF] In the same previous study, we observed an increased interaction between the CCK2R and 1,4-Tz compound [ 175 Lu]Lu-3 likely as the result of an additional cationp interaction between the triazole ring and Arg 356 of the receptor. [START_REF] Grob | Triazolo-Peptidomimetics: Novel Radiolabeled Minigastrin Analogs for Improved Tumor Targeting[END_REF] Similarly, the 1,5-Tz isomer [ 175 Lu]Lu-6 showed an improved affinity toward the CCK2R in comparison to the reference compound [ 175 Lu]Lu-1 (10 nM vs 15 nM). Even though the improvement of affinity of [ 175 Lu]Lu-6 was not as pronounced as in the case of the 1,4-Tz analog [ 175 Lu]Lu-3, the presence of an aromatic ring between Tyr 12 and Gly 13 indeed seems beneficial for the receptor interaction. Despite showing good cell internalization, [ 175 Lu]Lu-7 displayed a decreased affinity toward the CCK2R in comparison with its all-amide counterpart and its 1,4-Tz isomer (65 nM vs 15 nM and 21 nM, respectively). We conclude that the position between Ala 11 -Tyr [START_REF] Aucagne | Towards the Simplification of Protein Synthesis: Iterative Solid-Supported Ligations with Concomitant Purifications[END_REF] does not play a significant role in receptor binding due to its tolerance towards modifications.

Next, we studied the proteolytic degradation of compounds [ 177 Lu]Lu-5-7 in vitro to determine if the use of 1,5-Tz benefits the stability of this regulatory peptide. The metabolic half-lives (t1/2) of the conjugates were determined by incubation of the radiolabeled peptides in human blood plasma followed by quantification of their metabolites by g-HPLC over 24 h (Table 1).

The position of the amide-to-triazole switch utilizing 1,5-Tz had distinct effects on the metabolic stability of the resulting peptidomimetic conjugates in comparison to the previously studied 1,4-Tz-substituted minigastrins. [START_REF] Valverde | 1,2,3-Triazoles as Amide Bond Mimics: Triazole Scan Yields Protease-Resistant Peptidomimetics for Tumor Targeting[END_REF][START_REF] Valverde | Radiolabeled antagonistic bombesin peptidomimetics for tumor targeting[END_REF][START_REF] Mascarin | 1,2,3-Triazole Stabilized Neurotensin-Based Radiopeptidomimetics for Improved Tumor Targeting[END_REF] Peptidomimetics [ 177 Lu]Lu-5 and [ 177 Lu]Lu-7

showed an increased stability in comparison to reference compound [ 177 Lu]Lu-1, whereas the t1/2 was decreased for compound [ 177 Lu]Lu-6 (2.2 h vs. 3.9 h). It appears that, while the insertion of a Tz heterocycle can improve the affinity of the resulting peptidomimetic towards its receptor, it can also increase its susceptibility towards degradation by proteases. The degradation kinetics of 1,4-versus 1,5-Tz peptidomimetics of MG11 differed not only depending on the position of the heterocycle in the sequence but also on the substitution pattern of the triazole (Figure 5). Degradation of compound [ 177 Lu]Lu-5 was significantly slower than its 1,4-Tz counterpart [ 177 Lu]Lu-2 (t1/2 = 3.8 h vs. 91.1 h). The behavior was reversed in the case of conjugate [ 177 Lu]Lu-7, which was metabolized faster than 1,4-Tz compound [ 177 Lu]Lu-4 (t1/2 = 51.4 h vs. 6.7 h). Finally, degradation kinetics of compounds [ 177 Lu]Lu-3 and [ 177 Lu]Lu-6 were comparable (Figure 5). These results lead us to the conclusion that the contribution of differently substituted triazoles to the stability of a peptide may vary and does not follow a general pattern.

In vivo studies with mice bearing CCK2R-positive tumor xenografts were performed to investigate whether peptidomimetics [ 177 Lu]Lu-6 and 7 retained the biodistribution profile of reference compound [ 177 Lu]Lu-1 (see SI for details). Compound 5 was discarded for in vivo studies due to its reduced cell internalization and low affinity toward the CCK2R. Uptake of radioactivity in organs is expressed as % of injected dose per gram of organ or tissue (% ID/g). Female CD1 nu/nu mice were xenografted with A431-CCK2R cells and randomly assigned to groups of 4 mice. The results of the experiment are summarized in Figure 6 (see SI for complete data sets). At 4 h post-injection, all compounds showed typical biodistribution profiles of radiolabeled peptides with a fast clearance from the blood and unspecific uptake in the kidneys caused by renal excretion (Figure 6). Specific receptor-mediated uptake was demonstrated by co-injection with a 6000-fold excess of minigastrin (blocking), which led to a significant decrease of the uptake in the receptor-positive tumor and organs, e.g., the stomach (see SI). Compound [ 177 Lu]Lu-6 showed the highest tumor uptake of the investigated 1,5-Tz peptidomimetics (3.9 % ID/g). Similarly, a more than two-fold increased uptake in the tumor was also previously observed for the analogous 1,4-Tz isomer [ 177 Lu]Lu-3. This is consistent with their improved CCK2R affinity and cell internalization in comparison with reference compound [ 177 Lu]Lu-1 in vitro (Table 1). Tumor uptake of compound [ 177 Lu]Lu-7 was inferior to its 1,4-Tz counterpart [ 177 Lu]Lu-4 (1.7 % ID/g vs. 0.8 % ID/g,). The difference in the tumor uptake of compounds 7 and 4 could be due to a decreased receptor affinity of compound 7 in comparison with compounds 1 or 4 as well as to the superior proteolytic resistance of compound [ 177 Lu]Lu-4.

We herein report the use of 1,5-Tz heterocycles as metabolically stable mimics of cis-amide bonds in biologically active linear peptides. 1,5-Tz were placed in the backbone of the radiolabeled peptide MG11 in the vicinity of a proposed b-turn. The obtained peptidomimetics were compared side-by-side with the all-amide bond reference compound [Nle 15 ]MG11 as well as with previously studied 1,4-Tz analogs, which adopt transconformations for the probed amide bonds.

In comparison to the reference compound, some of the novel 1,5-Tz peptidomimetics showed improved biological properties such as enhanced plasma stability and affinity towards the CCKR2 in vitro or increased tumor uptake in vivo, thus providing another example for the successful application of 1,5-Tz as bioisosteres of amide bonds in bioactive peptides. Furthermore, the direct comparison of 1,4-Tz versus 1,5-Tz containing peptidomimetics revealed that a general prediction about the effect of the substitution pattern of the heterocycle on the biological properties of a peptide requires further investigations.

We encourage medicinal and peptide chemists who study 1,4-Tz-based peptidomimetics to consider also the use of 1,5-Tz as both heterocycles show high potential to enhance the pharmacological properties of bioactive peptides and thus, expand the peptidomimetic toolbox for drug discovery.

  and the SI). As anticipated from the internalization assays, compound [ 175 Lu]Lu-5 with the lowest cell internalization had the lowest affinity for its target (IC50 ~450 nM). Compounds [ 175 Lu]Lu-6 and [ 175 Lu]Lu-7 showed IC50 values of 11 nM and 65 nM, respectively, with only conjugate [ 175 Lu]Lu-6 approaching the affinity of the reference compound [ 175 Lu]Lu-1 (15 nM). Unlike compound [ 177 Lu]Lu-2, with a 1,4-Tz between Gly 13 -Trp 14 mimicking a trans-amide, compound [ 175 Lu]Lu-5, with a 1,5-Tz that mimics a cis-amide bond, did not show high nanomolar affinity toward the receptor. This result could indicate that the backbone of the peptide requires a linear (all-trans)

  ABBREVIATIONS

Figure 1 .

 1 Figure 1. Mimetics of amide bonds based on N-heterocycles.

Figure 2 .

 2 Figure 2. Comparison between trans-amide and cis-amide bonds and 1,4-(A) and 1,5-disubstituted 1,2,3-triazoles (B).

Figure 3 .

 3 Figure 3. Structure of [Nle 15 ]MG11 with numbering of the amino acid residues.

Figure 4 .

 4 Figure 4. Specific internalization kinetics of conjugates [ 177 Lu]Lu-1 to 7 in A431-CCK2R cells at 37 °C. Data of previously reported compounds [ 177 Lu]Lu-1-4 29 are included for comparison (n = 3 in triplicates).

9 [ 8 [ 6 [ 4 [ 1 [ 2 [ 7 a

 9864127 ]Lu-1 a [ 177 Lu]Lu-DOTA-DGlu-Ala-Tyr-Gly-Trp-Nle-Asp-Phe-177 Lu]Lu-2 a [ 177 Lu]Lu-DOTA-DGlu-Ala-Tyr-GlyΨ[1,4-Tz]Trp-Nle-Asp-Phe-177 Lu]Lu-3 a [ 177 Lu]Lu-DOTA-DGlu-Ala-TyrΨ[1,4-Tz]Gly-Trp-Nle-Asp-Phe-177 Lu]Lu-4 a [ 177 Lu]Lu-DOTA-DGlu-AlaΨ[1,4-Tz]Tyr-Gly-Trp-Nle-Asp-Phe-NH2 29.6 ±2.7 20.9 (17.0 -25.7) 51.177 Lu]Lu-5 [ 177 Lu]Lu-DOTA-DGlu-Ala-Tyr-GlyΨ[1,5-Tz]Trp-Nle-Asp-Phe-NH2 1.5 ± 0.9 447.7 (352.8 -568.1) 91.177 Lu]Lu-6 [ 177 Lu]Lu-DOTA-DGlu-Ala-TyrΨ[1,5-Tz]Gly-Trp-Nle-Asp-Phe-NH2 50.9 ± 4.1 10.5 (8.1 -13.7) 2.177 Lu]Lu-7 [ 177 Lu]Lu-DOTA-DGlu-AlaΨ[1,5-Tz]Tyr-Gly-Trp-Nle-Asp-Phe-NH2 30.3 ± 1.4 64.8 (51.7 -81.3) 6.Data of reference compound [ 177 Lu]Lu-1 and 1,4-Tz analogs [ 177 Lu]Lu-2-4 are reproduced for comparison. 29 b Competition experiments were performed with non-radioactive analogs of the compounds 1 to 7 labeled with 175 Lu.

Figure 5 .

 5 Figure 5. Degradation kinetics of 177 Lu-labeled peptide conjugates after incubation in human blood plasma for 24 h at 37 °C. Data points show mean ± SD. n = 2-3, for values of t½ see Table 1.

Figure 6 .

 6 Figure 6. Comparison of uptakes of compounds [ 177 Lu]Lu-1, -3, -4, -6, and -7 in selected organs of mice xenografted with A431-CCK2R positive cells at 4 h p.i. (n = 4). Data of compounds [ 177 Lu]Lu-1, -3, and -4 are reproduced for comparison. Superscripted pos indicates receptor-positive tissue. Data points show mean ± SD. n = 4, for values of all collected tissues and tumor-totissue ratios see the SI.

Table 1 .

 1

Summary of in vitro results of compounds [ 177 Lu]Lu-1-7.

Compound Sequence Internalization at 4 h (% a.

ACKNOWLEDGEMENTS

This work is part of the project "Pharmacoimagerie et Agents Théranostiques" supported by the Université de Bourgogne and Conseil Régional de Bourgogne through the Plan d'Action Régional pour l'Innovation (PARI), the Région Bourgogne Franche-Comté through the ANER program, and the European Union through the PO FEDER-FSE Bourgogne 2014/2020 programs. GDR CNRS "Agents d'Imagerie Moléculaire" 2037 is thanked for its interest in this research. We thank the "Plateforme d'Analyse Chimique et de Synthèse Moléculaire de l'Université de Bourgogne" (PACSMUB, http://www.wpcm.fr) for access to spectroscopy instrumentation. I.E.V. thanks Prof. Anthony Romieu and Dr. Adrien Normand for scientific discussions, Dr. Jerôme Bayardon for scientific discussions and technical support with chiral HPLC, and Dr. Quentin Bonnin and Marie-José Penouilh for HR-MS (Université de Bourgogne, Dijon, France). DOTA-tris(tBu)ester was a generous gift from Chematech (Dijon, France). A431-CCK2R cells were a kind gift of Dr. Luigi Aloj (University of Cambridge, United Kingdom). We thank Stefan Imobersteg (PSI, Villigen, Switzerland) for animal care.

FIGURES