Ecosystems are governed by dynamic processes such as competition for resources, reproduction and dispersal. These shape their biodiversity and how the system responds to change. Current approaches to modelling ecosystems, especially plants, focus on either describing fine-scale processes for individual species or broad-scale patterns for limited groups of plant functional types. Digitisation of herbarium and other plant records has unlocked a wealth of information that can be used to drive models of plant communities and make predictions for their future under different scenarios of climate change. The advent of increased computational capacity and fast, high level programming languages allows for simulation of such landscapes at unprecedented scales. Here, we demonstrate a tool for Ecosystem Simulation through Integrated Species Trait-Environment Modelling (EcoSISTEM), which models plant species across multiple ecosystem sizes, from patches and small islands to regions and entire continents. These simulated ecosystems support the ability to generate many different types of habitat, as well as reproducing different disturbance scenarios such as climate change, habitat loss and invasion. EcoSISTEM also reproduces examples of real-world species distributions by integrating plant occurrence records and global climate reconstructions to simulate plant species throughout the continent of Africa for the past century. EcoSISTEM allows us to flexibly explore the dynamics of tens of thousands of species interacting across a continent. The code parallelises efficiently across multiple nodes on high performance computing platforms, and has been scaled up to run on over 1000 cores. It allows us to study the impact of changes to climate, resources and habitat and investigate real-life mechanisms surrounding climate change and biodiversity loss. Following the recent Intergovernmental Panel on Climate Change (IPCC) report detailing how human-induced climate change has unequivocally altered our planet, we are now more aware than ever that deep changes need to be made to prevent future climate catastrophes (IPCC, 2021). Biodiversity loss, in particular, is expected to increase globally in response to loss of suitable habitats, in turn due to extreme weather events, sea level rise and ever-increasing temperatures. Plant ranges are limited by their dispersal capacities, often more so than animals, and for some species that exist in a particular climate niche, there may be no suitable alternative reachable habitat (Cunze et al., 2013). The International Union for the Conservation of Nature (IUCN) (SRLI) for Plants already estimates that at least one quarter of plants are threatened with extinction, while there are many thousands more species for which we do not have sufficient data to make an assessment (Brummitt et al., 2015). Great efforts have been undertaken over the past few decades to digitise herbaria and natural history records into online databases such as the (GBIF) (GBIF, 2019;Soltis & Soltis, 2016). Such data can be coupled with modelling efforts in order to fill these gaps. Africa is estimated to host approximately 45,000 plant species in total (Klopper et al., 2006) and has been the target of recent efforts to unite multiple datasets of vascular plant records (Dauby et al., 2016), although there are many areas estimated to have a higher species richness than the collection efforts reflect (Sosef et al., 2017). The continent is characterised by extremes of biodiversity, from the vast, species-poor Sahara desert, to the rain forests of the Congo Basin, the second largest rain forest region in the world storing more carbon per land area than the Amazon (Sullivan et al., 2017), to the Cape floristic region, the most speciesrich temperate area (Linder, 2014;Schnitzler et al., 2011). However, due to deforestation and destructive land-use practices, the tropical forests throughout Africa are now estimated to be net sources of carbon (Baccini et al., 2017). It is now more important than ever that we correctly characterise the biodiversity of these large-scale regions, so we can use this as a baseline for future measurements of loss. Advances in computing and big data approaches and the advent of fast, high-level programming languages such as Julia (Bezanson et al., 2017) allow for more sophisticated modelling to be carried out much more easily than could ever have been achieved before (Poisot et al., 2018). Although the need for large scale, species-focused modelling frameworks is urgent, those with greater geographical and taxonomic scope lack important biological processes such as demographics, dispersal and competition (Urban et al., 2016). Calls have already been made for more general models of biodiversity, highlighting the need for a framework similar to the IPCC to compare large-scale climate models, alongside monitoring of climate variables (Purves et al., 2013). Comparison with the IPCC is apposite, as climate change and land-use change are together the biggest threats to worldwide biodiversity. Dynamic models like the 'Madingley' model (Harfoot et al., 2014) already show success in replicating biodiversity patterns at broad functional group scale, despite reducing all terrestrial plants to a single group. Traditionally there have been three main approaches to modelling biological responses to climate change . Species-level approaches determine statistical correlations between current climate and species distributions to project future scenarios of environmental change (Keith et al., 2008;Thuiller et al., 2008). Due to their static nature, species distribution or environmental niche models (SDM/ENM) fail to account for important ecological processes affecting distributional change, such as competition, dispersal and evolution (Akcakaya et al., 2006;Hannah, 2014;Pearson & Dawson, 2003;Thuiller et al., 2008). Although SDMs are often developed species by species, this approach can be extended to global scope (e.g. Brum-mitt et al., 2021), but this would require modelling thousands of individual species to reveal community-level compositional change (Hannah, 2014). Furthermore, SDMs are often based upon incomplete species records and their interactions with climate or habitat and, without explicitly including biological mechanisms, predict new environments poorly. Many models do include whole-community dynamics at varying geographic and taxonomic scopes. Dynamic global vegetation models (DGVMs) simulate broad species classes as plant functional types (PFTs), with explicit mechanisms for processes including photosynthesis, growth and carbon cycling. DGVMs can be coupled with Earth system models incorporating feedbacks between environment and vegetation for improved climate predictions.
For example, the Joint UK Land and Environment Simulator (JULES) formerly classified plants into only five functional groups: broad-leaved and needle-leaved trees, C3 and C4 grasses, and shrubs (Clark et al., 2011). Coupled DGVMs are adept at modelling direct and indirect effects of climate change and disturbance on different vegetation types. However, although improvements include more detailed functional groups in JULES (Harper et al., 2016), these coarse PFTs limit the use of DGVMs in conservation planning, or differential within-group responses to changing environments (Thuiller et al., 2008). Equally, such models are primarily forced by detailed historical reconstructions and future climate predictions and often lack direct inputs of anthropogenic disturbance, such as habitat loss or land-use change (Quillet et al., 2010). At a much finer taxonomic and spatial scale there is the gap model, which simulates colonisation of forest gaps from falling trees. Although dynamic, including processes such as growth and competition, gap models act at the individual tree level to predict spatially-explicit changes to forest composition (Hannah, 2014;Thuiller et al., 2008). Gap models are often of commercially important tree species in small areas, for which demographic and dispersal data are readily available. They also historically focus on trees of temperate regions, although have been extended in recent years to grasslands, shrubs, and tropical trees (Fischer et al., 2016;Moorcroft et al., 2001). Scaling up to regional or continental levels is extremely computationally challenging, notwithstanding progress in simulating tree dynamics for the Amazon rain forest using parallel processing and high performance computing (Fischer et al., 2016;Shugart et al., 2018). Indeed gap models of larger areas often rely on aggregating species back into functional groups, similar to their DGVM counterparts (Shugart et al., 2018). No single approach currently provides all functionality necessary to model biodiversity change at large scales -a composite approach is needed. Such a model needs to sit in the overlap between DGVMs and traditional single-species models, operate at large scales, for many individual species, and include such biological processes as niche preferences, competition and dispersal . Here we present exactly such a composite dynamic ecosystem model -EcoSISTEM -and demonstrate its properties through idealised simulations of plant ecosystems across a range of scales and habitat types, from local patches with constant climates, to regional-or continental-scale areas with varied climates. We describe several example scenarios and provide information on the computational efficiency of the code. Finally, we demonstrate the capacity of the model to reconstruct real-world plant ecosystems across entire continents, drawing together plant records from GBIF and global climate reconstructions for the past century. Each characteristic can be found on an axis: geographic scale spanning from small to large, model approach spanning from mechanistic to correlative, and grouping of individuals from individual species to functional groups. Previous approaches, like dynamic global vegetation models (DGVMs), species distribution models (SDMs) and forest gap models tend to exist in the overlap of at least two of these properties, whereas EcoSISTEM unites all three, to operate over large scales, for individual species undergoing dynamic processes. The simulator is built around important, related facets of plant community composition, including trade-offs between resource availability and species' requirements, between the environmental variables and niche preferences. EcoSISTEM is a flexible and extensible framework used to simulate plant dynamics. The 2-D landscape is discretised across a regular grid, and accommodates an arbitrary number of plant species within an ecosystem of any size (we have tested it across Africa with over 50,000 species -the number of plant species on that continent). Within each grid cell we track the simulated change in abundance of each species over time. Resource (e.g. sunlight (kW) and water (mm)) availability and environmental data (e.g. temperature (  C)) are inputs into the model that affect the growth and survival of the plant populations. At each time step (e.g. monthly) the resource and environment is updated, and species' niche preferences, e.g. temperature, and resource requirements are traded off against availability in that grid cell. The match of species to environment and resource needs determines reproduction and death rates. At the end of each time step, after reproduction, seed dispersal takes place. An overview of the model is presented in First, we assume that plant dispersal only takes place via seed dispersal, i.e. dispersal follows immediately after seed production, as is observed in most plant species. Here we use a Gaussian dispersal kernel, which concentrates offspring close to the parent plant; however, the model can readily incorporate alternative choices of dispersal kernel. Second, reproduction is modelled as a Poisson process of germination rate per individual (or reproducing unit) and current population size, assuming that time until next seed production is exponentially distributed. This approximates the seed production processes, with an average number of germinations per adult plant per month determined by whether a plant is exposed to environmental conditions suitable for seed production, and whether there are sufficient resources to do so. Therefore, even though reproduction could occur at every time step in the simulation, environmental conditions and resource availability limit reproduction to windows of suitable conditions and show seasonal reproductive patterns in seasonally-fluctuating climates. Similarly, the probability of death also depends on environmental niche suitability and resource availability, and we model death as a Bernouilli process (Caswell, 2001). These distributions are standard in demographic modelling (Chisholm et al., 2014;Melbourne, 2012). The population dynamics of each plant species in the system are described by equation 1. The spatial domain is divided into N grid cells. The population at the next time step, t + t, at grid location X is a result of offspring dispersing and adults surviving from the previous time step t. The time step t can be altered, but is set at 1 month in simulations. The number of individuals of species s at time t in location X is pop s (t, X), Seeds s (t, X) is the number of seeds produced by species s and Deaths s (t, X) is the corresponding number of deaths.
Finally, k s (X, Y ) is the probability that a seed produced by a plant of species s at location Y disperses to X to grow. For each plant species, the number of individuals dying and the number of seeds produced at each time step is drawn from a Poisson and a binomial distribution, respectively, with the expected number of seeds per adult plant and probability of death during a time period of t as R b,s and 1  exp(P d,s ), respectively (Equations 2 and 3). where R b,s and P d,s depend on expected rates of seed production and death rates, b s and d s , in a suitable environment with adequate resources, moderated by the resource availability and environmental suitability via B e,s and D e,s : In the expressions for B e,s and D e,s (Equations 6 and 7, respectively), the last terms describe the effects of resource competition on each species, which depends on the ratios of resource i, K i (t, X), to the total of resource i used by all individuals of all species in the grid cell, E i (t, X) (= S s=1 s,i x pop s (t, X), where s denotes the species and s,i is the i-th resource requirement of species s). Here we assume there is a total of S species and R resources. Resource K i (X, t) is determined from resource (e.g. sunlight, water, etc.) data input into the model and is not affected by population abundance. There are no advantages when more resource is available than individuals to consume it. Species requiring high levels of resources reproduce less often, but also live longer, as observed in many species including plants (Reich, 2001) , and reflected in the first terms in Equations 6 and 7. This process is also modulated by fitness of a species in its environment (second terms in Equations 6 and 7), with those more mismatched facing decreased ability to reproduce and higher mortality through adjusted resource requirements,  i,s . For R resources and 1 environmental variable, then combined with resource competition, the above dynamics are captured in Equations 6 and 7. Resource competition , (7) where the average resource requirement of resource i for all species is  i (= 1 S S s=1 s,i ), and s,i  i denotes the relative resource i requirements of species s. The extent to which the species environmental niche ( Ts ) matches the experienced environment (T ) is captured in  i,x in Equation 8. The relative resource i requirement of species s, i,s  i is adjusted by the fit of the species niche preference to the current environment, f (T ; Ts ) (Equation 9), to penalise reproduction and increase mortality for species outside their natural habitats. f is modelled here as a Gaussian curve so that the penalty for being far away from the niche preference ( Ts ) is normally distributed (De Blasio et al., 2015), examples for which are given in Figure 3B. However, this can easily be modified to consider alternative distributions. where and T is the current environment in grid cell X, Ts is the niche mean for species s and  s the niche width (standard deviation) of the species. This rewards species staying within their niche more if the niche is narrow, thus giving an advantage to specialists over generalists. Only seeds produced in the current time step can disperse. We use a Gaussian dispersal kernel with standard deviation  s . The probability of a seed of species s dispersing from anywhere in the source grid cell, V , to destination grid cell, W , is k s (V, W ) as described by Equation 10: where (x V , y V ), (x W , y W ) are the grid square indices of the top right corner of the destination and source grid square, respectively. Space is scaled by grid cell width, , so that in the scaled space grid cells have unit width and height, as illustrated in Figure 4, and the mean dispersal distance of the species is  s = s  . The use of Gaussian dispersal kernels is contentious (Bullock et al., 2017), but it remains the simplest option in the absence of direct information about dispersal, such as dispersal modes and seed size. We illustrate results with the Gaussian kernel, but alternative dispersal kernels can be easily implemented in EcoSISTEM. All the code and underlying functionality for EcoSISTEM was written in Julia v1.5.3 (Besancon et al., 2021;Bezanson et al., 2017). The package is designed to interact with EcoJulia architecture and packages such as Diversity and Phylo (Reeve, 2021a,b), enabling the user to calculate diversity metrics for ecosystems and incorporate phylogenetic relationships between species. The code has been extensively tested and refined for speed and efficiency. It is multithreaded to run efficiently on multi-core workstations, and can be distributed across multiple nodes on high performance computing (HPC) platforms using the Message-Passing Interface (MPI) for inter-process communication (2015). Testing by the Numerical Analysis Group (NAG, Oxford, UK) on HPC found an essentially perfect (100%) global efficiency of the parallelisation in tests of up to 96 processes (Figure S1), and the code has been run on over 1000 cores. For the simulation of Africa, we used the European Centre for Medium-Range Weather Forecasting (ECMWF) climate reanalysis datasets (Dee et al., 2011;Laloyaux et al., 2018) to form the climate part of the ecosystem, combining the coarser 120km CERA-20C data with finer scale ERA-interim data to run simulations at an 80 km resolution. We used minimum monthly temperature (  C) and rainfall (mm) estimations for the full timespan, 1901-2018, to describe the climate niche. We also used the same combined ECMWF reanalysis datasets for the resources in the ecosystem. Here, they are surface solar radiation (kW/m 2 ) and soil water volume (m 3 water/m 2 ) at layer 1 (up to 7cm depth).
We filtered out those species with only a single recording point in the underlying GBIF data as they were too data deficient, leaving a total of 39,158 species recorded at least twice in Africa. The climate preferences were extracted from GBIF and ECMWF reanalysis datasets at the time and location of the GBIF record, corrected by sampling effort at that spatial location (C. Harris et al., 2021). We took the mean and standard deviation of the temperature and rainfall preferences as the niche mean and niche width for each species. There is little available information on requirements for solar energy and water for most plant species, so we took the environmental conditions they were exposed to as a proxy, assuming that they used all available resource in the areas they were found. Thus, we used the extracted information on resource needs for soil water volume and surface solar radiation, corrected for effort as in the niche preferences. Finally, we used randomised reproduction, mortality and dispersal parameters for this example as this information is unavailable for most plant species. To validate the model, simulations of small-scale ecosystem patches and islands were conducted. These were equal in size and seeded with 100 species / 100km 2 , but boundary conditions differed: islands had hard boundaries that species could not disperse beyond, whereas patches were toroidal and species could disperse in all directions. We set up each example with constant background temperature and ran each simulation for 10 years. We then expanded the case study to a grid with 80km squares and a total area similar to that of a continent. Using this continental backdrop, we investigated how quickly species could colonise large areas, whether they were previously populated or not. In order to initially seed the Africa simulations, we began by filling the ecosystems with the overall summed counts of each GBIF record according to the locations at which they had been previously found. These counts were adjusted by the availability of resource in each grid square, so that the ecosystem began with a number of individuals close to the equilibrium. Then we ran the simulations backwards in time from 2018 to 1901, adding in each species record at each location at the date it was recorded. This enabled us to simulate a baseline for plant diversity in 1901 from which to start the simulations, which could not be gained from a burn in forward through time. The number of individuals of each species added at each time point was scaled by the inverse of collection effort. This correction ensured that both under-and over-collected areas were seeded with similar numbers of plants despite dramatically different sampling efforts. Once the simulations reached 1901, they were cycled through this same year 10 times in order to ensure they were fully equilibrated. From this starting point, we ran the ecosystems ran through the full monthly climate history through to 2018. We performed 10 stochastic realisations, which was chosen after considering the computational expense of the simulations. We confirmed that abundance depends on available resources , here on an island ecosystem with two resources, water and sunlight, each on a gradient West to East and South to North, respectively. All species were seeded with the same resource requirements and vital rates. Abundance increased in cells with more water and sunlight, with some edge effects. As expected, abundances were invariant to the grid resolution , and ecosystems with greater areas could also support more individuals . We also tested ecosystems seeded with varying numbers of species, and similar numbers survive until the simulation ends in each case . Figure 6 shows overall abundance on an island with two species at opposite sides of the island after ten years of simulation: species with higher average dispersal distances moved faster and further into the unpopulated centre. We tested whether individual species' niche preferences and widths functioned as expected, on a small patch with species given different niche preferences and all other parameters kept equal. When all species had equal niche widths, species with niche preferences nearer the 25  C optimum were more abundant ; when all species had a preference for 25  C and varying niche widths, generalists with broader niches were less abundant than specialists with narrow niches . This is expected given the decline in species' match to their optimum environment as niche width increases in Figure 3B. However, if temperature was then increased by 1  C, species with narrowest niches went extinct, and generalists became more abundant , with a preference for species with a niche width of around 1  C, as predicted in Figure 3D. EcoSISTEM was designed to scale to much larger areas, supporting many more species. To illustrate this, we simulated 50,000 plant species at an 80km grid scale over a circular "continent" 8,000km across, with a constant environment of 25  C. When all species have equal fitness in the habitat, all 50,000 co-exist over long time scales of over 100 years . We also explored the selective advantage of invasive species: introducing an invasive with a narrow niche width tuned to the environment into an continent-sized landscape with an existing generalist, the invasive out-competes the generalist and spreads across the continent. Here, we define selective advantage as the difference in niche width between the invasive and generalist. The larger the invasive's selective advantage, the faster it will invade and colonise the landscape . When we simulate 50,000 generalist species coexisting across a landscape and introduce a single competitive invasive, the invasive can out-compete the local established populations . Figure 8D shows the Shannon entropy of communities within individual grid squares, using the Diversity package (Reeve, 2021a). The invasion centre-point, where the invasive dominates, is the most uneven species distribution, with a single species dominating whereas outside, every species coexists. The strength of competitive interaction can be adjusted so that invasives and generalists coexist in the system.
Finally, we simulated the competition, reproduction and dispersal of all of the plant species recorded on the African continent over the 20th and early 21st century, based upon digitised plant records and global climate reconstructions. As a result of the sparsity of GBIF records for plants in Africa many plant species had very low numbers of records, and competition between species during the burnin stage resulted in the survival of only approximately 33,000 species. Preliminary analyses seeding a threshold number of individuals from each species into the environment according to their climate niches suggests that the system can support the full complement of species, but that is not Instead we studied the 33,000 surviving species in our simulation study rather than the original 39,158 species. A comparison of their abundance between the beginning and end of the simulation suggests that they are relatively stable throughout the course of the simulation, with some stochastic gains and losses for rarer species saw similar losses (on average 1,968 (sd = +/-1452) species went extinct). Overall, simulated species richness decreased; however, local species richness was concentrated around well-known biodiversity hotspots in the central band of Africa, as well as the South . There were also some species rich regions that surrounded what are in reality urban areas, such as Lagos in Nigeria. Digitised taxonomic records from herbaria and natural history collections around the globe offer a unique opportunity, not just to monitor our current and historical plant biodiversity (T. Harris, Mulligan, & Brummitt, 2021), but also to make predictions for the future. In order for this to be feasible we must have models that operate over large geographic scales, at the species-level, that can incorporate a variety of mechanisms essential to plant growth, survival and reproduction. EcoSISTEM is such a tool, capable of simulating dynamics of tens of thousands of plant species at continental or even planetary scales. The initial case studies, 'Island, patch and continental dynamics', show that this model behaves as expected in various scenarios from small ecosystem patches to entire continents. On smaller-scale island and patch habitats, the simulator comfortably represents all desired properties from our (albeit non-exhaustive) list, and the code is tested automatically against these behaviours with an integrated test suite during development. Scaled up to whole continents, the simulator can sustain 50,000 species and demonstrates that an invasive can colonise such a habitat through competitive advantage. The main simulation study, 'Plant biodiversity of Africa', demonstrates the size and scale of modelling that can be achieved with EcoSISTEM. More complex environments and scenarios of environmental and biotic change can be readily simulated, together with the ability to incorporate climatic data layers, for example from WorldClim or ECMWF. Using such information, we were able to recreate patterns of species richness across the African continent, tracking the course of climate change across the past century. Here, we play out competitive, invasive dynamics across the African continent, comfortably supporting 33,000 species. The strength of the competitive interaction term can be tuned to simulate different propensities for co-existence. All species in the simulation were parameterised from analyses of online natural history records, but were also partly synthetic due to the absence of easily accessible data on reproduction strategies, dispersal and longevity. However, the model can be easily enhanced with any such data from real-life species as it becomes available. Despite the partly synthetic nature of our simulations, however, the model captured many features of estimated species distributions in Africa. Natural history collections offer a unique opportunity to further parameterise these life history traits, especially with the advent of digitised taxonomic descriptions, specimen records and photographs (T. Harris, Ottaviani, et al., 2021). Approaches to extract such information are currently underway, adapting natural language processing techniques to taxonomic descriptions. Future versions of the code will also seek to include important positive interactions between species, such as facilitation (Cavieres et al., 2014). As there are little data on diversity at this geographic scale, we compared patterns in species richness against those mapped globally ( Barthlott et al., 2007;Brummitt et al., 2021;Kier et al., 2005;Mutke & Barthlott, 2005) and for Africa (Dauby et al., 2016;T. Harris, Mulligan, & Brummitt, 2021;T. Harris, Ottaviani, et al., 2021;Linder, 2014). The species richness exhibited by the models broadly corresponds with well-known biodiversity hotspots. In Africa, the Cape floristic region showed as particularly diverse (Schnitzler et al., 2011(Schnitzler et al., , 2012)), as well as peaks in previously described highly rich regions of rain forest in Central Africa (Linder, 2014;Mutke & Barthlott, 2005). That the models broadly agree with global species richness patterns supports theories surrounding species richness gradients and climate, since these are the elements currently captured most faithfully by the model. There is thought to be a strong correlation between species diversity and climate when there is one limiting factor involved, such as energy or water availability (Mutke & Barthlott, 2005). The loss of biodiversity, about 5% of species, observed in simulations over the course of the century is impossible to verify given the absence of cross-sectional data that might inform us about current or historic African biodiversity at anything other than a very small geographic (and often taxonomic) scale. It is also difficult to set a baseline of climate stability, given that we cannot untangle the Earth's climate history from anthropogenic changes. Here, we chose to repeat the first year of climate reanalysis, 1901, as a proxy of pre-industrial climates. However, there are also evidently more factors involved in the current distribution of plants than climate. For example, urban areas and those undergoing land-use changes, such as those used for agriculture or forestry, are not accounted for in this model, but have a strong impact on the distribution and richness of plant species (Gerstner et al., 2014;Kier et al., 2009;Phillips et al., 2017). Indeed some regions, such as the area surrounding Lagos -the second most populous city in Africa -showed a richness of species much higher in our simulations than is observed in reality.
The simulation perhaps suggests that the area, with its tropical forest climate, might be very biodiverse were it not for the urban sprawl. More generally, Africa has seen widespread loss of native habitats due to deforestation and other destructive practices (Baccini et al., 2017;Gomes et al., 2019;Marengo et al., 2018;Perrings & Halkos, 2015;ter Steege et al., 2015). Elevation also plays a major component in species richness patterns (Mutke & Barthlott, 2005) and, although changes in climate are already accounted for in these models, there are other potential drivers of elevational effects, including boundary conditions (Colwell et al., 2004). There are other additional complexities related to the modelling of different types of biomes, particularly species-rich rain forests, where canopy cover plays an important role in regulating population dynamics and access to limiting resources, such as light (Meyer et al., 2019). EcoSISTEM is an open-source tool, freely available for use and extension via GitHub (C. Harris & Reeve, 2021), with several current applications beyond this continental biodiversity work, including simulating UK plant biodiversity under various climate and land-use change scenarios. Significant logistical challenges remain when modelling biodiversity at large scale, including parameterising this mechanistic, dynamic model for thousands or tens of thou-sands of species, extreme biases in available data on existing plant distributions, as well as computational logistics of running at such a large scale. This work constitutes the first attempt at a truly species-level, dynamic, continental-scale model of plant biodiversity parameterised with real plant records and climate reconstructions, and tested on the African continent. Although much work remains, we already see that it can capture broad patterns in plant biodiversity for these areas, and it identifies areas that may have undergone diversity change throughout the course of the 20 th century. Though we are not yet able to say the models are predictive, future directions involve fitting of the model to smaller regions with more complete data, and we are already able to run the model through projections of future climate change. The largest comparable plant vegetation models are DGVMs, such as JULES, which work with up to 9 plant functional types (Harper et al., 2016). Though global in scope, a lack of species-level information hinders their capacity to describe and predict trends in plant biodiversity (Mokany et al., 2016;Scheiter et al., 2013). Other attempts to model communities of plants concentrate on individuals, such as forest gap models, or use statistical methods to model individual species separately and together, as in SDMs and joint SDMs, respectively. Scaling up to thousands of species across large areas comprises a considerable computational demand for all of these approaches (Shugart et al., 2020;Wilkinson et al., 2019). In conclusion, here we present a simulator for plant biodiversity that includes demographic and competition processes and can be adapted for both size and type of habitat, as well as simulations of habitat and environmental change. Although the model is presented here with very incomplete information on the plant species simulated, it is nonetheless already parameterised on data from biological records and climate reconstructions. Using such data, leveraged from digital natural history collections, we show that changes in climate over the past century already have perceptible effects on the local plant biodiversity of Africa. This modelling approach also offers a route to exploit the value of natural history collections for novel information, including species traits and life history strategies. The source code for the analysis will be publicly available under an open source license on GitHub, and the release of the code associated with the paper will be provisioned with a DOI Table S1: Parameter details for case study 1: island and patch dynamics. Single values apply to all species, s, or grid cells, whereas others may vary within a given range. The source code for the analysis will be publicly available under an open source license on GitHub, and the release of the code associated with the paper will be provisioned with a DOI Table S1: Parameter details for case study 1: island and patch dynamics. Single values apply to all species, s, or grid cells, whereas others may vary within a given range.