Widespread decline in Central European plant diversity across six decades

Based on plant occurrence data covering all parts of Germany, we investigated changes in the distribution of 2136 plant species between 1960 and 2017. We analyzed 29 million occurrence records over an area of ~350,000 km2 on a 5 × 5 km grid using temporal and spatiotemporal models and accounting for sampling bias. Since the 1960s, more than 70% of investigated plant species showed declines in nationwide occurrence. Archaeophytes (species introduced before 1492) most strongly declined but also native plant species experienced severe declines. In contrast, neophytes (species introduced after 1492) increased in their nationwide occurrence but not homogeneously throughout the country. Our analysis suggests that the strongest declines in native species already happened in the 1960s–1980s, a time frame in which often few data exist. Increases in neophytic species were strongest in the 1990s and 2010s. Overall, the increase in neophytes did not compensate for the loss of other species, resulting in a decrease in mean grid cell species richness of −1.9% per decade. The decline in plant biodiversity is a widespread phenomenon occurring in different habitats and geographic regions. It is likely that this decline has major repercussions on ecosystem functioning and overall biodiversity, potentially with cascading effects across trophic levels. The approach used in this study is transferable to other large‐scale trend analyses using heterogeneous occurrence data.


| INTRODUC TI ON
Biodiversity loss is one of the greatest global challenges. Approximately 1 million species are threatened by extinction (Díaz et al., 2019). While the widespread decline in fauna (Dirzo et al., 2014) is discussed prominently in the scientific community and the general public, in recent times especially for terrestrial insects (e.g., Eisenhauer et al., 2019;Hallmann et al., 2017;Powney et al., 2015;van Klink et al., 2020), large-scale changes on the distribution of plants are less widely recognized. There are few examples of assessments of temporal trends in plant diversity over larger regions, such as whole countries (e.g., Finderup Nielsen et al., 2019;Rich & Woodruff, 1996;Walker & Preston, 2006). A deeper understanding of biodiversity change in plants is essential for predicting ecosystem-wide changes, including effects on ecosystem functioning (Hejda & de Bello, 2013) and the provision of ecosystem services relied on by humans (Guo et al., 2010).
A widespread problem in many datasets is limited spatial sampling, and often under-sampling of locations where change might have been the greatest. Moreover, many studies are based on short biodiversity time series that postdate time periods when important environmental changes have occurred, leading to underestimates of biodiversity change (Gonzalez et al., 2016). In addition, it has been noted that trends in species richness may also differ between native and non-native species (Cardinale et al., 2018). As pointed out by Sax and Gaines (2003), net gains or losses of native species in large study regions can be balanced or even overcompensated by the introduction of non-native species. For example, in a study on the flora of Denmark, Finderup Nielsen et al. (2019) reported an increased naturalization and subsequent spread in exotic species and increases in widespread native species, but reported declines in rare natives over the last 140 years. In total, this led to a net gain in species richness in the studied regions while species communities became more homogeneous across the study regions. A study by Winter et al. (2009) showed similar results on the European scale. A focus only on total species richness may lead to the conclusion that biodiversity in a region remains stable, whereas it is actually experiencing turnover in species composition. This turnover may drastically alter the floristic composition in a study region-and potentially affect ecosystem functioning and stability (Naeem et al., 2012).
Analysis of floristic changes at large scales is challenged by the lack of repeated surveys of large regions (Walker & Preston, 2006; but see Switzerland and the UK, Hintermann et al., 2000, Wood et al., 2017. Typically, species-level changes at large scales, such as national or global levels, are assessed via red lists as well as during structured monitoring projects. Data for the more common, non-threatened, or non-iconic taxa are mostly recorded through individual atlas projects that do not aim at repeated recording (but see Blockeel et al., 2014). Consequently, the different kinds of data on plant occurrences for large regions, for example, countries or continents may come with different protocols and methods for data collection or study focus. This often leads to heterogeneous data quality and structures, hampering the analysis of such data for temporal trends. Therefore, the integration of spatially and taxonomically comprehensive data across different sources of information-often including citizen science data-has become increasingly relevant for the assessment of biodiversity change across larger spatiotemporal scales Zipkin et al., 2019).
In Germany, the majority of knowledge on changes in plant biodiversity comes from plant community resurvey studies of (semi)permanent relevés. Conclusions are mixed, showing strong (Wesche et al., 2012), little (e.g., Jensen et al., 2012), or no significant changes in local species richness (e.g., Bruelheide, 1998;Diekmann et al., 2014;Litza & Diekmann, 2017), depending on the investigated habitats, landscapes, species groups, plot sizes, or the temporal extent of the study. Most studies report that declines in species richness and changes in community composition are strongest in agricultural landscape (e.g., Meyer et al., 2013;Meyer, Bergmeier, et al., 2015;. However, plot-level analyses may not allow for a straightforward extrapolation to larger scales due to several biases, for example, in the spatial (i.e., habitat) representativeness of the plots or local differences in disturbances or management (Cardinale et al., 2018;Sperle & Bruelheide, 2020). Ignoring the local, and potentially spatially biased, small-scale patterns may lead to erroneous conclusions on large-scale net changes in biodiversity (Cardinale et al., 2018;Gonzalez et al., 2016;Sax & Gaines, 2003). Moreover, to understand changes in largescale biodiversity, it is crucial to evaluate as many species as possible.
Germany, organized in 16 federal states, offers a good example to demonstrate how heterogeneous datasets on plant biodiversity can be combined and analyzed to detect changes in large-scale species distributions and, subsequently, species richness. Each federal state has carried out at least one floristic atlas survey (Bergmeier, 1992) followed by other survey projects, but protocols, time frames, and predefined species lists at least for the latter varied in and across states. Although the atlas surveys have not been systematically repeated, several local and often taxonomically less comprehensive surveys exist across Germany. The combination of these data has led to publication of the German Atlas of Vascular Plants (Netzwerk Phytodiversität Deutschland & Bundesamt für Naturschutz, 2013).
However, as the underlying data often show spatial and temporal inconsistencies in the sampling coverage of subregions within the complete study region, detection of changes in biodiversity over time can be very challenging (Hill, 2012;Pescott, Powney, et al., 2019). Relevé data on the plot level are available from research institutes, universities, online databases (e.g., German Vegetation Reference Jandt & Bruelheide, 2012;or veget web.de, Jansen et al., 2015), and private collections. These data mostly have accurate geographic information, are restricted in size (one to a few square meters) and often lack temporal replication and are thus often not suitable for large-scale temporal trend analyses on their own (Chytrý et al., 2014). Compiling and integrating different datasets (atlas, relevé and observational data from private observations, excursions, museum records, mobile apps, and from spatially referenced legacy collections) in a common analysis-and thus making use of the merits of each of these data types-may allow to quantify long-term changes of plant species distributions on large spatial and temporal scales, potentially also at a fine spatial grain. Meanwhile, modern statistical tools allow to incorporate different data types from different sources in robust analyses while accounting for their heterogeneity .
In the present study, we compiled an extensive dataset of the spatiotemporal occurrence of vascular plants in Germany. The dataset was collated from multiple sources on plant occurrence records and vegetation surveys in Germany, varying in taxonomic extent and sampling effort. After accounting for incomplete species recording across space and time, we (a) assess spatiotemporal changes in the occurrence of 2136 species on a 5 × 5 km grid cell basis; and (b) assess the balance between winners and losers. We (c) explore the temporal dynamics of these changes based on floristic status (natives vs. non-natives). Using spatiotemporal models, we (d) assess changes in mean grid cell species richness across the whole nation, accounting for spatiotemporal dependence in the data. Moreover, we (e) explore the spatial heterogeneity in the patterns of changes in grid cell species richness over the last six decades and identify the hotspots of biodiversity turnover.

| MATERIAL S AND ME THODS
We developed a workflow to harmonize the available data into a common format, taxonomically and spatially (Supporting Information; Figure S1). We accounted for the potential effects of imperfect detection, using the Frescalo algorithm (FREquency SCAling using Local Occupancy; Hill, 2012), which is provided in the R package "sparta" (August et al., 2015;v. 0.1.48). This algorithm has been specifically developed for repeated, large-scale surveys, such as atlases, and computes species occurrence probabilities (OPs) across periods of data availability on a defined grid size. Pescott, Humphrey, et al. (2019, p. 264ff) explain how slightly different widths of these periods should have no strong influence on the output of the Frescalo algorithm, given that these periods are selected with care (see Supporting Information as well as Hill, 2012or Pescott, Humphrey, et al., 2019 for some details on criteria). However, the current version (v. 0.1.48) of the "sparta" package does not account for possible temporal dependencies in the data. To this end, we further demonstrate the use of new approaches to account for temporal (and also spatiotemporal) autocorrelation in the data, thus enabling more reliable analyses of changes in species occurrences and richness across space and time.

| Data compilation and taxonomic harmonization
We compiled an extensive dataset of nearly 29 million occurrence records in Germany between 1960 and 2017 from 23 different data sources (Supporting Information; Table S1). The full dataset com- html)-is a compilation of occurrence records gathered from, inter alia, the Floristic Atlases of Western and Eastern Germany (Benkert et al., 1996;Haeupler & Schönfelder, 1989) and floristic data extracted from 74 mapping projects. We extended this dataset to the year 2017 by integrating data from more recent habitat mapping projects of federal states, vegetation relevés provided in two major German databases, GVRD (Jandt & Bruelheide, 2012) and vegetweb 2.0 (Jansen et al., 2015) and from universities and private collections. For all datasets, observations were georeferenced on a grid cell level (a quadrant of German ordinance maps, "MTBQ," approximately 5 × 5 km). In all spatially explicit analyses, central coordinates of each grid cell (in UTM, EPSG 4326) were used as sampling locations (Zuur & Ieno, 2017).
were aggregated to the species or, if necessary, to the aggregate level. For simplicity, we will refer to these as "species" in the following. Taxonomic harmonization was achieved using the R package "vegdata" (Jansen & Dengler, 2010). This resulted in data on the occurrence of 2976 vascular plant taxa, equaling 77% of all German vascular plant taxa (4305 taxa, including subspecies and varieties or 3868 taxa when raised to the taxonomic level used in this study). For analyzing trends, the dataset was binned into three periods (1960-1987, 1988-1996, and 1997-2017), each of them with similar number of total records and covering all 12,024 German grid cells. The temporal extent of these periods was determined by the need to find periods of similar coverage (spatial and taxonomically) of the whole nation (see Supporting Information; Technical details). We excluded species that were recorded in only one of the three periods or had fewer than 23 records in total (i.e., the 10% quantile of species frequency distribution data). Most of the excluded species were below this minimum threshold of records (308 species), and most of these were rated as either "not endangered" (91 species) or "very rare" (48 species) in the most recent German Red List of Endangered Plants (Metzing et al., 2018).
In general, most plant species listed in one of the threat categories of the German Red List (aside from species rated as "extinct or lost") were included in the present study (Table S2). Species were characterized according to their floristic status either as natives, archaeophytes (non-natives introduced before 1492) or as neophytes (non-natives introduced after 1492), using information available from the database BiolFlor (Kühn et al., 2004) and FloraWeb (Bundesamt für Naturschutz; http://www.flora web. de/). Species with an unknown floristic status were excluded (see species-level details in Supporting Information; Table S2). This left us with a total of 2234 species for analysis, equaling 58% of all vascular plant species known from Germany (harmonized to the taxonomic level of the present analysis).

| Correction for false absences
The majority of the data originated from grid-based occurrenceonly records (approximately 95%, Supporting Information; Table  S1). In cases where occurrence records do not originate from a project focusing on the complete floristic inventory or do not use complete checklists, false absences (i.e., not reporting a species that was present, but was either not detected or detected but not reported) may be an issue (Pescott et al., 2018). In addition, atlas projects, which aim at taxonomic completeness may not be finished in a federal state completely within one of the defined study periods, leading to taxonomic or spatial gaps in the data of a single study period. To correct for this so-called reporting bias, we used the Frescalo algorithm (Hill, 2012; see also Bijlsma, 2013;Blockeel et al., 2014;Fox et al., 2014;White et al., 2019) available in the R package "sparta" (August et al., 2015;v. 0.1.48). Briefly, the Frescalo algorithm calculates the OP of a species not detected or reported in a focal grid cell, based on the frequency of this species in the local neighborhood (here: 100 grid cells) of this cell while accounting for the ecological similarity of the neighborhood (Hill, 2012). Ecological similarity of the neighboring grid cells was calculated based on a set of 76 variables, comprising climatic, topographic, and edaphic measures. A detailed description of the specifications and considerations that are necessary when preparing data for an analysis using the Frescalo algorithm is given in Electronic appendix (Technical details; but also see Bijlsma, 2013;Hill, 2012). Correction for false absences increased the dataset to more than 41 million entries of OPs per species and grid cell between 1960 and 2017.

| Calculation of species-specific occurrence across space and time
Maps of the spatial distribution of occurrence estimates of a given species across the study region at a given period are not a direct outcome of the Frescalo algorithm as available in "sparta" (v. 0.1.48). However, the spatial distribution of the probability of a species being present at the focal grid cell of neighborhood in a certain period can be readily calculated from the available output using Equation 1 (cf. Bijlsma, 2013): Where OP jt is the occurrence probability of a species in the focal grid cell j of the respective neighborhood at time t; s jt is the sampling intensity (a measure of sampling completeness calculated by the Frescalo algorithm) for grid cell j at time step t, f j is the estimated frequency of the respective species in neighborhood j after rescaling, and tfact t is the time factor (the estimated relative frequency of the respective species; a parameter derived by the Frescalo algorithm) at time period t. OP jt was calculated for each species, separately. All variables are given in the Frescalo output file provided in R (for the respective R-Script see Supporting Information).
To account for uncertainty in the Frescalo estimates, calculations of species-specific occurrence probability OP jt were based on 1000 realizations of µ t (sampled from a species-specific normal distribution with mean µ t and σ t , cf. Equation 1). For each species and each realization, we calculated the nationwide occurrence of a species as the sum of all OPs of a species (SOP Spec ) across Germany for each period according to Equation 2: Where Spec represents species under consideration, j is the index of grid cell, t is the index of time period (i.e., 1 = 1960-1987; 2 = 1988-1996; 3 = 1997-2017), and OP Specjt is the occurrence probability of species in grid cell j at time t (Equation 1).

| Calculation of grid cell species richness
We summed up the OPs of all species within a grid cell as an estimate of species richness (SOP Grid ), while acknowledging that it is not species richness per se, since our analysis does not include the very rare species (see Sections 2.1, 2.2 and 4). Hence, our species richness values are most likely underestimating actual grid cell species richness. However, SOP Grid was found to be significantly correlated (r = 0.39, p < 0.001; Supporting Information; Figure S2) to species numbers in grid cells of the FlorKart dataset (cf. Section 2.1) that  identified to be well-sampled. Therefore, changes in SOP Grid can be interpreted as meaningfully representing relative changes between grid cells and time periods.
Where l is the index of the species and OP l t is the OP of this species l in a grid cell at time period t (as derived from Equation 1).

| Changes in species-specific occurrence
Changes in species-specific occurrences over time were evaluated at the species level, using a Bayesian log-linear mixed effects model and including a random walk of order 1 ("rw1") to account for temporal dependency. These specifications ensured that changes are bounded to −100% at the lower end, but not at the upper end.
Model 1: Estimation of changes in species-specific occupancy over time according to a random walk component of order 1. Here, µ t−1 is the estimated occurrence of the species for the preceding time period.
Of the total of 2234, 2206 species were found to have significant changes (for a definition of "significant," see below). A critical inspection showed that 102 species exhibited extreme trends (i.e., above or below the 95% quantile range of change). These extreme changes (2) were discussed in depth with taxon experts, and those considered to be unrealistic trends (70 species) were omitted from further analyses (cf. Supporting Information; Table S3).

| Changes in mean grid cell species richness
To analyze changes in nationwide mean grid cell species richness, we ran spatiotemporally explicit models, based on gamma distributions (residuals were clumped in time and space for the log-normal case). To account for temporal dependencies, we included a temporal correlation structure with a random walk of order 1. Moreover, we accounted for temporal autocorrelation in the spatial structure is termed "significant" hereafter, although within a Bayesian framework, the appropriateness of this word is debatable. Figure 1 shows the relative differences in species occurrences between the first  and the last (1997-2017) study period.

F I G U R E 1
Mean changes in occurrences of decreasing species were −28.3% (±19.8%, n = 1526) while changes in occurrences of the increasing species were +71.6% (±113.7%, n = 610). None of the studied species showed a decline of −100%, that is, went extinct in Germany. A detailed list of the winners and losers in occupancy for all three time steps is given in Table S3 and in an interactive plot (https://shiny.idiv. de/de25g eka/Winne rsLos ers/).
The species with strongest decline in mean OP (−99.8%) was Anagallis tenella L., an endangered native species of nutrient-poor fens and transition mires in Germany. In contrast, the neophyte Senecio inaequidens DC., showed the strongest increase (+696%; Figure 2).

| Changes in occurrences of floristic status groups
We found significant changes in mean occurrence among all grid cells in Germany for all floristic status groups (Table 1). Across all species, we found a total decline of −10.8% over the last six decades, with strongest losses occurring between the first  and the second (1988)(1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996) study period. While natives and archaeophytes decreased over the whole study period, neophytes showed an increase in mean occurrence. Native species and archaeophytes showed most change between the first and second periods. Among the floristic status groups, archaeophytes showed the strongest decrease across the whole study period (−21.6%). While losses were strongest between the first and second period, changes between the second (1988−1996) and third period  were marginally insignificant, that is, their 95% CIs overlapped with 0, their 90% CIs did not. There was a F I G U R E 2 Occurrence probability estimates for the three study periods on a 5 × 5 km grid. Top: the endangered native species Anagallis tenella L. Bottom: the neophyte Senecio inaequidens DC. Gray areas are outside of the range supported by the Frescalo algorithm TA B L E 1 Changes in mean occurrence as well as relative differences (%; incl lower and upper 95%CI) across the three study periods as predicted for the average species. Letters indicate significant differences in mean values based on the effect of the temporal component in Model 1 tendency of mean occurrence in archaeophytes to further decrease also in the last observation period. In contrast, neophytes showed increases in occurrence of +23.3% over the whole study period, with strongest increases between the second and third period. While the increase in neophytes compensated for the losses in natives and archaeophytes from the second to the third period, overall, the increase in neophytes did not level-off the decreases across the full study period (Table 1).

| Changes in mean grid cell species richness
We detected changes in the estimates of mean grid cell species richness in all floristic status groups as well as across all species ( Figure 3). While natives and archaeophytes showed consistent declines in mean grid cell species richness over the whole study period, neophyte species richness showed only slight increases between periods two and three (i.e., 1988-1996 vs. 1997-2017).
Losses were strongest for archaeophytes (−19.3%). As for changes in mean nationwide occurrence, increases in mean grid cell neophyte richness did not counterbalance the decreases in the other floristic status groups, causing a net decrease of −1.9% in mean grid cell species richness per decade across all species. Overall, the species richness trends were similar to mean species occurrence trends.

| Spatial patterns of species richness change
While the relative changes in grid cell species richness of archaeophytes and neophytes showed strong spatial heterogeneity, changes in native species richness were more uniform, with consistent declines in all grid cells (Figure 4; Supporting Information; Table S4). Archaeophytes were also consistently declining but the magnitude of the decline was more variable, reaching the largest declines out of all floristic status groups in some grid cells.
Neophyte changes were most spatially variable and included the regions of both decrease (especially in the northeast) and increase (especially in some southern regions). As for changes in occurrence, changes in species richness of archaeophytes and natives were strongest from the first to the second period, whereas for neophytes they were strongest between the second and third.
The spatial patterns of changes across all species closely followed those of natives (the most speciose floristic status group).
Therefore, we will omit maps on the spatial patterns across all species from the following. For species richness estimates and absolute changes, see Electronic appendix ( Figures S3 and S4 Native species richness also declined in the far east of Germany (Region 5) and the south (Regions 9 and 10). Regions of the strongest declines in archaeophytes were found in the southwestern parts of Germany (Region 6, but also 7 and 9) and in the northwest

| DISCUSS ION
Based on the collation of the largest databases on plant occurrence records in Germany to date, and correcting for reporting bias, we found significant declines in the German-wide occurrence of 71% of all investigated vascular plant species over the last six decades. (9) Alpine foothills, and (10) Lake Constance Region; Dotted lines demark the state boundaries. For interactive maps, see https://shiny.idiv. de/de25g eka/PCHM/ A much lower proportion (29%) of species showed increases in nationwide OPs. The increases in neophyte occurrence and species richness did not compensate for nationwide losses of other species, which led to a significant decrease of approximately −2% in mean species occurrence as well as mean grid cell species richness in Germany over the last six decades; temporal and spatiotemporal dynamics differed between floristic status groups. We provide evidence that the majority of plant species in Germany shows a decline in occurrence and that decreases in plant species diversity are widespread across most regions throughout Germany.
The present study overcomes some of the main critiques on large-scale studies of biodiversity change summarized by Cardinale et al. (2018). The spatial representation of our data is not restricted to certain plots or locations in Germany that may cause a spatial bias toward certain regions or habitat types. In fact, our data cover all 5 × 5 km German grid cells. Moreover, since our analysis does not treat datasets independently but as an amalgam of different sources, underlying biases in the potential drivers of certain datasets (e.g., some vegetation relevés may have originated from success control of restoration projects) are reduced by combining these data with, for example, grid mapping for atlases that do not include such biases. Moreover, our study differentiates between changes in natives vs. non-natives and we account for spatiotemporal dependencies, allowing for spatiotemporally explicit analyses.

| Species-specific changes in occurrence
With 55% of all vascular plant species in Germany (2136 out of 3868 taxa, data assigned to the harmonized taxonomy; the complete flora of Germany comprises a total of 4305 taxa, including subspecies and varieties; Metzing et al., 2018), our study covers a major part of the German flora. Most of the species omitted from this study occurred only in a single time period or were very rare, with less than 23 observations (cf . Table S2). Thus, our analysis included most of the rare, moderately common and common species, except some Oenothera, Taraxacum, and Rubus species groups, due to unstable taxonomical concepts. The study did not comprise very rare plant species, those that have gone extinct before 1987 or those that entered the German flora after 1996. A study by Lennon et al. (2004) demonstrated that, although rare species can constitute a major part of the species pool in a region, spatial structures of species richness are typically dominated by the more common species. Therefore, the detected temporal and spatiotemporal patterns can be assumed to be reliable estimates for change in species richness. Moreover, our approach is rather conservative in terms of species presence: a species that is absent from a grid cell in our original dataset, but detected (even only once) in the neighborhood of the respective grid cell causes the OP in the respective grid cell to be greater than zero, instead of it being rated as absent.
However, although we did everything to ameliorate biases due to differences in local recording effort in our data and accounted for spatiotemporal biases, we cannot rule out the possibility that our data still contain some artifacts. Indeed, we find some visible patterns of grid cell species richness that are well known to German botanic specialists. It is for example not clear whether the low numbers of archaeophytes in the federal state of Brandenburg (northeastern Germany in Figure 4) are due to systematic low sampling effort for archaeophytes in this area, or whether these low numbers reflect the reality. However, since our analysis focuses on relative changes among grid cells and time steps and includes a spatiotemporal component in the statistical models which correct for spatially structured systematic biases, such biases should not have large effects on the overall findings. They must, however, be kept in mind and should always be critically evaluated, especially if specific species-level changes are of interest. Nonetheless, our work represents a major advance in the investigation of German plant diversity, and complements existing German structured monitoring schemes that mainly focus on rare, threatened, or invasive species (Mitschke et al., 2005). However, we emphasize that close inspection of individual species trends and discussions with experts is crucial before making inferences. Moreover, the species-specific results provided in our analyses should be interpreted considering the spatial scale of the grid cell level. On the grid cell level, a species can only be rated as absent after the last individual in a grid cell has gone extinct (Chase et al., 2019). Thus, the results of our analyses refer to occurrence, and not abundance. For example, the species with the strongest decline of −99.9% (A. tenella), indicating a nearextinction, has been recognized as "threatened with extinction" (RL 1) in the German Red List of 1996 (Korneck et al., 1996) but has been changed to "endangered" (RL 2) in 2018 (Metzing et al., 2018).
While there are only three remnant occurrences of this plant species (shown in Figure 2) that have been decreasing in size in the last decades, populations within these remnants have stabilized due to nature protection measures (Metzing et al., 2018;Raabe et al., 2012). By contrast, the species with the most extreme increase, S.
inaequidens, has been recognized to be expanding since the early 1960s (Heger & Böhmer, 2005), mainly along the railway and road network. Our approach can identify such large-scale changes, but it cannot replace specialized, local-scale investigations such as population-based studies, for example, for red list assessments.

An increased number of endangered plant species in Germany has recently been reported in the German Red List of Endangered
Plants (Metzing et al., 2018). Our results report an even higher total number of species in decline. However, methods between the analysis approach demonstrated here and those of the red lists differ. Therefore, results cannot be directly compared. Many of species found to be in decline in our study are common species. This is congruent to the England vascular plant Red List (Stroh et al., 2014) which also used Frescalo and which led to the change in the threat level of many common species. Likewise, a study in north-east Germany reported that approximately 60% of the 355 studied species were declining and that moderately common species declined strongest . Similarly, in northwestern Germany, Bruelheide et al. (2020) reported significant declines for a large number of plant species, mainly herbs.

| Changes in occurrences of floristic status groups and their spatial patterns
The general loss in species richness across all species was dominated by the declines among native species, which is the most speciose group. We were able to show that the patterns of change do not only vary according to floristic status, but also in spatial and temporal patterns. An investigation of the causal relationships between large-scale measures and potential drivers was not in the scope of the present study. However, we can discuss the spatial patterns, especially the hotspots of change, with respect to knowledge from local, small-scale studies. For example, in the German coastline regions (also Region 1 in Figure 5), a number of coastal macrophytes, for example, × Calammophila baltica (Schrad.) Brand. and Leymus arenarius (L.) Hochst., two coastal grass species, were predicted to decline in climate envelope models for the German coastal regions due to climate warming (Metzing, 2010). Likewise, Kastler andMichaelis (1999) andEggert et al. (2006) reported that Zostera marina L. and Z. noltei Hornem., two submerse seagrasses, are declining due to increased sea temperatures, salinity, and eutrophication in the Wadden Sea. A study on the scale of the northeastern federal state of Mecklenburg-Western Pomerania by Jansen et al. (2020) found also other typical coastal species such as Salsola kali L. or Triglochin maritimum L. declining between 1980 and 2000, probably due to hampered coastal dynamics and reduced grazing of coastal grasslands. All mentioned species show declines in occurrence also in our analysis, indicating that the findings of the more local studies hold also on larger scales (Table S3).
The decline in the occurrence of archaeophytes shown in our study has also been reported from small-scale studies, and was often explained by an increase in land-use intensification (Baessler & Klotz, 2006;Comin & Poldini, 2009;Knapp et al., 2010;Leuschner et al., 2013;Meyer, Bergmeier, et al., 2015;Meyer, Wesch, et al., 2015). A study by Baessler and Klotz (2006) (Knapp et al., 2010; and human trade and transport (Rejmánek et al., 2005). In addition, several studies have shown that estuaries may especially be prone to the establishment of non-native species for various reasons (e.g., Wolff, 1999). For example, as demonstrated for the Elbe estuary, rivers may sediment excess nutrients there and the brackish water habitats often show the greatest "indigenous species minimum," so that more alien species can potentially establish (Nehring, 2006). In support for this, we detected strong increases in neophyte species richness around the city of Hamburg and the Elbe estuary (Region 2), which has the biggest German harbor, with strong international trade.
The area with strongest increases in neophyte species richness in the present study (Region 8, Figure 5) has been reported as a region of strong increases in the establishment and expansion of neophytes from the climatically mild Danube plains in a local-scale study by Sompek et al. (2017). The authors claim that this increased expansion is due to an increase in habitat suitability caused by climate change. Many neophytes are known to rapidly colonize newly available habitats, using mild valley refuges as a starting point for expansion (Rejmánek et al., 2005). Our spatial maps of changes in grid cell species richness confirm this spread and also show that the expansion of neophytes is widespread. Similarly, a study in the nature reserve in the northern part of the upper Rhine valley (close to Region 6) conducted by Vor and Schmidt (2008) reported an increase in neophyte species richness compared to atlas data from 1993 (Lang & Wolff, 1993). Based on our maps, we can show that this increase in neophytes is more widespread than the upper Rhine valley alone.
While acknowledging that our study cannot give empirical evidence for causal relationships, we conclude that the spatiotemporal patterns of change in the national plant biodiversity are highly variable, which is evidence for a complex interplay of drivers of biodiversity change. As demonstrated by the comparatively low level of spatial congruence in the detected hotspots of species richness and the lack of correlation between the hotspots of species richness change among the different floristic status groups, these factors act locally, and affect different species in different ways. Nevertheless, the fact that on the national level, net plant species richness is declining so pre-dominantly and apparently irreversibly is alarming.
While our study is focused on Germany, we have no reason to believe that these changes are only limited to this country. Plants, as primary producers, play pivotal roles in ecosystems and changes in their biodiversity may cascade throughout the food web and influence ecosystem functioning across trophic levels (Emmerson et al., 2016;Schuldt et al., 2018). For example, changes in the floristic composition of habitats in north Germany have been shown to result in lower potential nectar availability, with probably negative effects on pollinating insects . It is therefore possible that the detected large-scale changes in plant biodiversity is connected to recent insect declines (Hallmann et al., 2017;Seibold et al., 2019).

| CON CLUS ION
Declines in vascular plant biodiversity are widespread in Germany and apparent in more than 70% of plant species studied here. This includes approximately 40% of all moderately common to common vascular plant species in Germany. Urgent action is needed to halt this biodiversity loss. Our approach demonstrates how existing large datasets can be combined and used for reliable trend analysis.
Existing data should also be collated from other states in Europe and globally. The data integration and analysis approach used in this study is comprehensive and robust to different methodological biases. It can be applied to other large-scale research using heterogeneous occurrence data, given it can be harmonized to a common fine-scale grid. This makes our approach valuable for other projects, such as the growing Living Atlas community (Brenton et al., 2018) and may also help to inform and create new, more collaborative monitoring schemes that integrate knowledge and data from different actors in nature protection (e.g., governmental, academical, and volunteer; Kühl et al., 2020). Such new schemes should also include long-term monitoring of common species (see also Pescott et al., 2015). Clearly this is an ambitious endeavor, which can only be accomplished through joint efforts of a variety of stakeholders and should be underpinned with the financial and legislative power of (inter)national institutions. However, such approaches must not question the need for monitoring projects, which are still necessary and have the potential to identify the drivers of biodiversity change.

ACK N OWLED G EM ENTS
We are grateful to the many individuals who were involved in gathering the vast amount of data combined in our dataset and to those who helped to mobilize the data, especially to those who contributed to large data collections like FlorKart and the vegetation databases. Further, we are grateful to C. Storm, T. Heinken, T.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are openly available.
Species-specific occurrence estimates over time can be found in