Lévy Flights and Self-Similar Exploratory Behaviour of Termite Workers: Beyond Model Fitting

Animal movements have been related to optimal foraging strategies where self-similar trajectories are central. Most of the experimental studies done so far have focused mainly on fitting statistical models to data in order to test for movement patterns described by power-laws. Here we show by analyzing over half a million movement displacements that isolated termite workers actually exhibit a range of very interesting dynamical properties –including Lévy flights– in their exploratory behaviour. Going beyond the current trend of statistical model fitting alone, our study analyses anomalous diffusion and structure functions to estimate values of the scaling exponents describing displacement statistics. We evince the fractal nature of the movement patterns and show how the scaling exponents describing termite space exploration intriguingly comply with mathematical relations found in the physics of transport phenomena. By doing this, we rescue a rich variety of physical and biological phenomenology that can be potentially important and meaningful for the study of complex animal behavior and, in particular, for the study of how patterns of exploratory behaviour of individual social insects may impact not only their feeding demands but also nestmate encounter patterns and, hence, their dynamics at the social scale.


Introduction
Subterranean termites are tropical social insects having an economic impact in the billions of dollars all over the world each year [1]. For this reason, their patterns of movement have fomented both practical and theoretical studies in the past. Most of these studies, however, have concentrated in displacements at colony level and inside tunnels leaving the study at the individual level largely overlooked. Also, while animal movement in the context of resource acquisition has received a fair deal of attention [2][3][4][5], its links to individual interaction and -ultimately-to sociality still need better focus. Providing that interindividual interactions lie at the very heart of sociality [6] and that interactivity depends on movement and space [7], effective space exploration must play prominent role in social behaviour in general, and in termites in particular [8].
A common misconception about termite spatial orientation is that physical restrictions imposed by tunnel walls would guide effective space sweeping. While this could be true for termites foraging within tunnels underground, once inside their nest the scenario changes completely. There, galleries merge and split in unpredictable ways usually forming an entangled set of paths and chambers whose layout does not seem to help orientation in most cases. Termites could, then, be guided by specific chemical cues and vibratory signals to effectively find nestmates. But would that be entirely dependent on external clues? Or, as sustained by a long forgotten hypothesis [9], would termites orient themselves concatenating external stimuli with an optimized search strategy, as with other animals [10]?
Here we present evidence in favour of such a hypothesis, describing exploratory spatial behaviour in isolated termite workers kept in large containers, free from the constrained movements they experience within tunnels. In this way we were able to assess individual free exploratory behaviour in clueless environments and away from social interactions. We conclude that their searching patterns are compatible with scale-free strategies based on a fractal exploration of space and that these are key to the efficient flow of information between nestmates, thereby providing expressive hints on how self-organization underlies social cohesion [11].
Our study includes the analysis of anomalous diffusion where the mean squared displacement (MSD) is estimated as a base to identify superdiffusive aspects of termites movements. Within the framework of anomalous diffusion theory [12,13], the scaling exponent of the MSD is related mathematically to a Lévy probability distribution P(l)!l {m of the steps lengths l, with 1vmƒ3. We evaluate these relationships and compare them versus the results that can be extracted from a Kolmogorov structure functions analysis. We progress into exhibiting that a termite walking is a 1=f noise process with power-laws in the correlation function indicating long-memory that is compatible with the fractal structures revealed by an Iterated Function System algorithm. Finally, we use maximum likelyhood estimation (MLE) to show that waiting times with power-law scaling exponent values are also present in the termite exploratory behaviour with scaling exponents values compatible with those predicted by the physics of transport phenomena.

Experimental set-up
Cornitermes cumulans (Kollar) (Blattaria: Isoptera: Termitidae: Nasutitermitinae) workers were collected from a wild colony at the campus of the Federal University of Viçosa, Minas Gerais, Brazil (20u 459 S, 42u 529 W), on 14 Dec 2010. Cornitermes spp. are Neotropical termites occurring in several habitats, including forests, ''cerrados'' (Brazilian savannas), and man-modified habitats. No specific permissions were required for the locations or activities reported in this manuscript. The field studies did not involve endangered or protected species. The study did not involved human participants, specimens or tissue samples, or vertebrate animals, embryos or tissues.

Movement recordings
Experimental arena consisted of a glass Petri dish (Ø = 205 mm) upside down over a sand-blasted flat glass. To allow evaluation of free-walking behaviour, no obstacle, nestmates or food was present in the arenas. Within the arena, a single termite worker was inserted and its movement was tracked from above, continuously for 5-6 hours, with a closed-circuit video camera (Panasonic WV-BP334 B&W CCD) equipped with a macro lens (Fujinon YV562.7 R4B-SA2L). To allow for full contrast between termites and the background, thereby permitting appropriate video recording, termites were painted with nontoxic water soluble dye [14]. We also performed experiments in containers of diameters 140 mm, 90 mm and 55 mm and found that the size did not significantly altered the results here reported. In total, half a million worth of termite displacements were recorded. We concentrate in reporting the data from the largest containers and choose four representative examples for our analysis and discussion.
Termite trajectories were captured and digitized at a sample rate of one point every 0.5 s with an automatic video-tracking software (EthovisionXT version 8.5.614, Noldus Information Technology). The paths were converted as a series of Cartesian coordinates coupled with respective time records, which allowed posterior calculations of spatial-temporal displacements of termites and their associated parameters, as shown in Fig. 1 Recordings included all three displacing behaviours presented by termites in the arena: (i) walking along the border of the arena, (ii) free moving, and (iii) static waiting-times (see Fig. 2).

Path definitions
Paths were defined as the sequence of line segments walked by the termite individual every 0.5 seconds inside the Petri dish (Fig. 3). Waiting times are the time elapsed when the distance recorded was zero. Both quantities, distances and waiting times, are heavily dependent on the sampling rate. Under-sampling occurs if the time interval between two samples misses turning points and oversampling happens when two video shots are taken from the same line segment. Avoiding both, under-and oversampling, is a matter of prime importance. However, it is not always clear how to do it. From Fig. 4 it is obvious that different sampling rates would result in different conclusions if an histogram is the only choice made in order to measure potential scaling exponents, coming from hypothetical power-laws. Current trend of model fitting may even led toward selecting a model that excels the statistical criteria but has no phenomenological meaning. Here, we circumvent this, completing a series of alternate data analysis to evince key aspects of the termite walking behaviour while, at the same time, its dynamical richness is exhibited.

MSD and anomalous difussion
A widely-used method to evaluate the diffusive properties of mobile objects is the mean-square displacement (MSD) SDr r t (t) 2 T, where Dr r(t)~½r r(t){r r(tzt) is the object displacement at time t with a time delay t. Normal diffusion, as in the case of Brownian motion, occurs when the MSD increases linearly so that SDr r t (t) 2 T~2Dt, where D is the diffusion coefficient of the object. Diffusive regimes where the MSD is not linear in time but proportional to t a with a=1 are known as anomalous diffusion.  Furthermore, when a is the diffusion exponent and av1, the process is called subdiffusive, and if aw1, the process is called superdiffusive.
Cornitermes cumulans diffusion in the arena can be estimated and characterized when the value of the MSD is calculated and when, for the sake of clarity, it is compared against a simulated Brownian particle of a known a~1. It is easy to see, in Fig. 5, that while the Brownian particle has the expected scaling exponent of a~1, all termite's scaling exponents are aw1 and so their movements are properly superdiffusive. Termite step lengths and waiting times statistics ( Fig. 6) have broad distributions corresponding to superdiffusive (Lévy) process where the mean squared displacement grows faster than it does in the normal diffusion of a Brownian particle [15]. Phenomenologically, termite workers exhibit anomalous diffusion in their walking patterns because these individuals remain in motion without changing direction for a time that follows a scale-free distribution. These less-frequent long distance displacements are clearly visible in our experimental arenas (Fig. 2).
Physical systems where superdiffusion and Lévy fligths were first identified together and theoretically related to each other, have   been studied in great detail up to now. These systems include the transport of passive particles in non-linear fluids where fluid motion may be retarded at the border of regions having different velocity vectors. Such regions cause ''sticking'' that result in longtailed distributions, both in the distances traveled and the waitingtimes [16]. For all practical purposes both superdiffusion and Lévy fligths are two sides of the same coin, in such a way that a measure of one is related to a measure on the other. Superdiffusion and Lévy fligths scaling exponents (a and m) are know to be theoretically related by a simple expression [12,16]: that can be used to estimate Lévy fligths exponents from a MSD measure, as detailed below.

Structure functions and scaling
A well know method to characterize potential free-scale behaviour is due to Andrey Kolmogorov in the context of turbulent flows [17][18][19]. A generalization of the Kolmogorov method is known as structure functions of order q, being expressed as SDDl t D q T&t f(q) , where l is the traveled distance between two (x, y) points in the space, Dl is the distance increments, t is a time lag that does the role of spatially coarse-graining the data set at different scales, so that scale-invariant observables could be spotted through the scaling of the exponent function f(q) and, finally, the angled brackets denote an average operator. Since the 1941 Kolmogorov seminal work on turbulence, the q-order structure functions have been used extensively in physics [20][21][22][23][24][25], its use in biology for studying scale-free movement patterns in animals is relatively new and mostly applied to the analysis of invertebrates such as copepods [26] and fruit flies [27].
In order to study the free-scale patterns of termite movements under the Kolmogorov framework, we considered the trajectories yield by the walking of representative solitary individuals confined in a Petri dish arena where the values of the re-scaling lag t was estimated from Then the averages over all increments were calculated and elevated to the proper q power, so that the q-structure function is calculated: in such a way that log(SDDl t D q T)&log(t f(q) )&f(q)log(t). Notice that the last expression is a linear equation with slope f(q) that can be easily evaluated. The last step needed under this procedure is to realize that f(q)~dq and that the value of the slope d is related to the Lévy exponent m as d~(m{1) {1 . It is also important to remark that if f(q)~1 2 q, the process is normal-diffusive but if f(q)= 1 2 q then the process exhibits anomalous diffusion, with f(q)w 1 2 q being the superdiffusive regime (a condition fulfilled by all termites examined here). Equation 2 is used to estimate the structure functions SDDl t D q T of our time series for q-orders~f1,2::,8g. Results are plotted in Fig. 6 (left column) where linear scaling is seen on the log-log plot for values of tv10, therefore a power-law scaling is evident. The values of these slops are, in turn, plotted in Fig. 6 (right column) These plots evince the linear scaling behaviour of f(q) and so the value of the slopes d~(m{1) {1 , where m is the Lévy exponent.
Lévy flight exponents in the termite superdiffusive walking estimated by the two previous methods are compared in Table 1. Notice that the mean values for the m exponents gives Sm 1 T+sd = 2.2060.11 (from MSD) and Sm 2 T+sd = 2.1660.04 (from q-functions) that are remarkable similar.

Distance temporal fluctuations
The self-similar features of the distances traveled by a termite walking can be visualized by means of a Fast Fourier Transform that converts time fluctuations into a frequency fluctuations domain. This is useful for detecting temporal features such as periodicity or broad frequency spectra. The analysis is also very illustrative because it helps detecting self-similarity when the power spectrum scales as a power-law that, in turn, may hint the presence of colored noise, self-organization or critical dynamics. In addition, and for comparison purposes, it is possible to artificially generate a power-law time-series using well know algorithms in order to exhibit the resemblance of both power spectra. This procedure adds further support in characterizing the termite walking pattern as a fractal process in the traveled distance distributions.
For generating the artificial time-series, a transformation method between random probability distributions is preferred [28]. The method can generate a random power-law-distributed real number x~x min (1{r) ({1=(m{1)) in the range x min ƒxv? with exponent m, when feed with a random real number r uniformly distributed in the range 0ƒrv1. After the Fourier transform was applied to a 4096 long segment of the termite walking, its power-law scaling P(f )!f {b was evident (Fig. 7A), scaling with an exponent b~1:3. This information was used to generate a power-law with the same scaling exponent (Fig. 7B). The resemblance between both spectra is remarkable, adding support for considering the termite walking pattern as a power-law relaxation process.

Fractal signatures and IFS
Lévy flights have spatial self-similarity due to their power-law nature. Because of this, time-series containing Lévy-distributed random variables can be visualized to evince their subtle structural details, a pattern shared with time-series containing correlated coloured noise that can be generated artificially. An example of this is given in Fig. 8 where a fractal recurrence method known as Iterated Function System (IFS), has been used to reveal the detailed self-similarity found in a termite walking.
IFS fractal reconstruction is a method widely used to visualize biological data such as in DNA sequences [29], insect population dynamics [30] or the temporal activity of social insects [6,31].
Here the method is used for the first time to visualize the selfsimilar features of an animal walking while performing exploratory behaviour. Fig. 8 depicts one of the termite time-series. It reveals in full detail the self-similarity of it, on a set of fractal triangular structures. A similar picture arises when the IFS algorithm is used to visualize a time series of an artificial power-law distribution of a known FFT scaling exponent b~1:3. The foregoing result is sound evidence in favor of considering the termite walking pattern as a self-similar process with power-law distribution of step lengths.
Step correlations Successive steps in a Cornitermes cumulans termite motion are correlated to a certain extent, as evident from the fact that these resemble correlated pink noise (Fig. 8). In a termite walking, successive steps are no independent from each other with temporal dependencies that decrease along the trajectory. One method to identify and characterize such correlating distances is the correlation function (CF): where N is the time-series length, t is a time increment, s 2 is the variance and h h is the mean of the steps l. Correlation function is useful for finding temporal patterns in the time series such as periodicity or independence. In the case of white noise DC(t)D~0, so no correlation is present and the successive values are The values of exponents a and m 1 are related by theoretical results [16], and are also compatible with the results of a Lévy exponent m 2 via a Kolmogorov structure functions analysis. doi:10.1371/journal.pone.0111183.t001 independent. We learn from Fig. 9 that this is not the case of termites. On the other hand, it is common that CF decays exponentially with rapid rates as DC(t)D!e {t=ts , in which case the time-series have correlations only among short distances t s . Nontrivial correlation over long distances usually decay asymptotically as power-laws DC(t)D!t {c , with 0vcv1.
Cornitermes cumulans walking patterns have been found here to have long-range correlations decaying as power-laws. This means that a step length chosen by the termite at a particular time, will influence the length of future steps in a process akin to longmemory. For years physicists have been studying systems that have power-law decaying CFs and have linked them to either systems which are far from equilibrium or very close to a critical point (e.g., a phase transition). We do not know the details of the locomotive physiological mechanism at play that generates a correlated termite walk but we can compare it to an artificially generated time series in order to gain more insights. Consider Fig. 9 where the same auto-correlated time-series explored before (Fig. 8) is again analyzed. We now know that this artificial timeseries has a power-law decaying CF but has originated by a multiplicative growth process. There are differences worth mentioning. Notice that the termite CF varies softly in contiguous values of t while the artificial signal has a CF varying wildly in amplitude along of all its decay, evidencing its stochastic origin.

Waiting times
Foraging patterns of mobile individuals often exhibit resting periods known as waiting times. In an ecological context, it has been speculated that these may result from the encounter of places of interest that need time for being explored and exploited so that observed statistical properties of waiting times will reveal the distribution of resource patch sizes [10]. In a behavioural approach waiting times may be originated in the fact that a mobile organism actually need periods of resting between activities. A number of foraging studies on the statistical properties of waiting times report that these behave as random-like burstiness that appear to be scale-free and so well-approximated by power-laws, meaning that short waiting periods are frequent while long are rare [13]. It has also been recently speculated that the powerlaw scaling of waiting times is a behavioural ''rule of thumb'' evolved to optimize move-wait decisions in unpredictable environments [32]. Remarkably these waiting times resemble the scaling laws found in the study of superdiffusive passive trackers moving in non-linear flows and where these waiting times are a well-known indivisible aspect of anomalous diffusion [16]. In fact, let's return to equation 1 where a relationship between the values of the MSD and Lévy scaling exponents was established. As a matter of fact, equation 1 holds best if lw2 [13,33]. Since the equation actually holds in our case because the estimated values of m 1 and m 2 are consistent (Table 1), we expect that the termite waiting times will be described as power-laws w(t)!t {l with lw2:0: We show in Fig. 10 that Cornitermes cumulans worker termites indeed exhibit power-law waiting times when performing exploratory behaviour in navigational clueless environment, where no food or nestmates are present. Waiting times were calculated from the termite walking time series when the isolated termite workers did not changed position. Before walking again a termite may stop and wait for as long as 0.5 seconds to as much as 200 seconds. The log-log plot of the waiting time data, for the example shown, resulted in a negative power-law function with an exponent l~2:64; other measured values were 2.28, 2.52 and 2.96, with SlT+sd = 2.6060.28.

Turning angle distribution
Inside the entangled network of tunnels, termites follow well defined routes for traveling between different places of interest, from the foraging areas, to the royal chamber, etc. Here we are interested in exploring to what extent there is persistence in the direction of consecutive steps in termites that would reveal, for example, if termites move exploring space in preferential angles or not. This is because it has been speculated that they may be doing movements with angles of about 40-60 degrees as an optimal searching strategy that may be related to the branching angle of termite tunnels [34]. When exploring the distribution of turning (B) An artificially correlated time-series generated with a relaxation return map x tz1~g x t z(1{g)E t where E t is a normally distributed random variable with zero mean and unit variance and g is real valued parameter whose value determines the colour of the resulting time-series scaling. Colour in this context means the classification of a noise mode 1=f b , where f is the frequency in a Fourier transformed space. b is the scaling exponent and when b~0, the process is uncorrelated white noise (C), b~1 is correlated pink noise (D) and b~2 is a correlated brown noise (E). Termite walking (A) lies in between a pink (D) and brown noise (E) scaling, being compatible with the fact that the termite scaling exponent in the FFT is b~1:3. For details on how the IFS algorithm operates, see [30]. doi:10.1371/journal.pone.0111183.g008 Figure 9. Long-range correlations. Termite walking exhibit powerlaw decaying long-range correlations as measured by a correlation function along the walk time-series (blue). An artificial correlated timeseries, as explained in Fig. 8 was used also to compare a correlated decaying process (red). The black line is a power-law with a scaling exponent c~{0:3. doi:10.1371/journal.pone.0111183.g009 angles between successive steps, it is found here that Cornitermes cumulans workers in open spaces do not show evidence of prewired preferential angles. There is a directional trend to move forwards but the resulting probability distribution is a uniform bell shaped curve centered in 0 degrees, as seen in Fig. 11. The absence of preferential angles of around 40-60 degrees may imply that the geometry of the tunnels may be either the outcome of social interactions or a response of termites to terrain clues, or a combination of both.

Discussion
The study of animal movements is of prime importance for understanding ecological and behavioural traits of individual displacements, needed to the efficient use of space. These  movements in a social context are extremely important because these regulate the rate of interactions that are the basis for building self-organization global patterns of behaviour [6] which are, ultimately, essential to cooperation phenomena and social stasis. Social individuals exchange information related to both the status of individuals and the social group and this information flows between individuals who move and encounter each other.
Our results show that free termite displacement, in space and time, present self-similar patterns highly consistent with anomalous diffusion. In social animals interacting nestmates are valuable resources to look for, and that is achieved by moving about in an efficient manner. Isolated termite displacement patterns revealed in this study, add support for the hypothesis that animals adopt Lévy flights because these confer an advantage in terms of higher fitness resulting from greater efficiency in finding nestmates to whom socially interact. The current trend that highlights statistical model fitting has perhaps put to much attention in fitting methodologies, which are important of course, but other qualitative and quantitative aspects of the phenomena should not remain overlooked, specially the biological phenomenology. Superdiffusion, Lévy flights and waiting times are quantitatively interconnected in ways that must be explored further. It is intriguing that passive trackers in fluid dynamics do behave qualitatively and quantitatively the same way that termites do. Superdiffusion is a movement pattern observed in a wide range of organisms thought to be related to optimal search strategies. Here we show that termites perform superdiffusive displacements, suggesting that this movement pattern is important for developing social interactions upon which self-organization at the colony level is built. Our results suggest that scale-free superdiffusion may have played an important role in the evolution of societies and cooperation because it may confer an advantage in terms of fitness to individuals that profit from efficient social encounter rates, that translates into efficient information transfer.
What is remarkable in our study, is the fact that assayed termites exhibited superdiffusion and Lévy flights even in the absence of external stimuli or clues: the arenas contained no food nor other nestmates. It seems, therefore, that such a movement pattern is indeed hardwired in termite ''instincts''. Or recalling a long forgotten idea in termitology, termites orient themselves concatenating external stimuli with some form of internal momentum [9]. The important connotation arising from such a result is that termite diffusion can proceed independently of a reactive phasis, as expected for a walker with limited (if any) cognitive abilities [35]. More specifically, blind termites seem to be equipped with strategies which secure them to find their targets (food or nestmates) even in environments where these are scarce or criptically located.

Supporting Information
Datasets S1 Four examples of termite displacement data in x,y format. (ZIP)