Quantifying spatial and temporal patterns of prey searching is of primary importance for understanding animals' critical habitat and foraging specialization. In patchy environments, animals forage by exhibiting movement patterns consisting of area-restricted searching (ARS) at various scales. Here, we present a new method, the fractal landscape method, which describes the peaks and valleys of fractal dimension along the animal path. We describe and test the method on simulated tracks, and quantify the effect of track inaccuracies. We show that the ARS zones correspond to the peaks from this fractal landscape and that the method is near error-free when analyzing high-resolution tracks, such as those obtained using the Global Positioning System (GPS). When we used tracks of lower resolution, such as those obtained with the Argos system, 9.6–16.3% of ARS were not identified, and 1–25% of the ARS were found erroneously. The later type of error can be partially flagged and corrected. In addition, track inaccuracies erroneously increased the measured ARS size by a factor of 1.2 to 2.2. Regardless, the majority of the times and locations were correctly flagged as being in or out of ARS (from 83.8 to 89.5% depending on track quality). The method provides a significant new tool for studies of animals' foraging behavior and habitat selection, because it provides a method to precisely quantify each ARS separately, which is not possible with existing methods.

The movement patterns of foraging animals are affected by the availability and distribution of resources, as well as the animals' ability to forage coupled with the level of physiological need(Bell, 1991). Ultimately,efficient searching behavior results in higher energy gains relative to the time and/or energy expended (Bell,1991; Pyke et al.,1977). Therefore, greater understandings of the factors that shape foraging behavior are crucial to evaluate ecological specialization and adaptation.

The use of electronic devices to track animal movement has provided detailed information on the movement patterns and prey-searching behavior of a wide variety of organisms (Austin et al.,2004; Block et al.,2002; Costa, 1993; Hays et al., 2004). As top predators forage in a patchy, hierarchical environment(Fauchald, 1999), they display alterations in traveling speed and direction between different behavioral phases in the track. These phases are scale dependent(Fauchald et al., 2000), and result in hierarchical area-restricted searching (ARS) pattern of movement. This pattern has been shown to be a response to patchy resource distribution(Benhamou, 1992; Fauchald and Tveraa, 2003; Lode, 2000). To investigate foraging behavior in a patchy environment it is necessary to independently and accurately describe each ARS zone along a track.

Several techniques have been developed to measure the scales at which ARS occurs (Robinson et al.,2007). One of the simplest methods is to measure changes in transit speed along a track (LeBoeuf et al., 2000). The assumption here is that intensive searching is associated with reduced speed; other approaches measure the changes in sinuosity and/or angularity of the track(Erlandsson and Kostylev, 1995; Laidre et al., 2004). However,none of these methods provide an index of the scale or domain of the region of ARS. Several approaches have been developed to measure the region or spatial scale of the ARS. First-passage time (FPT) provides a measure of the time an animal takes to cross a circle of a certain radius that is moved along the track (Fauchald and Tveraa,2003). The circle radius corresponding to the highest variance in FPT is used to define the animal's operational spatial scale. One type of fractal analysis uses a similar approach, but uses fractal dimension (D) (or fractal dimension estimator) instead of time(Nams, 1996). A segment of a given length is moved along the track, and D is calculated for each segment. As D generally increases with track convolutions, the segment lengths corresponding to highest average D, and/or highest variance in D, are used as a cut-off to identify the operational spatial scales(Nams, 1996; Nams, 2005). These two methods are similar in principle, but produce different results. FPT analysis identifies a radius (Pinaud and Weimerskirch, 2005) whereas the above fractal analysis identifies a segment length (Fritz et al.,2003; Nams, 2005). Intense searching behavior is often a combination of longer time spent and higher track convolutions in an ARS zone, which in theory should be identified by both methods. However, these methods have fundamental differences. One focuses on time whereas the other focuses on space coverage, both in an exclusive way. Therefore, the two methods are not really measuring the same thing (Robinson et al., 2007)and the use of one or the other has to be driven by underlying questions and/or theory about ARS behavior. Both methods may accurately detect the major ARS zones but they may differ in describing them, and they have a major flaw:they both identify scales that are roughly constant throughout the track(Fauchald and Tveraa, 2003; Fritz et al., 2003; Pinaud and Weimerskirch,2005). In order to finely describe ARS behavior, it is important to identify the scale and position of each ARS zone individually.

The present study describes and tests a new method for quantifying ARS behavior in tracking data, and tests the effect of the track quality on the efficiency of the method.

Theory and analytical strategy

Imagine you lost your watch in a field. To search for it, you zigzag here and there, but you do not move through the field along a straight path. Unconsciously, your behavior will aim at covering the space you are searching in, not just turning more. This is an important consideration because track analyses often involve turning rate and sinuosity indexes(Erlandsson and Kostylev, 1995; Laidre et al., 2004), but none focused strictly on space coverage. The fractal dimension of a set of two-dimensional points can be seen as a measure of its propensity to cover the plane, being a value of one for no plane coverage (a straight line or a circle, for example) and two for full coverage of some area in the plane (a Hilbert curve, for example). The fractal dimension does not measure track sinuosity per se, and, if sinuosity is the question of interest, some other sinuosity indexes should be used instead(Benhamou, 2004). When analyzing searching behavior it makes sense to quantify how well an animal covers the region or plane of interest both spatially and temporally: for example, the greatest chance of recovering your watch requires that you maximize the amount of area searched within the plane and this will require more time. If your watch is not precious, you will give up earlier and perhaps not cover the plane as well. Since plane coverage and time are two fundamental ingredients for searching behavior, we argue that methods of track analysis should quantify both components of searching behavior: time and space coverage.

Based on these considerations, we develop a new method for detecting and quantifying searching behavior focusing on space coverage and time. The underlying idea is that animals should increase their local plane coverage(i.e. fractal dimension of the track) when they exhibit ARS behavior. After testing the accuracy of the method on simulated tracks of known characteristics, we tested the effect of track inaccuracy on the process by altering the perfect simulated tracks in a way that resembles Argos track qualities, and reanalyzed these tracks to determine the effect of track quality.

Argos track qualities vary considerably with the type of animal tracked(Tremblay et al., 2006). Consequently, differences in track qualities can affect apparent plane coverage (and sinuosity) and/or spatial distribution of presence (i.e. time),and it is therefore important to quantify these effects. We therefore simulated albatross-like and elephant seal-like tracks. These two species were chosen because they represent two extremes in animal speed and in Argos track qualities (Tremblay et al.,2006).

Note that we focused here on the description of small-scale ARS, where animals are presumably searching for prey. At a larger scale, the animal is probably searching for suitable areas that include small-scale ARS(Pinaud and Weimerskirch,2005).

Track simulations

Tracks were simulated by concatenating portions of tracks determined by correlated random walks (CRW) with different degrees of correlation between successive steps. The degree of the correlation between successive angles was controlled by the standard deviation (s.d.) used to generate the circular normal distribution from which the turning angles were randomly selected. When s.d. is small, the path is straighter, and when s.d. is large, the track approaches an uncorrelated random walk (i.e. Brownian motion). The process essentially follows the simulation process described by Fauchald and Tveraa,with some adjustments (Fauchald and Tveraa, 2003).

Adjustments include the following:

  1. Speed was not considered constant between steps; instead it was selected from truncated normal distributions of known mean and s.d., based on real data(Table 1).

  2. Tracks were built by concatenating three kinds of phases to form a type of three-state Markov model for the track: (i) very straight, directed movements at the beginning and at the end of the track to simulate large-scale dispersal of the animal towards a hypothetical foraging zone. These portions had a lower turning rate and higher speed; (ii) ARS movements, with a higher turning rate and slower speeds; and (iii) between straight and ARS portions and between two ARS portions, movements of intermediate turning rates and intermediate speed were used. This aims to simulate animal searching at larger spatial scales.

  3. The time spent in each ARS movement was controlled, and the ARS movement was restricted to a circle of known size. The circle's size did not necessarily correspond to the ARS size because we did not force the ARS movement to fill the circle. Therefore, the circle dimensions represent the maximum possible size of the ARS. Circles were filled only when the time spent in them was sufficient. If not, the ARS sizes were necessarily smaller than the corresponding circle and were not known precisely. Although this could be seen as a major flaw, we argue that forcing the ARS to describe a circle of know dimension is much more problematic. In such a case we would generate what the method is aimed at describing (not what the animal actually does), leading us to artificially enhance the apparent performance of the method. Circles are used as an approximation of the ARS size and shape, and were determined using a statistical method (detailed below).

Table 1.

Parameters used to simulate animal tracks containing area-restricted search (ARS) patterns

ParameterCharacteristicCalculation for Albatross-like tracksCalculation for Elephant seal-like tracks
Interval between steps Constant 5 min 5 min 
Average animal speed (km h-1Random (truncated normal) 30±20 (0.1-100) 3.5±2.3 (0.1-7) 
Speed during directed movements Average speed + random change Random change = 10-30% (uniform) Random change = 10-30% (uniform) 
Sinuosity during directed movements Random (circular normal). Mean=0 Random s.d. = 0.005-0.010 (uniform) Random s.d. = 0.005-0.010 (uniform) 
Speed during intermediate movements Random (truncated normal) Same as average speed (above) Same as average speed (above) 
Sinuosity during intermediate movements Random (circular normal). Mean=0 Random s.d. = 0.04-0.06 (uniform) Random s.d. = 0.04-0.06 (uniform) 
Speed during ARS movements Average speed - random change Random change = 10-30% (uniform) Random change = 10-30% (uniform) 
Sinuosity during ARS movements Random (circular normal). Mean=0 Random s.d. = 0.5-0.7 (uniform) Random s.d. = 0.5-0.7 (uniform) 
Number of ARS in the track Random (uniform) 1-20 1-20 
Maximum size of ARS (radius) Random (uniform) 1-60 km 1-60 km 
Time spent in ARS Random (uniform) 3-56 h 10-480 h 
ParameterCharacteristicCalculation for Albatross-like tracksCalculation for Elephant seal-like tracks
Interval between steps Constant 5 min 5 min 
Average animal speed (km h-1Random (truncated normal) 30±20 (0.1-100) 3.5±2.3 (0.1-7) 
Speed during directed movements Average speed + random change Random change = 10-30% (uniform) Random change = 10-30% (uniform) 
Sinuosity during directed movements Random (circular normal). Mean=0 Random s.d. = 0.005-0.010 (uniform) Random s.d. = 0.005-0.010 (uniform) 
Speed during intermediate movements Random (truncated normal) Same as average speed (above) Same as average speed (above) 
Sinuosity during intermediate movements Random (circular normal). Mean=0 Random s.d. = 0.04-0.06 (uniform) Random s.d. = 0.04-0.06 (uniform) 
Speed during ARS movements Average speed - random change Random change = 10-30% (uniform) Random change = 10-30% (uniform) 
Sinuosity during ARS movements Random (circular normal). Mean=0 Random s.d. = 0.5-0.7 (uniform) Random s.d. = 0.5-0.7 (uniform) 
Number of ARS in the track Random (uniform) 1-20 1-20 
Maximum size of ARS (radius) Random (uniform) 1-60 km 1-60 km 
Time spent in ARS Random (uniform) 3-56 h 10-480 h 

When randomly selected values from a normal distribution were used, the mean and s.d. of the distribution is given. Minimum and maximum values are given in parentheses. The terms `normal' and `uniform' refer to the type of distribution.

The parameters were chosen for simulating each type of track and for each portion according to Table 1. An example of a simulated track is given in Fig. 1A. The time between each step was 5 min.

Alteration of tracks

Intact simulated tracks were considered to be near error-free and similar to Global Positioning System (GPS) tracks. These tracks were then altered spatially and temporally, in order to introduce the same kind of error typically obtained within albatross and elephant seal Argos tracks.

The number of locations per day used to sub-sample each track was randomly selected from a truncated normal distribution centered on 17.770 and 5.125 with a s.d. of 4.450 and 2.826 in albatross and elephant seal tracks,respectively. Minima were set to six and three locations per day and maxima were set to 26 and 16 locations per day, for albatross and elephant seals,respectively. These values were based on data collected for the TOPP program during 2002–2005(http://www.toppcensus.org/).

Similarly, we introduced spatial error to each selected point(Fig. 1B). For this, a quality class was randomly attributed to each selected point, in order to match mean proportions of the various quality classes in real deployments(Table 2). An error was then randomly selected from a normal distribution centered on the mean error and with the same s.d. as the errors for the corresponding class(Table 2), as obtained during static tests at the laboratory [similar results to those found in tests by Vincent et al. (Vincent et al.,2002)].

Table 2.

Parameters used to spatially alter the simulated tracks

Argos classMean accuracy (km)s.d. accuracy (km)Albatross tracks (%)Elephant seal tracks (%)
0.178 0.386 1.150 2.149 
0.372 0.433 3.620 3.799 
0.770 0.958 11.054 5.893 
2.539 6.549 47.159 9.028 
0.720 2.523 16.274 23.006 
5.692 41.414 17.601 48.210 
205.943 382.994 3.142 7.914 
Argos classMean accuracy (km)s.d. accuracy (km)Albatross tracks (%)Elephant seal tracks (%)
0.178 0.386 1.150 2.149 
0.372 0.433 3.620 3.799 
0.770 0.958 11.054 5.893 
2.539 6.549 47.159 9.028 
0.720 2.523 16.274 23.006 
5.692 41.414 17.601 48.210 
205.943 382.994 3.142 7.914 

Accuracy data were obtained using our own static tests (Y.T., A.J.R. and D.P.C., unpublished). Percentages were calculated using 24 564 hits from 153 deployments on albatrosses and 78 565 hits from 106 deployments on elephant seals.

The resulting tracks were then considered to be similar to Argos tracks,and were therefore filtered and further interpolated along a Bezier curve(μ=0.2), following Tremblay et al.(Tremblay et al., 2006)(Fig. 1C).

The filtering process involved a speed filter set at 80 and 10 km h–1 for albatross and elephant seal tracks, respectively. An angle filter set at 170° was also used in both types of track, to reject small-scale spikes not removed by the speed filter. In addition, when locations were closer than 10 min, the location with the least quality was removed irrespective of its quality class. This allows for the removal of induced sinuosity, created by location proximity, and avoids taking into account locations that are calculated using common hits from the transmitters(in real data). All quality classes were allowed to be removed.

Fig. 1.

Example of a simulated track of an elephant seal. (A) A Global Positioning System (GPS)-like track containing two area-restricted searching (ARS) zones.(B) This shows the track after the introduction of spatial error and temporal sub-sampling, as obtained in real deployments using the Argos tracking technique. (C) The same track after the filtering and the interpolation process (see Materials and methods for details).

Fig. 1.

Example of a simulated track of an elephant seal. (A) A Global Positioning System (GPS)-like track containing two area-restricted searching (ARS) zones.(B) This shows the track after the introduction of spatial error and temporal sub-sampling, as obtained in real deployments using the Argos tracking technique. (C) The same track after the filtering and the interpolation process (see Materials and methods for details).

Calculating fractal dimension

A fractal dimension (D) can be calculated from a set of points using several distinct methods that differ in accuracy, in the sensitivity to the number of points used and computing time required(Esteller et al., 1999; Kallimanis et al., 2002; Nams, 2006). It is therefore essential to carefully select the appropriate algorithm for a given application (Esteller et al.,1999; Jelinek et al.,1998). In this study, fractal D was calculated using the information dimension (Halley et al.,2004) as it is less prone to errors in finite datasets (A. J. Roberts: http://www.sci.usq.edu.au/staff/robertsa/soft.html). We also tested another algorithm implementing the box-counting method to calculate fractal D (Giorgilli et al.,1986; Liebovitch and Toth,1989). The results were noticeably less accurate (data not shown).

The accuracy of our calculation was assessed using five mathematical curves of known fractal dimension (Fig. 2): a Hilbert curve composed of 4096 points (D=2), a Sierpinski triangle of 2409 points (D=1.58), a Koch's snowflake of 3075 points (D=1.26),a circle of 2000 points (D=1) and a straight line of 2000 points (D=1). For these curves there is no difference between the information dimension and any other fractal dimension. The information dimensions D calculated on these sets were 2.01, 1.56, 1.29, 1.05 and 1.05, respectively, corresponding to a mean error of 3%.

To evaluate the effect of the number of input points on fractal D estimates, each curve was randomly sub-sampled by a different number of points and the fractal D was calculated on each subset. Since each dataset was altered, its apparent fractal dimension was also modified(Fig. 2), but the relative order of the determined fractal D of the sets was apparent with at least 50 points, and unequivocal when using 100 points or more. In our case study, 100 points represented retention of 2.4–20% of the total number of data points. This shows that the algorithm detected the differences in relative complexity between the sets even when using a relatively small proportion of the data. The estimate of fractal D stabilizes and the s.d. is acceptable even when a minimum of ∼100 points are used(Fig. 2). The stability of the fractal D estimates, with respect to the number of points, reflects the sensitivity of the shape of the curve to sub-sampling. As an illustration,consider the decrease in estimated fractal D of the circle in Fig. 2: the shape of a circle is indeed affected by the number of points that defines it, being here a polyhedron of number_of_points –1 number of edges. The higher the numbers of points, the closer to the real shape of the circle and, therefore,the more the estimated fractal D tends to 1. The fact that fractal D is <1 for the straight line data set is normal, because, in this case, it is a measure of line filling (D≤1) rather than plane filling (1≤D≤2).

Fig. 2.

Test of the accuracy and the sensitivity to sub-sampling of the algorithm used to calculate fractal dimension. Tests were run on five datasets of known theoretical fractal dimensions. The thick line is a smoothed representation(moving average) of the thin line.

Fig. 2.

Test of the accuracy and the sensitivity to sub-sampling of the algorithm used to calculate fractal dimension. Tests were run on five datasets of known theoretical fractal dimensions. The thick line is a smoothed representation(moving average) of the thin line.

A fractal dimension does not characterize a dataset at some specific length scale; instead it characterizes a self-similar pattern that occurs in the data across many length scales. Thus, basic fractal dimension estimates come from fitting a straight line to data on a log-log plot where one axis covers a range of length scales. As part of a least-squares fit, the algorithm not only determines the two parameters of the best fractal similarity, but also the two length scales between which the fractal similarity holds. For the data in this work, we found that fractal self-similarity holds over an order of magnitude in scale. Such a decade of length scales is acceptable to claim fractality(see Feder, 1988). In our application, it is necessary to note that consistent estimates of fractal dimension (i.e. higher values for higher plane coverage) are more important than strictly accurate estimates. Therefore, if the self-similarity in our tracks is different to the mathematical curves, the accuracy of fractal dimension may change slightly, but this is not detrimental because consistency is maintained (see below).

Implementation of the fractal landscape method

The fractal landscape method consisted of four distinct steps. Step 1: we determined a particular segment length (i.e. a distance) for the track at which the largest most common ARS occurred. Since speed is reduced in ARS, we determined the segment length as follows. First, the rank of speed between successive steps was computed. Second, the portion of the consecutive lowest third of ranks was determined. Third, the distance traveled during these portions of relative low speed was calculated. Finally, we used the longest of these distances as the segment length. This segment length is used as the window size for the next step.

Step 2: we calculated fractal D along the track. For each point, D was calculated using the locations belonging to the segment centered on each point(the window size). When the number of points was more than 200, we randomly selected 200 points for the given window. This considerably reduced the computing time while preserving reasonable accuracy(Fig. 2), ensuring that sufficient data were used in each D estimate. The plot of D against time describes fractal peaks and valleys, i.e. the fractal landscape(Fig. 3A).

Step 3: we automatically detected a threshold value that distinguished relevant fractal peaks from fractal valleys. In order to distinguish a minor oscillation from a real fractal peak, we found the threshold value over which one oscillation in fractal D could be considered a peak. For this, we counted the number of peaks above each potential threshold value. When the potential threshold value was small and within the baseline, the number of peaks above it increased dramatically (Fig. 3B). This point was chosen as the threshold cut-off.

Fig. 3.

Graphic output of the fractal landscape analysis. (A) The fractal landscape: seven peaks in fractal dimension are distinguishable. (B) This illustrates the way the threshold for determining the peaks is determined. The blue line is a smoothing (moving average) of the calculated number of peaks over the threshold (gray line). (C) This shows the corresponding track, with the seven zones corresponding to the seven fractal peaks highlighted and numbered. The color of the circles delimiting the area-restricted search zones is scaled based on the area over the threshold of the corresponding fractal peak (green shade in A). This is considered as an index of searching intensity.

Fig. 3.

Graphic output of the fractal landscape analysis. (A) The fractal landscape: seven peaks in fractal dimension are distinguishable. (B) This illustrates the way the threshold for determining the peaks is determined. The blue line is a smoothing (moving average) of the calculated number of peaks over the threshold (gray line). (C) This shows the corresponding track, with the seven zones corresponding to the seven fractal peaks highlighted and numbered. The color of the circles delimiting the area-restricted search zones is scaled based on the area over the threshold of the corresponding fractal peak (green shade in A). This is considered as an index of searching intensity.

Step 4: based on this threshold, we then extracted fractal peak characteristics, which were used to visualize and quantify ARS behavior. For each peak we extracted its duration, area above the threshold value, mean fractal D, and the circle position and size where each fractal peak was overlaid on the mapped track (Fig. 3C). The zone corresponding to each fractal peak is determined by a circle centered on its mean coordinates. Standard deviations of the coordinates were calculated for the x and y axes. The radius of the circle was then calculated as three times the smaller of the two s.d. This reduced the effects of outliers, and kept the majority of ARS coordinates within the circle. High-intensity-searching behavior can be seen as the most efficient method of covering the localized circle (high fractal D) over a long time interval and is therefore defined by a high area of fractal peak. This area was thus calculated and used as an index of searching intensity, and color coded for visualization in Fig. 3C.

Analysis of intact simulated tracks

As expected, fractal dimensions increase when points from the track segments include points from an ARS. A plot of fractal D against time clearly reveals peaks and valleys, i.e. a fractal landscape(Fig. 3A). For each fractal peak, we extracted the corresponding time and position data from the track. The positions of the peaks corresponded very closely with the most tortuous portions of the track (black sections of track in Fig. 3C). Overall, 98.7% of the 456 simulated ARS zones were correctly identified(Fig. 4A). We examined the few ARS that were not correctly identified, and found that they were all concatenated with another ARS because the time (and distance) between them was very short. Consequently, 100% of the ARS found by the method overlapped with at least one simulated ARS. As expected, sizes of the ARS were found to be mostly inferior to the size of the circle in which the simulated ARS was forced to be (Fig. 5), and visual inspection shows that indeed, the calculated circle size matched the actual size of the ARS (Fig. 5). Approximately 98% (38.8+59) and 99% (37.5+62) of the track locations (i.e. time) were correctly identified as being part of an ARS or not, in albatross-like and elephant seal-like tracks, respectively(Table 3). This shows that the fractal landscape method is near error-free when used with near-perfect tracks(e.g. using GPS).

Fig. 4.

Performance of the fractal landscape method in the detection of area-restricted searching (ARS) in intact simulated tracks (Global Positioning System (GPS)-like tracks (A) and in corrected tracks (Argos-like tracks) (B). Open circles, albatross; blue triangles, elephant seal-like tracks.

Fig. 4.

Performance of the fractal landscape method in the detection of area-restricted searching (ARS) in intact simulated tracks (Global Positioning System (GPS)-like tracks (A) and in corrected tracks (Argos-like tracks) (B). Open circles, albatross; blue triangles, elephant seal-like tracks.

Table 3.

Percentages of track locations in relation to initial status [in area-restricted search (ARS) versus not in ARS] and status found by the fractal landscape method

Results from the fractal landscape method
Intact tracks
Degraded tracks
Location statusIn ARSNot in ARSIn ARSNot in ARS
Albatross-like simulated tracks (N=10) In ARS 38.8±8.2 1.5±1.7 30.8±7.2 9.5±5.7 
 Not in ARS 0.7±0.6 59.0±7.1 1.0±0.7 58.7±6.3 
Elephant seal-like simulated tracks (N=10) In ARS 37.5±11.3 0.2±0.1 26.1±8.8 11.6±9.7 
 Not in ARS 0.3±0.5 62.0±11.4 4.7±5.5 57.7±11.9 
Results from the fractal landscape method
Intact tracks
Degraded tracks
Location statusIn ARSNot in ARSIn ARSNot in ARS
Albatross-like simulated tracks (N=10) In ARS 38.8±8.2 1.5±1.7 30.8±7.2 9.5±5.7 
 Not in ARS 0.7±0.6 59.0±7.1 1.0±0.7 58.7±6.3 
Elephant seal-like simulated tracks (N=10) In ARS 37.5±11.3 0.2±0.1 26.1±8.8 11.6±9.7 
 Not in ARS 0.3±0.5 62.0±11.4 4.7±5.5 57.7±11.9 

Note that the number of location is strictly equivalent to time, because track locations were equally spaced in time.

Fig. 5.

Zoomed portion of a simulated track showing the circles in which the area-restricted searching (ARS) was constrained to stay within (gray circles,see Materials and methods for details) and the circles corresponding to the output of the fractal landscape method (yellow circles). The yellow circles are visually very close to the spatial extent of the simulated ARS.

Fig. 5.

Zoomed portion of a simulated track showing the circles in which the area-restricted searching (ARS) was constrained to stay within (gray circles,see Materials and methods for details) and the circles corresponding to the output of the fractal landscape method (yellow circles). The yellow circles are visually very close to the spatial extent of the simulated ARS.

Analysis of corrected simulated tracks

The introduction of Argos-like errors into the tracks resulted in two different effects. Some ARS were lost (i.e. existent but not found), and some ARS were produced (i.e. found but non-existent). As a result, the number of ARS regions was more variable when Argos tracks were used(Fig. 4B).

In albatross-like degraded tracks, 9.6% of the track's ARS were lost whereas 90.4% appropriately identified areas of actual ARS. In elephant seal-like degraded tracks, these values were 16.3 and 83.7%, respectively. In addition, the number of produced ARS regions was 1 out of 105 and 26 out of 108 in albatross- and elephant seal-like tracks, respectively. As a result of the differences in track qualities between albatross- and elephant seal-like tracks, elephant seal tracks had more lost and produced ARS regions (i.e. more errors) than in albatross-like tracks. However, some of these errors can be easily flagged: produced ARS were generally short in time and intensity, with 65.6% of them belonging to the lowest quintile of ARS intensities (area under the peak, see Materials and methods). The deletion of all ARS regions where the search-intensity index falls within the lowest quintile can remove many of the wrong ARS regions, while losing only a few of the good but not intense ARS regions, in elephant seal Argos tracks.

Nevertheless, the proportion of locations (i.e. time) correctly identified as in or out of an ARS was still reasonably high(Table 3; 89.5% in albatrosses and 83.8% in elephant seals). Note that Argos inaccuracies caused a disproportionate error towards lost ARS rather than produced ARS in both kinds of tracks (Table 3). Logically,ARS are missing because of temporal degradation, whereas they are produced through spatial degradation introducing more spatial coverage. This indicates that the biggest problem with using Argos data is not the spatial inaccuracy but the relatively low sampling rate. In this respect, recently developed methods for filtering Argos data without deleting points(Jonsen et al., 2005) should considerably improve the situation, since it is common to delete a third to a half of the raw Argos data with classical filters(Austin et al., 2003; McConnell et al., 1992).

However, differences in spatial accuracies between albatross-like and elephant seal-like corrected tracks had a visible effect: 4.7 times more locations were incorrectly attributed to belonging to an ARS in elephant seal-compared with albatross-like tracks. This was the largest difference in the results between the two kinds of tracks(Table 3).

Spatial accuracy also affected ARS size. Not surprisingly, ARS radii are larger in degraded tracks than in intact tracks, and this difference was higher in elephant seal-than in albatross-like tracks(Fig. 6). We determined that ARS radii in corrected tracks could be corrected by a factor of 1.2 and 2.2 in albatross- and elephant seal-like tracks, respectively. Dividing the radii by these factors allowed us to match the distribution of the ARS radii obtained with intact tracks (Fig. 6). Note that as we chose two extreme cases in terms of Argos track qualities, it is expected that correction factors applied to datasets of intermediate qualities should necessarily be intermediate as well.

Fig. 6.

Distribution of the sizes of the area-restricted searching (ARS) calculated by the fractal landscape method in albatross- (A) and elephant seal-like tracks (B). Kernel densities are shown instead of histograms for clarity. The blue, solid line is the distribution obtained with intact tracks [Global Positioning System (GPS)-like tracks] and the red, dotted line was obtained with corrected tracks (Argos-like tracks). The shift towards larger ARS is a consequence of the Argos spatial inaccuracies. A correction factor applied to these data improved results (black, broken line; see Results and discussion for details).

Fig. 6.

Distribution of the sizes of the area-restricted searching (ARS) calculated by the fractal landscape method in albatross- (A) and elephant seal-like tracks (B). Kernel densities are shown instead of histograms for clarity. The blue, solid line is the distribution obtained with intact tracks [Global Positioning System (GPS)-like tracks] and the red, dotted line was obtained with corrected tracks (Argos-like tracks). The shift towards larger ARS is a consequence of the Argos spatial inaccuracies. A correction factor applied to these data improved results (black, broken line; see Results and discussion for details).

Does the method work on real data?

We applied the method on real non-duty-cycled Argos tracks of 31 albatrosses and 20 elephant seals and the output was consistent with results from this study. In the Laysan albatross (Phoebastria immutabilis)track shown in Fig. 7, nine peaks were clearly identified in the fractal landscape, each peak corresponding to an ARS on the track. Interestingly, the ARS regions with the highest intensity index (highest area under the peak, see Materials and methods) were distributed along a west–east transect whereas the lowest indexes were found in ARS regions north and south of this line(Fig. 7). This line is situated approximately 1700–2000 km north of the breeding colony (Tern island,Hawaii), and corresponds with the position of the food-rich North Pacific Transition Zone between sub-Arctic and subtropical waters at this time period(Hyrenbach et al., 2002; Polovina et al., 2000). This result indicates that more intensive ARS behavior occurred when closer to the front, and shows the utility of this index in identifying biologically important regions.

How sensitive is the method to the window size?

We ran a sensitivity analysis to investigate how the segment length used for the estimation of fractal dimension affects the output of the method. We ran the analysis 60 times on one elephant seal-like track containing seven ARS(30 times for both the intact and the corrected track versions). For each iteration, the segment length was changed from 20–600 km by steps of 20 km. The segment length had little effect on the ARS detection ability of the method. The range over which 6–8 ARS were detected was between 40–500 and 100–600 km in intact and corrected tracks, respectively(Fig. 8). This indicates that the segment length can be fixed for a set of tracks from different animals,and does not need to be calculated for each individual track.

Technical problems encountered and how they were solved

When we implemented the method, the first problem was the large amount of computing time initially required. This was solved as described in the methodology, by sub-sampling the data if more than 200 points were available. The change in computing time was important (empirically), and a 6-month track interpolated every 5 min can be processed in just a couple of minutes with a standard laptop.

Later, we realized that degrading the track affects the speed calculation along the track, and this sometimes posed problems for the calculation of the window's size (step 1 of the method). This was solved by smoothing the speed data using a classic moving average filter. The window used for the filter was simply chosen by visually interpreting the plot of speed against time, in a way that maintains the major oscillations and removed high-frequency noise. The process was always obvious and straightforward to achieve and therefore we did not automate it.

Fig. 7.

Example of the fractal landscape analysis output on a real track of Laysan albatross (Phoebastria immutabilis). The inset figure shows the entire track. The color of the area-restricted searched (ARS) zones corresponds to a searching-intensity index (see Materials and methods, and Fig. 3). The ARS with the highest intensities are spread along a west-east transect corresponding to a food-rich frontal zone (see Results and discussion). The broken perimeter corresponds to the area in which the small-scale searching (the ARS) occurred. It is proposed as a possible (but non-exclusive) method to extrapolate the second (larger) scale at which the individual operated (see Results and discussion for details).

Fig. 7.

Example of the fractal landscape analysis output on a real track of Laysan albatross (Phoebastria immutabilis). The inset figure shows the entire track. The color of the area-restricted searched (ARS) zones corresponds to a searching-intensity index (see Materials and methods, and Fig. 3). The ARS with the highest intensities are spread along a west-east transect corresponding to a food-rich frontal zone (see Results and discussion). The broken perimeter corresponds to the area in which the small-scale searching (the ARS) occurred. It is proposed as a possible (but non-exclusive) method to extrapolate the second (larger) scale at which the individual operated (see Results and discussion for details).

A similar problem arose with the high-frequency noise in the consecutive fractal D calculations (quick oscillations in the fractal landscape), and in the determination of the threshold. The problem was solved in the same manner,by using a moving average smoothing.

Comments on the method

Concerns have been raised about the use of fractal geometry in biological sciences (Halley et al., 2004; Turchin, 1996) and, in particular, the use of fractal dimension (D) to estimate track sinuosity(Benhamou, 2004). Some criticism results from the confusion over what the fractal dimension, or a fractal dimension analysis, measures biologically(Jelinek et al., 1998). In the analysis of animal distribution there are several ways of applying fractal geometry concepts (Laidre et al.,2004; Nams, 1996; Russell et al., 1992; With, 1994). A common mistake is that D measures sinuosity (Benhamou,2004; Laidre et al.,2004; Nams, 1996). The track of an animal that turns 180° at every movement step will show very high sinuosity but minimum D. Similarly, a track describing a circle and a track describing a straight line both have a low D value (namely one) but will have different sinuosity. A high D value will only result when a track's convolutions lead to reasonably efficient coverage of an area in the plane. In our study, D should be considered as an area-filling index, and is therefore particularly suitable for the analysis of searching behavior. Fractal dimension is classically used in studies interested in plane or volume filling measurements (Chmel et al.,2005; Kim and Kim,2005; Phattaralerphong and Sinoquet, 2005; Uttieri et al., 2005). A recent study using fractal dimension in three-dimensions highlighted the subtle differences between degree of convolution and fractal dimension (Uttieri et al., 2005).

Fig. 8.

Number of area-restricted searching (ARS) detected using the fractal landscape method using different segment lengths for the calculation of the fractal dimensions along the track.

Fig. 8.

Number of area-restricted searching (ARS) detected using the fractal landscape method using different segment lengths for the calculation of the fractal dimensions along the track.

Calculating a fractal dimension for non-fractal objects is also problematic. In our case, the range of length scales over which our calculation was done was always >10, which indicates that our segments were near-fractal. Although this was true in simulated intact tracks, this was artificial in corrected tracks and resulted from the interpolation process. The interpolation process was, however, crucial in order to obtain accurate estimates of D.

The main advantage of the fractal landscape method over others is that it provides a means to find and precisely describe each ARS zone separately. The method is an application of fractal geometry rather than a fractal analysis per se. It permits measurement of the variance in the characteristics of ARS zones, which is very informative and not possible with the previously available methods. Our method also produced results that were spatially explicit and quantitatively straightforward to interpret.

The principle of the method was focused on the animal's efficiency in covering some area of the plane. The places where this efficiency was relatively higher in the track were considered to be part of an ARS. It is important to note that this is different to focusing on the time the animal spends in some zones, as this time can be spent immobile or in rectilinear movement. The relationship between these two aspects can be easily confused as they are potentially linked. However, the link goes only one way: high plane coverage efficiency involves more time, but the reverse is not true. The reverse statement can, however, appear true when too much random inaccuracy affects the track. Together with our results, this indicates that the method's interest increases with the tracks' quality because the ability to discriminate between the effects of time and plane coverage increases. With the increasing usage of GPS tags, this method could potentially help in understanding subtle aspects of foraging behavior, such as the temporal optimization of searching efficiency (or how to search more space in a minimum amount of time). With lower quality tracks, the same can be achieved, but only at larger scales, for which track inaccuracy is negligible.

Comments on the concept of scale

Our method provides an efficient approach to finding and describing regions of small-scale ARS but not in determining other, larger scales of movement. By contrast, the first-passage time method and the traditional fractal analysis method are more designed to identifying average scales at which animals modify their behavior. In both methods, this is done in different ways, using and outputting different measures. Nobody can tell which is more appropriate,because there is no consensus on what really is a scale or an ARS. Therefore,fundamental to this discussion is the question: what are the appropriate scales of measurement of an animal track, and in which units should they be measured?

Using a circle to approximate the size of an ARS region appears to be an acceptable measure because, by definition, it is of small scale and it generally looks circular. However, when an animal forages at larger scales, a circle is generally not an appropriate descriptor of the region of movement(Pinaud and Weimerskirch,2005) and thus does not accurately describe the animal's movement. The segment length used in classical fractal analysis appears more appropriate in this respect. However, animal tracking data are by nature two-dimensional and because filling the searched space as efficiently as possible is a key component of searching behavior, we believe that the most relevant metric of ARS is the amount of coverage of a plane (i.e. an area) (taking the earlier example of the lost watch, this makes intuitive sense). For ARS, the area of the circle could be used as an estimate of the smallest scale. For the larger scale, techniques used in estimates of territory size and home range could be used, such as the convex polygon method (example given in Fig. 7) or a more sophisticated kernel analysis that describes the distribution of ARS or track location. This would allow accurate description of areas that match the shape of the animal space-use pattern at all scales. Nonetheless, the fractal landscape method derived here provides a measure of the ARS of small-scale events with good precision, and with the added advantage of providing a measure of the time spent searching within that localized region. Depending on the questions of interest, different metrics could potentially be used to create a similar landscape and to find associated zones (speed, sinuosity index, etc). For example, as speed is often related to energy expense, a plot of speed against time could be used to determine the portion of lower versus higher speeds (i.e. potential energy expense). The principle is simple, intuitive and different metrics can be adapted to the question initially asked.

This research was part of the Tagging of Pacific Pelagics (TOPP) program,funded by the National Ocean Partnership Program (N00014-02-1-1012), the Office of Naval Research (N00014-00-1-0880 and N00014-03-1-0651) and the Alfred P. Sloan, the Packard and The Gordon and Betty Moore Foundations.

Austin, D., McMillan, J. I. and Bowen, W. D.(
2003
). A three-stage algorithm for filtering erroneous Argos satellite locations.
Mar. Mamm. Sci.
19
,
371
-383.
Austin, D., Bowen, W. D. and McMillan, J. I.(
2004
). Intraspecific variation in movement patterns: modeling individual behaviour in a large marine predator.
Oikos
105
,
15
-30.
Bell, W. J. (
1991
).
Searching Behaviour: The Behavioural Ecology of Finding Resources
. London:Chapman & Hall.
Benhamou, S. (
1992
). Efficiency of area-concentrated searching behavior in a continuous patchy environment.
J. Theor. Biol.
159
,
67
-81.
Benhamou, S. (
2004
). How to reliably estimate the tortuosity of an animal's path: straightness, sinuosity, or fractal dimension?
J. Theor. Biol.
229
,
209
-220.
Block, B. A., Costa, D. P., Boehlert, G. W. and Kochevar, R. E. (
2002
). Revealing pelagic habitat use: the tagging of Pacific pelagics program.
Oceanol. Acta
25
,
255
-266.
Chmel, A., Smirnov, V. N. and Astakhov, M. P.(
2005
). The Arctic sea-ice cover: fractal space-time domain.
Physica A
357
,
556
-564.
Costa, D. P. (
1993
). The secret life of marine mammals: novel tools for studying their behavior and biology at sea.
Oceanography
6
,
120
-128.
Erlandsson, J. and Kostylev, V. (
1995
). Trail following, speed and fractal dimension of movement in a marine prosobranch, Littorina littorea, during a mating and a nonmating season.
Mar. Biol.
122
,
87
-94.
Esteller, R., Vachtsevanos, G., Echauz, J. and Litt, B.(
1999
). A comparison of fractal dimension algorithms using synthetic and experimental data. In
Proceedings of the IEEE International Symposium on Circuits and Systems: Adaptive Digital Signal Processing
. Vol.
3
, pp.
199
-202. Orlando, FL, USA.
Fauchald, P. (
1999
). Foraging in a hierarchical patch system.
Am. Nat.
153
,
603
-613.
Fauchald, P. and Tveraa, T. (
2003
). Using first-passage time in the analysis of area-restricted search and habitat selection.
Ecology
84
,
282
-288.
Fauchald, P., Erikstad, K. E. and Skarsfjord, H.(
2000
). Scale-dependent predator-prey interactions: the hierarchical spatial distribution of seabirds and prey.
Ecology
81
,
773
-783.
Feder, J. (
1988
).
Fractals
. New York: Plenum Press.
Fritz, H., Said, S. and Weimerskirch, H.(
2003
). Scale-dependent hierarchical adjustments of movement patterns in a long-range foraging seabird.
Proc. R. Soc. Lond. B Biol. Sci.
270
,
1143
-1148.
Giorgilli, A., Casati, D., Sironi, L. and Galgani, L.(
1986
). An efficient procedure to compute fractal dimensions by box counting.
Phys. Lett. A
115
,
202
-206.
Halley, J. M., Hartley, S., Kallimanis, A. S., Kunin, W. E.,Lennon, J. J. and Sgardelis, S. P. (
2004
). Uses and abuses of fractal methodology in ecology.
Ecol. Lett.
7
,
254
-271.
Hays, G. C., Houghton, J. D. R. and Myers, A. E.(
2004
). Pan-Atlantic leatherback turtle movements.
Nature
429
,
522
.
Hyrenbach, K. D., Fernández, P. and Anderson, D. J.(
2002
). Oceanographic habitats of two sympatric North Pacific albatrosses during the breeding season.
Mar. Ecol. Prog. Ser.
233
,
283
-301.
Jelinek, H. F., Jones, C. L. and Warfel, M. D. (1998). Is there meaning in fractal analysis? Complexity Int. 6, http://www.complexity.org.au/ci/vol06/jelinek/jelinek.html.
Jonsen, I. D., Flenming, J. M. and Myers, R. A.(
2005
). Robust state-space modeling of animal movement data.
Ecology
86
,
2874
-2880.
Kallimanis, A. S., Sgardelis, S. P. and Halley, J. M.(
2002
). Accuracy of fractal dimension estimates for small samples of ecological distributions.
Landscape Ecol.
17
,
281
-297.
Kim, H. J. and Kim, W. H. (
2005
). Automatic detection of spiculated masses using fractal analysis in digital mammography.
Lect. Notes Comput. Sci
.
3691
,
256
-263.
La Boeuf, B. J., Crocker, D. E., Costa, D. P., Blackwell, S. B.,Webb, P. M. and Houser, D. S. (
2000
). Foraging ecology of northern elephant seals.
Ecol. Monog.
70
,
353
-382.
Laidre, K. L., Heide-Jorgensen, M. P., Logsdon, M. L., Hobbs, R. C., Dietz, R. and VanBlaricom, G. R. (
2004
). Fractal analysis of narwhal space use patterns.
Zoology
107
,
3
-11.
Liebovitch, L. S. and Toth, T. (
1989
). A fast algorithm to determine fractal dimensions by box counting.
Phys. Lett. A
141
,
386
-390.
Lode, T. (
2000
). Functional response and area-restricted search in a predator: seasonal exploitation of anurans by the European polecat, Mustela putorius.
Austral Ecol.
25
,
223
-231.
McConnell, B. J., Chambers, C. and Fedak, M. A.(
1992
). Foraging ecology of southern elephant seals in relation to the bathymetry and productivity of the Southern Ocean.
Antarct. Sci.
4
,
393
-398.
Nams, V. O. (
1996
). The VFractal: a new estimator for fractal dimension of animal movement paths.
Landscape Ecol.
11
,
289
-297.
Nams, V. O. (
2005
). Using animal movement paths to measure response to spatial scale.
Oecologia
143
,
179
-188.
Nams, V. O. (
2006
). Improving accuracy and precision in estimating fractal dimension of animal movement paths.
Acta Biotheor.
54
,
1
-11.
Phattaralerphong, J. and Sinoquet, H. (
2005
). A method for 3D reconstruction of tree crown volume from photographs: assessment with 3D-digitized plants.
Tree Physiol.
25
,
1229
-1242.
Pinaud, D. and Weimerskirch, H. (
2005
). Scale-dependent habitat use in a long-ranging central place predator.
J. Anim. Ecol.
74
,
852
-863.
Polovina, J. J., Kobayashi, D. R., Parker, D. M., Seki, M. P. and Balazs, G. H. (
2000
). Turtles on the edge: movement of loggerhead turtles (Caretta caretta) along oceanic fronts, spanning longline fishing grounds in the central North Pacific, 1997-1998.
Fish. Oceanogr.
9
,
71
-82.
Pyke, G. H., Pulliam, H. R. and Charnov, E. L.(
1977
). Optimal foraging: a selective review of theory and tests.
Q. Rev. Biol.
52
,
137
-154.
Robinson, P. W., Tremblay, Y., Antolos, M., Crocker, D. E.,Kuhn, C. E., Shaffer, S. A., Simmons, S. and Costa, D. P.(
2007
). A comparison of indirect measures of feeding behavior based on ARGOS tracking data.
Deep Sea Res. Part II Top. Stud. Oceanogr
. In Press.
Russell, R. W., Hunt, G. L., Coyle, K. O. and Cooney, R. T.(
1992
). Foraging in a fractal environment: Spatial patterns in a marine predator-prey system.
Landscape Ecol.
7
,
195
-209.
Tremblay, Y., Shaffer, S. A., Fowler, S. L., Kuhn, C. E.,McDonald, B. I., Weise, M. J., Bost, C. A., Weimerskirch, H., Crocker, D. E.,Goebel, M. E. et al. (
2006
). Interpolation of animal tracking data in a fluid environment.
J. Exp. Biol.
209
,
128
-140.
Turchin, P. (
1996
). Fractal analyses of animal movement: a critique.
Ecology
77
,
2086
-2090.
Uttieri, M., Zambianchi, E., Strickler, J. R. and Mazzocchi, M. G. (
2005
). Fractal characterization of three-dimensional zooplankton swimming trajectories.
Ecol. Modell.
185
,
51
-63.
Vincent, C., McConnell, B. J., Ridoux, V. and Fedak, M. A.(
2002
). Assessment of Argos location accuracy from satellite tags deployed on captive gray seals.
Mar. Mamm. Sci.
18
,
156
-166.
With, K. A. (
1994
). Using fractal analysis to assess how species perceive landscape structure.
Landscape Ecol.
9
,
25
-36.