Introduction

Conservation has traditionally responded to potential loses of biota by establishing legally protected areas (Pressey 1994; Rodrigues et al. 2004). Wittemyer et al. (2008) have shown that average human population growth rates on the borders of protected areas in Africa and Latin America were nearly double the average rural growth, suggesting that protected areas attracted human settlement. People perceive or obtain benefit from their proximity to such areas (de Sherbinin and Freudenberger 1998; Scholte 2003) but, there could be a concomitant threat to biodiversity within them. Many species are continuing to decrease within protected areas (Brashares et al. 2001; Newmark 2008) often due to the illegal wildlife harvesting for meat and trophies (Milner-Gulland et al. 2003). This is particularly true for African nature reserves where local species extinctions are directly linked to human population proximity, high reserve perimeter to area ratios, and bushmeat hunting (Brashares et al. 2001; Ogutu et al. 2009).

In the Serengeti ecosystem, Tanzania, there have been marked declines in black rhino (Diceros bicornis), elephant (Loxodonta africana) and African buffalo (Syncerus caffer) inside the protected area (Dublin et al. 1990b; Metzger et al. 2007; Sinclair et al. 2007). Declines in the numbers of large herbivores were attributed to cessation of anti-poaching activities during a period of economic decline. Analysis of the trends in the buffalo population over the whole area has suggested that population change was primarily due to illegal hunting, and that enforcement of wildlife laws reduced the illegal offtake (Hilborn et al. 2006) a conclusion also reached for other areas (Hilborn et al. 2006; Jachmann and Billiouw 1997; Keane et al. 2008; Leader-Williams and Milner-Gulland 1993). Using 50 years of buffalo census data, Hilborn et al. (2006) established that illegal hunting and enforcement activities could account for the overall trends in buffalo population yet examination of the buffalo total counts indicated variation in the buffalo population recovery; some areas have almost completely recovered from the population low of 1994 and other areas have failed to recover. Therefore, the main purpose of this paper is to analyse the possible causes of these spatial differences.

Buffalo are known to be targeted by illegal hunters (Sinclair 1977). Park rangers who actively search for snares and signs of illegal hunting have identified buffalo carcasses in the field (Hilborn personal observation) and buffalo meat appears in villagers bushmeat diets (Ndibalema and Songorwa 2007). Illegal hunting remains a large threat to conservation efforts in the Serengeti (Holmern et al. 2007; Kaltenborn et al. 2005; Loibooki et al. 2002) and, therefore, we determined whether illegal hunting was a contributing factor to the spatial differences in buffalo recovery.

Many factors can contribute to variation in animal population change including disease, food supply, drought, and natural predation. Disease has been well monitored in the buffalo and we have no indication that the observed changes in population have been caused by epizootics. Rinderpest virus originally caused major declines in buffalo numbers after 1890 but the virus has not caused declines since the 1960s (Dobson 1995; Dublin et al. 1990a; Rossiter et al. 1983; Sinclair et al. 2008), and indeed it is now globally extinct (Normille 2008). Bovine tuberculosis (Myobacterium bovis), although prevalent in South Africa (Cross et al. 2009), has not been found in Serengeti buffalo (Cleaveland et al. 2008; Sinclair 1977). Drought can be a major controlling factor and drought induced mortality occurred in 1993 causing approximately 40% mortality in the buffalo population. This mortality was equally distributed across the ecosystem and therefore cannot be responsible for the spatial patterns in recovery (Dublin et al.1990a; A. Sinclair unpublished data). While it is possible that other factors may contribute to the spatial variation of buffalo recovery, the major controlling factors are likely to be food supply, natural predation and illegal hunting. We analyzed the impacts of these three factors—hunting, food supply and natural predation—using a spatial analysis to separate out their effects. Thus, human population density and rate of increase, which we show are related to hunting within the reserve (Campbell and Hofer 1995; Hofer et al. 2000), are greatest in the west and northwest. In contrast, food limitation, which is a function of rainfall (Sinclair 1977), is most severe in the east and south, while predation is evenly spread over the buffalo range. The greatest food supply is in the north where rainfall is highest (Fig. 1). The next highest food levels are in the west, while the lowest food supplies are in the east (Sinclair 1977). During the 1960s, prior to the population collapse, these northern areas supported the highest densities of buffalo recorded in Africa, and in general Serengeti buffalo are limited by food and not by predation (Sinclair 1977).

Fig. 1
figure 1

The Serengeti ecosystem in Tanzania, East Africa includes the Serengeti National Park as a protected area, and the game reserves and conservations areas. These are the Ikorongo Game Reserve, Grumeti Game Reserve, Maswa Game Reserve, Ngorongoro and Loliondo Conservation Areas, which surround Serengeti National Park and have restrictions on settlement within their borders. The Serengeti National Park is divided up into zones (north, far west, centre, far east and south). Rainfall isohyets, showing the highest rainfall in the northwest and the lowest in the southeast. Rainfall data collected at local rainfall stations across the Serengeti ecosystem has been interpolated to produce the isohyets

Full size image

Materials and methods

Study area

The Serengeti-Mara ecosystem is located east of Lake Victoria and northwest of the Ngorongoro highlands and the Rift Valley (Fig. 1) and is described elsewhere (Sinclair and Arcese 1995b; Sinclair et al. 2007; Sinclair and Norton-Griffiths 1979). Serengeti National Park is designated IUCN land category II and is managed for ecosystem protection and recreation. A network of game reserves and conservation areas are located to the west and east of Serengeti National Park (Fig. 1). This whole area is known as the Greater Serengeti Ecosystem. The east of the national park boundary is settled by Maasai pastoralists who rarely hunt for wild meat and their lifestyles tend to be consistent with conservation of wildlife (Polansky et al. 2008). In contrast, human settlements to the west of the park boundary do consume game meat regularly (Holmern et al. 2006; Loibooki et al. 2002; Nyahongo et al. 2005).

Buffalo total counts

Beginning in the early 1960s, buffalo populations were censused by aerial survey every few years. A detailed description of methods is given in Sinclair (1977). In 1970 all observations of buffalo (individuals and herds) in the Greater Serengeti Ecosystem were plotted on a map of the ecosystem. These observations were later incorporated into a GIS using the Universal Transverse Mercator (UTM) coordinates. From the 1992, 1998, 2000, 2003, and 2008 censuses similar data were obtained using global positioning system (GPS) technology. The buffalo population was close to its maximum in 1970 and this census was therefore used as the baseline with which we compared the following years. We determined the instantaneous rate of change in the buffalo population from 1970 to 2008 by zone. Zones within the park (Fig. 1) represent distinct geographical and ecological areas. Buffalo herds are relatively sedentary, confine themselves to a home range of less than 20 km in diameter, and so rarely cross over zone boundaries (Sinclair 1977). These zones were the north, far east, far west, center, south and short grass plains. Because buffalo do not use the short grass plains we did not include this area in our analysis. We summed buffalo numbers within each zone for each year that we had census data and compared these numbers with those in 1970 to show the relative change. A major drought in 1993 affected all zones and caused a 40% mortality (Sinclair et al. 2007, 2008).

Spatial population dynamics model

We used a spatially structured population dynamics model to determine the trends in buffalo abundance in the five different regions between 1965 and 2008 (Hilborn et al. 2006). We examined a range of possible influences on abundance. These factors included carrying capacity, which is a function of size of zone times rainfall (a surrogate for food supply, Sinclair and Arcese 1995a), lion predation, and hunting effort. The population \( \hat{N}_{a,y} \), the predicted number of buffalo in zone a year y is given by;

$$ \begin{aligned} \hat{N}_{a,y + 1}& = \left[ {\hat{N}_{a,y} + r\hat{N}_{a,y} \left( {1 - {\frac{{\hat{N}_{a,y} }}{{k_{a} }}}} \right)} \right]S_{y} - u_{a,y} \hat{N}_{a,y} - E_{a,y} \\ {\text{while}} \\ u_{a,y} &= qv_{a} P_{y} \hfill \\ {\text{and}} \\ E_{a,y}& = \hat{N}_{a,y} \left[ {1 - \exp \left( {zL_{y} } \right)} \right] \\ \end{aligned} $$
(1)

where r is an intrinsic rate of growth assumed to be the same in all zones. Carrying capacity for zone a is k a and S y is survival from drought in year y, assumed to be 1.0 for all years except 1993, the year of the drought. The exploitation rate from hunting in zone a and year y is u a , P y is the relative hunting effort in year y, v a is the relative hunting effort for zone a, and q is a scalar relating hunting effort and area specific vulnerability to the exploitation rate. E ay is the number of buffalo in zone a killed by lions in year y, L y is an index of the number of lions in buffalo habitat in year y, and z scales the lion abundance index to lion mortality rate.

We explored a range of nested models, in various configurations that either included or excluded hunting, lion predation, and rainfall. We estimated the parameters using census data for each of five zones assuming a lognormal likelihood

$$ L\left( {N_{a,y} |{\text{parameters}}} \right) = {\frac{1}{{\sigma \sqrt {2\pi } }}}\exp \left( { - {\frac{{\left[ {\ln \left( {N_{a,y} - \hat{N}_{a,y} } \right)} \right]^{2} }}{{2\sigma^{2} }}}} \right) $$
(2)

where N ay is the observed number of buffalo in zone a, year y, and σ the standard deviation of the lognormal observation process.

The relative hunting effort (P) is poachers arrested per number of patrols day−1 (see Hilborn et al. 2006. Figure 1b). The zone specific vulnerability parameters (v a ) were estimated relative to that in the north which was fixed at 1.0. The parameter q is the harvest rate per unit of hunting effort (P) in a zone with v = 1.

Food supply and rainfall

We also considered a range of hypotheses regarding carrying capacity. First, we assumed all zones had the same carrying capacity. Secondly, we assumed that carrying capacity in each zone (k a ) was proportional to the size of the zone and the rainfall. Thus,

$$ k_{a} = pA_{a} R_{a} $$
(3)

where A a is the area in square km of zone a, R a is the average dry season rainfall in zone a, and p is a scalar to relate the product of area and rainfall to the carrying capacity.

While rainfall was the primary determinant of the food supply in most of Serengeti (Fig. 1), the far east differed by lacking riverine grassland. In this zone rainfall was less suitable as a predictor of resources (Sinclair 1977). Hence, thirdly we estimated the carrying capacity for each zone independent of its size and rainfall.

Intrinsic rate of increase and lion predation

While we could, in theory, estimate the intrinsic rate of increase (r) from the spatial data using the likelihood in Eq. 2 we found that the estimates obtained in that fashion were much lower than the total population growth rate in the 1960s and 1970s. This is because the variability of the data by zone is much higher than the variability for the total population. We estimated the intrinsic rate of increase (r = 0.092) from the total census between 1965 and 1976. We fixed r in all scenarios except that when lion predation was added we needed to increase the r value to account for additional lion mortality

$$ r = 0.092 + \left[ {1 - \exp \left( {zL_{1958} } \right)} \right] $$
(4)

Each model scenario included different parameters estimated or fixed. In all cases we estimated the number of buffalo in each zone in the first year, and the parameter σ. Thus, the simplest model has 6 parameters, plus a single carrying capacity (k) for a total of 7 parameters. As r is fixed in all models it is not considered an estimated parameter.

Fine scale buffalo population rate of increase (1970–1998 and 2000–2008)

The spatial trend in buffalo population was examined by comparing two time periods; 1970–1992 and 2000–2008 by creating a fine resolution map of buffalo population change across the park. To do this we first constructed a buffalo density map. In the GIS we divided the Serengeti National Park into 5 × 5 km areas and all observations within each 25 km2 area were summed. These numbers were then transformed to density (animal km−2) within each 25 km2 area. In order to smooth across the 25 km2 cell boundary the whole park was subdivided into 1 km2 units. The 30 nearest neighbor 1 km2 cells were averaged for each 1 km2 cell using the neighborhood analysis tool in ArcGIS 9.2. This allowed us to reduce the heterogeneity created from the clumping effect of large herds in some grid cells adjacent to empty cells. The 30 km2 area was of a similar magnitude to the maximum home range of buffalo (Sinclair 1977). We calculated the instantaneous rate of population change per year (r) using the raster calculator tool in ArcGIS 9.2 spatial analyst. Instantaneous rate of population change is defined as:

$$ r = ln\left( {N_{t} /N_{0} } \right)/t $$
(5)

where N t is the population size at time t, N 0 is the population size at the start of the time period, and t is the number of years between the two. The r calculation was performed on each cell in the density map for the two time periods 1970–1992 and 2000–2008.

Relation between buffalo numbers and human densities

We calculated the distance of each buffalo observation to the nearest edge of the park where there was human settlement in 1970, 1992, 1998, 2000, 2003 and 2008. Using Pearson’s correlation coefficient we determined the spatial correlation between buffalo counts and distance to humans (see below).

Hunter population estimates

We used two years of human census data, 1978 and 2002 (Bureau of Statistics, Dar es Salaam) for the area west of the Serengeti National Park boundary to Lake Victoria. Census data were organized by local areas called wards (similar to US counties). The area (km2) of each ward was known and we converted the ward population to density (humans km−2). From the human density we calculated the hunter density. Hunter density is a proportion of human density, which changes with the distance from the protected area boundary. The equation is that given in Campbell and Hofer (1995). For each ward, we determined the instantaneous rate of change of hunter density between 1978 and 2002.

Results

Total population numbers of buffalo in the protected area

Figure 2 shows the changes in total numbers of buffalo in the Serengeti National Park since 1965. At that time the population was recovering from the impacts of the viral disease, rinderpest, and numbers subsequently increased to a peak of 74,237 in 1975 (Sinclair 1977). Shortly after that, in 1977, anti-hunting activities were severely restricted by an economic crisis in Tanzania (Hilborn et al. 2006; Sinclair and Arcese 1995b) and widespread hunting on this species (and others) followed (Dublin et al. 1990a). By 1992 anti-hunting efforts had returned but the population had been reduced to 36,119 animals, some 49% of the peak number. The sharp decline to 21,186 buffalo in 1994 reflects the effect of a severe drought amounting to an additional mortality of 42% of the remaining population. Since 1998 the population has slowly increased. The most recent census (2008) of buffalo recorded 28,524 individuals in Serengeti National Park. Because these are total counts there are no sampling errors associated with the data. However, bias errors have been calculated, accommodated by technique design and kept constant over the years (Mduma and Hopcraft 2008; Sinclair 1972).

Fig. 2
figure 2

Buffalo population trends for the Serengeti National Park. Data from Sinclair et al. (2007)

Full size image

Buffalo population trends by region

Figure 3a, b presents the distribution of buffalo herds in 1970 before the main hunting period, and in 2003 during the recovery phase. These show that northern and western parts of the protected area have lost herds while the center and east have developed larger herds. Figure 4 shows the proportional changes in buffalo population in each zone (see Fig. 1), relative to 1970, the year when we have complete spatial distribution of animals prior to the onset of hunting. By 1992 the north had lost 84% (±5% (95%CL)) the far west some 38% (±9%) and the center 29% (±7%) of their numbers. In contrast, the south had lost 23% (±10%) and the far east only 12% (±6%) of the population. Since the drought of 1993, the south has increased above the 1970 level (120% ±3%) and the far east is at 62% (±6%) of those levels. The far west and center, although just beginning to recover by 2008, are still only at 54% (±9%) and 40% (±8%) of their original numbers. The northern population has been unable to recover at all and remains at a mere 2% (±0.3%) of original numbers. In summary, the three regions with western borders had consistently lower recovery throughout the period 1970–2008 than the east and south.

Fig. 3
figure 3

The location of all buffalo herds in the park-wide censuses of (a) 1970, and (b) 2008 showing the loss of herds in the north and far west. Data from Sinclair (1977), S. A. R. Mduma unpublished. Dots represent the size of the buffalo herds at each location. The largest dot represents the largest herd at 2,800 individuals and the smallest a single buffalo

Full size image
Fig. 4
figure 4

Proportional changes in buffalo population in the five zones of the Serengeti at different times relative to the starting number in 1970. Ninety-five percent confidence intervals were calculated (largest was 0.12%) but were too small to show

Full size image

Spatial population dynamics model

Details of the model (Eq. 1) can be found in Table 1. In our basic model configuration we assumed that the carrying capacity of a zone was proportional to the area and rainfall (Eq. 3). The second model included the same hunting effort in each zone of the park with no lion predation and no drought. The third model included lion predation (Eq. 5) but no hunting effort and no drought effect. These first three models fitted the data poorly. In model 4 hunting differed in each zone but had no lion predation and the fit of the model improved greatly. Model 5 was similar to model 4 but included the mortality from the 1993 drought (Eq. 1) and again the fit of the model improved. In model 6 we allowed the carrying capacity in the far east to be different from that of other areas (for the reasons explained above that resources differed), and this provided another significant improvement in fit. Again building on model 6, in model 7 we included the impact of lion predation and this too provided an improvement. Thus, the model incorporating unequal hunting effort, survival rates resulting from drought, carrying capacity in the far east estimated separately, and lion predation provided the best fit to the census data (Fig. 4). Using the likelihood ratio model 7 would be the preferred model.

Table 1 Candidate models of buffalo population changes over the last 50 years in the five regions of the Serengeti
Full size table

Final model parameter estimates

The model that explained the most variation in population across the zones was model 7 (Fig. 5). Using this model we estimated that the north had the highest intensity of hunting with the exploitation rate in 1982 (the worst year for hunting) being 31%. The far west had the second highest hunting being 16% in 1982, while the far east was estimated to have zero hunting mortality (Table 2). In addition, the carrying capacity in the far east was not adequately estimated from area and rainfall, and so was estimated independently in model 7. Lion predation rate was estimated to be 10% (assumed constant in all areas), and the 1993 drought mortality was estimated to be 48%.

Fig. 5
figure 5

Observed abundance of African buffalo (dots) and model predictions (solid line) for the zones of the Serengeti and for the total population

Full size image
Table 2 Final ‘best’ model parameter estimates that predict population changes for the five different regions (L was 10% for the final model). Hunting was greatest in the North zone
Full size table

Fine-scale analysis of buffalo and human population changes

The fine scale spatial analysis produced a gradation in the rates of buffalo population increase (Fig. 6) during the hunting period (1970–1992). There were negative rates of increase in the northwest and positive rates of increase in the east and south. The far west was more complex but rates of increase were still lower there than in the east.

Fig. 6
figure 6

Fine scale spatial differences in the rate of population change 1970–1992 showing the greatest loss in the north and far west. Dark areas represent negative population increases and light areas represent higher values (r = –0.3 to +0.05)

Full size image

A similar pattern (Fig. 7a) is exhibited during the increase phase (1998–2008) with population decreases in the northwest and west and population increases in the east. In the increase phase, the areas of population decreases were more concentrated and restricted to the northwest and west of the park compared to the hunting phase. While there were areas in the western corridor that still exhibited population decreases the area south of Grumeti Game Reserve shows population increases compared to the hunting phase.

Fig. 7
figure 7

(a) Fine scale spatial differences in the rate of population change 2000–2008 showing the slowest increase in the north and far west. Dark areas represent negative population increases and light areas represent higher values (r = –0.9 to +0.48). (b) Instantaneous rate of population change of hunter population densities to the west of Serengeti National Park. Dark areas represent high population growth whereas light areas represent low population growth (r = –0.6 to +0.59). Location of fastest increase is adjacent to areas of slowest increase in buffalo seen in Fig. 7a

Full size image

This pattern of buffalo population growth is the converse of the human population growth adjacent to the protected area (Fig. 7b). Hunters living within 40 km of the protected area were estimated as 20,000 in 1973 and 36,000 in 2002. The instantaneous rate of increase was 0.03 per year, similar to the national average. Numbers within 10 km of the protected area increased from 13,000 to 24,000, a similar rate of change (3% per year). However, the population adjacent to the northwest border of the park has increased at a faster rate (4.4% per year) mostly through immigration (Campbell and Hofer 1995). In 1978 human population densities varied ranging from 0.1–954 people km−2 and in 2002 the range in densities spanned 15–2,840 people km−2.

Buffalo increase and distance to reserve boundary

In 1970, prior to the decrease phase there was no spatial relationship between buffalo occupancy and distance to the reserve boundary. During the hunting phase there was a positive relationship between distance and buffalo numbers (1992, r = 0.15, P-value < 0.001, and 1998 r = 0.23, P-value = 0.03). When enforcement was increased (2000–2008) there was no relationship between distance and buffalo.

Discussion

Our results explain the spatial variation in buffalo population recovery across the protected area and elaborate on the work of Hilborn et al. (2006), which confined itself to the time trends of the whole buffalo population. Buffalo population changes are best explained largely by hunting but model fit was improved with the addition of predation mortality. Food supply was only a factor in areas where hunting was least, namely the east and south. In addition, our spatial results are consistent with the trends in the elephant population (Sinclair et al. 2007, 2008). Elephants behave more as one cohesive population, which moves away from disturbed areas and finds sanctuary in more peaceful areas, so that densities are a reflection of movements.

Buffalo on the other hand are extremely philopatric and remain within their home range irrespective of the disturbance (Sinclair 1977). Therefore, buffalo numbers reflect the local dynamics of an area. The buffalo populations in close proximity to areas with higher rates of human settlement had low or negative intrinsic rates of population increases and, therefore, were either slowest to recover or failing to recover at all. Areas of slow buffalo recovery are consistent with the previous analysis by Campbell and Hofer (1995), which identified similar areas of high human exploitation on a different suite of species, namely the resident antelopes. Hunter populations reside along the western and north-western boundaries of the protected area, and incursions are made into the park from the west. Human populations have increased considerably in the past 30 years (Fig. 7b) and buffalo numbers in the north and the far west reflect these changes in human populations along the boundaries. In contrast, the strong gradient in food supply, which is determined by the rainfall gradient (Fig. 1), is opposite to the trends in population recovery, i.e. areas of high food supply are those with the least recovery.

In conclusion, the increase in human populations along the western boundaries of the Serengeti ecosystem has led to negative consequences within the protected area on wildlife populations, as indicated by trends in the buffalo population. This result is consistent with the predicted impacts from increases of human settlement around protected areas elsewhere in Africa (Harcourt et al. 2001; Wittemyer et al. 2008). Hilborn et al. (2006) suggested that these negative consequences are mitigated by increases in enforcement of wildlife laws by protected area authorities.