*1.3. Trampling*

Another explanation for the increased abundance of woody vegetation is the effect of trampling. Trampling from the high frequency and density of pastoral farming causes significant declines in cyanobacterial soil crust [32,33]. Savannas are characterized by low soil nutrient content [34–36], although many areas have biological soil crusts that increase soil surface stability, thereby reducing nutrient loss by erosion and atmospheric nitrogen fixation [33]. Studies have found that the soil crust is greatly influenced by this pastoral trampling within 2 to 8 km of boreholes [37], and that *Acacia* (new *Senegalia* and *Vachellia* classifications) species are often found in higher abundances within areas closer to boreholes, due to their low palatability and the positive species-specific association between canopy and soil crust development [38]. Boreholes are narrow shafts drilled into the ground in order to extract water and are the primary source of water for livestock farmers in southern Africa. Furthermore, cattle rarely stray more than 13–18 km from these water sources in Africa [39], meaning areas closer to boreholes may have increased woody vegetation cover.

## *1.4. Fire*

Fire is a factor that restricts woody vegetation diversity and abundance, preventing the formation of canopies [40–42] as well as removing seedlings and subsequently preventing the establishment of new trees [43]. Furthermore, for certain species fire can also kill the larger trees [44,45]. Seymour and Huyser [45] found that infrequent fires were enough to kill established *Vachellia erioloba* trees, which are an important keystone species in the region, meaning an increase in fire frequency could have implications on biodiversity. In unmanaged areas, the build-up of large quantities of grass biomass in the understory results in high-intensity fires that are capable of destroying juvenile trees [46]. For example, Sankaran et al. [11] studied the effect of fire return intervals on the percentage of woody cover in African savannas and found that a shorter return interval reduced established woody cover, which kept the community in a juvenile state by 'top-killing' seedlings. In managed landscapes, fires are not as frequent or intense enough to have a discernible impact on mature trees [40], and a common feature of savannas is the reduction of fires due to mitigation strategies [47]. However, Joubert et al. [48] note that fire is crucial to disrupt transition from grassy savanna to thicket, and that managers who prevent fires at this stage are likely to experience bush thickening in the future.

#### *1.5. Research Gap and Questions*

Variation in species characteristics is fundamental to understanding biogeographic patterns [49]. One reason for the possible lack of conclusive evidence explaining the main drivers of di fferent woody vegetation patterns in previous research is the variation in how vegetation has been measured (e.g., single species, multiple species, richness, percent woody cover), as well as the di fferences in spatial scales of the previous studies (ranging from garden experiments to coarse continental extents). Assessing diversity as total species richness does not always adequately characterize the way in which species di ffer from each other, and it is these di fferences in traits, which often indicate that species respond in di fferent ways to changes in the environment [50,51]. Alternatively, studying only one species in isolation could lead to species-specific results that are not generalizable to the larger system or to other species. Several mechanisms (outlined above) have been invoked to explain the processes responsible for woody vegetation composition; however, these are often investigated separately at scales not best suited to land-managers, thereby impeding the evaluation of their relative importance.

Subsequently, this study focuses on the vegetation composition of the Botswana Kalahari, with the aim to investigate the relative influence of the environmental drivers of woody vegetation at a regional scale. By classifying species into morphological groups based on shared physiological traits, the drivers of woody vegetation richness and abundance can be interpreted more meaningfully at a regional scale that is more appropriate for landscape managemen<sup>t</sup> decisions. This study will explore three main questions: (1) what is the woody vegetation composition of the Kalahari in western Botswana? (2) What are the environmental drivers of woody species richness? and (3) what are the environmental drivers of woody species abundance?

#### **2. Materials and Methods**

#### *2.1. Study Area*

We conducted our research in western Botswana between 2009 and 2011 (Figure 1). We created a 950 km transect following the observed rainfall gradient along the western part of the Kalahari. This transect ran from Shakawe in the northwest of the country to Bokspits in the southwest of the country. Rainfall along the transect decreases from the north to south, ranging from a MAP of 550 mm to 350 mm [52]. Along this transect, we identified 15 regions (Figure 1) where we conducted multiple vegetation surveys. We selected regions on their accessibility and a minimum distance of 75 km to the previous region.

**Figure 1.** Location of the 15 regions along the Kalahari Transect where fieldwork was undertaken.

#### *2.2. Data Collection*

Vegetation was surveyed using the line interception transect (LIT) method. Within each region, we fixed six 100 m transects radially from a center point. For the dry season, the direction of the first transect was determined by a random number (between 0 and 360), and the further two transects were offset by 120 degrees. Transects of the wet season were spaced exactly between dry season transects, resulting in an offset of 60 degrees from the very first transect laid. Transects were placed 200 m from the center point to avoid over-sampling a small area. See Krebs [53] for a further description of the LIT methodology. We recorded all woody vegetation that was taller than 25 cm following the nomenclature provided by Palgrave [54], whereby average height, distance covered over the transect line, and distance and direction of the stem(s) were documented. Species richness and abundance were recorded at all transects, and species identity were recorded at all sites, with the exception of the wet season transects at Sites 1, 3, 4, and 5 due to uncertain species identification resulting from missing leaves. The results of the vegetation survey meant we had data from 78 transects for use in the statistical analysis.

Species were categorized into five morphological groups based on the classification guidelines outlined by Meyer et al. [55]. Morphological group I consisted of species characterized by bipinnate leaf structures and growth form ranging from multi-stemmed shrub like appearance to single-stemmed trees. Morphological group II included broad leaf species forming dense canopy structures where the majority of the growth form is either multi-stemmed (generally less than five stems) or single-stemmed. Morphological group III contained multi-stemmed broad leaf shrubs with closed canopies, seldom exceeding 2 m in height. In contrast, morphological group IV contained shrub species characterized by open canopies. Morphological Group V included relatively short shrub species (<1.5 m) with small, open canopies (<0.5 m in diameter). We also obtained data on precipitation, fire frequency, cattle density, and borehole locations that represent the possible drivers of diversity and abundance of woody vegetation (Table 1).


**Table 1.** Description of the environmental drivers used to explore the diversity and abundance of woody vegetation in western Botswana.

#### *2.3. Data Analysis*

We performed regression analysis in order to explore the environmental drivers of woody species richness and species abundances. Environmental variables were checked for multicollinearity using variance inflation factor, then standardized using z-scores in order to compare their relative influence on the ecological indicators. We performed all regression analyses using R 3.3.0. [61]. We selected regression analyses based on a preliminary evaluation of the data and their error distribution. Histogram exploration identified a mixture of Poisson and censored Gaussian distributions. We subsequently used a combination of generalized linear models with Poisson error distributions and Tobit regression models to analyze our data. For data that had a Poisson distribution, a Generalized Linear Model procedure with a Poisson error distribution and a log link function was used:

$$\log(y) = \beta\_0 + \beta\_1 X\_1 + \dots + \beta\_n X\_n \tag{1}$$

where *y* is the abundances, *Xn* is the *n*th predictor, and β*n* is the Poisson regression coefficient.

A censored Gaussian distribution represents a dataset that has a normal error distribution, but has some limit, either from below or above. Ecological data is often collected with a large proportion of the observations just above zero, while data cannot extend below zero or above certain thresholds (e.g., percentage cover). Tobit regression overcomes this bias and has been shown to perform better than ordinary least squares (OLS) (e.g., [62]) and is widely used in criminology (e.g., [63]) and land use change research (e.g., [64]). Species richness of woody vegetation is censored at zero (i.e., there cannot be a species richness of −1), and so any parameter estimates obtained by conventional OLS would be biased. Developed by Tobin [65], the Tobit regression model fits a set of parameters to where the dependent variable is left-censored at zero:

$$y\_i^\* = x\_i \beta + \varepsilon\_i \tag{2a}$$

$$y\_i = \begin{cases} \ 0 \text{ if } y\_i^\* \le 0 \\ \ y\_i^\* \text{ if } y\_i^\* > 0 \end{cases} \tag{2b}$$

where the subscript *i* = 1, 2, 3 ... *n*, indicates the observation, *y*∗*i* is an unobservable variable, *xi* is a vector of explanatory variables, β is a vector of unknown parameters, and ε*i* is the error term. To estimate the censored regression models, we used the censReg [66] and MaxLik [67] packages. Final models were selected based on Akaike Information Criterion (AIC) using both forwards and backwards stepwise selection of models. We investigated third and fourth order interactions, but these did not improve the final models. Therefore, our final models only include main effects and second order interactions.
