Animals need to navigate between resources such as water, food and shelter, and how they achieve this is likely to vary with species. Here, using high-accuracy GPS data, we studied repeated journeys made by wild plains zebra (Equus quagga) through a naturally vegetated environment to explore whether they consistently follow the same route through the area or whether they use a range of routes to reach their goal. We used a model to distinguish and quantify these two possibilities and show that our observations are consistent with the use of multiple routes. Our model performs better than assuming a uniform angular distribution of trajectories. The typical separation of the routes was found to be small (1.96 m), while the scale at which neighbouring trajectories are informative to direction of travel was found to be large (with a confidence interval of 1.19–26.4 m). Our observations are consistent with the hypothesis that zebra are able to navigate without having to return to previously used routes, instead using numerous different routes of similar trajectories.

In natural ecosystems, resources are typically widely distributed, meaning that animals must move through the landscape to access them (Wiens, 1976). To do this while minimising use of time and energy, animals need to be able to navigate between key locations.

They could achieve this using a number of strategies. The first is through innate navigational abilities, for instance through knowledge of their location, destination and the use of magnetic or sun compasses. This has been reported for a variety of birds and aquatic species but not for most terrestrial mammals (Lohmann et al., 2007; Geva-Sagiv et al., 2015).

A second approach (beaconing) relies on moving in towards a directly perceptible cue. It ranges from using a single cue, situated close to the ultimate goal, to the use of a series of perceptual cues to allow them to locate the ultimate goal relative to their current position. Beaconing has been demonstrated in many species (Zollner and Lima, 1999) and can use a range of senses (Croney et al., 2003); however, it is limited to goals within an animal's immediate perceptual range and can therefore only be used for long-distance navigation in species with large perceptual ranges such as birds and bats (Williams et al., 1966).

It is hypothesised that a third approach, mental maps, allows animals to navigate to destinations outside of their immediate perceptual range by memorising the relative location of a series of landscape features and then using this neural representation to navigate through a large-scale environment. The most complex form of a mental map is a cognitive map (O'Keefe and Nadel, 1978) where it is theorised that an animal mentally stores information about all important features in its landscape, as well as the relative angle and distance between those features, similar to an actual physical map. Although a popular concept, the cognitive cost and mental capacity necessary for such a map may be prohibitive (Bennett, 1996) and at present there is no evidence of its use in any animal taxon (Bennett, 1996; Janson, 2000; Byrne, 2000).

A suggested alternative that is supported by a number of biologists is the memorisation of a limited number of specific routes or navigational cues between different sites within an animal's range (Byrne, 1979; Milton, 1980; Poucet, 1993). Such a map could be formed of routes, multiple successive beacons, or a combination of routes and beacons (Dumont and Petit, 1998). The use of routes, the repeated use of the same course or way between a starting point and a destination, has been recorded in numerous species, including ungulates living in savannah, forest and mountain environments (Agnew, 1966; Kays et al., 2011; Newmark and Rickart, 2012). Route use has also been widely recorded in arboreal primates; for instance, wild ranging howler monkeys repeatedly used the same routes that either were close to areas of high resource density or were elevated (de Guinea et al., 2019), whilst spider and woolly monkeys repeatedly used the same routes, many of which were associated with distinct topographical features, whilst travelling through their home range (Di Fiore and Suarez, 2007).

Repeated use of the same route can lead to the creation of physical tracks; the track network in African landscapes is extensive, with the same tracks often being used for many decades (Haynes, 2006), including by savannah-dwelling ungulates such as zebra (Smuts, 1976). It has been hypothesised that route use facilitates navigation by reducing the complexity of a landscape through limiting movement decisions at each step and therefore reducing the amount of information an animal needs to retain to successfully navigate between two points (Newmark and Rickart, 2012). It also reduces movement costs by compacting substrate and reducing superstrate barriers (Shepard et al., 2013).

Day to day navigational strategies utilised by terrestrial animals have been studied in laboratory (Morris et al., 1982) and engineered settings (Beecham, 2001; Mueller and Fagan, 2008; Ulanovsky and Moss, 2008) but rarely in natural landscapes and then mostly at much finer scales; for instance, in ants, where heading sensing and awareness of distance covered are key (Collett et al., 2013). Other work has focused on flying, rather than terrestrial, species such as birds and bats (Tsoar et al., 2011). Many theoretical studies have considered the importance of perceptual cues and cognitive ability on animal movement patterns (Beecham, 2001; Mueller and Fagan, 2008; Hirvonen et al., 1999); however, these capabilities have not been fully evaluated and are challenging to quantify in field settings.

In this study, we used a combination of high-rate and high-accuracy GPS movement data and mathematical modelling to analyse the movement of zebra navigating over long distances between resources to differentiate between potential navigation strategies. Zebra exist in stable family groups, harems, consisting of a stallion, two to five adult females and their offspring (Skinner and Smithers, 1990). Whilst harems occasionally join to form larger groups (herds), movement decisions are made at a harem level (Smuts, 1976) so tracking one individual is equivalent to tracking the whole harem. Zebra, like many herd-dwelling herbivores, travel one behind another in a linear fashion when moving longer distances (Smuts, 1976); it can therefore be expected that all zebra of the same harem follow the same linear route. The zebra studied here spend the dry season in the Moremi Game Reserve in the Okavango Delta, Botswana. When in their dry season range, the zebra drink at a single water source but use disparate grazing areas. They move between these distinct grazing and drinking areas (around 10 km apart; see Figs 1 and 3) every 2–3 days. The regularity and frequency of these journeys provide a high number of repeated movements, making this an excellent system to investigate the navigational strategies utilised by wild ranging herbivores moving distances greater than their immediate perceptual range.

We used a variant of the technique of function approximation by radial basis function networks (Broomhead and Lowe, 1988) to investigate how routes (for this paper, defined as a sequence of occupied spatial positions, i.e. GPS fixes) recently used by a zebra influence its future movements when returning to a known location.

We characterised how the zebra move in terms of two possible navigation strategies, one where zebra consistently follow the same route through an area and an alternative strategy where there is minimal preference for particular routes and multiple routes are utilised. The first scenario would indicate reliance on, and movement between, a limited number of learnt landscape features, the second an ability to navigate to the distant destination through a more advanced knowledge of location and orientation or through the learning of numerous different routes. Our model also allows us to determine how predictive other routes are: is the heading used by a zebra predicted by a contemporaneous route when such a route is nearby?

We hypothesised that zebra utilise multiple routes to reach their destination, but that routes are highly predicted by other nearby routes, as a consequence of zebra using knowledge gained from previous journeys to the destination to aid the efficiency of subsequent journeys.

Subjects

Nine female zebra, Equus quagga Boddaert 1785, each from a different harem, were fitted with GPS collars of our own design (Wilson et al., 2013, 2018). The zebra selected were members of a migratory sub-population that move between the Moremi Game Reserve, part of the Okavango Delta, and the Makgadikagdi National Park, a large salt pan and grassland system, in northern Botswana (Bartlam-Brooks et al., 2011). These two protected areas fall within one of the largest continuous tracts of land available to wildlife in southern Africa. Numerous species, both herbivores and carnivores, utilise movement corridors within this greater system as part of either seasonal migratory or dispersal movements.

Zebra live in stable harems, composed of an adult male, up to five adult females and their juveniles. Unfortunately, because of harems temporarily mixing during darting, the precise harem structure for collared individuals could not be recorded but as most foals are born in the wet season in this zebra population (Bartlam-Brooks et al., 2011), the age structure of harems was broadly similar with no young foals. Adult females were selected to reduce the risk of collar damage from intra-sexual fighting; because zebra live in harems that include different age and sex cohorts, collaring only adult females does not introduce bias due to sex- and age-specific movement patterns.

Animal handling

Zebra were darted from a stationary vehicle by an experienced wildlife veterinarian. Zebra were sedated with 7 mg Etorphine, 1667 i.u. Hyalase and 80 mg Azaperone. During sedation, collars were fitted snugly at the top of the zebra's neck and basic size measurements were made (neck circumference, hindleg and foreleg length). Dart sites were infused with intra-mammary antibiotics to guard against infection. Immobilisation was then reversed with 80 mg Naltrexone. All zebra recovered rapidly and were observed rejoining their harems, and none showed any lasting effects. All collars were fitted with mechanical drop-offs (Sirtrack, Hawkes Bay, New Zealand) to ensure collars fell off after 18 months.

All animal handling procedures were carried out under ethical approval from the Ethics and Welfare Committee of the Royal Veterinary College, London (URN 2013 1233) and via Botswana Department of Wildlife and National Parks Research Permits held by A.M.W. (EWT 8/36/4 plus additions) and Botswana Veterinary Registration held by A.M.W. Darting permits were provided by the Department of Wildlife and National Parks (DWNP) based on research permit EWT 8/36/4 XXIV (193) prior to all darting.

Study area

The Moremi Game reserve is located in the Okavango Delta in northern Botswana, between 22.0°–23.5°E and 18.5°–20.0°S. A satellite image of the study region can be found in Fig. 1A. The Okavango Delta is a large landlocked alluvial fan, covering 22,000 km2 (McCarthy and Ellery, 1998). It is fed by the Okavango River System, which originates in the Angolan highlands. Rainwater falling in the catchment basin creates a flood surge that peaks in the Delta in June, the middle of the dry season. The Delta therefore experiences two moisture regimes: the annual rains and the annual flood (McCarthy and Ellery, 1998). The system supports a large variety of wildlife species; 32 species of large mammals and 72 species of small mammals have been recorded (Ramberg et al., 2006). High habitat heterogeneity, a shortened limiting season due to the dual moisture regime and the dynamic vegetation succession resulting from the variable flooding pattern may all contribute to the Delta's ability to support higher biodiversity than the region's poor quality soil and low annual rainfall would predict (Bonyongo and Harris, 2007; Ramberg et al., 2006).

Fig. 1.

Zebra movements recorded between grazelands and water, and model predictions. (A) Trajectories (red) superimposed on a Universal Transverse Mercator projection of north Botswana (Google Maps). The bottom left corner of this map is at 23°29′E, 19°00′S. (B) An aerial 3D photogrammetric survey reconstruction of one such track in a typical area of the survey region imaged with 3D photogrammetry. (C) Example depiction of the effect on a trajectory on our model predictions with different values for σ and λ (see Materials and Methods). The blue line is a trajectory with an observation in the top left of each plot (one of the circles in from the top left) to the bottom right. These trajectories are identical in all four sub-panels. The brightness of the circles indicates the probability of a trajectory at that angle at that location. The angle in the circle is the angle of movement: red is high probability, blue is low probability, grey is a value comparable to the uniform distribution. In each panel, there is a black and a red line: the black line is of length λ while the red line is of length σ. σ and λ were selected by hand for illustrative purposes. σ controls the scale of the basis. The top row of C has a short range basis and as a result the angular distribution tends towards a uniform distribution rapidly as we move away from the blue trajectory. In the bottom row, σ is large and the distribution is non-uniform over much larger ranges. λ controls how rapidly the trajectories transition from parallel movement to perpendicular return to a track. On the left, this scale is short and our model predicts a rapid return to a trajectory. On the right panel, the scale is long and the model predicts movement parallel to the blue trajectory everywhere. For zebra, the parameters are similar to those displayed in the bottom right panel. (D) Position traces of the zebra, and the distribution of predicted angles of movement for the zebra (key as for C).

Fig. 1.

Zebra movements recorded between grazelands and water, and model predictions. (A) Trajectories (red) superimposed on a Universal Transverse Mercator projection of north Botswana (Google Maps). The bottom left corner of this map is at 23°29′E, 19°00′S. (B) An aerial 3D photogrammetric survey reconstruction of one such track in a typical area of the survey region imaged with 3D photogrammetry. (C) Example depiction of the effect on a trajectory on our model predictions with different values for σ and λ (see Materials and Methods). The blue line is a trajectory with an observation in the top left of each plot (one of the circles in from the top left) to the bottom right. These trajectories are identical in all four sub-panels. The brightness of the circles indicates the probability of a trajectory at that angle at that location. The angle in the circle is the angle of movement: red is high probability, blue is low probability, grey is a value comparable to the uniform distribution. In each panel, there is a black and a red line: the black line is of length λ while the red line is of length σ. σ and λ were selected by hand for illustrative purposes. σ controls the scale of the basis. The top row of C has a short range basis and as a result the angular distribution tends towards a uniform distribution rapidly as we move away from the blue trajectory. In the bottom row, σ is large and the distribution is non-uniform over much larger ranges. λ controls how rapidly the trajectories transition from parallel movement to perpendicular return to a track. On the left, this scale is short and our model predicts a rapid return to a trajectory. On the right panel, the scale is long and the model predicts movement parallel to the blue trajectory everywhere. For zebra, the parameters are similar to those displayed in the bottom right panel. (D) Position traces of the zebra, and the distribution of predicted angles of movement for the zebra (key as for C).

Data collection

Data were collected between September and November 2015, a period of the dry season when all the zebra were exclusively drinking in one location.

Collar design

RVC collars used solar cells and a rechargeable battery in order to maintain a sufficient collar life span at a 5 min sample rate. All collars were constructed in-house. The collar mass was 930 g plus a 50 g drop off; this equates to only 0.3% of an adult female zebra's body mass (320 kg) (Skinner and Smithers, 1990). The collar circuit was based around a low power MSP430 16-bit micro-controller (Texas Instruments Inc., Dallas, TX, USA), running custom software written in the ‘C’ programming language developed using an integrated development system from IAR Systems.

The microcontroller contains several internal peripheral blocks, including an 8-channel 12-bit analog-to-digital converter (ADC), four serial communications modules, plus various timers, general-purpose digital input and output lines, and other support modules. A connected 2 GB micro-SD flash memory card (Sandisk Corp., Milpitas, CA, USA) provided data storage. GPS position was obtained from a NEO-6T or NEO-M8N GPS module (u-Blox AG); other satellite constellations were not used. In addition to internally computed position and velocity, the module is able to generate raw pseudo-range, phase and Doppler data for the signal from each satellite, enabling detailed GPS performance evaluation, and use of customised differential techniques for increased accuracy.

The resulting GPS tracks are displayed in Fig. 1A. A 2.4 GHz chirp spread spectrum communication module (Nanotron Technologies GmbH, Berlin, Germany) communicating at 1 Mbit per second was used to download data and upload software configurations and firmware updates. A conventional wildlife tracking transmitter in the 149 MHz band (Sirtrack) was used for long-range animal tracking using conventional direction-finding techniques from the ground or the air. Power for the collars was provided by two 13 A h lithium thionyl chloride primary battery (Saft Groupe SA, Bagnolet, France) and a 900 mA h lithium-polymer rechargeable battery (Active Robots, Radstock, UK), charged by a solar cell array consisting of 10 monocrystalline silicon solar cells (Ixys, Milpitas, CA, USA). The switch of the collars electrical load from one battery to the other occurred depending on battery state (voltage over time), monitored by the microcontroller.

The collar provides GPS position using a NEO-6T or NEO-M8N GPS module (u-Blox AG). The NAV-POSLLH data message includes hAcc, a module-derived 2D (i.e. horizontal) accuracy estimate. This corresponds to a 1-sigma value when six or more satellites are tracked and is calculated from the number of satellites tracked, satellite signal to noise data (effectively signal strength), geometric distribution of tracked satellites (DOP value) and the residuals on the pseudorange values. Over the whole dataset, this had a median value of 3.93 m with quartiles of 2.48 m and 6.77 m, so slightly above the data sheet CEP of 2.5 m (NEO-6T) and 2.0 m (Neo-M8N), respectively. As typical for GPS, the position error was highly non-gaussian with a heavy tail. The most extreme observations were removed along with those outside of the study area (see below).

Model approach

The model attempts to predict the heading, estimated from successive GPS fixes, of a zebra as a function of location by assuming that movement in a particular direction in the vicinity of a prior observation increases the probability of later observing a similar heading. The input to the model is each GPS fix along with the heading vector to the subsequent GPS fix. The model initially assumes that the heading vector has a uniform distribution. Then, at every pair of successive GPS observations, a radial basis function is added which increases the likelihood of movement in this direction. As such, the model requires that the observed tracks must be reasonably approximated by linear segments, and the spatial and temporal sample frequency must be reasonably high. The model can be adjusted to change the scale of the radial basis term, which determines over what distance a successive pair of GPS observations are informative, and two parameters which allow for different navigation strategies.

Two model parameters were used to characterise movement: σ and λ. The first, σ, is a measure of how the distance from a previous route predicts future animal movement. The second model parameter, λ, tells us how animals move relative to previous routes. Consider a scenario where we record a zebra's route, and then, a few months later, the zebra is very close to the same location. It is reasonable to expect that accurate predictions about the new route can be derived from the previous route. If on another day, the animal is kilometres away from any previous route, it is reasonable to expect that few, if any, of the previous routes will tell us much about the animal's heading. σ is the distance (m) where the model transitions from one regime (previous tracks are very useful for predicting new tracks) to the other (previous tracks are less useful for predicting tracks). Observations that are separated by distances much greater than σ yield little information about each other's orientation, while observations that are separated by only a fraction of σ are informative of one another's orientation.

For small values of λ, animals are predicted to move strongly towards (converge towards) previous routes. For larger values of λ, the animals are predicted to move parallel to existing routes. Thus, λ controls which strategy zebra use in our model: small values of λ correspond to behaviour consistent with usage of a small number of routes; larger values are consistent with moving to the final goal but ignoring nearby routes that were used previously, i.e. using a more complex mental map of the environment. This model parameter is the length scale (m) where we transition between these two extremes. Observations that are separated by distances much greater than λ will, if they are informative of one another, tend to be perpendicular, while observations that are separated by only a fraction of λ, if informative of one another's orientation, will tend to be parallel.

If animals have preferred routes that they gravitate towards, then we expect λ to be small, in this case compared with σ. Values of λ much larger than σ predict the use of new routes that will be parallel but not co-incident with existing routes, and values between these two extremes predict angles between these two extremes (somewhat attracted towards a previous route, but not directly at it).

Data analysis

We collected GPS data with samples every 5 min from nine zebra over a period of 499 days (1106 zebra days total, minimum 39 days per harem, median 60 days per harem). We filtered the data by removing erroneous positions, those with a horizontal position accuracy value greater than 15 m, and then extracted GPS position (Wilson et al., 2013); subsequent positions were differentiated to obtain trajectory segments – the heading used was the orientation of this vector. Because of this filtering and occasional missed satellite fixes, some small gaps in the data were therefore present. A total of 29 GPS fixes (around 0.37%) were separated by more than two sample intervals. These samples were differentiated as above. No interpolation was performed as only the direction contributed to the fit; this could be obtained even from these more widely separated samples (albeit with less local reliability), and interpolation would have given these outliers undue weight as a result of duplication.

The spatial trajectories were manually sifted through looking for regions where the zebra were travelling between known grazelands and water sources. The start of zebra movements was identified by a speed of greater than 1 m s−1 combined with a sustained directional movement in at least three successive GPS points. The end of movements was identified as either arrival at a water source or when sustained movements slowed to a speed of less than 1 m s−1 combined with high tortuosity, indicating zebras had arrived at a grazing or resting patch. This process delivered 217 journeys, 15–34 per zebra, to and from the grazing areas to the water source across a 15× 15 km area of interest.

We used a set of basis terms (we use the term ‘basis’ informally here as in ‘radial basis functions’; technically, the terms used do not constitute a basis, but rather are a collection of functions whose finite linear spans are dense in the function space of interest), centred on the observations of the training data, to construct a likelihood model for the observed trajectories in the testing data. This can be viewed as an application of the technique of function approximation by radial basis function networks (Broomhead and Lowe, 1988). For every two subsequent observations in the training dataset (where the zebra was observed to move between two locations x0 and x1; see Fig. 2), we have a basis function which calculates the weight, W, associated with moving in a direction θ, at location x, given that we have two successive GPS observations at x0 and x1, with λ and σ parameters as previously discussed:
formula
(1)
Fig. 2.

Plot of successive GPS observations showing the co-ordinates and terms used to construct the radial basis function-like terms used to model the probability distribution of the directionality of the navigating harems. Black crosses are observed positions of a harem, at 5 min intervals, joined by blue lines. The purple line indicates a contour which is equidistant from the closest point on the middle pair of observations in this sequence, x0 and x1. These locations have the same value of the distance, d, which functions as the ‘radius’ in our terms (Eqn 3). The angle of the movement, θ, is that made against the x-axis (horizontal black line), and the line joining the pair of observations.

Fig. 2.

Plot of successive GPS observations showing the co-ordinates and terms used to construct the radial basis function-like terms used to model the probability distribution of the directionality of the navigating harems. Black crosses are observed positions of a harem, at 5 min intervals, joined by blue lines. The purple line indicates a contour which is equidistant from the closest point on the middle pair of observations in this sequence, x0 and x1. These locations have the same value of the distance, d, which functions as the ‘radius’ in our terms (Eqn 3). The angle of the movement, θ, is that made against the x-axis (horizontal black line), and the line joining the pair of observations.

The spatial weight of basis, Wx, is given by:
formula
(2)
where d is given by:
formula
(3)
with:
formula
(4)
and:
formula
(5)
where μ is the projection of the point x onto the line between x0 and x1 such that at x=x0, μ is zero, and at x=x1, μ is 1; u is the distance from x0 to x1.
One can think of d as the shortest distance between the location of our test observation and the line between the two points in our training observation. The weight falls off as a negative squared exponent as this increases. The angular weight is given by:
formula
(6)
where:
formula
(7)
if τ(x0, x1, x, θ, λ)≤b, and:
formula
if τ(x0, x1, x,θ, λ)>b, where:
formula
(8)
where tan−1 is the arc tangent of the vector components with appropriate sign, and:
formula
(9)
formula
(10)
formula
(11)
formula
(12)
and
formula
(13)
where p is the distance of the point from the line between x0 and x1 in the natural units of this line; q is a factor whose sign changes as we cross the line; α is a squared exponential which decays as you move away from the line with a constant offset, the rate of decay being determined by λ; b and c are angular terms which together control whether the weight is directed towards the line between x0 and x1 and ensure that the cosine has a periodic argument in angle as one moves x around the line between x0 and x1; and Rπ/2 is the quarter circle rotation matrix:
formula
(14)
The infinitesimal probability, P(θ|x,λ,σ), to observe a zebra at point x moving at angle θ is given by:
formula
(15)
where we sum over the training data. We can reinterpret this equation as a likelihood, L(λ,σ,θ,x}, and estimate the parameters λ and σ by the process of maximum likelihood; as dθ will be common to all terms, we can drop this factor, giving:
formula
(16)
For a set of test observations enumerated by a, xa, θa, we aim to solve:
formula
(17)
and use bootstrapping to estimate errors on these parameters as reported in the Results section.
This model was compared with a uniform model where the likelihood Lu is given by:
formula
(18)

This likelihood Lu along with the likelihood L were used to calculate the Akaike information criterion (AIC) for the purposes of model comparison. Details of this procedure are given below.

Model fitting and performance

The zebra trajectories were split into testing and training datasets. Our model was fitted to the training data to produce an angular probability field (in the geometric sense) for the motion of the zebra at any point close to an observed point in space. We bootstrapped the fitting procedure for our model by splitting our data into testing and training data sets (150 trajectories were sampled with repetition for the training data, and 67 were sampled with repetition for the testing datasets). This was done 20,000 times.

The model parameters were estimated by maximising the likelihood of the test data using the training data by a simplex minimiser (Nelder and Mead, 1965). The 95% confidence intervals for the parameters λ and σ were estimated, along with the ratio λ/σ. If this ratio is much greater than 1, then the animals are rapidly returning to fixed tracks, if it is much less than 1, then they are following different but similarly orientated tracks as the move (as illustrated in Fig. 1C). We note that this ratio is only capable of capturing the global strategy for isolated tracks. Where tracks intersect, there will be some interference, as we do not know from a single position estimate along which track the animal is moving when they intersect. Crossing tracks could not have contributed significantly to our parameter estimates, however, as crossing tracks correspond, as far as the model is concerned, to returning in a highly directed manner to a track that has been deviated from, and the opposite effect was observed.

For λ and σ, we require a characteristic scale to compare these parameters to. To estimate this scale, we compute the minimum distance between a trajectory and its nearest neighbour and calculate the median of this distance; this typical distance between our trajectories was approximately 1.96 m. If λ and σ are comparable to this scale (a ratio of greater than one-third being reasonable as these are scales for a squared exponential decay), then the majority of our trajectories will be at least partially predicted by our model.

The quality of the fit was also assessed by evaluating an AIC (Akaike, 1974). This was estimated by using the mean parameter values for λ and σ from the bootstrapping process. To estimate our information criterion, we took our 217 trajectories, removed one trajectory and estimated the likelihood of that trajectory given the other trajectories.

This model, fitting procedures and statistics were implemented in the python programming language (Python Software Foundation) using the scipy (Jones et al., 2001) extensions. The source code for these programs is available from GitHub (https://github.com/swilshin/trajbasis).

An example of this fit is provided in Fig. 1D, where the bootstrapped mean values of the model parameters were used and basis functions placed on every example trajectory.

We included 217 journeys from the grazing area towards the waterhole from 9 zebra in the analysis, totalling 931 km of routes (see Fig. 1A and Fig. 3; GPS data for these journeys can be found in Table S1). The median journey length was 3.96 km and typical tortuosity was 1.14, indicating that the trajectories used were fairly straight (a combined histogram, kernel density plot and box plot of the path lengths and tortuousities is included in Fig. 4). The median minimal route separation across all route pairings, calculated from the distance of closest approach for each pair of routes and taking the median of these values, was 1.96 m. This median separation remained small, 2.39 m, even if a 1 km area around the point of convergence was removed from the observations.

Fig. 3.

Plot of the routes between grazing areas and the water source at different spatial scales. Each harem is shown in a separate colour. (A,B) Journeys starting at the water source. (C,D) Return trips. The red square in A and C shown on an enlarged scale in B and D, respectively.

Fig. 3.

Plot of the routes between grazing areas and the water source at different spatial scales. Each harem is shown in a separate colour. (A,B) Journeys starting at the water source. (C,D) Return trips. The red square in A and C shown on an enlarged scale in B and D, respectively.

Fig. 4.

Length and tortuosity of recorded zebra movements between their grazelands and water. (A) Top: plot of the frequency of tortuosities (blue is a density-normalised histogram and black kernel density plot with Gaussian kernel with bandwidth estimated by Scott's rule). Bottom: boxplot of tortuosities of the zebra movements, operationalised by the arc-chord ratio; that is, the ratio of the length of the curve, measured by taking the sum of the distances between successive GPS observations, and the distance between the start and end of the curve, measured by taking the distance between the last and first GPS observation. This ratio cannot be less than one, and indicates how indirect the route taken was. (B) Histogram and boxplot of the frequency of path length, which is the sum of the distances between successive GPS observations. Subjects generally took a reasonably direct route between water and grazelands, and vice versa, but many highly indirect routes were taken; this is reflected in Fig. 1.

Fig. 4.

Length and tortuosity of recorded zebra movements between their grazelands and water. (A) Top: plot of the frequency of tortuosities (blue is a density-normalised histogram and black kernel density plot with Gaussian kernel with bandwidth estimated by Scott's rule). Bottom: boxplot of tortuosities of the zebra movements, operationalised by the arc-chord ratio; that is, the ratio of the length of the curve, measured by taking the sum of the distances between successive GPS observations, and the distance between the start and end of the curve, measured by taking the distance between the last and first GPS observation. This ratio cannot be less than one, and indicates how indirect the route taken was. (B) Histogram and boxplot of the frequency of path length, which is the sum of the distances between successive GPS observations. Subjects generally took a reasonably direct route between water and grazelands, and vice versa, but many highly indirect routes were taken; this is reflected in Fig. 1.

The bootstrapped 95% confidence intervals on σ, λ and the ratio σ/λ were 1.19–26.4 m, 68.4–1.02×108 m and 1.20×10−7–0.162, respectively. These results are consistent with σ (a measure of how large an area around a previous route we can use to predict future animal movement) being comparable to the characteristic scale of the distance between neighbouring trajectories (1.96 m). Therefore, a route that is separated from a second route by this characteristic scale can be used to infer how the zebra will move along this second route. As the model includes a squared exponential even at a location several multiples of σ from a route, we can estimate how a zebra will move based on that route. The ratio implies that the zebra do not gravitate towards particular paths but rather follow a number of tracks as they travel rather than preferring any particular route. This is visually illustrated in Fig. 1C,D, where the blue line represents a zebra route under four different model scenarios. The colour intensity in each circle around the line indicates travel direction probability, red being a high probability of travelling in that direction and blue being a low probability. The direction of movement is more likely to be towards the route when λ is low, as seen in the left hand panels of Fig. 1C, whilst the range at which routes can influence the direction of new routes increases with σ, as seen in the bottom two panels of Fig. 1C. Our results are most similar to those in the bottom right panel, with a high probability of moving parallel to the existing track, as shown for true routes in Fig. 1D.

The AIC for our model was −26,000, while the AIC of a uniform model was −32,000. Our model is therefore vastly superior to a uniform model of the zebra trajectories. The same is implied by more conservative methods of model comparison such as the Schwarz Bayesian information criterion (also −26,000 and −32,000 for the two models to three significant digits) (Schwarz, 1978).

The most probable heading direction determined from the model matched the observed trajectories more reliably that those from the uniform model; the average absolute error was 17.6 deg, while the average absolute error from the uniform model was 45 deg, showing that majority of zebra routes went in the same direction as previously used routes.

We hypothesised that zebra utilise multiple routes to reach their destination, but that routes are highly predicted by other nearby routes, as a result of zebra using knowledge gained from previous journeys to the destination to aid the efficiency of subsequent journeys. Our findings supported this hypothesis: zebra did not repeatedly use the same routes when moving through their environment; instead, they used a series of routes that shared directional properties with previous routes. Even when routes were close together, less than 2 m at their closest points, they did not converge onto each other.

Many animals, including humans, prefer to repeatedly use the same routes when travelling through their environment. Route use can reduce energetic costs as cost of transport (cost to move 1 kg of body mass a distance of 1 m) varies substantially with substrate: in humans, walking cost increased 2.5 times when walking on sand versus solid ground (Lejeune et al., 1998), and the repeated use of the same route creates physical trails with a denser substrate that is free of vegetation (Shepard et al., 2013). Route use that results in the creation of well-defined trails can also simplify navigation, reducing landscape complexity by reducing the number of navigational decisions from a step time scale to a junction time scale (Newmark and Rickart, 2012). The creation of physical trails through repeated use should be particularly beneficial to animals in habitats with loose substrates or those that are hard to travel over, such as sand or snow. However, route use could also increase journey length, especially if an animal only uses a small number of known routes to move between different locations, somewhat like a motorway network, and is not able to make novel shortcuts. So, whilst potentially energetically and navigationally beneficial, the strategy that the zebra in this study utilise of not following a small number of specific routes does have potential advantages: time and energy are not wasted walking extra distances to join a known route, navigation is more resilient to environmental perturbations (for instance, the loss of a landmark through environmental change) and prey species (such as zebra) are less predictable in their movement patterns, making it harder for an ambush hunting strategy to be effective.

The area that the zebra in this study moved through is covered in large numbers of game trails of diverse orientations (see Fig. 1 for an aerial 3D photogrammetric survey of such trails), yet all zebra routes were highly directed towards their end point. The level of directedness suggests that zebra have a good knowledge about their spatial environment and are not misled by trails in a different orientation. The method utilised by zebra to achieve these highly directed yet variable tracks is unknown. However, the lack of local topography, with an overall height variation of less than 2 m across the 15 km square (McCarthy and Ellery, 1998), combined with distance travelled, eliminates the use of one single visual beacon close to the destination. An olfactory or auditory beacon would be less impacted by lack of topography and as such may allow for beaconing from a greater distance. The use of such stimuli in spatial navigation is well documented in rodents (Lavenex and Schenk, 1998) and in small-scale experiments with domestic ungulates (Edwards et al., 1997) but has not been documented across a large spatial scale in wild-ranging ungulates.

We note that as σ is bounded from below the confidence interval, it can never include zero. As such, standard hypothesis testing cannot be applied to determine whether this parameters is substantially different from zero. However, we also note that λ is extremely large, especially when compared with σ. The proposed model, when scored on an AIC, is vastly superior to that of the uniform model with only these two parameters. This implies not only that the zebras are not converging on specific routes but also that the model where λ approaches infinity is a reasonable proxy for our model. As such, the large drop in the AIC is mostly due to the effect of the σ term. The confidence interval for σ is 1.19–26.4 m, which suggests that we should be able to look at the observed movements of the zebra and identify multiple, similarly orientated tracks separated by approximately this distance to a small multiple of this distance (the squared exponential still has a weight of around 2% at 3 times the characteristic scale), around 10–30 m. We note that these movements, if they converge, do so very slowly over large distances. This is exactly what we would expect to see given the parameter values observed for the model. It should be mentioned that the heading is derived from fixes 5 min apart, so it is the overall chosen direction not instantaneous heading at that point which would be disrupted (e.g. walking round a bush).

It is possible the zebra achieve this directional movement by knowing the relative location of various features or other cues within their landscape and use this knowledge to select a trail that leads to the desired destination. It is also possible they use some form of innate navigation, such as the position of the sun, to identify a trail that is orientated in the correct direction. Either strategy would require a zebra to make a navigational decision at each trail intersection to ensure they remain on an efficient route and have knowledge of position heading along with a spatial map.

Perturbation studies where zebra were intentionally driven off known tracks (emulating, for instance, a predator or group of animals) might reveal how this highly directed movement was achieved. Such a study would also serve to test the predictions of this model about how the zebra should move in response to such a perturbation when in transit. Specifically, zebra should, if moved sufficiently far from their original route/trail, continue along another similarly orientated route or trail, rather than return to their original track.

We inferred the zebra’s navigation strategy by constructing a model of the direction in which zebra move using a method similar to radial basis function networks (Broomhead and Lowe, 1988). Such an approach has several advantages. First, it places two navigation hypotheses in direct opposition, permitting us to determine which is a better explanation of the observed behaviour. Second, this model could easily be modified with additional parameters, which would allow for wider application and permit the model to capture more variability. For example, in a multi-species study, the model could include a factor for species, or if spatial inhomogeneity were suspected, perhaps as a result of variable terrain, then the model could include a terrain factor to account for this. The model does, however, have several important limitations. It models the heading of the animal but not the speed, and it requires that the proposed navigation hypothesis be describable in spatially extended terms. This would make incorporating navigation strategies such as beaconing more difficult (though not necessarily impossible, because beacons could be estimated and spatially located).

This model could be further refined to permit other factors, beyond spatial location, that may affect animal movement. One obvious extension would be to include a state variable for the level of hydration of the subject or when it last drank. The model predicts with roughly equal probability a trajectory at each location with the subject moving either towards the grazing lands or towards the water source. If the level of hydration of the zebra were known, these two could be disambiguated.

We thank J. Lowe and S. Amos for their contribution to the development of the technology and field assistants M. J. Nkape, N. Terry and M. Classe for downloading data from the collars.

Author contributions

Conceptualization: H.B.-B., S.D.W., T.Y.H., A.M.W.; Methodology: H.B.-B., S.D.W., T.Y.H., A.M.W.; Software: S.D.W.; Formal analysis: H.B.-B., S.D.W., T.Y.H., S.H.; Investigation: H.B.-B., S.D.W., T.Y.H., E.B., A.M.W.; Resources: A.M.W.; Data curation: H.B.-B., S.D.W., E.B.; Writing - original draft: H.B.-B., S.D.W.; Writing - review & editing: H.B.-B., S.D.W., T.Y.H., S.H., E.B., A.M.W.; Supervision: S.H., A.M.W.; Project administration: A.M.W.; Funding acquisition: A.M.W.

Funding

We extend our thanks to the EPSRC (EP/H013016/1), BBSRC (BB/J018007/1) and ERC (323041) for funding.

Data availability

Code for model, fitting procedures and statistics is available from GitHub: https://github.com/swilshin/trajbasis.

Agnew
,
A. D. Q.
(
1966
).
The use of game trails as a possible measure of habitat utilisation by larger mammals
.
Afr. J. Ecol.
4
,
38
-
46
.
Akaike
,
H.
(
1974
).
A new look at the statistical model identification
.
IEEE Trans. Autom. Control
19
,
716
-
723
.
Bartlam-Brooks
,
H. L. A.
,
Bonyongo
,
M. C.
and
Harris
,
S.
(
2011
).
Will reconnecting ecosystems allow long-distance mammal migrations to resume? A case study of a zebra Equus burchelli migration in botswana
.
Oryx
45
,
210
-
216
.
Beecham
,
J. A.
(
2001
).
Towards a cognitive niche: divergent foraging strategies resulting from limited cognitive ability of foraging herbivores in a spatially complex environment
.
Biosystems
61
,
55
-
68
.
Bennett
,
A. T.
(
1996
).
Do animals have cognitive maps?
J. Exp. Biol.
199
,
219
-
224
.
Bonyongo
,
M. C.
and
Harris
,
S.
(
2007
).
Grazers species-packing in the okavango delta, botswana
.
Afr. J. Ecol.
45
,
527
-
534
.
Broomhead
,
D. S.
and
Lowe
,
D.
(
1988
).
Multivariate functional interpolation and adaptive networks
.
Complex Syst.
2
,
321
-
355
.
Byrne
,
R. W.
(
1979
).
Memory for urban geography
.
Q J. Exp. Psychol. (Colchester)
31
,
147
-
154
.
Byrne
,
R. W.
(
2000
).
Evolution of primate cognition
.
Cogn. Sci.
24
,
543
-
570
.
Collett
,
M.
,
Chittka
,
L.
and
Collett
,
T. S.
(
2013
).
Spatial memory in insect navigation
.
Curr. Biol.
23
,
R789
-
R800
.
Croney
,
C. C.
,
Adams
,
K. M.
,
Washington
,
C. G.
and
Stricklin
,
W. R.
(
2003
).
A note on visual, olfactory and spatial cue use in foraging behavior of pigs: indirectly assessing cognitive abilities
.
Appl. Anim. Behav. Sci.
83
,
303
-
308
.
de Guinea
,
M.
,
Estrada
,
A.
,
Nekaris
,
K. A.-I.
and
Van Belle
,
S.
(
2019
).
Arboreal route navigation in a neotropical mammal: energetic implications associated with tree monitoring and landscape attributes
.
Mov. Ecol.
7
,
39
.
Di Fiore
,
A.
and
Suarez
,
S. A.
(
2007
).
Route-based travel and shared routes in sympatric spider and woolly monkeys: cognitive and evolutionary implications
.
Anim. Cogn.
10
,
317
-
329
.
Dumont
,
B.
and
Petit
,
M.
(
1998
).
Spatial memory of sheep at pasture
.
Appl. Anim. Behav. Sci.
60
,
43
-
53
.
Edwards
,
G. R.
,
Newman
,
J. A.
,
Parsons
,
A. J.
and
Krebs
,
J. R.
(
1997
).
Use of cues by grazing animals to locate food patches: an example with sheep
.
Appl. Anim. Behav. Sci.
51
,
59
-
68
.
Geva-Sagiv
,
M.
,
Las
,
L.
,
Yovel
,
Y.
and
Ulanovsky
,
N.
(
2015
).
Spatial cognition in bats and rats: from sensory acquisition to multiscale maps and navigation
.
Nat. Rev. Neurosci.
16
,
94
-
108
.
Haynes
,
G.
(
2006
).
Mammoth landscapes: good country for hunter-gatherers
.
Quatern. Int.
142-143
,
20
-
29
.
Hirvonen
,
H.
,
Ranta
,
E.
,
Rita
,
H.
and
Peuhkuri
,
N.
(
1999
).
Significance of memory properties in prey choice decisions
.
Ecol. Model.
115
,
177
-
189
.
Janson
,
C.
(
2000
).
Spatial movement strategies: theory, evidence, and challenges
. In
On the Move: How and Why Animals Travel in Groups
(ed.
S.
Boinski
and
P. A.
Garber
),
pp. 165-203
.
Chicago
:
University of Chicago Press
.
Jones
,
E.
,
Oliphant
,
T.
and
Peterson
,
P.
(
2001
).
SciPy: Open Source Scientifica Tools for Python
.
Kays
,
R.
,
Tilak
,
S.
,
Kranstauber
,
B.
,
Jansen
,
P. A.
,
Carbone
,
C.
,
Rowcliffe
,
M.
,
Fountain
,
T.
,
Eggert
,
J.
and
He
,
Z.
(
2011
).
Monitoring wild animal communities with arrays of motion sensitive camera traps
.
Int. J. Res. Rev. Wireless Sensor Netw.
1
,
19
-
29
.
Lavenex
,
P.
and
Schenk
,
F.
(
1998
).
Olfactory traces and spatial learning in rats
.
Anim. Behav.
56
,
1129
-
1136
.
Lejeune
,
T. M.
,
Willems
,
P. A.
and
Heglund
,
N. C.
(
1998
).
Mechanics and energetics of human locomotion on sand
.
J. Exp. Biol.
201
,
2071
-
2080
.
Lohmann
,
K. J.
,
Lohmann
,
C. M. F.
and
Putman
,
N. F.
(
2007
).
Magnetic maps in animals: nature's GPS
.
J. Exp. Biol.
210
,
3697
-
3705
.
Mccarthy
,
T. S.
and
Ellery
,
W. N.
(
1998
).
The okavango delta
.
Trans. R. Soc. S. Afr.
53
,
157
-
182
.
Milton
,
K.
(
1980
).
The Foraging Strategy of Howler Monkeys: a Study in Primate Economics
.
Columbia University Press
.
Morris
,
R. G. M.
,
Garrud
,
P.
,
Rawlins
,
J. N. P.
and
O'keefe
,
J.
(
1982
).
Place navigation impaired in rats with hippocampal lesions
.
Nature
297
,
681
-
683
.
Mueller
,
T.
and
Fagan
,
W. F.
(
2008
).
Search and navigation in dynamic environments – from individual behaviors to population distributions
.
Oikos
117
,
654
-
664
.
Nelder
,
J. A.
and
Mead
,
R.
(
1965
).
A simplex method for function minimization
.
Comput. J.
7
,
308
-
313
.
Newmark
,
W. D.
and
Rickart
,
E. A.
(
2012
).
High-use movement pathways and habitat selection by ungulates
.
Mamm. Biol.
77
,
293
-
298
.
O'keefe
,
J.
and
Nadel
,
L.
(
1978
).
The Hippocampus as a Cognitive Map
.
Oxford
:
University Press
.
Poucet
,
B.
(
1993
).
Spatial cognitive maps in animals: new hypotheses on their structure and neural mechanisms
.
Psychol. Rev.
100
,
163
-
182
.
Ramberg
,
L.
,
Hancock
,
P.
,
Lindholm
,
M.
,
Meyer
,
T.
,
Ringrose
,
S.
,
Sliva
,
J.
,
Van As
,
J.
and
Vander Post
,
C.
(
2006
).
Species diversity of the okavango delta, botswana
.
Aquat. Sci.
68
,
310
-
337
.
Schwarz
,
G.
(
1978
).
Estimating the dimension of a model
.
Ann. Stat.
6
,
461
-
464
.
Shepard
,
E. L. C.
,
Wilson
,
R. P.
,
Rees
,
W. G.
,
Grundy
,
E.
,
Lambertucci
,
S. A.
and
Vosper
,
S. B.
(
2013
).
Energy landscapes shape animal movement ecology
.
Am. Nat.
182
,
298
-
312
.
Skinner
,
J. D.
and
Smithers
,
R. H. N.
(
1990
).
The Mammals of the Southern African Subregion
.
Pretoria
:
University Of Pretoria
.
Smuts
,
G. L.
(
1976
).
Population characteristics of Burchell's zebra (Equus burchelli antiquorum, H. Smith, 1841) in the kruger national park
.
S. Afr. J. Wildl. Res.
6
,
99
-
112
.
Tsoar
,
A.
,
Nathan
,
R.
,
Bartan
,
Y.
,
Vyssotski
,
A.
,
Dell'Omo
,
G.
and
Ulanovsky
,
N.
(
2011
).
Large-scale navigational map in a mammal
.
Proc. Natl Acad. Sci. USA
108
,
E718
-
E724
.
Ulanovsky
,
N.
and
Moss
,
C. F.
(
2008
).
What the bat's voice tells the bat's brain
.
Proc. Natl Acad. Sci. USA
105
,
8491
-
8498
.
Wiens
,
J. A.
(
1976
).
Population responses to patchy environments
.
Annu. Rev. Ecol. Syst.
7
,
81
-
120
.
Williams
,
T. C.
,
Williams
,
J. M.
and
Griffin
,
D. R.
(
1966
).
The homing ability of the neotropical bat Phyllostomus hastatus, with evidence for visual orientation
.
Anim. Behav.
14
,
468
-
473
.
Wilson
,
A. M.
,
Lowe
,
J. C.
,
Roskilly
,
K.
,
Hudson
,
P. E.
,
Golabek
,
K. A.
and
Mcnutt
,
J. W.
(
2013
).
Locomotion dynamics of hunting in wild cheetahs
.
Nature
498
,
185
-
189
.
Wilson
,
A. M.
,
Hubel
,
T. Y.
,
Wilshin
,
S. D.
,
Lowe
,
J. C.
,
Lorenc
,
M.
,
Dewhirst
,
O. P.
,
Bartlam-Brooks
,
H. L. A.
,
Diack
,
R.
,
Bennitt
,
E.
,
Golabek
,
K. A.
, et al. 
(
2018
).
Biomechanics of predator–prey arms race in lion, zebra, cheetah and impala
.
Nature
554
,
183
-
188
.
Zollner
,
P. A.
and
Lima
,
S. L.
(
1999
).
Search strategies for landscape-level interpatch movements
.
Ecology
80
,
1019
-
1030
.

Competing interests

The authors declare no competing or financial interests.

Supplementary information