A random walk description of individual animal movement accounting for periods of rest.

Animals do not move all the time but alternate the period of actual movement (foraging) with periods of rest (e.g. eating or sleeping). Although the existence of rest times is widely acknowledged in the literature and has even become a focus of increased attention recently, the theoretical approaches to describe animal movement by calculating the dispersal kernel and/or the mean squared displacement (MSD) rarely take rests into account. In this study, we aim to bridge this gap. We consider a composite stochastic process where the periods of active dispersal or 'bouts' (described by a certain baseline probability density function (pdf) of animal dispersal) alternate with periods of immobility. For this process, we derive a general equation that determines the pdf of this composite movement. The equation is analysed in detail in two special but important cases such as the standard Brownian motion described by a Gaussian kernel and the Levy flight described by a Cauchy distribution. For the Brownian motion, we show that in the large-time asymptotics the effect of rests results in a rescaling of the diffusion coefficient. The movement occurs as a subdiffusive transition between the two diffusive asymptotics. Interestingly, the Levy flight case shows similar properties, which indicates a certain universality of our findings.

Introduction

Dispersal is a common property of all ecological species [1,2]. Understanding its properties is important for many reasons, in particular because dispersal is known to be a primary factor controlling the rate and extent of species’ spatial spread [3-5]. While dispersal of plants occurs mostly due to the seeds’ transport by environmental flows (such as, depending on the species traits, wind or water flows), dispersal of animals mostly happens due to their self-motion. In this paper, we focus on the dispersal of animals.

Dispersal is usually described by the dispersal kernel[1] which is essentially the probability density to find the animal at a given location (say x) after a certain time (say t) of its movement away from the point of its release. A comprehensive approach to dispersal is therefore hardly possible without a good understanding of individual animal movement, and indeed the latter has been a focus of considerable attention and controversy over the last two decades [7-12].

A traditional way to consider the Brownian motion as an appropriate model of individual movement [13,14] (with the dispersal kernel having ‘thin’ exponential or even Gaussian tail corresponding to the ordinary diffusion) was eventually challenged by experimental observations that movement of some animals apparently follow a different pattern where the dispersal kernel has a ‘fat’ tail with a lower rate of decay described by a power-law [11,15-17]. The latter is usually associated with the Levy walk [12] or, more generically, with anomalous diffusion [8]. The rate of decay at the dispersal tail can have huge implications as it is believed to ultimately define the rate of species spread [3]. Theoretical approaches were developed to show that the foraging following the Levy-type movement pattern can optimize the search [18] and hence may arise as a result of evolution [15]. The issue as a whole remains highly controversial though. For some animal species that originally were claimed to perform Levy walks, later analysis re-classified them to the ‘thin-tailed’ Brownian walkers; see [7,11] for albatrosses, [15,19] for mussels, [20,21] for Tenebrio beetles.

Whatever the specific movement pattern is, a straightforward interpretation of a theoretical probability density function (pdf) as an ecologically meaningful dispersal kernel may be misleading because of the existence of different time scales. Standard approaches of statistical mechanics to calculate the propagator (dating back to the original work by Einstein [22]) usually consider a particle performing either uninterrupted movement or movement occasionally interrupted by short stops. Importantly, the duration of these short stops is on the same time scale as the duration of the short steps (e.g. between turning points in the so-called velocity model [6]) along the movement path. The combination of short steps and short stops defines the movement to which we, for convenience, will refer as continuing or uninterrupted. Also in biology, theoretical studies on animal movement aiming to calculate the dispersal kernel and/or the rate of animal spread (e.g. the mean squared displacement (MSD) as a function of time) usually consider an animal performing continuing non-stop movement [14,23-27]. However, this approach is only valid when animal movement is considered over a relatively short time. Indeed, if considered on a larger time scale, the animal would normally alternate long periods of continuing movement (‘bouts’) with long periods of rest [28,29]. It remains unclear how the properties of the dispersal kernel calculated for the uninterrupted movement can change as a result of upscaling when the animal movement is considered over a longer time (or at a larger spatial scale) that includes periods of rest. The goal of this paper is to bridge this gap. We are particularly concerned with the following questions:

— how the presence of rest times can affect the baseline dispersal kernel of uninterrupted movement; in particular, whether a thin-tailed kernel may become fat-tailed and vice versa;

— how the rate of spread (e.g. the dependence of the MSD on time) can be modified.

This paper is organized as follows. In §2, we introduce a model to describe the composite movement consisting of periods of uninterrupted movement and periods of rest and provide a semi-intuitive analysis of its properties. In §3, we perform a rigorous analytical study of the composite movement to obtain a general equation for the corresponding pdf. The mathematical analysis is illustrated and complemented with simulation results. In §4, our mathematical analysis of the movement is extended to address the rate of convergence of the pdf to its large-time asymptotics. Finally, in §5 we discuss the implications of our findings.

Composite movement: semi-intuitive analysis

We consider an idealized situation of animal foraging in a uniform and stationary environment. Although the assumption of uniformity is likely to be inadequate if the movement is considered in a natural environment on large spatial and temporal scales, arguably it should work well at small and intermediate scales (e.g. imaging a bumble-bee or a small rodent foraging in a large meadow or crop field), and it certainly works well in the case of a laboratory environment. We neglect the size of the animal compared with the typical travel distances, so that its location in space is given by a single point, e.g. =(x,y).

Let be the pdf of the animal position at a time t≥0 in the baseline case of ‘uninterrupted’ movement (i.e. without long periods of rest, see §1), where D is the distribution parameter; in particular, in the case of diffusion, D is the diffusion coefficient.

Let us now consider the situation where the movement stops suddenly at t=τ1, so that the animal is at rest for any t>τ1. Compared with the uninterrupted movement, this is a different stochastic process and it is described by its own pdf which is related to ρ as follows: as time after t=τ1 does not contribute to the pdf, because the position of the animal does not change for t>τ1. In the case that τ1 is a fixed parameter, equations (2.1) provide the full description of the movement (assuming that the baseline pdf is known). In the case that τ1 is a random variable, equations (2.1) give the description only for one particular realization of the process and should be averaged over all possible values of τ1 (in the manner shown in the next section).

Consider now a sequence of several bouts, τ1,…,τ, and the corresponding rests separating them, ξ1,…,ξ. Consider the moment t corresponding to the end of the last bout in this sequence, that is Using the same arguments as above, it is readily seen that the pdf of the animal position at time t is where T=τ1+⋯+τ is the total duration of the actual movement, i.e. total time excluding periods of rest. Moreover,

As an instructive example, let us consider the case where the baseline case of uninterrupted movement is described by a Gaussian dispersal kernel. We then observe that where Therefore, for the stochastic process defined by a given sequence of bouts and rests, the pdf of the animal position at time t coincides with the pdf of the uninterrupted movement up to the rescaling of the diffusion coefficient.

In any realistic animal movement, the durations of bouts and rests are of course not fixed. Instead, they can be regarded as random variables and described accordingly by their probability density functions, say ϕ(τ,λ) and ψ(ξ,ω), where λ−1 and ω−1 are the average duration of bouts and rests, respectively. Correspondingly, expression (2.5) is not the actual probability density function of this composite movement but the conditional one, say , i.e. the pdf conditioned by the given sequence of bouts and rests. This conditional pdf then must be averaged over all possible sequences {τ} and {ξ} in order to obtain the averaged or expected pdf, which we denote as . It is this expected pdf that gives a global description of the movement process taking into account details of the given stochastic process across different time scales.

For a small or intermediate t, the averaging is likely to turn the simple rescaling relation (2.6) into something more complicated. However, when time becomes large, i.e. when t≫λ−1+ω−1, the movement includes many bouts and rests of various duration. Intuitively, this is equivalent to averaging over different stochastic realizations. One therefore can expect that (2.6) is true asymptotically, i.e. for sufficiently large time. In the next section, these heuristic arguments will be confirmed by a rigorous analysis.

Rigorous analysis and results

Our analysis is based on the following two non-restrictive assumptions:

— The bout and rest durations are statistically independent random variables fully determined by their probability density functions, i.e. by ϕ and ψ, respectively;

— Bouts and rests do not present either dependence or influence on other variables of the movement, in particular on the functional form of the baseline pdf .

It does not seem possible to perform the analysis in the case of unspecified pdfs of bouts and rests duration, therefore we have to decide about appropriate functions ϕ and ψ. Based on substantial experimental evidence [29-31], a suitable choice for the probability distribution of bouts is the exponential distribution

The situation with the duration of rests is less straightforward. A standard assumption often made in the theoretical literature (sometimes implicitly, see section 2.3.2 in [31] for a discussion of this issue) is that the distribution of the rest times is exponential too Recently, there has been growing empirical evidence that, at least in some situations, the distribution of rest times can significantly deviate from the exponential, in particular to exhibit fat tails, cf. [29]. For the sake of analytical tractability, in this section we stick to the standard choice (3.2); possible modification of our results in the case that ψ has a fatter tail (e.g. being described by a power law) will be briefly discussed in §5.

Once the statistical properties of bouts and rests are thus defined, our approach consists of computing the probability of having k bouts within the observation time t, and then gathering the contributions of all possible number of bouts, which goes from one to infinity, to compose the expected quantity . To proper illustrate the procedure, let us first consider in detail the contribution from a realization with only one bout. If we assume the composite movement to always start with a bout, then at the observation time t there are two situations that must be considered: one in which a bout started at t=0 is still going and the other in which the bout period has stopped at some time t′

From these two simple cases, we can extend the formalism onto the general case in which an arbitrary number k of bout periods are obtained from the statistical ensemble. Again there are two slightly different situations, so let us first analyse the case in which the kth bout period is still ongoing at the observation time t. Prior to t, there must be exactly k−1 bout periods and k−1 rest periods, with a total duration no larger than t. But as the functional contribution in t does not depend on the specific order in which the periods occur, but rather on the total amount of time spent in bouts, the periods can be joined into two separate intervals with durations τ and ξ, corresponding to the total duration of bout and rest activities, respectively. As each period is obtained from the sum of k−1 exponentially distributed random variables, the resulting description of the total duration is given by a Gamma distribution. Hence, the probabilities of observing a total duration τ of bouts and ξ of rests obtained from k distinct periods are given, respectively, as

The same approach may be used to obtain the expected value of the function, which consists of assuming the total duration of the k−1 periods of bout and rest to be smaller than observation time, τ+ξunspecified fraction of the last one, but this issue is easily resolved by realizing that the whole bout period must be equal to t−ξ. So the full contribution of the case in which the kth bout period is still ongoing at t may be written as The case in which the kth rest period is still ongoing is a very similar, with the only differences being that there must be k bout periods prior to t, and that the functional contribution may be written directly as dependent on the total bout duration τ. Hence this contribution may be written as

Once the method for computing the contribution from a general ensemble of k bout periods has been described, the only thing left to do is to sum all contributions to obtain the expected value . As the case for k=1 is already accounted for by equations (3.3) and (3.4), the contribution of equations (3.6) and (3.7) should be within a summation that goes from k=2 up to infinity

Now, if we recall the series representation of the modified Bessel function of the first kind, then the summations in equation (3.8) easily turn into well-known functions and When all terms are combined, from (3.8) we obtain an expression that allows for a compact integral representation

It is possible to further simplify equation (3.11) if we consider a change in the integration variables in order to obtain a single integral representation, which is done by defining ξ=t(1−y) and τ=ty. We may also consider the following rescaled parameters which reduce the dimensionality of the parameter space. As a result, the expected value of the dispersal process may be written as where . Equation (3.13) gives the exact behaviour of the expected function at some observation rescaled time ζ. However, because of its complexity, the spatial behaviour of the distribution is not readily seen. In order to reveal this behaviour, we have to use some approximations and asymptotic analysis. In the following analysis, it is convenient to consider three time scales: (i) small, (ii) intermediate, and (iii) large.

(i) Small time asymptotics. When considering the movement at a very early stage of dispersal, ζ≪1 (i.e. t≪1/λ), the system is predominantly driven by the first bout. Therefore, at the corresponding observation time the process is described by the pdf of non-interrupted movement .

(ii) Intermediate time asymptotics. For a larger observation time, i.e. ζ∼1/λ, the activity switches between bout and rest period become statistically relevant, and their influence on the expected value should be properly addressed. To perform such analysis, it will require an approach that approximates the integrals without obscuring (or locally approximating) the spatial dependence, and at the same time it should be valid for any choice of . One way of accomplishing such feature is to approximate the kernel of the integral around the region where it exhibits its maximum value, in a way that the integrand would become a power series multiplied by the original dispersal function. By considering only the kernel of the integral in equation (3.13), i.e. then the objective is to obtain a power series in ζ expanded around the maximum value of the function defined on the interval [0,1], but the precise location of the maximum depends on both ζ and α. For ζ sufficiently small, the kernel presents its maximum value either at y=0 for , or y=1 for . In the intermediary range there is a transition in the behaviour, as the kernel function presents a concavity with the maximum, located inside the interval (0,1), being determined by ∂K(α,ζ,y)|=0. In the general case, it is always possible to determine, at least numerically, the precise position y(α,ζ) of the maximum, so (3.13)–(3.14) turn into the following: where (auxiliary figure 1 shows the precise location of the maximum in the initial dynamics regime, and the several values of α illustrating this transition behaviour indicate the neighbourhood of y where the power series should be expanded around). Equation (3.15) is generic and is valid for any dispersal process. However, it is only possible to obtain an explicit expression for when an explicit expression for the baseline pdf of the non-interrupted movement is available. In §§3.1 and 3.2 below, we consider equation (3.15) for two qualitatively different processes, i.e. for the Brownian motion when is the normal distribution and for the Levy flight when is the Cauchy distribution.

Figure 1.

Position of the maximum of the kernel function (3.14) for several values of α. On the right, from bottom to top, α=0.1,0.2,0.3,0.4,0.5,0.6,0.75,1.0,1.5,2.5 and 5.0. For small values of α, the maximum is located at , as long as ζ is also sufficiently small. As α increases, the position continuously shifts towards .

(iii) Large-time asymptotics. At sufficiently large times, ζ≫1, when a relatively large number of bout and rest periods have already taken place, we can use asymptotic methods to approximate equation (3.13). If we consider the asymptotic expansion of the modified Bessel functions, then we may consider only the first term of the expansion in order to write In this presented form, the integral (3.17) may be associated with the general expression which can be analysed by using the Laplace method to estimate the leading asymptotics for large ζ. As the function f e tends to zero for and (f alone diverges at y=1, but the limit guarantees the convergence), and S(y) has a unique maximum at the leading order of the integral at large observation times is obtained from

With the aid of the aforementioned approximation methods, we are now ready to further analyse the properties of the dispersal kernel by considering two specific representations for . Having assumed that the environment is uniform and isotropic, we will restrict our analysis to the one-dimensional case. We will consider two qualitatively different dispersal processes, i.e. the Brownian motion (‘Gaussian walk’), representing a thin-tailed distribution, and a Levy flight (‘Cauchy walk’), with its power law decay representing a fat-tailed distribution.

Gaussian diffusion

A one-dimensional Gaussian diffusion process is defined by the following dispersal function: where D is the diffusion coefficient. If we start the analysis by considering first the behaviour of the dispersal function 〈ρG(x,βζ)〉 in the small time regime, for ζ sufficiently small it is expected that the main contribution would come from the damped Gaussian term e−ρG(x,βζ), as the integral term in equation (3.13) is associated with multiple bout sequences. In order to address its contribution, we must look at values of ζ where the stochasticity of the bout and rest sequences start to become relevant, which typically occurs at . The dispersal function in this regime is given by equation (3.15), with the integral being written as where Γ(k,z) is the upper incomplete Gamma function. If we look first at the asymptotic spatial regime, then with the aid of the expansion for large real arguments (), it is possible to conclude that, in this intermediary time regime, the original Gaussian process is always asymptotically dominant because the leading order of every term defined by equation (3.22) presents a faster decay for large spatial values due to the presence of an x−2 factor. In a more formal way, the complete sum of the series expansion (3.15) for the Gaussian dispersal behaves asymptotically as which decays faster for large x than e− present in the original process. As a result, even though for small times the original Gaussian diffusion accounts for most of the dispersal curve, the exponential dumping diminishes its contribution as the system evolves, which leads to a restructure of the dispersal function that happens at small spatial scales and tends to spread to large ones as the time grows. For any fixed time, the bout and rest sequences contribute to a reshaping of the dispersal function for small x while keeping the tail of the original Gaussian process for large x: but, the larger the time ζ we look at the function, then we need to look at bigger values of x in order to have its tail described by the original process. So the picture we get at these intermediary times is one in which a slower spreading dispersal function 〈ρG(x,βζ)〉, slower because of the presence of rest periods, is associated with fatter tails coming from the original faster dispersal ρG(x,βζ), as shown in figure 2.

Figure 2.

Rescaled dispersal functions for the Gaussian process at small time regime, for β=1 and α=2. From left to right, ζ=0.5,1,2,4. With the scaling factor S(ζ)=〈ρG(0,ζ)〉, blue (thick) curves represent the rescaled dispersal , and black (thin) curves the rescaled damped Gaussian S−1 e−ρG(x,ζ). The tail of the distribution is always determined by the tail of the original process, but as the system evolves its contribution becomes less significant, mostly restricted to large distances.

To understand what happens to the dispersal function at large times, we must resort to the Laplace method (3.20) to estimate the leading asymptotics. In this regime, the damped Gaussian term should not contribute any more, and the dispersal function may be approximated only by the asymptotic contribution of integral (3.17), resulting in From equation (3.25), it is straightforward to conclude that asymptotically we obtain the following rescaling of the diffusion constant, which coincides with the results of the semi-intuitive analysis in §2. Correspondingly, it leads to a renormalized diffusion process being described by

As the time evolution continuously transforms the dispersal function from the original ρG(x,Dt) process into a slower dispersing process defined by ρG(x,D*t), because D>D*, it also becomes evident why we observe a dispersal function with a somewhat fatter tail: whenever we combine the functional forms of two Gaussian processes characterized by different diffusivities, the result is not exactly a Gaussian process but rather a dispersal process characterized by the slower diffusion at small spatial scales, with fatter tails (from the faster diffusion) present at large spatial scales. Figure 3 shows the dispersal function with fatter than Gaussian tails, and how the convergence takes place for large times, thus confirming our conclusions.

Figure 3.

Rescaled dispersal functions for the Gaussian process at large-time regimes, for β=1 and α=2. From left to right, ζ=4,8,16,32. Blue (thick) curves represent the rescaled dispersal , while black (thin) ones represent the rescaled asymptotic diffusive behaviour ρG(0,β*ζ)−1ρG(x,β*ζ), as defined by equation (3.27), where β*=βα(1+α)−1. As the dispersal function converges to the asymptotic diffusion from above (at sufficiently large spatial values), for finite times the core of the function may be approximated by a Gaussian diffusion, but it still presents a faster than normal decay at the tails.

In the case that the baseline process (non-interrupted movement) is a Brownian motion, so that the pdf 〈ρ(x,⋅)〉 is thin-tailed, a convenient alternative way to describe dispersal is to consider the MSD, i.e. the second moment of the pdf. It is well known that a non-interrupted diffusive process (i.e. bouts only, no rests) results in the MSD depending linearly on time, 〈x2〉∼Dt, where D is the diffusion coefficient. Now we are going to consider the MSD for the composite movement combining alternating bouts and rests in order to reveal how this linear dependence can be modified.

From the definition of the MSD for the dispersal process, we obtain Based on our above analysis of small-time and large-time asymptotics, it is readily seen that where β*=βα/(1+α).

Comparing the small-time and large-time asymptotics in (3.29), intuitively, one can expect that, at intermediate time, the movement slows down. Unfortunately, for the intermediate asymptotics ζ∼1, analytical investigation does not seem to be possible because of the complexity of the corresponding equation (3.15). The dynamics on the intermediate time scale is therefore studied through Monte-Carlo simulations. The results are shown in figure 4. It appears that, on the intermediate time scale which can span over a few orders of magnitude, the movement shows clear subdiffusive properties so that 〈x2〉∼t, where the exponent γ can be anywhere between 0 and 1, being an increasing function of α (figure 5). The exponent tends to zero as , because the system tends to stay in rest for most of the time, and as , because of the opposite situation where the animal keeps moving almost all the time as the proportion of bout duration becomes much greater than the proportion of rests. Interestingly, if this result is considered together with the properties of the dispersal kernel, the movement has somewhat counter-intuitive properties as the subdiffusive spread occurs for the kernel with the tail fatter than Gaussian.

Figure 4.

Mean squared displacement of the dispersal process, (a) in the original variables, (b) in rescaled variables. In both panels, curves are shown for (from bottom to top) α=10−3,10−2,10−1 and 100. The dashed line shows a power law proportional to ζ−1. When the expected value of the rest periods is much larger than the bout ones (α≪1), the system exhibits long periods of almost no movement, which accounts for the 〈x2〉≈const. If α is increased until the expected periods become comparable (α∼1), then the MSD presents a power law type of transient behaviour of the form 〈x2〉∼ζ, with γ∈(0,1), which is characteristic from a subdiffusive type of spread.

Figure 5.

Exponent of the transient subdiffusive spread 〈x2〉∼ζ as a function of the ratio α of expected periods of movement and rest. The exponent γ is calculated as the maximum slope of 〈x2〉(t), cf. figure 4.

Levy flights

As an example of a fat-tailed ‘superdiffusing’ process, we are going to consider the Cauchy walk that is often regarded as a standard model of Levy flights [12]. In the one-dimensional case, the corresponding dispersal kernel is given by the following expression: To describe the behaviour of the dispersal process in intermediary times, we must consider the expansion (3.15), with the integral being written as where 2F1(a,b;c;z) is the hypergeometric function, as defined by and (q) is the Pochhammer symbol representing the rising factorial, Differently from the Gaussian process, where we were able to show that the presence of bout and rest periods tend to affect the dispersal kernel in a way that propagates from the origin (bulk) towards the tails, the effect on the Cauchy dispersal appears in a global way, because the perturbations induced are of the same order as the dispersal process itself, as may be observed from the spatially asymptotic contribution of the integral term, i.e. Even though we are observing a process where the dispersal kernel is being restructured across all spatial scales simultaneously (which is probably related to the fact that the dispersal function is also scale-free), the dispersal behaviour at large times shows the same general property as the Gaussian process, as asymptotically we also recover a renormalized Cauchy type of dispersal defined by where the scaling relation for the distribution parameter, that is is exactly the same as in the case of the Brownian motion.

Universality of rescaling

From the analysis of the two above-mentioned examples of dispersal processes, one cannot leave unnoted the emerging behaviour at large-time regimes, as both dispersal functions recover their initial functional forms in a renormalized time scale. In fact, if the observation time ζ is sufficiently large, it should be expected that the fluctuations arising from the random sampling tends to disappear, as the system would settle into a regime where it splits its whole activity into well-defined portions of rest and bout periods. If this statement is indeed true, then the rescaling feature should be dependent only on the statistical nature of the activity periods sampling, thus being observable with any dispersal function at hand. Although we cannot formally prove this conjecture, we can show some evidence strongly supporting the statement if we look at the expected duration of the total amount of bout period that takes place within a given observation time ζ, which may be written as For ζ sufficiently large, having applied the Laplace method, one can show that this function tends to a linear time evolution that, when written in the form enables us to interpret the asymptotic dynamics as a constant fraction of the movement spent on bout activity. Hence, any type of dispersal function that has its proper behaviour altered by the statistical distribution of bout and rest periods should be described by its original function with this rescaled time, and if the degree of spatial dispersal of the function at a given time ζ is characterized by a single constant entering the functional form only as a time multiplier, then the whole change in behaviour should be characterized by the renormalization process independent of the specifics of the function.

Convergence rates

As was shown above, the dispersal process 〈ρG(x,βζ)〉 converges to the asymptotic diffusion ρG(x,β*ζ). However, details of the convergence have remained unclear. One way of approaching this problem is to define a measure of how different the two functions are and look for general properties of its time evolution. With that in mind, we are going to define a global convergence function (or global error function) as which gives a measure of the proximity of the two functional forms at a given time ζ, and enable us to determine via an asymptotic analysis the typical time scale of the convergence. But there is a scaling issue with this function, because the divergence appearing as does not allow one to define a proper (natural) scale for comparison. To circumvent this problem, we may consider a scaling function defined by the integrated squared difference between two arbitrary Gaussian processes, i.e. As this scaling function presents the same divergence when two distinct diffusivities (β1 and β2) are considered, we may look at the ratio between the global convergence function R(α,β,ζ) and the maximum possible difference between 〈ρ(x,βζ)〉 and ρ(x,β*ζ), as defined by Λ(β,β*,ζ), which gives a natural way of measuring the time evolution of global functional distances because the ratio is always defined on the interval (0,1). If we perform a numerical evaluation (figure 6), then it is possible to observe that intermediate times are characterized by a non-uniform evolution of the ratio, but asymptotically the convergence presents a well-defined power law type of convergence. As our main interest lies in determining the typical behaviour on this asymptotic limit, we must look for a way to determine higher order contributions coming from integrals of the type (3.18), but this involves a re-examination of the Laplace method.

Figure 6.

Ratio between the global convergence (R) and maximum difference (Λ) functions for the Gaussian dispersal. From left to right, α=10, with n= 1,0,−1,−2,−3. Continuous curves (blue) are the numerical evaluation of the ratio function, while dashed (black) lines represent the asymptotic behaviour as defined by equation (4.7).

When considering integrals of the type (3.18), we typically consider a series expansion around the maximum y0 of the function S(y), followed by a Gaussian type of approximation of the integral, which is only justified if the kernel is sufficiently sharp around its neighbourhood. Because we are interested in determining the first few orders of contribution in the asymptotic limit , let us consider the following expansion: where the notation used is defined by Now, if the derivatives do not show an increase in the order in ζ when compared with the leading order term f0, then the next leading order will come from terms whose expansion in ζz satisfy δ2=2(δ1+1), because each integration of z2 generates a ζ−( multiplying factor. When taking this into account, it is possible to conclude that the second leading order will come from contributions of the following approximation: If we come back now to the specific case of the Gaussian dispersal, and remember to take into account the second term of the asymptotic Bessel expansion (3.16), then it is possible to compute the asymptotic behaviour of the convergence function as thus leading to the asymptotic ratio which matches perfectly the numerical results (figure 6). Our result shows that the convergence time asymptotically decays proportionally to ζ−2, so in order to obtain a typical time scale of convergence one needs only to define an arbitrarily small value ϵ of comparison for the ratio and invert equation (4.7) to determine the time ζ(ϵ) required for the system to reach that confidence level.

In the case of the Cauchy walk, as this dispersal also shows a convergence between two functional forms, it is possible to use the same convergence function (4.1) to determine the typical time scale of convergence for the Cauchy dispersal. The numerical evaluation of the ratio function is shown in figure 7, and its asymptotic behaviour may be obtained, after some calculations, as which clearly shows the same type of time scaling ζ−2 as the Gaussian dispersal.

Figure 7.

Ratio between the global convergence (R) and maximum difference (Λ) functions for the Cauchy type of dispersal. From left to right, α=10, with n=1,0,−1,−2,−3. Continuous curves (blue) are the numerical evaluation of the ratio function, while dashed (black) lines represent the asymptotic behaviour as defined by equation (4.8).

Discussion and concluding remarks

In this paper, we considered the properties of the composite random walk consisting of alternating periods of movement (bouts) and rest. Having used the time scale separation argument, i.e. that the periods of rest between the bouts occur on a larger time scale compared with possible short stops of the baseline ‘uninterrupted’ random movement, we performed a rigorous mathematical analysis of this composite walk. Namely, we calculated the ‘expected’ probability density function of the position of the walker which provides a global description of the composite movement across different time scales and over different stochastic realizations of bouts and rests. This global pdf 〈ρ〉 appears to be related to the probability density function of the corresponding uninterrupted movement. In a general case, the pdf of the composite movement is given by the integral relation (3.11). We then focused on two specific cases where the pdf of the uninterrupted movement is (i) a normal distribution quantified by the diffusion coefficient and (ii) the Cauchy distribution quantified by its parameter. For each of these two cases, we obtained the explicit expressions for the probability distribution of the composite walk in the large-time asymptotics, see equations (3.27) and (3.35), respectively.

With regard to the questions itemized in the introduction, we therefore have obtained the following answers:

— Regardless of the properties of the uninterrupted movement, i.e. whether it is the Brownian motion described by (3.21) or the Levy walk described by (3.30), the evolution of the pdf in time follows the same scenario resulting, in the large-time asymptotics, in the rescaling of the distribution parameter in the same way, cf. equations (3.26) and (3.36). The functional form of the composite walk’s pdf therefore coincides with the pdf of the uninterrupted movement subject to the ‘universal’ rescaling of the distribution parameter; see (3.39). For the intermediate time scale, the pdf of the composite movement is given by the integral expression (3.11) which, as we have shown by numerical simulations, predicts a somewhat fatter tail compared with the pdf of the uninterrupted movement (figure 3).

— In the case that the ‘baseline’ uninterrupted movement is described by the Brownian motion, on the intermediate time scale the MSD shows clear subdiffusive properties, i.e. 〈x2〉∼t with 0<γ<1, where the value of the exponent γ depends on the ratio of the mean bout duration and the mean rest time.

Note that the analysis done in §3.3 suggests that the rescaling (3.39) should remain the same regardless the properties of the dispersal kernel. Regarding the convergence rate in the large-time asymptotics (see §4), although we cannot prove universality in a general case, the fact that the rate of convergence is the same (i.e. ∼ζ−2) in the two very different cases such as the Brownian motion and the Cauchy walk suggests that this is a general result too.

The subdiffusive movement on the intermediate time scale is a somewhat surprising result. We mention here that subdiffusion is widely observed in physical and chemical systems [31], where it usually appears as a result of particle trapping by certain spatial structures of the system (referred to as ‘traps’). By contrast, in our study the animal moves in a perfectly uniform environment without traps or refuges. The periods of immobility (i.e. rest) are attributed to purely behavioural reasons. Indeed, the existence of behaviour is a factor that can make movement of animals significantly different from movement of inanimate particles, cf. [26].

An interesting question is how the composite movement can be modified if the distribution of the rest times is not exponential but is described by a fat-tailed distribution (e.g. a power law). Although we cannot currently provide a quantitative study into this issue, the results that we obtained in this paper make it possible to figure out what the movement properties are going to be. A fatter-tailed distribution of the rest times would obviously result in further slowing down the movement on the intermediate scale, hence making the emerging subdiffusive behaviour more prominent. In the large-time asymptotics, however, one should expect that the functional form of the pdf of the uninterrupted movement is restored, subject to the rescaling of the corresponding distribution parameter as given by the universal relation (3.39).

We point out that our findings (such as the general expression (3.11) for the dispersal kernel and the universal rescaling (3.39) valid in the large-time asymptotics) have immediate ecological applications. Provided the statistical properties of the bouts and rests are known—which can be obtained from a small-scale experimental study—and the pattern of uninterrupted individual animal movement is established, it can be upscaled to the large-time ecological scale to give the dispersal kernel. Once the dispersal kernel is thus obtained, it can be used to construct kernel-based mathematical models describing the spatio-temporal dynamics of a given population over a broad range of spatial and temporal scales, e.g. see [4].

A question may arise here about the generality of the above message. Indeed, our results were obtained under some simplifying assumptions. For instance, our model of composite movement did not take into account the heterogeneity of the environment where the animal lives and moves; meanwhile, it is known that the complexity of the environment can affect the movement pattern [15,32]. However, here we argue that, unless the environment has a large-scale structure, e.g. as given by the presence of mountains or lakes, a small-scale heterogeneity is likely to affect only the choice of the movement step, which is already described by the pdf.

In conclusion, we mention that our study obviously leaves a few open problems. Firstly, we assumed that the animal uses only one movement mode, e.g. as described by a single probability distribution. In reality, however, animals often have more than one mode [19,33,34] and switching between the modes may depend on a variety of factors, in particular on the spatial and temporal scales of the movement [35]. Whether the existence of different movement modes may or may not affect the upscaling relation (4.1) is not immediately clear. Secondly, here we have considered only the individual movement, thus neglecting any interaction between the animals; however, the latter is known to affect their movement [36], for instance, through formation of groups such as schools, packs, herds, flocks, etc. [37]. Moreover, some recent results show that there may exist a more intricate collective behaviour such as density-dependent cut-off when the movement type changes at a certain spatial scale because of the interactions with conspecific animals [38]. Although in the first instance we can argue, similarly to the above, that the interactions with the other animals should already be taken into account by a properly chosen pdf of the uninterrupted movement, this issue requires a more careful consideration. All these issues should become a focus of the future study.