Water availability, bedrock, disturbance by herbivores, and climate determine plant diversity in South-African savanna

To identify factors that drive plant species richness in South-African savanna and explore their relative importance, we sampled plant communities across habitats differing in water availability, disturbance, and bedrock, using the Kruger National Park as a model system. We made plant inventories in 60 plots of 50 × 50 m, located in three distinct habitats: (i) at perennial rivers, (ii) at seasonal rivers with water available only during the rainy season, and (iii) on crests, at least ~ 5 km away from any water source. We predicted that large herbivores would utilise seasonal rivers’ habitats less intensely than those along perennial rivers where water is available throughout the year, including dry periods. Plots on granite harboured more herbaceous and shrub species than plots on basalt. The dry crests were poorer in herb species than both seasonal and perennial rivers. Seasonal rivers harboured the highest numbers of shrub species, in accordance with the prediction of the highest species richness at relatively low levels of disturbance and low stress from the lack of water. The crests, exposed to relatively low pressure from grazing but stressed by the lack of water, are important from the conservation perspective because they harbour typical, sometimes rare savanna species, and so are seasonal rivers whose shrub richness is stimulated and maintained by the combination of moderate disturbance imposed by herbivores and position in the middle of the water availability gradient. To capture the complexity of determinants of species richness in KNP, we complemented the analysis of the above local factors by exploring large-scale factors related to climate, vegetation productivity, the character of dominant vegetation, and landscape features. The strongest factor was temperature; areas with the highest temperatures reveal lower species richness. Our results also suggest that Colophospermum mopane, a dominant woody species in the north of KNP is not the ultimate cause of the lower plant diversity in this part of the park.

Main effect of bedrock. The results of statistical analyses are summarized in Table 2. Plots on granite harboured more herbs and shrubs than those on basalt (Fig. 3A, D; mean ± SD: 80.9 ± 16.8 and 12.3 ± 6.9 vs. 62.8 ± 18.2 and 8.2 ± 3.7; p = 0.001 and p = 0.0756, respectively). Similarly, plots on granite had a higher Shannon diversity H' of herbs than basalt plots (2.12 ± 0.60 vs. 1.65 ± 0.65; p = 0.0113, respectively; note that due to similar results of analyses of species richness and Shannon H' only the former are presented graphically). Plots on granite had more indicator species of herbs and shrubs compared to plots on basalt (37.2 ± 7.9 and 8.8 ± 3.8 vs. 28.5 ± 8.4 and 6.4 ± 2.7; p = 0.001 and p = 0.0697, respectively). However, neither herb and shrub cover (Fig. 3B, E), nor the proportions of indicator species significantly differed between bedrocks (Fig. 3C, F).
On granites, vegetation on crests and at seasonal rivers contained more indicator species of herbs compared to perennial rivers ( Fig. 5C; 38.8 ± 7.9 and 38.5 ± 5.2 vs. 34.2 ± 9.2, respectively), but on basalts, plots at perennial and seasonal rivers harboured more herb indicator species than plots on crests (29.1 ± 4.5 and 32.4 ± 9.8 vs. 24.1 ± 7.9, respectively; p = 0.0512 for the overall models on the effects of bedrock, habitat type and their interactions). The proportions of indicator shrub species in different habitats did not significantly differ between the two bedrocks (Fig. 5F).
Effects of bedrock and habitat on grass versus herb covers. The relative cover of grasses differed among bedrocks and habitats (p = 0.0721, and p = 0.0003, respectively), with higher values recorded on basalt (74%, range 7-99%) than granite (60%, range 11-97%) (Fig. 6). Concerning the type of habitat, the largest relative cover of grasses was recorded on crests (83%, range 48-99%), differing significantly from perennial . Richness, cover and proportion of indicator species of herbs and shrubs on granitic and basaltic bedrock. Bars show means for each type of bedrock (n = 30) and error bars are associated with standard deviation (SD). The bold letters above bars indicate significant differences between bedrocks (including marginally significant, 0.05 < p < 0.1, in normal font). www.nature.com/scientificreports/ rivers (50%, range 7-93%, p = 0.0001) and marginally significantly from seasonal rivers (68%, range 11-98%, p = 0.055). Perennial rivers also marginally significantly differed from seasonal rivers in their relative grass cover (p = 0.056).
Effects of bedrock and habitat on species composition. The most common species on particular bedrocks and in the habitats sampled are listed in Table 1. The multivariate ordination of data showed that the composition of herb communities differed between granites and basalts, but the difference was only marginally significant (p = 0.094); the species composition of shrubs did not differ between the two bedrocks (p = 0.13).
In contrast, the composition of both herbs and shrubs significantly differed among habitats (p = 0.002 and p = 0.002, respectively). For herbs, the crests differed from both seasonal and perennial rivers, and perennial rivers differed from seasonal rivers. For shrubs, crests significantly differed from seasonal and perennial rivers, but seasonal rivers did not differ significantly from perennial rivers. In the ordination plot, the first axis was associated with the presence of a river (as opposed to dry crest), the second with its type, separating seasonal and perennial rivers. Herb species such as Acanthospermum hispidum, Alternanthera pungens, and Amaranthus praetermissus were typically associated with perennial river plots, Abutilon fruticosum, Corchorus asplenifolius, and Panicum maximum with seasonal rivers, and Brachiaria nigropedata, Phyllanthus incurvus, and Pogonarthria squarrosa were typically found on crests (Fig. 7A). Shrubs Lippia javanica, Searsia gueinzii, and S. pentheri were typically found at perennial rivers, Hippocratea crenata and Ximenia americana at seasonal rivers, and Combretum zeyheri, Ormocarpum trichocarpum, and Terminalia sericea on crests (Fig. 7B).
Of the total 540 species recorded, 48.0% occurred in all three habitats, and 25.2% were specific to one of them. From the life-history perspective, similar proportions of forbs, grasses and woody species shared all three habitats (46.1%, 51.5% and 51.1%, respectively). Still, seasonal rivers harboured the highest number of woody species that occurred exclusively in that habitat (14, i.e. 10.1% of the total), as compared to perennial rivers (12 , and error bars show the associated standard deviation (SD). The letters above bars indicate significant differences among the three types of plots (including marginally significant, 0.05 < p < 0.1). Asterisks indicate the significance of the overall model when pairwise contrasts were not significant. www.nature.com/scientificreports/ species, i.e. 8.6%) and crest (6 species, i.e. 4.3% of all woody species). In contrast, perennial rivers were the richest in forbs recorded exclusively in that habitat (31, corresponding to 11.0% of all forbs) ( Supplementary Fig. 1). For both herbs and shrubs, a significant interaction between bedrock and habitat was detected (p = 0.002 and p = 0.008, respectively), indicating that species associated with different habitats reflecting water availability also differed between the two bedrocks. For example, Jatropha zeyheri, Panicum coloratum, and Zanthoxylum humile are associated with crests on basalt, while Commiphora africana, Cucumis hirsutus and Pogonarthria squarrosa

Effect of large-scale factors. The ordination model (RDA) for plant species richness and cover with all
large-scale predictors (see Supplementary Table 1), and bedrock and habitat type as covariables yielded nonsignificant results on canonical axes, either with or without covariables representing the spatial effects (Table 3). However, the forward selection identified mean temperature, variation in Enhanced Vegetation Index, and the distance from dirt roads as the most important predictors. Accordingly, the direct gradient ordination model (RDA) with these predictors was highly significant (p = 0.002) and remained so with the spatial effects accounted for (p = 0.008; Fig. 8).
The effect of temperature-related variables on plant species richness and covers was highly significant (p = 0.008) and remained so with the spatial effects included (p = 0.006). The effect of distances was non-significant but turned marginally significant after including the spatial effects (p = 0.086; Table 3).
In the linear regression models with residuals from the linear mixed-effect models (LMM) as responses, the total plant richness and herb richness were negatively and significantly related to all four measures used to characterize the long-term temperature, except for the total plant richness being marginally significantly related to temperature minima and maxima (p = 0.057 and 0.062, respectively). All significant relationships turned to be marginally significant after including the spatial effects into the LMM model (Table 4). Besides temperature, there were marginally significant effects of other variables on grass and herb covers in the LMM models; grass cover decreased with the distance of the plot from a gravel road (p = 0.092 and p = 0.08 with and without spatial effects, respectively). With spatial effects included, herb cover was positively and marginally significantly related to the variation in the Enhanced Vegetation Index (p = 0.082; Table 4). It should be also noted that the amount of variance in the data explained by the large-scale predictors is generally low, even in case of predictors with Table 3. Results of the ordination models (direct gradient analysis-RDA) exploring and testing the relations between large-scale predictors and plant species richness/cover in the Kruger National Park. In all ordination models, total plant richness, richness of herbs, richness of shrubs, cover of grasses and cover of shrubs were used as the response variables. Spatial effects are represented by the three main principle coordinates (PCO1, PCO2, PCO3). Significant effects are displayed in bold, marginally significant effects in italics. The significance values refer to individual predictors in case of forward selection, and to the significance of all canonical axes in the remaining ordination models. See Supplementary Table 1 for detailed description of variables. 1 Summary effect of rainfall (rainSum, rainMean, rainSD). 2 Summary effect of distance from potential sources of propagules-roads, rivers, rest camps and KNP boundary (DistBnd, distCamp, distTar, distDir, distRiv, distStrm; see Supplementary Table 1). 3 Summary effect of temperature and its variation over time (tempMean, tempSD, tempMin, tempMax). 4 Summary effect of surface water occurrence (waterSum, waterMean, waterSD). 5 Summary effect of Enhanced Vegetation Index (eviSum, eviMean, eviSD). 6 Summary effect of fire (FireSum, FireMean, FireSD).   Table 4. The results of univariate regression models testing the relationships between the large-scale predictors, identified as significant by the ordination models (Table 3), and measures of plant species richness and cover in plots in the Kruger National Park. For each response variable, values from models not considering spatial effects as well as those with significant principal coordinates (PCO1, PCO2, PCO3) included to account for spatial effects are presented. Significant effects are displayed in bold, marginally significant in italics. See Supplementary Table 1  Effect of Colophospermum mopane cover. The marginal effects, i.e. not considering the confounding factors, of the mopane cover and north-south gradient on total, herb, and shrub species richness were significant ( Table 5). The partial effects of the mopane cover (after excluding the effect of the north-south gradient) on all three measures of plant richness were not significant, but the partial effect of the north-south gradient (after excluding the effect of mopane cover) was significant, with p = 0.001, p = 0.004, and p = 0.011 for the total, herb and shrub species richness, respectively. After excluding the effects of bedrock, habitat, and triplet as a grouping variable, the partial effects of mopane cover were not significant for total species richness, but marginally significant for herbs and shrubs (p = 0.095 and p = 0.059, respectively). The partial effects of the north-south gradient (excluding the effects of bedrock, habitat, and triplet as a grouping variable) were marginally significant for total species richness and significant for herbs (p = 0.077 and p = 0.049, respectively), and non-significant for shrubs.

Discussion
Determinants of plant species richness: interaction of bedrock and habitat. Our study is based on the sampling of the complete floristic composition, both woody plants and herbs, of plant communities in those Kruger National Park habitats that are assumed to underlie the major gradient of variation in plant species richness and composition. The robustness of the data is illustrated by the fact that with the 540 taxa recorded in our plots, we captured over a quarter of the total plant diversity of KNP (1903 taxa in Eckhardt et al. 35 ; most recent number: ~ 2000 taxa following G. Zambatis et al., unpublished data), even though the total area of the plots represented less than 1/100,000 of the total area of KNP. Savannas are inherently heterogeneous ecosystems driven by complex interactions between physical landscape properties, rainfall, herbivores and fire to produce dynamic landscape patterns 5,36,37 . Correspondingly, we found that bedrock and water availability affect plant species richness, diversity, and composition of savanna communities in KNP. We found that the vegetation on granites is richer in species and more diverse than that on basalts, which have a more varied physical landscape template 37 but more fertile soils 38,39 . We suggest that the low number of species on basalts occur because of less favourable conditions, in terms of water supply and mechanical stress, for most plant species. Different soil types develop on each bedrock-clayey soils are typically formed by the weathering of basalts, while sandy soils are formed on granites. At low rainfall, clayey soils are relatively drier for plants than sandy ones because their permanent wilting point occurs at higher soil water content than at sandy soils 40,41 . Further, the reduced infiltration on clays results in quick runoff during heavy rains and the inability to store water in soil 9 ; shallow water penetration is more favourable for grasses, whose root zones are located in the upper soil layers 42 . In addition, basalts weather into more fertile soils than granites 38,39 , which contributes to the dominance of C4 grasses 17 , and in turn, may result in low plant diversity on this bedrock. Low water requirements of C4 grasses enable them to coexist with trees, the strongest competitors for light from the plant kingdom, whose competitive ability is reduced under limited water supply 15,43,44 . Lastly, soils differ in relation to geomorphology especially on the granitic bedrock; these catenas include coarse sandy soils on crests and clayey soils at the footslopes and river bottoms. Soil sequences primarily depend on terrain modulation and character Table 5. The results of variance partitioning between the effects of mopane (Colophospermum mopane) cover and the north-south (N-S) gradient in the Kruger National Park on plant species richness. The marginal effects represent the coarse effects, without accounting for the confounded factors, while partial effects represent the net effects after accounting for the confounded factors. Significant effects are in bold, marginally significant in italics, the direction of the effect is indicated by the plus/minus symbol. www.nature.com/scientificreports/ of watercourse network; they are more developed in hilly landscapes with big rivers 17 . In our study, we did not collect soil samples, thus we were unable to classify soils and assess their effect on plants. Most of our plots were located either on flat valley bottoms close to rivers or on flat crests. Therefore, we believe that soil variation in our plots could be roughly described as a combination of habitat type and bedrock (i.e. coarse sandy soils found typically on granite crests and fine clayey soils on basaltic crests). Nevertheless, the character of soil on plots near perennial and seasonal rivers may have been influenced, besides the bedrock, also by the alluvial effects. Besides, the effect of bedrock interacts with that of water availability; this is most strongly manifest by crests on basalts that host the lowest numbers of herbs, and their species poverty can be attributed to the combined effect of mechanical stress, drought, and competition of dominant grasses 24 . Desiccation of basaltic soils, leading to large cracks, may also damage the roots mechanically 45 . Such mechanical stress may favour the dominance of grasses, which have many fibrous roots, and vertical cracks damage only part of them, while forbs with a single taproot may be damaged more seriously 9 . The strong dominance of grasses on basalts also reduces space and resource capacity for other species, thereby decreasing richness and diversity. In our basalt plots, grasses made up on average 74% of the herb layer cover, markedly more than on granite, where the relative cover of grasses reached 60%. Accordingly, the basaltic crests also hosted the lowest numbers of indicator species typical for given vegetation units (as defined in Mucina and Rutherford 9 ).
Species' cover showed a pattern different from that of species richness and diversity: herbs, including grasses, reached the highest cover on basaltic crests, where their species richness was lowest of all habitats. The high grass cover on basalts may be explained by their fast recovery after heavy grazing due to a high amount of nutrients in soil on this bedrock (in contrast to granites that are poor in nutrients and the recovery is slower 35,46 ) and C4 metabolism of grasses that enables them to photosynthesize more effectively under dry and hot conditions than in C3 forbs 47 . This may suggest that the diversity of herbs is limited by the dominance of a few grass species rather than by environmental stress, but most likely, these two factors act in concert-mechanical and physiological stress acts as an environmental filter, selecting few grass species that reach dominance and further reduce the richness and diversity of other species.

Structure of plant communities: differences in species composition.
Only marginally significant differences in the species composition of herbs were detected between the two types of bedrocks. This result is somewhat unexpected, as acidophilous floras usually differ substantially from basiphilous floras worldwide 48 . However, the results from KNP suggest that the two bedrocks share most species; the vegetation on basalts represents an impoverished version of that on granite, without a major qualitative difference-the basaltic community is "nested" within the granitic community ( Supplementary Fig. 1).
On the contrary, compositional differences between the three habitats were significant for herbs as well as for shrubs. In general, plots near perennial rivers are rich in species but host many ruderal weeds or naturalized aliens, reflecting that such sites are both disturbed and rich in nutrients 34 . Contrary to that, plots on crests are less diverse but host the largest share of species considered important indicators of vegetation units typical for KNP 9 .
The significant interaction between bedrock and water availability shows that the response of both herbs and shrubs to water availability is specific for granites and basalts. The most apparent aspect of this difference is the remarkable dominance of grasses (Bothriochloa radicans, Panicum coloratum, Themeda triandra) on basaltic crests compared to those on granites.
The importance of seasonal rivers: disturbance by herbivores and water availability. The effect of herbivores, elephants in particular, on vegetation has been thoroughly studied 49 , yet not in the perspective of the two factors addressed in our study or directly exploring the role of seasonal rivers. With high population densities, elephant herds influence the savanna ecosystem not only by consuming large amounts of plant tissue (mainly leaves) but also by damaging or uprooting grown trees 18,50 , hence changing woody savanna into more open grass dominated states 51,52 . While some plant species cope with these disturbances well, other important tree canopy species regenerate poorly 53 . Other species regenerate more easily but cannot cope with bark-stripping and thus continued elephant use keeps them from attaining the size they would reach in the absence of elephants 54 . Bark-stripping also makes the trees more vulnerable to fungal infection, insect attacks, and fire 55,56 . Elephants also affect trees by pulling out seedlings or eating their apical meristems, disabling them to grow to a size at which they can survive fires 28,57 . The impact of other large herbivores can be even stronger. For example, Scogings et al. 58 reported that species richness and density of woody plants increased five years after exclusion of all large herbivores, but not after the exclusion of elephants alone (see also Barnes 57 ). The joint impact of all herbivore species on vegetation influenced not only species richness but also vegetation structure, as documented by Asner et al. 19 who found that 3-D structure of woody vegetation differed significantly between areas protected from and accessible to herbivores, with up to 11-fold greater woody canopy cover in the former areas.
The vegetation on crests, where water availability is most limiting of the three habitats we sampled, hosts a lower diversity of herbs than seasonal or perennial rivers. Shrubs represent a life history that is most exposed to herbivores' impacts, and as such, they followed a different pattern in our study-they reached ~ 26-27% greater species richness near seasonal rivers than in the other two habitats. While the difference in shrub species richness between seasonal rivers vs. perennial rivers and crests was only marginally significant (which can be attributed to a low sample size given by the field research logistics), it seems to hold for both bedrock types. Seasonal river plots on granites harboured on average 31% more shrub species than perennial rivers and 24% more than crests, with corresponding figures on basalt being 18% and 32%, respectively (Fig. 5D).
This finding supports our prediction; we hypothesized that shrubs near perennial rivers might be damaged by animals or riparian disturbances 59 , while shrubs on crests may suffer from low water availability, mechanical root damage and the influence of fire 60 . Therefore, seasonal rivers where the relatively low disturbance combined www.nature.com/scientificreports/ with a low stress from the lack of water seem to be the habitat that supports the greatest richness and diversity of shrubs, while the total shrub cover remains at a comparable level to the other two habitats (~ 35%). This conclusion is further supported by the fact that seasonal rivers harboured the greatest number of woody species that occurred only in this habitat (Supplementary Material File 2). Rivers are an important factor driving elephant distribution and densities, as these animals generally concentrate closer to major rivers 61 . Smit and Ferreira 29 further suggested that based on aerial census data collected annually during the mid-dry season of 1985-2007 across KNP, rivers of different sizes can be used as a proxy for elephant densities in the dry period of the year. They showed that since 1999 when water provision in KNP was reduced and culling ceased, elephant densities at seasonal rivers (both large and intermediate) were lower than at perennial rivers and higher than at sites far from rivers. Assuming a correlation between elephant density and elephant impact, these authors report a disproportionally greater increase in impact around large rather than small rivers, including seasonal, and crests 29 . Together with our data on occurrences of all herbivores recorded by camera traps, these findings support our hypothesis of the greatest shrub species richness at seasonal rivers due to relatively lower disturbance, interacting with water availability gradient.
The signal we observed, pointing to the role of seasonal rivers as supporting high species richness of shrubs, a key plant life form of the savanna ecosystem, is consistent across landsystems and supported by data on the herbivores' habitat affinities. From the management perspective, the role of seasonal rivers in maintaining high diversity of woody vegetation should be taken into account when predicting the risk of local extirpation of specific woody species 29,62 , especially in light of increasing elephant densities closer to rivers 61 .
Too hot to handle: the effect of large-scale factors. The only strong signal from the large-scale analyses was related to temperature. The variation in Enhanced Vegetation Index, representing the fluctuation in vegetation phenology over time and space 63 , and distance from the gravel roads, which can be interpreted as opportunities for rapid propagule dispersal, were also identified as important in determining plant species richness, and cover. However, their effects were only marginally significant and not consistent, depending on whether or not the spatial effects were included in the models. Other variables such as the number and intensity of fires, rainfall characteristics, or the presence of surface waters did not show significant effects in our models. This does not mean they are not important, rather the opposite is true as evidenced by decades of botanical research in South-African savannas 5,11,12 . Within our system, however, the results strongly suggest that habitat, bedrock, and herbivores are the key determinants of local plants species richness, cover, and community composition. We argue that large-scale factors, such as temperature, further influence the effects of these factors thereby explaining part of the residual variation. The non-significant effects of some park-wide large-scale factors can be further related to complex spatiotemporal interactions. For example, the possible effects of the presence of surface water, related to soil moisture as a critical factor for plants, could be overwhelmed by the water availability in nearby rivers. Similarly, the effect of rainfall and possibly fires will also manifest in the variation in vegetation productivity. Or, the effect of fire on plant richness could be masked by the direct effect of bedrock-this would be indicated by the number of fires on granites in our system being marginally significantly greater than on basalts (p = 0.099). In addition, we faced a statistical limitation associated with the sampling design-the power of all tests was limited due to rather low number of independent replicates (20 triplets, 60 plots altogether). On the other hand, the plots were spread across the whole area of KNP and therefore can be assumed to be representative for the effects of large-scale predictors, such as climate or fire regime.
The temperature in KNP is negatively related to the species richness of all target groups of plants, i.e., shrubs, herbs, and all species. Here, it needs to be emphasized that our models allowed us to identify the pure effects of temperature, not confounded by other variables and geography of KNP. This shows that high temperatures are associated with low plant diversity in a situation when other important drivers (bedrock, water availability, geography) are constant, suggesting that it is the heat stress what limits the plant diversity. Indeed, plant diversity in warm regions with a limited water availability decreases with increasing temperature 64,65 . This pattern, however, is only poorly documented in savannas. For example, Knight et al. 66 found a strong negative relationship between annual solar radiation and plant richness in southern-African savannas.
As expected, the effect of the mopane cover on plant characteristics was heavily confounded with the north-south gradient in KNP. However, the partial effects of the mopane cover (after accounting for the north-south gradient) were not significant and were featured by low portions of explained variability. On the contrary, the effects of the north-south gradient were significant even after accounting for mopane. Although the south of KNP hosts richer flora, and mopane dominates the vegetation only in the north, we suggest that the high cover of mopane does not reduce the floristic richness in the north of the park, at least not so dramatically as previously thought 33 . Besides the statistical analyses we present here, we have also recorded several plots with a high cover of mopane rich in herb flora. For example, plots within a triplet near Phalaborwa (PHA2) hosted 75, 88 and 87 herb species at perennial river, seasonal river and on crest, respectively, and the covers of mopane at these plots ranges between ~ 40-70%.

Species indicating the representativeness of savanna vegetation.
Interestingly, crests showed the highest proportion of herb and shrub species considered as indicators of typical vegetation (as defined in Mucina and Rutherford 9 ), contrary to plots at perennial rivers, which harboured the least indicator species typical for those habitats. This result shows that even though the plots around perennial rivers were species-rich, a part of their diversity was due to the occurrence of widespread weeds (Acalypha indica, Amaranthus praetermissus, Pupalia lappacea) or invasive aliens (Acanthospermum hispidum, Alternanthera pungens, Gomphrena celosioides 34 ). The affinity of weedy species to perennial rivers is also the reason why this habitat hosted the highest numbers of forbs occurring only there. www.nature.com/scientificreports/ To conclude, plant species richness, diversity, composition, and cover all respond to bedrock and water availability in South African savannas, as well as to their interactions. For the woody life histories, the effects of these factors are fine-tuned by herbivory. The plant community characteristics, however, differ in their response to the main factors. For example, crests exhibit the lowest species richness and diversity, but the proportion of 'important species' (sensu Mucina and Rutherford 9 ) indicating typical, well-developed vegetation units is the highest of all habitats. Similarly, herbs, including grasses, reach the highest cover on basaltic crests, suggesting that their lowest diversity in this habitat may be associated with the high biomass of a few dominants.
These results indicate that the typical vegetation of different habitats is associated with manifold characteristics. The crests, exposed to relatively low pressure from grazing but stressed by the lack of water, are important from the conservation perspective because of the presence of typical, sometimes rare savanna species, and so are seasonal rivers whose shrub richness is increased and maintained by the combination of low levels of disturbance and position in the middle of the water availability gradient. Sampling and field data. MOSAIK's primary objective is to sample plant and animal (mammal, bird, and insect) communities in habitats across the four main landsystems in KNP. To account for differences in water availability, we used plots arranged in spatially defined triplets, each comprising three habitats: (i) near a perennial river or another permanent water source, such as artificial water points or dams, (ii) near a seasonal river with a lack of water during dry periods, and (iii) on dry crests, at least 5 km from any source of water (Fig. 1). Each triplet was located on either granitic (Skukuza and Phalaborwa landsystems) or basaltic (Satara and Letaba landsystems) bedrock 68 (Fig. 1). Therefore, we captured the effects of two factors: type of bedrock, defined at the level of triplets, and habitat representing water availability, defined at the level of individual plots within each triplet. The plots within each triplet were selected within a distance of ~ 7-13 km between them. In total, we sampled vegetation in 60 plots, 50 × 50 m in size, each bedrock represented by 30 plots, and each habitat by 20 plots across the entire KNP (Fig. 2). The plots were selected with a focus on shrubby savanna; those by perennial rivers were located behind the riparian gallery forest, and although some trees occasionally occurred in plant inventories, they were not included in analyses here.
Plants were sampled during two rainy seasons (when the vegetation is in optimal phenological stage allowing sampling), 16 January to 4 February 2019 and 17 January to 3 February 2020; we sampled 33 and 27 plots, respectively, stratified so as to proportionally represent all landsystems and habitats in each sampling period. The January rainfall for KNP in 2019 and 2020 was 88 mm and 97 mm, respectively (averaged across all KNP stations). All vascular plant species (for simplicity, intraspecific taxa such as subspecies were included under the term 'species') were recorded in each 2500 m 2 plot. Their abundance was visually estimated using the Braun-Blanquet cover-abundance seven-grade scale 84 . To quantify the occurrence of species in plots, the Braun-Blanquet scores were transformed to percentage values as follows: 5 = 87.5%, 4 = 62.5%, 3 = 37.5%, 2 = 15%, 1 = 2.5%, + = 1.0%, r = 0.02% 85 . Abundance within the herb, shrub, and tree layers was estimated separately. The time spent by three botanists to sample a plot ranged from 1 to 7 h, with an average of 2:15 ± 1:01 h (mean ± S.D.); the same three botanists sampled all sites. www.nature.com/scientificreports/ Plant species were classified according to life form (herbs, shrubs) and in each plot, the following characteristics were calculated separately for each of these life forms and used as response variables in statistical analyses: (i) species richness expressed as the total number of species in a plot; (ii) species diversity expressed as the Shannon index H' , calculated using species percentage cover as importance values 86 ; (iii) total cover of the herb and shrub communities (%) estimated in the field 84 . To assess the extent to which the vegetation in plots represents typical natural vegetation as defined for South Africa, we identified species that are labelled as 'important taxa' for the respective vegetation units and indicate their phytosociological assignment 9 . Species labelled as such in any of the 13 vegetation units captured by our sampling are further termed as 'indicator species, ' and (iv) their proportion was used as an additional characteristic analysed.
To estimate the numbers of herbivores in the three habitats and relate the extent of disturbance to vegetation with botanical data, we used camera traps located in the same plots (n = 60), one in each (Bushnell Essential E3 Camera Trap with low glow IR flash). The cameras were serviced once in ~ 3 months, starting in June 2018. The monitoring was part of the MOSAIK project, aimed at recording the diversity and activity of mammals. Here we use animal records over 140 days, including both the dry season (August-October 2018) and the rainy season (December 2018-February 2019). For each species, we calculated the number of animal records taken by the camera in a given site. The overall herbivore load per plot was then expressed as the sum of all the herbivorespecies' records, i.e. reflecting their abundances. We considered large herbivores and megaherbivores, both browsers, grazers and mixed feeders, where an impact on vegetation is expected (18 species in total): buffalo, bushbuck, common duiker, elephant, giraffe, grysbok, hippo, impala, kudu, nyala, black rhino, white rhino, sable antelope, steenbok, tsessebe, waterbuck, wildebeest, and zebra.

Large-scale variables.
To analyse the effect of large-scale factors on plant species richness and cover, we collected large-scale variables related to the fire history (number of fires and their frequency), vegetation productivity (Enhanced Vegetation Index, EVI), rainfall (sum and mean), temperature (mean, minimum, maximum) and surface-water occurrence density (sum, mean). All these variables were measured in a 4-km 2 grid cell surrounding the sample plot and over the long-term (2000-2019). For all variables, we also included the measure of their variation over time (see Supplementary Table 1 for detailed description and data sources). In addition, we also measured the distance of the plot from rivers, roads, restcamps, and Kruger boundary, to account for opportunities of propagule dispersal, both naturally and by human action (Supplementary Table 1). Data analysis. The differences in herb and shrub species richness, Shannon diversity H' , herb and shrub total covers, proportions of indicator species, and herbivore numbers were analysed using LMM regression models within the nlme package of R software 91 . The identity of a triplet was arranged as a random factor (grouping variable), and the spatial distances between individual triplets were modelled as a continuous function, based on the GPS coordinates. As a result, the univariate LMM model on species richness, diversity, or total cover was: lme(response variable ~ bedrock*habitat, random = ~ 1|triplet identity, cor = corGaus(form = ~ gpsn + gps e)). The differences between habitats (perennial rivers, seasonal rivers, crests), if significant, were subsequently tested by the Tukey posthoc pairwise comparison of estimated marginal means, using the emmeans package in R 92 . The data on total covers and proportions of indicator species were arcsine transformed. The data on species richness and Shannon diversity H' were square-rooted if Shapiro-Wilk test revealed significant deviations from normality.
The differences in species composition between the two bedrock types (basalt vs. granite) and the three habitats reflecting water availability (perennial rivers, seasonal rivers, crests) were tested using multivariate ordination methods and Monte-Carlo permutation tests in Canoco 5 software 93 . A hierarchical split-plot arrangement was used to account for the pseudoreplication structure in the data: individual triplets were defined as whole plots and the permutations at the whole-plot level tested the main effect of bedrock. In contrast, permutations within the whole plots tested the main effect of habitat (water availability), a factor defined at the split-plot level. Finally, a model testing the bedrock × habitat interaction was created by permuting the whole-plots as well as split-plots.
The spatial effects, given by the location of individual triplets (whole-plots) were accounted for by the method of PCNM (Principal Coordinates of Neighbour Matrices 9 ), which uses the Euclidian spatial distances to create a matrix of spatial vectors, independently at different scales 94,95 . The scores of the first three spatial PCoA vectors were then used as covariables in the ordination models that accounted for the spatial effects.
To explore how the large-scale factors are correlated in affecting plant species richness and cover, we used the direct gradient ordination analyses (RDA). Total plant richness, herb and shrub species richness, and covers of herbs, shrubs and grasses (in %) were response variables; the large-scale factors listed in Supplementary Table 1 were predictors (Table 1). Here, we also analysed the total vascular plant species richness, which is traditionally the main variable of interest in large-scale analyses, i.e. not only shrub and herb as in analyses of local factors where we were more interested in vegetation structure defined by vegetation layers; for the same reasons the Shannon diversity H' was only tested in local-scale analyses, not in the large-scale one. The plots were permuted in a split-plot scheme to account for the arrangement of individual plots into triplets. Individual triplets were defined as whole plots, while the split-plots were represented by individual vegetation plots within triplets. The manual forward selection of environmental variables was used to identify the large-scale predictors with the strongest explanatory powers. Further, the bedrock (basalt, granite) and habitat (perennial river, seasonal river, crest) were set as covariables because we were interested in the variability beyond that explained by these local factors. The method of PCNM was used to account for spatial effects, beyond those given by the arrangement of individual plots into triplets. As a result, the first three principal coordinates (PCO1, PCO2, PCO3) were included in the ordination models as covariables representing spatial effects. www.nature.com/scientificreports/ To test the significance of the relationships between large-scale factors and plant species richness/cover, indicated by the ordination methods (Fig. 8), we used linear regression models. First, the effects of principal coordinates (PCO1, PCO2, PCO3) on a particular richness/cover measure were tested by a linear regression model. Then, we ran an LMM model with bedrock, habitat, bedrock × habitat interaction, principal coordinates with significant effects, triplet as a random effect, and richness/cover as a response variable. Finally, the effects of large-scale factors were tested in a regression model, with residuals from the corresponding LMM model as a response variable. By doing this, we accounted for (i) the effect of bedrock, habitat and their interaction, (ii) the effect of spatial autocorrelation caused by the arrangements of individual plots in triplets, and (iii) the autocorrelation given by the distribution of individual plots as well as triplets in the whole target area of KNP.
The effect of the Colophospermum mopane cover (estimated in the sampling plots in the field) on species richness (total, herbs, shrubs) was tested by using a linear regression model. Variance partitioning was applied to separate the effects of the mopane cover from those of the north-south gradient. Further, the partial effects of mopane and the north-south gradient were tested after accounting for the effects of bedrock (granite, basalt), habitat (perennial river, seasonal river, crest), and their interaction, similarly to the effects of the large-scale predictors. The effects of the mopane cover and north-south gradient were tested on residuals from the LMM model with bedrock and habitat as predictors, and triplet as a random effect.