- Email: [email protected]
Optimal cooperative searching using purely repulsive interactions
Journal of Theoretical Biology 361 (2014) 159–164
Contents lists available at ScienceDirect
Journal of Theoretical Biology journal homepage: www.elsevier.com/locate/yjtbi
Optimal cooperative searching using purely repulsive interactions Noriyuki P. Tani, Alan Blatt, David A. Quint n, Ajay Gopinathan Department of Physics, University of California Merced, United States
H I G H L I G H T S
We explore how repulsive interactions could enhance cooperativity in groups of diffusing foragers. We find an optimal repulsive range, where the time to find a resource/target is minimal. Optimal repulsion strikes a balance between minimizing redundancy and maximizing total area searched. Our results are insensitive to the exact form of the repulsive interaction and scale in a simple manner with forager density.
art ic l e i nf o
a b s t r a c t
Article history: Received 7 November 2013 Received in revised form 24 July 2014 Accepted 25 July 2014 Available online 2 August 2014
Foraging, either solitarily or collectively, is a necessary behavior for survival that is demonstrated by many organisms. Foraging can be collectively optimized by utilizing communication between the organisms. Examples of such communication range from high level strategic foraging by animal groups to rudimentary signaling among unicellular organisms. Here we systematically study the simplest form of communication via long range repulsive interactions between multiple diffusing Brownian searchers on a one-dimensional lattice. We show that the mean first passage time for any one of them to reach a fixed target depends non-monotonically on the range of the interaction and can be optimized for a repulsive range that is comparable to the average spacing between searchers. Our results suggest that even the most rudimentary form of collective searching does in fact lower the search time for the foragers suggesting robust mechanisms for search optimization in cellular communities. & 2014 Elsevier Ltd. All rights reserved.
Keywords: Collective behavior Foraging
1. Introduction Understanding the process of searching or foraging in living systems has been of great interest in many disciplines, such as biology, physics, computer science, and robotics (Viswanathan et al., 2011). The mechanisms by which different organisms forage for food can be quite varied, for example bears and wolves use their sense of smell in order to acquire food (Mattson, 2005), while bats and dolphins use echolocation to locate their food (Au and Snyder, 1980; Schnitzler and Kalko, 2001). Some animals have the ability to search for food individually; however, many other organisms must work in tandem in order to efficiently find food, such as ants and fish (Ioannou et al., 2012; Jackson and Ratnieks, 2006). Studying collective foraging patterns in nature can reveal basic algorithmic features that can be directly compared with artificial searching algorithms used in computer science and robotics (Shoghian and Kouzehgar, 2012). This type of analysis can help animal behavioral scientists and n
Corresponding author. Tel.: þ 1 209 228 4048. E-mail addresses: [email protected] (D.A. Quint), [email protected] (A. Gopinathan). http://dx.doi.org/10.1016/j.jtbi.2014.07.027 0022-5193/& 2014 Elsevier Ltd. All rights reserved.
computer scientists understand how these algorithms evolved over time and became robust over the wide range of environmental scenarios. One application of these searching process is currently used in robotics, where robots can utilize collective searching motifs that help them to navigate unexplored terrain and also assist in search and rescue efforts (Saeedi et al., 2009; Ko and Lau, 2009; Reich and Sklar, 2006; Hoff et al., 2010). Cooperation among both living systems and artificial ones strives for the same goals, such as minimizing the search time (i.e minimize energetic cost) while maximizing the search space. Collective foragers or searchers, found in nature, display a high degree of coordination and communication within the collective as compared to a single searcher on its own. In fact, movement at the individual organism level within a collective is strongly correlated with the information that is being derived from their surrounding neighbors. There are several biological systems that generically display foraging behavior at the individual level; however the type of foraging strategies that these various systems utilize span a wide range of spatial and temporal scales, such as sub-diffusive, diffusive and super-diffusive (Bartumeus et al., 2002; Golé et al., 2011; Seuront and Stanley, 2013).
160
N.P. Tani et al. / Journal of Theoretical Biology 361 (2014) 159–164
A major theme arose from the suggestion that Lèvy walks can be seen in the foraging ants because their fractal and superdiffusive properties are likely to be advantageous when searching (Shlesinger and Klafter, 1986). Analytical and experimental studies (Schuster and Levendowsky, 1996; Levendowsky et al., 1997; Viswanathan and et al., 1999; Bartumeus et al., 2002; Humphries et al., 2012; Raichlena and et al., 2013) subsequently showed that Levy walks are indeed advantageous when searching for randomly distributed resources, and are typically more effective than other simple search patterns, such as Brownian walks and ballistic (straight-line) movements. Natural selection may therefore have led to adaptations for Levy searching. This prediction is now supported by a wide variety of empirical studies (Viswanathan et al., 2011). However, it has been recently shown that such strategies are only an optimal modality of foraging for certain environmental conditions and can also depend on the species under study (Seuront and Stanley, 2013). In the case of “blind searchers” (e.g. basking sharks, jellyfish predators, leatherback turtles, and southern elephant seals), if one considers discrete Lèvy flights and not continuous walks, the optimal search strategy depends on the location of the intended target (food), and for targets that are in close proximity to the searcher Brownian motion is the optimal search strategy (Palyulina et al., 2013). While Lèvy strategies have enhanced efficiency for individual foraging, a wide range of organisms ranging from single cells to mammals perform movements that are predominantly modeled by correlated or persistent random walks, which, at large spatiotemporal scales, become Gaussian or Brownian random walks (Turchin, 1998). One question that arises in this context is whether such Brownian searchers can collectively improve their efficiency by sharing information via interacting in some manner. An example model system, which spans both types of foraging behavior (solitary or collective), is the eukaryotic cell Dictyostelium discoideum (Dicty). Solitary Dicty can search their environment in a persistent Brownian random walk when food sources are high (Shenderov and Sheetz, 1997; Selmeczi et al., 2005). On the contrary, Dicty can form colonies that can search their environment collectively by utilizing chemotactic signaling (Keating and Bonner, 1977; Konijn et al., 1968; Pan et al., 1972). These chemotatic chemical interactions, which tend to correlate the motions of individual Dicty, allow for collectively foraging that can optimize foraging efforts, such that the entire colony benefits (Gelbart, 2010). At much larger length scales and with higher quality of information, Mongolian gazelles perform movements that are well approximated by simple Brownian motion but enhance their collective efficiency in finding sparse patches of vegetation by calling to each other when they find food (MartinezGarcia et al., 2013). Intriguingly the study finds that there is an optimal intermediate range of communication that maximizes efficiency. Furthermore, this enhancement via communication extends qualitatively to Lèvy foragers as well but with significantly less quantitative improvement (Martinez-Garcia et al., 2014). While high-quality information sharing can lead to higher group efficiency, an intriguing question is whether there exist rudimentary or minimal interactions during the search process that can speed up the collective search for resources. For this to happen, it would be desirable for the interactions to minimize redundancy in searching while still maintaining sufficient exploration. We hypothesized that a simple and easily implementable interaction with this property that would prove advantageous for Brownian searchers, such as the ones mentioned above, would be repulsion or mutual avoidance. In this study, we examine how Brownian foragers cooperate with each other via purely repulsive interactions, while searching for a single target (e.g. food). We ask how the search time depends on the density of searchers, the range of repulsion and the functional form of the repulsion.
We address these questions by simulating multiple Brownian particles that search for a fixed target on a closed one-dimensional lattice. In the simplest case we first study this system without any interaction between two searchers and then compare this with the more complex system of two interacting Brownian particles by measuring the average mean first passage time (MFPT) to the target (Cepa and Lepingle, 2001; Sokolov et al., 2005). We also studied systems with three searchers and different forms of the interaction potential. We found that interactions among the searchers affected the search time, and an optimal repulsive range for foraging was found. This suggests that in order to optimize collective foraging, organism should interact such that they minimize redundant search patterns and maximize the search area in their environment. In Section 1, we discuss our model and the dynamics of our simulation. In Section 2 we present our results; in Section 3.1, we compare the MFPT of three different systems; one searcher, two searchers without interactions two, searchers with interactions as well as varying the form of the interaction potential; in Section 3.2, we present, by dimensional analysis, the relationship between the optimal repulsive strength and the lattice size; In Section 3.3, we present the relationship between the average encounter time and the lattice size and discuss the extension of these results to the case with three interacting searchers. In Section 3, we discuss the implications of our results.
2. Model and simulation We study a discrete system consisting of two interacting Brownian searchers (random walkers) that move along a onedimensional periodic lattice with N sites (Fig. 1). Initialization of both searchers and the target is selected from a uniform random distribution, such that the domain of the distribution corresponds to the lattice size, N. The dynamics of this model are such that the bare diffusion constant for both searchers when they are not interacting is D¼
a2 ¼ 1; T
ð1Þ
where the lattice spacing is a¼ 1. Repulsive interactions are considered only between the two mobile searchers. Specifically we use an inverse power of the distance between the two searchers similar in form to an electrostatic potential between two like charges. The form of the potential in general is V ¼ α=r γ and for our simulations γ ¼ 1. The range of the potential is set by the parameter α and in our simulations it is given in terms of the lattice spacing a, which can range from 0 to 2N. The distance between the two searchers is
Fig. 1. The pictorial representation of our simulation model. S1 and S2 are the positions of the two searchers (red circles), the green square is the target and the black dotted circles are the repulsive boundary set by the value of α. r1 and r2 are the distances between the two searchers, and N=2π is the effective radius of periodic system. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)
N.P. Tani et al. / Journal of Theoretical Biology 361 (2014) 159–164
given as
12
r 1 ¼ jS1 S2 j;
ð2Þ
r2 ¼ N r1 :
ð3Þ
8
Here the j j represents the absolute value and the subscripts refer to the relative distance between the searchers on either side of the periodic boundaries (Fig. 1). The distance to the target relative to either of the searchers is only used to end the simulation when one of the searchers “finds” the target. Once the target is found we record the time (i.e. the number of time steps) it took for that specific realization of the simulation to end. Searchers' positional updates were found by evaluating the energy difference between the current system configuration and a randomly chosen proposed new configuration of the system (Si ðt þ 1Þ ¼ Si ðtÞ 7 a). In this way the system has no memory beyond the previous time step. The energies of either configuration are given as
6
Ei ¼ Ef ¼
α r 1;i
α r 1;f
þ
161
α α r 2;f
4 2 0
200
400
600
800 1000 1200 1400 1600 1800
Fig. 2. Semi-log plot of the frequency of search times for the single searcher (red solid line), two non-interacting searchers (green solid line) and two interacting searchers (blue solid line). The slope of the fitted dashed lines is the MFPT for each 1 system and for the single searcher is τ1;sim ¼ 0:002, for the two non-interacting 1 1 searchers τ2;sim ¼ 0:004 and for the two-interacting searchers is τ2i;sim ¼ 0:0005. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)
ð4Þ
r 2;i
þ
10
;
ð5Þ
where i refers to the previous configuration before the searcher's positional update and f refers to the new proposed configuration. Once these energies are calculated we then calculate the probability for the change in the configuration of the system using a Boltzmann distribution: P ¼ expð ðEf Ei ÞÞ
ð6Þ
Employing the METROPOLIS Monte Carlo (MMC) method the value of P is then compared to a random number s, which is drawn from a uniform distribution between 0 and 1. If s o P, then the proposed new configuration of the system is accepted, which corresponds to a lower energy of the system. This procedure is carried out until the target is found. In light of the fact that will be measuring the MFPT of our system, the choice of MMC over the usage of kinetic Monte Carlo (kMC) is worth a brief discussion. Although kMC explicitly measures the dynamical transition rates of our simulation, which are related to the real passage times of the system, in this case, time as calculated via kMC turns out to be equivalent to the actual computer simulation time, measured in MC simulation iteration steps, for sufficiently long simulation times . All statistical quantities that are presented here were computed by averaging over many initial conditions of the system for a fixed lattice size.
3. Results A single particle diffusing on a closed d-dimensional space has been studied in great detail. Here we focus on a random walk of two particles on a 1d closed lattice, or a ring where the two particles interact via a repulsive interaction.
their repulsive energies are of order N/2. In Fig. 2 we plot the logarithm of the frequency of target encounters versus the time that was taken to reach the target (first passage times) for all three cases and find that indeed the repulsive random walkers find the target faster (pink dashed line). In the first two cases, for the single searcher and two non-interacting searchers, we should expect that the mean first passage time to the target's position, averaged over all initial starting positions of the target, to be proportional to the square of the system size for large lattice sizes. More precisely we should expect for finite lattices that the mean number of steps taken before the target is found is Montroll (1969) 〈n〉 ¼ NðN þ 1Þ=6: In our simulations, for the case when the searchers are noninteracting, we have that the number of searcher steps taken equals the number of time steps in our simulation before the target is found, hence 〈τ 〉 ¼ 〈n〉. In Fig. 2 we find an excellent agreement with our numerical simulation for both the single random walker and two non-interacting walkers, which yield a theoretical MFPT for a system size of N ¼ 50 as 〈τ〉1 1 ¼ 0:002 and 〈τ〉2 1 ¼ 0:004 (black dashed lines). We find in general that the addition of a repulsive interaction lowers the mean first passage time in comparison to the noninteracting cases. We do expect that as we approach the limit where α goes to zero we should recover the mean first passage time for the two non-interacting searchers. However, in the limit that 2α=N becomes comparable to the thermal noise in the system (kb T ¼ 1), we should expect that either searcher exhibits slower effective diffusive motion (number of steps2 =total time o 1). This implies that the mean first passage time to the target should increase, thus we should expect to find an optimal repulsive strength where the searchers discover the target in the least amount of time. 3.2. Results II – repulsive strength and lattice size
3.1. Results I – distribution of first passage times To compare our results with previous studies for non-interacting Brownian searchers on a 1d ring, we calculate numerically the distribution of first passage times for a single searcher, two noninteracting searchers and two interacting searchers to reach the target over all random realizations of the target and searcher starting positions. Intuitively, one should expect that the time for repulsive searchers, on average, requires less time than the noninteracting cases. Since the repulsive interaction encourages the two searchers to avoid covering the same locations on the lattice, when
To test the assertion that there exists an optimal repulsive strength, we performed simulations over the entire range of repulsive strengths such that the repulsive energy goes from α ¼ ½0; 2N. Fig. 3a shows ensemble averaged data for the mean first passage time as a function of the repulsive strength for two different lattice sizes (see inset of Fig. 3a). In each case we find that there is indeed a minimum mean first passage time that occurs at an optimal repulsive strength near α N=2. Increasing the repulsion leads to an increase in the mean first passage time and a decrease in the effective diffusion rate for
162
N.P. Tani et al. / Journal of Theoretical Biology 361 (2014) 159–164
Fig. 3. (a) The MFPT over the range of repulsive strengths α ¼ ½0; 2N for a lattice size of N ¼ 50. The red marker indicates approximately where the minimum MFPT occurs at a point α⋆ ((a) The MPFT verses repulsive strength for a system size of N ¼ 150.) b() Changing the form of the repulsive interaction (Eqs. (4) and (5)) potential γ ¼ ½1=2 red; 1 blue; 2 black. The horizontal axis has been re-scaled so that each curve collapses indicating the results in (a) are independent of the exponent of the potential. All three curves are for a lattice size of N ¼ 50. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)
either of the two searchers. The results in Fig. 3a are also independent of the form of the repulsive interaction as seen in Fig. 3b, where we have varied the exponent of the interaction potential to range from 1/2 to 2. The horizontal axis from Fig. 3a has been re-normalized by a factor rmin, which is the range when the searchers begin to sense each other and is given by the following equation: 2 31=ðγ þ 1Þ 6 ð4γαδÞ 7 7 r min ðα; γ Þ 6 4 4γαδ 5 1 γ þ1 N
ð7Þ
This minimum interaction distance in Eq. (7) is derived from considering the case when the two searchers are within a range such that their interaction energy (Eqs. (4) and (5)) is of the same order as the fluctuations in their motion (O 1), then expanding to first order in their displacement away from this minimum range (i.e. the lattice spacing δ ¼ 1). After re-scaling the repulsive minima in Fig. 3b we find that r min ðα⋆ ; γ Þ const. It is interesting to note that the data collapse under rescaling the repulsive range by rmin implies that the finding time to the target depends only on this length scale, i.e. the form of the potential is not as important as the range of the interaction. At low values of the repulsion strength we find that the trajectories for each searcher overlap frequently (Fig. 4a). As the repulsion strength is increased we find that the amount of overlap between the two searchers' trajectories becomes significantly less (Fig. 4b) near the optimum shown in Fig. 3a. Increasing the repulsive interaction past the optimum we find that any overlap of the two trajectories is suppressed (Fig. 4c). In this regime, the energy associated with such a large value of α makes any over lap unfavorable during the simulation. To understand this, we can examine the allowed update modes to the system (there are four in total) as the repulsion is increased. We should expect that allowed movements of the two searchers become restricted to particular set of modes allowed by the partition function. There are two translational modes, where the change in the distance between the two searchers is unchanged Δr i ¼ 0, hence the energy of the system is unchanged. In contrast there exist two modes where the energy of the system changes, such that the searchers move apart or move closer with respect to each other, Δr i a0. For values of repulsion near or above the optimum, modes where Δr i a 0 are unfavorable. In Fig. 4d–f we show the frequency of searcher separations using the relative coordinate, r1 modulo N, for a given repulsive
strength (α ¼ ½0; N=2; N; respectively). For the case with no repulsion we find that all possible particle separations are possible, hence the flat distribution in Fig. 4d. As we turn up the repulsive interaction (α ¼ N=2 (Fig. 4e) and α ¼ N (Fig. 4f)) we find that the distribution of allowed particle separations is limited to only a set that is maximal at r 1 ¼ 0. This observation suggests that for repulsive strengths that are near or above the optimal the searchers are more likely to be found at diametrically opposed locations on the lattice with very little over lap of their search area, which is consistent with Fig. 4a and b. 3.3. Results III – average encounter time and lattice size To shed light on the reason for the optimal value of the repulsive interaction, one can use dimensional analysis to determine the location optimal value of α as a function of the system size. There are two important length scales in our system, one of which is the lattice size N and the other is determined by the diffusion constant D (length2/time). We posit that the critical (optimal) value of the repulsive energy is a function of both these parameters such that α⋆ ¼ f ðN; DÞ. Since the repulsive energy is dimensionless (Eqs. (4) and (5)), α is defined to have units of length. The bare diffusion of the two searchers is defined in terms of the underlying length and time scales of our simulation and is unity. Therefore we expect the form of the scaling function to be linear in the lattice size
α⋆ ¼ bN; where b is a unit-less constant. In Fig. 5a we plot the location of the optimal value of α versus the system size N. Fitting the data using the linear fitting function suggested by our dimensional scaling analysis, we find an excellent agreement between the data and the scaling function (red dashed line) proposed by our dimensional analysis. The fitting free parameter b is found to be about 0.6 which is very close to what we expect intuitively since at the value of α ¼ N=2 is where the two searchers begin to interact at an energy scale where the probability for accepting/rejecting a proposed move is of order e 1 . This optimal repulsive strength indicates that the interaction between the searchers are such that both do not diffuse over the same space (Fig. 4b) as doing so would be energetically unfavorable. To further solidify that our intuition is correct in regard to the scaling of the two searcher case, we also looked at the three searcher system (Fig. 5a, green squares). We find that indeed the system picks the minimum between all three searchers
N.P. Tani et al. / Journal of Theoretical Biology 361 (2014) 159–164
163
Fig. 4. (a) The frequency of lattice position occupation for both searchers (red-S1) and (blue-S2) for α ¼ 0 displays a high value of overlap. (b) At the optimal value of the repulsion strength α⋆ N=2 we find that the overlap between the two searchers is minimal while maintaining a good search coverage over the lattice. (c) Past the optimal repulsive strength we find that while the overlap is minimal, the search coverage of the lattice is poor because the searchers are more localized. (d–f) The frequency of particle separation as a function of the absolute value of the relative coordinate jS1 S2 j. For the non-interaction case all searcher separations are equally likely, consistent with the overlap shown in (a). In contrast, at the optimal and large repulsion values (as in b–c) we find that the most likely separation is when S1 S2 , which is consistent with the optimal repulsion value of α N=2. All data shown was for a lattice size of N¼ 200. For (a–c) simulations were run for 105 time steps. For (d–f) simulations were run for 106 time steps. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)
to be roughly N/3, thus for a lattice size of N ¼150 we find that α 50, which is the result for the two searcher case for N ¼100. These results suggest that the minimal searcher number that needs to be considered is two. In addition to the existence of an optimal repulsive interaction, which lowers the average mean first passage time to the target, we also expect there will be an average (minimal) encounter time to the target, τ . This characteristic time is related to the intrinsic length and time scales that define our system and should also be related to other relevant quantities, like the bare diffusion constant, D. Using dimensional analysis, we can write down a scaling ansatz similar to the one for the optimal repulsive strength α⋆ . Since we are interested in a quantity that has units of time, the only two parameters that the mean encounter time can depend on is the lattice size and the bare diffusion rate, T ¼ f ðN; DÞ. Given that the diffusion constant has units 2 of ½Length nTime 1 we expect that the dependence of τ to scale as the square of the system size. Therefore we expect that
τ ¼ an
N2 : D
In Fig. 5b we find the average mean first encounter time to the target as a function of the lattice size squared (blue points). Fitting this data with the scaling function derived above we found an excellent agreement between the data and the best fit line (red dashed line).
4. Discussion Foraging strategies can vary among species of animals from bacteria to vertebrates. Foraging strategies could in principle vary across a single species, which could allow for better target tracking as well as coping with changing environmental conditions
(Bartumeus et al., 2002; Seuront and Stanley, 2013). Aside from solitary foraging there also exist animal groups that forage collectively. The same principles of optimization of searching strategies will also apply to these systems as well, but will be conditional upon the collective migration aspect of foraging, which will in general depend on the type of signaling that is used between members of the collective. In summary, we have shown that for multiple Brownian searchers foraging for a single target that the mean first passage time (MFPT) to the target is minimized when there is a mutually repulsive interaction between them and furthermore that there is an optimal range for this repulsion that depends on the density of searchers. These results suggest that the mutually repulsive interaction introduces a cooperative effect that allows the two searchers to optimize both their search time and the search area. The cooperativity in this system is established by the mutual repulsion that each particle exerts on the other, as their motion tends to become correlated and redundancy is minimized when the repulsive range is large enough. However if the repulsive range is too large, the searchers influence each other too strongly and effective search area of each searcher becomes more localized. The balance dictates the optimum range which is roughly the average separation of searchers. In a sense, communication between particles is mediated by the repulsive interaction and the flow of information about their relative positions on the lattice only occurs when their repulsive zones are overlapping and the optimal point is when they are just overlapping. Although we have focused on Brownian dynamics in our model, it should be noted that in many biological systems it has been found that foraging at the individual level can be optimal for other step-size distributions (i.e. Lèvy flights). However, a wide range of organisms ranging from single cells, like Dicty and
164
N.P. Tani et al. / Journal of Theoretical Biology 361 (2014) 159–164
poor environment will remain high (Gelbart, 2010). Similar considerations apply to the behavior of mongolian gazelles that call to each other when they find food (Martinez-Garcia et al., 2013). While the context and mode of interaction are different, it is interesting to note that the study also finds that there is an optimal intermediate range of communication that maximizes efficiency. This raises an interesting evolutionary question in that not only do organisms have to utilize some intercommunication signal to optimize their resource gathering but that interaction must all be optimally tuned in order for collective foraging to always promote a high fitness value in a wide range of environmental conditions.
Acknowledgments The authors would like to thank both Anatoly Kolomeisky and K.C. Huang for their insightful comments and discussions. This work was partially supported by National Science Foundation Grant EF-1038697 and a James S. McDonnell Foundation Award. Undergraduate support was funded by the Undergraduate Research and Mentoring (URM) Program Sponsored by the National Science Foundation NSF Grant DBI-1040962.
References
Fig. 5. (a) Plot of optimal repulsive range αn versus lattice size N for both two searchers (blue circles) and three searchers (green squares). Blue and green bars represent the error and red line is the best fit line α⋆ ¼ 0:5755 7 0:05nN weighted by the error for each data point shown. Error bars were calculated using data points that were within 2% of the minimum (Fig. 3) for each lattice size. (b) Graph of N2 and τ . Blue bars represent the error for each data point and the red line is the best fit line using to our scaling ansatz, which gives τ ¼ 0:0507nN 2 . (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)
bacteria, to mammals, like gazelles, perform movements that are predominantly modeled by correlated or persistent random walks, which, at large spatiotemporal scales, become Gaussian or Brownian random walks (Turchin, 1998). Our results indicate that such Brownian searchers can derive a significant advantage in searching by implementing an optimal range of repulsion. Optimal spacing of foragers via interactions through chemical signaling has been noted in real systems of slime molds and flour beetles (Naylor, 1959, 1961, 1965; Bonner, 1963). This type of behavior in biological systems has been thought to be needed for maximal utilization of food sources as well as searching for good habitats (Shorey, 1976; Martinez-Garcia et al., 2013) and our results show that these interactions can be optimally tuned. Minimizing the search time by each searcher translates to real systems as a minimization of energy consumed by each searcher when foraging for food. This is certainly an important consideration in many biological systems. For example, when food supplies are low single cell Dicty begins to produce chemo-repellent in places that it has searched for food (Keating and Bonner, 1977; Konijn et al., 1968; Pan et al., 1972). This has the effect of informing other single Dicty cells in the population to avoid these areas as a way to conserve energy. This strategy allows the entire population to search for resources collectively as opposed to when food supplies are abundant, where each cell is free to search randomly, although not pure diffusive as in the case for Dicty. From this point of view it is easy to see why having an antagonist interaction among foragers is evolutionarily advantageous. By minimizing energy consumption Dicty can ensure that its fitness in a nutrient
Au, W.W.L., Snyder, K.J., 1980. J. Acoust. Soc. Am. 68, 1077. Bartumeus, F., Catalan, J., Fulco, U.L., Lyra, M.L., Viswanathan, G.M., 2002. Phys. Rev. Lett. 88, 097901. Bonner, J.T., 1963. Sci. Am. 209, 84. Cepa, E., Lepingle, D., 2001. ESAIM: Probab. Stat. 5, 203. Gelbart, M., 2010. Foraging Strategies for Dictyostelium discoideum. Princeton University, Princeton, NJ, USA. Golé, L., Rivière, C., Hayakawa, Y., Rieu, J.-P., 2011. PLoS ONE 6, e26901. Hoff, N.R., Sagoff, A., Wood, R.J., Nagpal, R., 2010. In: IEEE International Conference on Robotics and Biomimetics, vol. 123. Humphries, N.E., et al., 2012. Proc. Natl. Acad. Sci. 109, 7169–7174. Ioannou, C.C., Guttal, V., Couzin, I., 2012. Science 337, 1212. Jackson, D.E., Ratnieks, F.L., 2006. Curr. Biol. 16, R570. Keating, M.T., Bonner, J.T., 1977. J. Bacteriol. 130 (1), 144. Ko, A., Lau, H.Y., 2009. Int. J. Adv. Sci. Technol. 3, 69. Konijn, T., Barkley, D., Chang, Y., Bonner, J., 1968. Am. Natural. 102, 225. Levendowsky, M., White, B.S., Schuster, F., 1997. Acta Protozool. 36, 237. Martinez-Garcia, R., et al., 2014. Phys. Rev. E 89, 032718. Martinez-Garcia, R., et al., 2013. Phys. Rev. Lett. 110, 248106. Mattson, D.J., 2005. USGS. Montroll, E.W., 1969. J. Math. Phys. 10, 753. Naylor, A.F., 1959. Ecology 40, 453. Naylor, A.F., 1961. Ecology 42, 231. Naylor, A.F., 1965. Ecology 46, 341. Palyulina, V.V., Chechkinb, A.V., Metzler, R., 2013. Proc. Natl. Acad. Sci. 111, 2931. Pan, P., Hall, E., Bonner, J., 1972. Nat.: New Biol. 237 (75), 181. Raichlena, D.A., et al., 2013. Proc. Natl. Acad. Sci. 111, 728. Reich, J., Sklar, E., 2006. In: Proceedings of the IEEE International Workshop on Safety, Security and Rescue Robotics. Saeedi, P., Sorensen, S.-A., Hailes, S., 2009. In: Safety, Security and Rescue Robotics (SSRR), IEEE International Workshop. Schnitzler, H.U., Kalko, E.K.V., 2001. BioScience 51 (7), 557. Schuster, F., Levendowsky, M., 1996. J. Euk. Microbiol. 43, 150. Selmeczi, D., Mosler, S., Hagedorn, P.H., Larsen, N.B., Flyvbjerg, H., 2005. Biophys. J. 89, 912. Seuront, L., Stanley, H.E., 2013. Proc. Natl. Acad. Sci. 111, 2206. Shenderov, A., Sheetz, M., 1997. Biophys. J. 72, 2382. Shlesinger, M.F., Klafter, J., 1986. On Growth and Form: Lévy Walks versus Lévy Flights. Springer, Netherlands. Shoghian, S., Kouzehgar, M., 2012. World Acad. Sci. Eng. Technol. 72, 447. Shorey, H.H., 1976. Animal Communication by Pheromones. Academic Press, Waltham, Massachusetts, USA. Sokolov, I.M., Metzler, R., Pant, K., William, M.C., 2005. Phys. Rev. E 72, 041102. Turchin, P., 1998. Quantitative Analysis of Movement: Measuring and Modeling Population Redistribution in Plants and Animals. Sinauer Associates, Sunderland, MA, USA. Viswanathan, G.M., et al., 1999. Nature 401, 911. Viswanathan, G.M., et al., 2011. The Physics of Foraging: An Introduction to Random Searches and Biological Encounters, 1st ed. Cambridge University Press, Cambridge, England, p. 2011.
Contents lists available at ScienceDirect
Journal of Theoretical Biology journal homepage: www.elsevier.com/locate/yjtbi
Optimal cooperative searching using purely repulsive interactions Noriyuki P. Tani, Alan Blatt, David A. Quint n, Ajay Gopinathan Department of Physics, University of California Merced, United States
H I G H L I G H T S
We explore how repulsive interactions could enhance cooperativity in groups of diffusing foragers. We find an optimal repulsive range, where the time to find a resource/target is minimal. Optimal repulsion strikes a balance between minimizing redundancy and maximizing total area searched. Our results are insensitive to the exact form of the repulsive interaction and scale in a simple manner with forager density.
art ic l e i nf o
a b s t r a c t
Article history: Received 7 November 2013 Received in revised form 24 July 2014 Accepted 25 July 2014 Available online 2 August 2014
Foraging, either solitarily or collectively, is a necessary behavior for survival that is demonstrated by many organisms. Foraging can be collectively optimized by utilizing communication between the organisms. Examples of such communication range from high level strategic foraging by animal groups to rudimentary signaling among unicellular organisms. Here we systematically study the simplest form of communication via long range repulsive interactions between multiple diffusing Brownian searchers on a one-dimensional lattice. We show that the mean first passage time for any one of them to reach a fixed target depends non-monotonically on the range of the interaction and can be optimized for a repulsive range that is comparable to the average spacing between searchers. Our results suggest that even the most rudimentary form of collective searching does in fact lower the search time for the foragers suggesting robust mechanisms for search optimization in cellular communities. & 2014 Elsevier Ltd. All rights reserved.
Keywords: Collective behavior Foraging
1. Introduction Understanding the process of searching or foraging in living systems has been of great interest in many disciplines, such as biology, physics, computer science, and robotics (Viswanathan et al., 2011). The mechanisms by which different organisms forage for food can be quite varied, for example bears and wolves use their sense of smell in order to acquire food (Mattson, 2005), while bats and dolphins use echolocation to locate their food (Au and Snyder, 1980; Schnitzler and Kalko, 2001). Some animals have the ability to search for food individually; however, many other organisms must work in tandem in order to efficiently find food, such as ants and fish (Ioannou et al., 2012; Jackson and Ratnieks, 2006). Studying collective foraging patterns in nature can reveal basic algorithmic features that can be directly compared with artificial searching algorithms used in computer science and robotics (Shoghian and Kouzehgar, 2012). This type of analysis can help animal behavioral scientists and n
Corresponding author. Tel.: þ 1 209 228 4048. E-mail addresses: [email protected] (D.A. Quint), [email protected] (A. Gopinathan). http://dx.doi.org/10.1016/j.jtbi.2014.07.027 0022-5193/& 2014 Elsevier Ltd. All rights reserved.
computer scientists understand how these algorithms evolved over time and became robust over the wide range of environmental scenarios. One application of these searching process is currently used in robotics, where robots can utilize collective searching motifs that help them to navigate unexplored terrain and also assist in search and rescue efforts (Saeedi et al., 2009; Ko and Lau, 2009; Reich and Sklar, 2006; Hoff et al., 2010). Cooperation among both living systems and artificial ones strives for the same goals, such as minimizing the search time (i.e minimize energetic cost) while maximizing the search space. Collective foragers or searchers, found in nature, display a high degree of coordination and communication within the collective as compared to a single searcher on its own. In fact, movement at the individual organism level within a collective is strongly correlated with the information that is being derived from their surrounding neighbors. There are several biological systems that generically display foraging behavior at the individual level; however the type of foraging strategies that these various systems utilize span a wide range of spatial and temporal scales, such as sub-diffusive, diffusive and super-diffusive (Bartumeus et al., 2002; Golé et al., 2011; Seuront and Stanley, 2013).
160
N.P. Tani et al. / Journal of Theoretical Biology 361 (2014) 159–164
A major theme arose from the suggestion that Lèvy walks can be seen in the foraging ants because their fractal and superdiffusive properties are likely to be advantageous when searching (Shlesinger and Klafter, 1986). Analytical and experimental studies (Schuster and Levendowsky, 1996; Levendowsky et al., 1997; Viswanathan and et al., 1999; Bartumeus et al., 2002; Humphries et al., 2012; Raichlena and et al., 2013) subsequently showed that Levy walks are indeed advantageous when searching for randomly distributed resources, and are typically more effective than other simple search patterns, such as Brownian walks and ballistic (straight-line) movements. Natural selection may therefore have led to adaptations for Levy searching. This prediction is now supported by a wide variety of empirical studies (Viswanathan et al., 2011). However, it has been recently shown that such strategies are only an optimal modality of foraging for certain environmental conditions and can also depend on the species under study (Seuront and Stanley, 2013). In the case of “blind searchers” (e.g. basking sharks, jellyfish predators, leatherback turtles, and southern elephant seals), if one considers discrete Lèvy flights and not continuous walks, the optimal search strategy depends on the location of the intended target (food), and for targets that are in close proximity to the searcher Brownian motion is the optimal search strategy (Palyulina et al., 2013). While Lèvy strategies have enhanced efficiency for individual foraging, a wide range of organisms ranging from single cells to mammals perform movements that are predominantly modeled by correlated or persistent random walks, which, at large spatiotemporal scales, become Gaussian or Brownian random walks (Turchin, 1998). One question that arises in this context is whether such Brownian searchers can collectively improve their efficiency by sharing information via interacting in some manner. An example model system, which spans both types of foraging behavior (solitary or collective), is the eukaryotic cell Dictyostelium discoideum (Dicty). Solitary Dicty can search their environment in a persistent Brownian random walk when food sources are high (Shenderov and Sheetz, 1997; Selmeczi et al., 2005). On the contrary, Dicty can form colonies that can search their environment collectively by utilizing chemotactic signaling (Keating and Bonner, 1977; Konijn et al., 1968; Pan et al., 1972). These chemotatic chemical interactions, which tend to correlate the motions of individual Dicty, allow for collectively foraging that can optimize foraging efforts, such that the entire colony benefits (Gelbart, 2010). At much larger length scales and with higher quality of information, Mongolian gazelles perform movements that are well approximated by simple Brownian motion but enhance their collective efficiency in finding sparse patches of vegetation by calling to each other when they find food (MartinezGarcia et al., 2013). Intriguingly the study finds that there is an optimal intermediate range of communication that maximizes efficiency. Furthermore, this enhancement via communication extends qualitatively to Lèvy foragers as well but with significantly less quantitative improvement (Martinez-Garcia et al., 2014). While high-quality information sharing can lead to higher group efficiency, an intriguing question is whether there exist rudimentary or minimal interactions during the search process that can speed up the collective search for resources. For this to happen, it would be desirable for the interactions to minimize redundancy in searching while still maintaining sufficient exploration. We hypothesized that a simple and easily implementable interaction with this property that would prove advantageous for Brownian searchers, such as the ones mentioned above, would be repulsion or mutual avoidance. In this study, we examine how Brownian foragers cooperate with each other via purely repulsive interactions, while searching for a single target (e.g. food). We ask how the search time depends on the density of searchers, the range of repulsion and the functional form of the repulsion.
We address these questions by simulating multiple Brownian particles that search for a fixed target on a closed one-dimensional lattice. In the simplest case we first study this system without any interaction between two searchers and then compare this with the more complex system of two interacting Brownian particles by measuring the average mean first passage time (MFPT) to the target (Cepa and Lepingle, 2001; Sokolov et al., 2005). We also studied systems with three searchers and different forms of the interaction potential. We found that interactions among the searchers affected the search time, and an optimal repulsive range for foraging was found. This suggests that in order to optimize collective foraging, organism should interact such that they minimize redundant search patterns and maximize the search area in their environment. In Section 1, we discuss our model and the dynamics of our simulation. In Section 2 we present our results; in Section 3.1, we compare the MFPT of three different systems; one searcher, two searchers without interactions two, searchers with interactions as well as varying the form of the interaction potential; in Section 3.2, we present, by dimensional analysis, the relationship between the optimal repulsive strength and the lattice size; In Section 3.3, we present the relationship between the average encounter time and the lattice size and discuss the extension of these results to the case with three interacting searchers. In Section 3, we discuss the implications of our results.
2. Model and simulation We study a discrete system consisting of two interacting Brownian searchers (random walkers) that move along a onedimensional periodic lattice with N sites (Fig. 1). Initialization of both searchers and the target is selected from a uniform random distribution, such that the domain of the distribution corresponds to the lattice size, N. The dynamics of this model are such that the bare diffusion constant for both searchers when they are not interacting is D¼
a2 ¼ 1; T
ð1Þ
where the lattice spacing is a¼ 1. Repulsive interactions are considered only between the two mobile searchers. Specifically we use an inverse power of the distance between the two searchers similar in form to an electrostatic potential between two like charges. The form of the potential in general is V ¼ α=r γ and for our simulations γ ¼ 1. The range of the potential is set by the parameter α and in our simulations it is given in terms of the lattice spacing a, which can range from 0 to 2N. The distance between the two searchers is
Fig. 1. The pictorial representation of our simulation model. S1 and S2 are the positions of the two searchers (red circles), the green square is the target and the black dotted circles are the repulsive boundary set by the value of α. r1 and r2 are the distances between the two searchers, and N=2π is the effective radius of periodic system. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)
N.P. Tani et al. / Journal of Theoretical Biology 361 (2014) 159–164
given as
12
r 1 ¼ jS1 S2 j;
ð2Þ
r2 ¼ N r1 :
ð3Þ
8
Here the j j represents the absolute value and the subscripts refer to the relative distance between the searchers on either side of the periodic boundaries (Fig. 1). The distance to the target relative to either of the searchers is only used to end the simulation when one of the searchers “finds” the target. Once the target is found we record the time (i.e. the number of time steps) it took for that specific realization of the simulation to end. Searchers' positional updates were found by evaluating the energy difference between the current system configuration and a randomly chosen proposed new configuration of the system (Si ðt þ 1Þ ¼ Si ðtÞ 7 a). In this way the system has no memory beyond the previous time step. The energies of either configuration are given as
6
Ei ¼ Ef ¼
α r 1;i
α r 1;f
þ
161
α α r 2;f
4 2 0
200
400
600
800 1000 1200 1400 1600 1800
Fig. 2. Semi-log plot of the frequency of search times for the single searcher (red solid line), two non-interacting searchers (green solid line) and two interacting searchers (blue solid line). The slope of the fitted dashed lines is the MFPT for each 1 system and for the single searcher is τ1;sim ¼ 0:002, for the two non-interacting 1 1 searchers τ2;sim ¼ 0:004 and for the two-interacting searchers is τ2i;sim ¼ 0:0005. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)
ð4Þ
r 2;i
þ
10
;
ð5Þ
where i refers to the previous configuration before the searcher's positional update and f refers to the new proposed configuration. Once these energies are calculated we then calculate the probability for the change in the configuration of the system using a Boltzmann distribution: P ¼ expð ðEf Ei ÞÞ
ð6Þ
Employing the METROPOLIS Monte Carlo (MMC) method the value of P is then compared to a random number s, which is drawn from a uniform distribution between 0 and 1. If s o P, then the proposed new configuration of the system is accepted, which corresponds to a lower energy of the system. This procedure is carried out until the target is found. In light of the fact that will be measuring the MFPT of our system, the choice of MMC over the usage of kinetic Monte Carlo (kMC) is worth a brief discussion. Although kMC explicitly measures the dynamical transition rates of our simulation, which are related to the real passage times of the system, in this case, time as calculated via kMC turns out to be equivalent to the actual computer simulation time, measured in MC simulation iteration steps, for sufficiently long simulation times . All statistical quantities that are presented here were computed by averaging over many initial conditions of the system for a fixed lattice size.
3. Results A single particle diffusing on a closed d-dimensional space has been studied in great detail. Here we focus on a random walk of two particles on a 1d closed lattice, or a ring where the two particles interact via a repulsive interaction.
their repulsive energies are of order N/2. In Fig. 2 we plot the logarithm of the frequency of target encounters versus the time that was taken to reach the target (first passage times) for all three cases and find that indeed the repulsive random walkers find the target faster (pink dashed line). In the first two cases, for the single searcher and two non-interacting searchers, we should expect that the mean first passage time to the target's position, averaged over all initial starting positions of the target, to be proportional to the square of the system size for large lattice sizes. More precisely we should expect for finite lattices that the mean number of steps taken before the target is found is Montroll (1969) 〈n〉 ¼ NðN þ 1Þ=6: In our simulations, for the case when the searchers are noninteracting, we have that the number of searcher steps taken equals the number of time steps in our simulation before the target is found, hence 〈τ 〉 ¼ 〈n〉. In Fig. 2 we find an excellent agreement with our numerical simulation for both the single random walker and two non-interacting walkers, which yield a theoretical MFPT for a system size of N ¼ 50 as 〈τ〉1 1 ¼ 0:002 and 〈τ〉2 1 ¼ 0:004 (black dashed lines). We find in general that the addition of a repulsive interaction lowers the mean first passage time in comparison to the noninteracting cases. We do expect that as we approach the limit where α goes to zero we should recover the mean first passage time for the two non-interacting searchers. However, in the limit that 2α=N becomes comparable to the thermal noise in the system (kb T ¼ 1), we should expect that either searcher exhibits slower effective diffusive motion (number of steps2 =total time o 1). This implies that the mean first passage time to the target should increase, thus we should expect to find an optimal repulsive strength where the searchers discover the target in the least amount of time. 3.2. Results II – repulsive strength and lattice size
3.1. Results I – distribution of first passage times To compare our results with previous studies for non-interacting Brownian searchers on a 1d ring, we calculate numerically the distribution of first passage times for a single searcher, two noninteracting searchers and two interacting searchers to reach the target over all random realizations of the target and searcher starting positions. Intuitively, one should expect that the time for repulsive searchers, on average, requires less time than the noninteracting cases. Since the repulsive interaction encourages the two searchers to avoid covering the same locations on the lattice, when
To test the assertion that there exists an optimal repulsive strength, we performed simulations over the entire range of repulsive strengths such that the repulsive energy goes from α ¼ ½0; 2N. Fig. 3a shows ensemble averaged data for the mean first passage time as a function of the repulsive strength for two different lattice sizes (see inset of Fig. 3a). In each case we find that there is indeed a minimum mean first passage time that occurs at an optimal repulsive strength near α N=2. Increasing the repulsion leads to an increase in the mean first passage time and a decrease in the effective diffusion rate for
162
N.P. Tani et al. / Journal of Theoretical Biology 361 (2014) 159–164
Fig. 3. (a) The MFPT over the range of repulsive strengths α ¼ ½0; 2N for a lattice size of N ¼ 50. The red marker indicates approximately where the minimum MFPT occurs at a point α⋆ ((a) The MPFT verses repulsive strength for a system size of N ¼ 150.) b() Changing the form of the repulsive interaction (Eqs. (4) and (5)) potential γ ¼ ½1=2 red; 1 blue; 2 black. The horizontal axis has been re-scaled so that each curve collapses indicating the results in (a) are independent of the exponent of the potential. All three curves are for a lattice size of N ¼ 50. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)
either of the two searchers. The results in Fig. 3a are also independent of the form of the repulsive interaction as seen in Fig. 3b, where we have varied the exponent of the interaction potential to range from 1/2 to 2. The horizontal axis from Fig. 3a has been re-normalized by a factor rmin, which is the range when the searchers begin to sense each other and is given by the following equation: 2 31=ðγ þ 1Þ 6 ð4γαδÞ 7 7 r min ðα; γ Þ 6 4 4γαδ 5 1 γ þ1 N
ð7Þ
This minimum interaction distance in Eq. (7) is derived from considering the case when the two searchers are within a range such that their interaction energy (Eqs. (4) and (5)) is of the same order as the fluctuations in their motion (O 1), then expanding to first order in their displacement away from this minimum range (i.e. the lattice spacing δ ¼ 1). After re-scaling the repulsive minima in Fig. 3b we find that r min ðα⋆ ; γ Þ const. It is interesting to note that the data collapse under rescaling the repulsive range by rmin implies that the finding time to the target depends only on this length scale, i.e. the form of the potential is not as important as the range of the interaction. At low values of the repulsion strength we find that the trajectories for each searcher overlap frequently (Fig. 4a). As the repulsion strength is increased we find that the amount of overlap between the two searchers' trajectories becomes significantly less (Fig. 4b) near the optimum shown in Fig. 3a. Increasing the repulsive interaction past the optimum we find that any overlap of the two trajectories is suppressed (Fig. 4c). In this regime, the energy associated with such a large value of α makes any over lap unfavorable during the simulation. To understand this, we can examine the allowed update modes to the system (there are four in total) as the repulsion is increased. We should expect that allowed movements of the two searchers become restricted to particular set of modes allowed by the partition function. There are two translational modes, where the change in the distance between the two searchers is unchanged Δr i ¼ 0, hence the energy of the system is unchanged. In contrast there exist two modes where the energy of the system changes, such that the searchers move apart or move closer with respect to each other, Δr i a0. For values of repulsion near or above the optimum, modes where Δr i a 0 are unfavorable. In Fig. 4d–f we show the frequency of searcher separations using the relative coordinate, r1 modulo N, for a given repulsive
strength (α ¼ ½0; N=2; N; respectively). For the case with no repulsion we find that all possible particle separations are possible, hence the flat distribution in Fig. 4d. As we turn up the repulsive interaction (α ¼ N=2 (Fig. 4e) and α ¼ N (Fig. 4f)) we find that the distribution of allowed particle separations is limited to only a set that is maximal at r 1 ¼ 0. This observation suggests that for repulsive strengths that are near or above the optimal the searchers are more likely to be found at diametrically opposed locations on the lattice with very little over lap of their search area, which is consistent with Fig. 4a and b. 3.3. Results III – average encounter time and lattice size To shed light on the reason for the optimal value of the repulsive interaction, one can use dimensional analysis to determine the location optimal value of α as a function of the system size. There are two important length scales in our system, one of which is the lattice size N and the other is determined by the diffusion constant D (length2/time). We posit that the critical (optimal) value of the repulsive energy is a function of both these parameters such that α⋆ ¼ f ðN; DÞ. Since the repulsive energy is dimensionless (Eqs. (4) and (5)), α is defined to have units of length. The bare diffusion of the two searchers is defined in terms of the underlying length and time scales of our simulation and is unity. Therefore we expect the form of the scaling function to be linear in the lattice size
α⋆ ¼ bN; where b is a unit-less constant. In Fig. 5a we plot the location of the optimal value of α versus the system size N. Fitting the data using the linear fitting function suggested by our dimensional scaling analysis, we find an excellent agreement between the data and the scaling function (red dashed line) proposed by our dimensional analysis. The fitting free parameter b is found to be about 0.6 which is very close to what we expect intuitively since at the value of α ¼ N=2 is where the two searchers begin to interact at an energy scale where the probability for accepting/rejecting a proposed move is of order e 1 . This optimal repulsive strength indicates that the interaction between the searchers are such that both do not diffuse over the same space (Fig. 4b) as doing so would be energetically unfavorable. To further solidify that our intuition is correct in regard to the scaling of the two searcher case, we also looked at the three searcher system (Fig. 5a, green squares). We find that indeed the system picks the minimum between all three searchers
N.P. Tani et al. / Journal of Theoretical Biology 361 (2014) 159–164
163
Fig. 4. (a) The frequency of lattice position occupation for both searchers (red-S1) and (blue-S2) for α ¼ 0 displays a high value of overlap. (b) At the optimal value of the repulsion strength α⋆ N=2 we find that the overlap between the two searchers is minimal while maintaining a good search coverage over the lattice. (c) Past the optimal repulsive strength we find that while the overlap is minimal, the search coverage of the lattice is poor because the searchers are more localized. (d–f) The frequency of particle separation as a function of the absolute value of the relative coordinate jS1 S2 j. For the non-interaction case all searcher separations are equally likely, consistent with the overlap shown in (a). In contrast, at the optimal and large repulsion values (as in b–c) we find that the most likely separation is when S1 S2 , which is consistent with the optimal repulsion value of α N=2. All data shown was for a lattice size of N¼ 200. For (a–c) simulations were run for 105 time steps. For (d–f) simulations were run for 106 time steps. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)
to be roughly N/3, thus for a lattice size of N ¼150 we find that α 50, which is the result for the two searcher case for N ¼100. These results suggest that the minimal searcher number that needs to be considered is two. In addition to the existence of an optimal repulsive interaction, which lowers the average mean first passage time to the target, we also expect there will be an average (minimal) encounter time to the target, τ . This characteristic time is related to the intrinsic length and time scales that define our system and should also be related to other relevant quantities, like the bare diffusion constant, D. Using dimensional analysis, we can write down a scaling ansatz similar to the one for the optimal repulsive strength α⋆ . Since we are interested in a quantity that has units of time, the only two parameters that the mean encounter time can depend on is the lattice size and the bare diffusion rate, T ¼ f ðN; DÞ. Given that the diffusion constant has units 2 of ½Length nTime 1 we expect that the dependence of τ to scale as the square of the system size. Therefore we expect that
τ ¼ an
N2 : D
In Fig. 5b we find the average mean first encounter time to the target as a function of the lattice size squared (blue points). Fitting this data with the scaling function derived above we found an excellent agreement between the data and the best fit line (red dashed line).
4. Discussion Foraging strategies can vary among species of animals from bacteria to vertebrates. Foraging strategies could in principle vary across a single species, which could allow for better target tracking as well as coping with changing environmental conditions
(Bartumeus et al., 2002; Seuront and Stanley, 2013). Aside from solitary foraging there also exist animal groups that forage collectively. The same principles of optimization of searching strategies will also apply to these systems as well, but will be conditional upon the collective migration aspect of foraging, which will in general depend on the type of signaling that is used between members of the collective. In summary, we have shown that for multiple Brownian searchers foraging for a single target that the mean first passage time (MFPT) to the target is minimized when there is a mutually repulsive interaction between them and furthermore that there is an optimal range for this repulsion that depends on the density of searchers. These results suggest that the mutually repulsive interaction introduces a cooperative effect that allows the two searchers to optimize both their search time and the search area. The cooperativity in this system is established by the mutual repulsion that each particle exerts on the other, as their motion tends to become correlated and redundancy is minimized when the repulsive range is large enough. However if the repulsive range is too large, the searchers influence each other too strongly and effective search area of each searcher becomes more localized. The balance dictates the optimum range which is roughly the average separation of searchers. In a sense, communication between particles is mediated by the repulsive interaction and the flow of information about their relative positions on the lattice only occurs when their repulsive zones are overlapping and the optimal point is when they are just overlapping. Although we have focused on Brownian dynamics in our model, it should be noted that in many biological systems it has been found that foraging at the individual level can be optimal for other step-size distributions (i.e. Lèvy flights). However, a wide range of organisms ranging from single cells, like Dicty and
164
N.P. Tani et al. / Journal of Theoretical Biology 361 (2014) 159–164
poor environment will remain high (Gelbart, 2010). Similar considerations apply to the behavior of mongolian gazelles that call to each other when they find food (Martinez-Garcia et al., 2013). While the context and mode of interaction are different, it is interesting to note that the study also finds that there is an optimal intermediate range of communication that maximizes efficiency. This raises an interesting evolutionary question in that not only do organisms have to utilize some intercommunication signal to optimize their resource gathering but that interaction must all be optimally tuned in order for collective foraging to always promote a high fitness value in a wide range of environmental conditions.
Acknowledgments The authors would like to thank both Anatoly Kolomeisky and K.C. Huang for their insightful comments and discussions. This work was partially supported by National Science Foundation Grant EF-1038697 and a James S. McDonnell Foundation Award. Undergraduate support was funded by the Undergraduate Research and Mentoring (URM) Program Sponsored by the National Science Foundation NSF Grant DBI-1040962.
References
Fig. 5. (a) Plot of optimal repulsive range αn versus lattice size N for both two searchers (blue circles) and three searchers (green squares). Blue and green bars represent the error and red line is the best fit line α⋆ ¼ 0:5755 7 0:05nN weighted by the error for each data point shown. Error bars were calculated using data points that were within 2% of the minimum (Fig. 3) for each lattice size. (b) Graph of N2 and τ . Blue bars represent the error for each data point and the red line is the best fit line using to our scaling ansatz, which gives τ ¼ 0:0507nN 2 . (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)
bacteria, to mammals, like gazelles, perform movements that are predominantly modeled by correlated or persistent random walks, which, at large spatiotemporal scales, become Gaussian or Brownian random walks (Turchin, 1998). Our results indicate that such Brownian searchers can derive a significant advantage in searching by implementing an optimal range of repulsion. Optimal spacing of foragers via interactions through chemical signaling has been noted in real systems of slime molds and flour beetles (Naylor, 1959, 1961, 1965; Bonner, 1963). This type of behavior in biological systems has been thought to be needed for maximal utilization of food sources as well as searching for good habitats (Shorey, 1976; Martinez-Garcia et al., 2013) and our results show that these interactions can be optimally tuned. Minimizing the search time by each searcher translates to real systems as a minimization of energy consumed by each searcher when foraging for food. This is certainly an important consideration in many biological systems. For example, when food supplies are low single cell Dicty begins to produce chemo-repellent in places that it has searched for food (Keating and Bonner, 1977; Konijn et al., 1968; Pan et al., 1972). This has the effect of informing other single Dicty cells in the population to avoid these areas as a way to conserve energy. This strategy allows the entire population to search for resources collectively as opposed to when food supplies are abundant, where each cell is free to search randomly, although not pure diffusive as in the case for Dicty. From this point of view it is easy to see why having an antagonist interaction among foragers is evolutionarily advantageous. By minimizing energy consumption Dicty can ensure that its fitness in a nutrient
Au, W.W.L., Snyder, K.J., 1980. J. Acoust. Soc. Am. 68, 1077. Bartumeus, F., Catalan, J., Fulco, U.L., Lyra, M.L., Viswanathan, G.M., 2002. Phys. Rev. Lett. 88, 097901. Bonner, J.T., 1963. Sci. Am. 209, 84. Cepa, E., Lepingle, D., 2001. ESAIM: Probab. Stat. 5, 203. Gelbart, M., 2010. Foraging Strategies for Dictyostelium discoideum. Princeton University, Princeton, NJ, USA. Golé, L., Rivière, C., Hayakawa, Y., Rieu, J.-P., 2011. PLoS ONE 6, e26901. Hoff, N.R., Sagoff, A., Wood, R.J., Nagpal, R., 2010. In: IEEE International Conference on Robotics and Biomimetics, vol. 123. Humphries, N.E., et al., 2012. Proc. Natl. Acad. Sci. 109, 7169–7174. Ioannou, C.C., Guttal, V., Couzin, I., 2012. Science 337, 1212. Jackson, D.E., Ratnieks, F.L., 2006. Curr. Biol. 16, R570. Keating, M.T., Bonner, J.T., 1977. J. Bacteriol. 130 (1), 144. Ko, A., Lau, H.Y., 2009. Int. J. Adv. Sci. Technol. 3, 69. Konijn, T., Barkley, D., Chang, Y., Bonner, J., 1968. Am. Natural. 102, 225. Levendowsky, M., White, B.S., Schuster, F., 1997. Acta Protozool. 36, 237. Martinez-Garcia, R., et al., 2014. Phys. Rev. E 89, 032718. Martinez-Garcia, R., et al., 2013. Phys. Rev. Lett. 110, 248106. Mattson, D.J., 2005. USGS. Montroll, E.W., 1969. J. Math. Phys. 10, 753. Naylor, A.F., 1959. Ecology 40, 453. Naylor, A.F., 1961. Ecology 42, 231. Naylor, A.F., 1965. Ecology 46, 341. Palyulina, V.V., Chechkinb, A.V., Metzler, R., 2013. Proc. Natl. Acad. Sci. 111, 2931. Pan, P., Hall, E., Bonner, J., 1972. Nat.: New Biol. 237 (75), 181. Raichlena, D.A., et al., 2013. Proc. Natl. Acad. Sci. 111, 728. Reich, J., Sklar, E., 2006. In: Proceedings of the IEEE International Workshop on Safety, Security and Rescue Robotics. Saeedi, P., Sorensen, S.-A., Hailes, S., 2009. In: Safety, Security and Rescue Robotics (SSRR), IEEE International Workshop. Schnitzler, H.U., Kalko, E.K.V., 2001. BioScience 51 (7), 557. Schuster, F., Levendowsky, M., 1996. J. Euk. Microbiol. 43, 150. Selmeczi, D., Mosler, S., Hagedorn, P.H., Larsen, N.B., Flyvbjerg, H., 2005. Biophys. J. 89, 912. Seuront, L., Stanley, H.E., 2013. Proc. Natl. Acad. Sci. 111, 2206. Shenderov, A., Sheetz, M., 1997. Biophys. J. 72, 2382. Shlesinger, M.F., Klafter, J., 1986. On Growth and Form: Lévy Walks versus Lévy Flights. Springer, Netherlands. Shoghian, S., Kouzehgar, M., 2012. World Acad. Sci. Eng. Technol. 72, 447. Shorey, H.H., 1976. Animal Communication by Pheromones. Academic Press, Waltham, Massachusetts, USA. Sokolov, I.M., Metzler, R., Pant, K., William, M.C., 2005. Phys. Rev. E 72, 041102. Turchin, P., 1998. Quantitative Analysis of Movement: Measuring and Modeling Population Redistribution in Plants and Animals. Sinauer Associates, Sunderland, MA, USA. Viswanathan, G.M., et al., 1999. Nature 401, 911. Viswanathan, G.M., et al., 2011. The Physics of Foraging: An Introduction to Random Searches and Biological Encounters, 1st ed. Cambridge University Press, Cambridge, England, p. 2011.