## Abstract

The Lévy-flight foraging hypothesis states that because Lévy flights can optimize search efficiencies, natural selection should have led to adaptations for Lévy flight foraging. Some of the strongest evidence for this hypothesis has come from telemetry data for sharks, bony fish, sea turtles and penguins. Here, I show that the programming for these Lévy movement patterns does not need to be very sophisticated or clever on the predator's part, as these movement patterns would arise naturally if the predators change their direction of travel only after encountering patches of relatively strong turbulence (a seemingly natural response to buffeting). This is established with the aid of kinematic simulations of three-dimensional turbulence. Lévy flights movement patterns are predicted to arise in all but the most quiescent of oceanic waters.

## 1. Introduction

Lévy flights, named after the French mathematician Paul Lévy, arose in a purely mathematical context in the first half of the past century [1]. They first entered the biological literature when Shlesinger & Klafter [2] proposed that they can be observed in the movement patterns of foraging ants. Lévy flights comprise clusters of short step-lengths with longer movements between them. This pattern is repeated across all scales with the resultant clusters creating fractal patterns. The hallmark of a Lévy flight is a step-length distribution with a heavy power-law tail, *p*_{l}(*l*)∼*l*^{−μ}, where 1<*μ*<3. The divergence of the second-moment is indicative of the absence of a characteristic scale and the associated formal divergence of the diffusion constant is indicative of ‘superdiffusion’. It has long been recognized that the fractal and superdiffusive properties of Lévy flights can be advantageous when searching for randomly distributed resources [3]. This theoretical realization has led to the ‘Lévy-flight foraging hypothesis’ [4]. This hypothesis states that because Lévy flights can optimize search efficiencies, natural selection should have led to adaptations for Lévy flight foraging. The evolution of optimal search patterns is predicted, because natural selection will favour individuals that are best able to find resources critical to survival. Initial evidence for Lévy flights in wandering albatrosses [5] together with the realization that these movement patterns can be advantageous in random search scenarios [3] led to an explosion of interest in Lévy flights as models of animal movement patterns [4,6]. Lévy flights became contentious when it was found that they had been wrongly attributed to the wandering albatross and to many other species through the employment of inappropriate statistical techniques and through misinterpretations of the data [7,8]. More recent studies have however provided seemingly compelling evidence that many organisms, including T cells, *Escherichia coli*, mussels, honeybees, extinct marine creatures and even human hunter–gatherers, have movement patterns that can be approximated by Lévy flights [9–15].

Some of the strongest evidence for the Lévy-flight foraging hypothesis has come from telemetry data for marine predators such sharks, bony fish, sea turtles and penguins [16,17]. Nonetheless, the key to prediction and understanding lies in the elucidation of mechanisms underlying the observed patterns [18]. ‘Without an understanding of mechanisms, one must evaluate the stress on each new system de novo, without any scientific basis for extrapolation; with such an understanding, one has the foundation for understanding’ [18]. This sentiment was recently echoed by Stumpf & Porter [19], who rightly noted that ‘a statistically sound power-law is no evidence of universality without a concrete underlying “generative mechanism” to support it’.

Here, I show that the programming for Lévy movement patterns in marine predators does not need to be very sophisticated or clever on the predator's part, as these movement patterns would arise naturally if the predators change their direction of travel only after encountering patches of relatively strong turbulent velocity fluctuations (a seemingly natural response to buffeting). This is done with the aid of kinematic simulations of turbulence [20]. Lévy movement patterns are also predicted to arise in marine predators that change their direction of travel on encountering large fluctuations in the pressure gradients (i.e. large turbulent accelerations).

The analysis stems from Kolmogorov's similarity theory of turbulence [21]. Kolmogorov similarity theory has provided a cornerstone for the interpretation of oceanic spectral measurements of turbulence and some of its predictions have been convincingly verified by observations [22]. According to this theory, the turbulent velocities, *u*(*t*), of the water in the immediate vicinity of a moving predator are characterized by 〈Δ*u*^{2}(Δ*t*)〉∝*ε*^{2/3}*V* ^{2/3}Δ*t*^{2/3}, where *t* is time, *ε* is the mean rate of dissipation of turbulent kinetic, *V* is the speed of the predator relative to the water and Δ*u*(Δ*t*)=*u*(*t*+Δ*t*)−*u*(*t*). Directly analogous scaling was obtained by Fung *et al.* [23] for inertial particles falling through turbulent flows. The scaling 〈Δ*u*^{2}(Δ*t*)〉∝Δ*t*^{2/3} can also be derived directly from the observed form of the energy spectrum without invoking Kolmogorov similarity theory. This is because the energy spectrum exhibits Kolmogorov −5/3 power-law scaling over an extended range of scales [22], and so its Fourier transform, i.e. 〈Δ*u*^{2}(Δ*t*)〉, will be proportional to Δ*t*^{2/3} over an associated range of scales [24]. It follows from 〈Δ*u*^{2}(Δ*t*)〉∝Δ*t*^{2/3} that the distribution of time-intervals between consecutive bouts of strong velocity fluctuations will have a power-law tail [25,26]. The characteristic power-law exponent is dependent on the underling dynamics. Predators that move with near constant speed and only change their direction of travel after encountering strong turbulence are therefore expected to have movement patterns that can be approximated by Lévy flights. For directly analogous reasons, predators that change their direction of travel on encountering large fluctuations in pressure gradients are also predicted to have movement patterns that can be approximated by Lévy flights. This is because Kolmogorov similarity theory predicts that the pressure gradients in the vicinity of a moving predator scale in the same way as turbulent velocity fluctuations, namely as 〈|Δ*p*(Δ*t*)|〉∝*ε*^{2/3}*V* ^{2/3}Δ*t*^{2/3}.

Here, I test these predictions with the aid of numerical simulations of turbulent flows. Turbulence can be calculated from first principles without approximation in direct numerical simulations, but these simulations become computationally prohibitive for Reynolds numbers, *Re*∼10^{5}. This is restrictive, because oceanic flows can have much higher Reynolds numbers; the gulf stream, for example, has *Re*=*uL*/*ν*∼10^{11}(*u*∼1 ms^{−1}, *L*∼10^{5} m and *ν*∼10^{−6} m^{2} s^{−1}). This restriction on Reynolds number can be alleviated by computing solutions to the dynamical equations only for the large scales and by modelling of smaller scales. These large-eddy simulations, although very promising, are also computationally very expensive. Here, instead, the velocity field of three-dimensional, incompressible homogeneous, isotropic turbulence is modelled as a sum of independent unsteady random Fourier modes that vary in space and time and over different realizations [20]. The modes are chosen so that spectra have the Kolmogorov form which pertains to turbulence in the ‘inertial subrange’ that is predicted to give rise Lévy flights, i.e. pertains to turbulent structures having sizes, *l*, in the range *L*≫*l*≫*l*_{η}, where *l*_{η}=*L*/*Re*^{3/4} is the scale on and below which dissipation occurs. The sizes of marine predators (sharks, bony fish, sea turtles and penguins) displaying Lévy flight movement characteristics fall into this range of scales. The structure of the simulated velocity field is found to be similar to that computed by direct numerical simulations with the same spectrum, but the flow fields do not necessarily satisfy the dynamical equations nor do they have all the known statistical distributions [20]. Nonetheless, particle tracking and pressure statistics agree with experimental measurements and this suggests that these detailed aspects of the flow are not so important when modelling many aspects of turbulent phenomena [20]. Notably, structure functions, 〈Δ*u*^{2}(Δ*t*)〉, characterizing fluid-velocities experienced by particles along their trajectories are consistent with Kolmogorov similarity theory and with observations [20]. Nonetheless, the absence of inertial subrange intermittency makes incremental changes in Eulerian and Lagrangian (tracer-particle) velocities resolutely Gaussian, in contrast with laboratory results which exhibit highly non-Gaussian statistics [20]. The same shortcoming can be expected for fluid velocities along the trajectories of heavy particles, and this may impact on model predictions for Lévy exponents [25,26].

## 2. Methods

In the numerical simulations, a model predator swims in a straightline through synthetic turbulence until encountering a patch of relatively strong turbulence when upon it randomly changes its direction of travel to avoid travelling directly through that patch. Distances travelled between consecutive turns were recorded and then compared with the expectations for Lévy walks and Brownian walks, a null model of the simulated movement patterns. These numerical simulations encapsulate key properties of turbulence, namely its incompressibility and its spectral properties which theory suggests will give rise Lévy flight movement patterns. Other aspects of turbulence such as the advection of smaller eddies by larger eddies (known as sweeping) that do not enter into the theory are not captured because they will incur prohibitively high computational costs. For modest Reynolds numbers, *Re*∼100, sweeping can be incorporated into computationally tractable kinematic simulations [20]. But the approach is not readily extended to higher Reynolds numbers where Kolmogorov scaling (the key ingredient in the theory predicting Levy flights) extends over several decades.

Here, following Fung *et al.* [20], velocity fields of synthetic three-dimensional, incompressible, homogeneous, isotropic turbulence are modelled as the sums of *N* unsteady random Fourier modes
*k*_{n} are random independent vectors that are uniformly distributed over a sphere with radius *k*_{n}≡|*k*_{n}|=*k*_{1}+(*k*_{η}−*k*_{1})((*n*−1)/(*N*−1))≡*k*_{1}+(*n*−1)Δ*k*, where *k*_{η} is the wavenumber corresponding to the Kolmogorov dissipation scale *l*_{η}=2*π*/*k*_{N}. The vectors *a*_{n} and *b*_{n} are independent random vectors with amplitudes chosen according to |*a*_{n}|^{2}=|*b*_{n}|^{2}=2*E*_{s}(*k*_{n})Δ*k*_{n}, where *k*_{n} and *E*(*k*) is the energy spectrum which here is taken to have the Kolmogorov form
*C*_{k}=1.5. The frequency-dependence is modelled by assuming that the energy at each wavenumber is spread over a range of frequencies around a characteristic frequency, *N* increases, the velocity field becomes closely Gaussian (as a consequence of the central limit theorem) and independent of the choice of distribution of the wavenumbers. The frequencies *ω* account for temporal fluctuations of the small-scale turbulence that are partly caused by the non-uniform advection by the large scales and partly by their dynamical interactions. The mean flow is neglected because it carries with it both the predator and the turbulence.

Reynolds numbers associated with the kinematic simulations are determined by *k*_{η}*L*∼*Re*^{3/4}, where *L*∼*σ*^{3}/*ε* is the lengthscale of the largest eddy and *Re*′≈4*Re*. The kinetic simulations were implemented with *N*=100, *k*_{1}=1, *ε*^{2/3}=0.1 and with *k*_{η}=100 or *k*_{η}=1000 in arbitrary length-time units, corresponding to Reynolds *Re*∼10^{3} and *Re*∼10^{5}. Going to higher Reynolds numbers is computationally prohibitive.

The choice of wavevectors, although random, determines the structure of the turbulent flow. For this reason, simulation data for each Reynolds number were obtained for 50 different realizations of the turbulent flows each of which was simulated for 10^{5} arbitrary time units.

For the sake of simplicity, it is assumed that predators swim at constant speed (speed *s*=2.5*σ*). Model predictions do not change significantly when the assumption of a constant swimming speed is relaxed, e.g. when a new speed is selected at random from an exponential distribution, *p*(*s*)=*e*^{−(s−2.5σ)} for *s*≥2.5*σ* otherwise *p*(*s*)=0, after each turn. Turbulence was taken to be weak when the local speed was less than 2.5*σ*, and strong when the local speed was greater than 2.5*σ*. Three scenarios were examined; movements confined to a fixed depth, one-dimensional vertical dive movements and freely roaming three-dimensional movements. In other kinematic simulations, the absolute value of the pressure gradient (i.e. the local acceleration of the turbulence) rather than speed was used to distinguish between weak and strong turbulence and so trigger changes in a predators' direction of travel. Pressure gradients were calculated using the prescription given in Fung *et al.* [20]. Note, however, that the pressure–frequency (*p*−*ω*) spectrum obtained from kinematic simulations follows the law *ω*^{−7/3} predicted by Kolmogorov's similarity theory, but that the constant of proportionality does not agree with theoretical expectations [20].

The Akaike information criterion [27] was used to test whether the simulation data provided more evidence for the step-length distributions having power-law
*et al.* [7], such comparisons are now routinely made in the analysis of empirical movement pattern data [9–15]. The power-law exponent, *μ*, and the exponential decay rate, λ, were determined using log-maximum-likelihood methods [28]. The start of the tail of the distributions (*a*≈10) was ascertained by visual inspection of the survival function (the complement of the cumulative distribution function). To construct the survival function, the simulation data for the step-lengths *l*_{i} were first ranked from largest to smallest {*i*=1…*n*}. The probability that a length is greater than or equal to *l*_{i} (the survival function) is then estimated as *i*/*n*.

## 3. Results

The results of the numerical simulations reveal that step-length distributions have exponential tails at low Reynolds number (*Re*∼10^{3}) and power-law tails at high Reynolds number (*Re*∼10^{5}; figures 1 and 2). Consequently, movement patterns at low Reynolds number are better approximated by Brownian walks, while those at high Reynolds number are better approximated by Lévy flights. The maximum-likelihood estimates for the Lévy exponents are *μ*=2.4. These findings are not dependent on: the dimensionality of the movement pattern; the speed of the predator, as expected on physical grounds [11]; whether turning is triggered by the presence of relatively strong turbulent velocity fluctuations or by the presence of strong pressure gradients (turbulent accelerations).

## 4. Discussion

The identification of Lévy flight movement patterns in marine predators—sharks, bony fish, sea turtles and penguins [15,16]—shifted the focus of Lévy flight research. The situation was succinctly encapsulated by Frederic Bartumeus, who said that ‘the question now isn't whether animals perform Lévy flight but when they do—and why’ [29] and by David Sims, who suggested rightly that ‘We’re interested now in understanding how these patterns arise. Not just the where and when, but also the how and the why. Those are the important questions' [30].

In this paper, a candidate generative mechanism for the occurrence of Lévy flight movement patterns in marine predators was identified with the aid of turbulence theory and numerical simulations. Predators that change their direction of travel only after encountering patches of relatively strong turbulence were predicted to have Lévy flight movement patterns when the Reynolds number of the turbulence *Re*≥10^{5}; a condition that is expected to be satisfied in all but the most quiescent of oceanic waters (figures 1 and 2). A region of the ocean with size *L*∼10^{4} m and having mean flow *u*∼0.1 ms^{−1} would, for example, have *Re*=10^{9}. When the turbulence is globally strong, then turning away from patches of relatively strongly turbulence could just be a natural response to buffeting. But when turbulence is globally very weak, e.g. when *R*∼10^{5}, then patches of relatively strong turbulence may be barely discernible and in these cases turning cannot be attributed to a buffeting-avoidance behaviour but may instead be selected for because it is advantageous.

The scaling exponent characterizing the Lévy flights in the numerical simulations, *μ*≈2.4, also characterizes the Lévy flight dive patterns of basking sharks and bigeye tuna (*Thunnus obesus*) [16]. The Lévy flight dive patterns of Atlantic cod (*Gadus morthua*) and leatherback turtles and Magellanic penguins are characterized by *μ*=2.0, 1.9 and 1.7, respectively [16]. Note, however, that for air-breathers (turtles, penguins) the change in chosen maximum depth between successive dives in single foraging trips was used as a proxy for searching to remove the anomalous effect on the step-length frequency distribution caused by the necessity of leaving and returning to the surface to breathe air [16].

The new finding suggests that Lévy flight movement patterns in marine predators could be a fortuitous accident rather than the product of natural selection acting on plastic neurological and/or physiological processes that can tune up for Lévy flights. The emergent movement patterns just happen to be near optimal when foraging. In the presence of highly abundant prey, these movement patterns will, if co-opted as search patterns, get arrested and as a result appear Brownian-like [11], in accordance with observations [17]. de Jager *et al.* [11] showed that any search pattern will resemble Brownian motion at high-resource densities, provided that movement is interrupted upon encounters. This is simply because the step-length distribution will be truncated just beyond the mean-free path length (the average straightline distance between adjacent resources). If the mean-free path length is sufficiently short then power-law tails, the hallmark of Lévy walks, will be entirely absent.

The predicted emergence of Lévy flight movement patterns in marine predators has resonance with the occurrence of Lévy flights in honeybees [10] which may be derived from the Weber–Fechner law in a bee's odometer (i.e. from errors in the estimation of distance) [31], with Lévy walks in T cells [12] which may be a by-product of their stick–slip locomotion [32], and Lévy walks in *E. coli* [9] which can be derived from noise in the chemotactic pathway [33]. Taken together, these findings call for a reappraisal of the Lévy flight foraging hypothesis. They suggest that this hypothesis can be replaced or supplemented by a broader, simpler, more general hypothesis. Namely that Lévy flights can emerge spontaneously and naturally from innate behaviours and banal, innocuous responses to the environment, and if advantageous then there can be selection *against* losing them.

To falsify the hypothesized mechanism, it will be necessary to show that an important class of marine predator has Lévy movements when swimming in still water or show that the changes of directions could not be turbulence induced. This would be very difficult to test with sharks, bony fish, sea turtles and penguins as most of these organisms cannot be observed across typical foraging scales under controlled (i.e. laboratory) non-turbulent conditions. But it would be possible with some of the minute invertebrates that populate the oceans as their movement patterns can be accurately recorded [34].

## Data accessibility

Computer codes can be obtained from the author.

## Funding statement

This research is supported by the Biotechnology and Biological Sciences Research Council (BBSRC).

- Received May 22, 2014.
- Accepted August 19, 2014.

- © 2014 The Author(s) Published by the Royal Society. All rights reserved.