Towards Swarm Calculus: Urn Models of Collective Decisions and Universal Properties of Swarm Performance
aa r X i v : . [ c s . N E ] M a r Towards Swarm Calculus:Urn Models of Collective Decisions and UniversalProperties of Swarm Performance
Heiko HamannDepartment of Computer ScienceUniversity of PaderbornZukunftsmeile 133102 Paderborn, Germany [email protected]
March 19, 2013
Abstract
Methods of general applicability are searched for in swarm intelligencewith the aim of gaining new insights about natural swarms and to developdesign methodologies for artificial swarms. An ideal solution could be a‘swarm calculus’ that allows to calculate key features of swarms such asexpected swarm performance and robustness based on only a few param-eters. To work towards this ideal, one needs to find methods and modelswith high degrees of generality. In this paper, we report two models thatmight be examples of exceptional generality. First, an abstract model ispresented that describes swarm performance depending on swarm densitybased on the dichotomy between cooperation and interference. Typicalswarm experiments are given as examples to show how the model fits toseveral different results. Second, we give an abstract model of collectivedecision making that is inspired by urn models. The effects of positivefeedback probability, that is increasing over time in a decision making sys-tem, are understood by the help of a parameter that controls the feedbackbased on the swarm’s current consensus. Several applicable methods, suchas the description as Markov process, calculation of splitting probabilities,mean first passage times, and measurements of positive feedback, are dis-cussed and applications to artificial and natural swarms are reported.
The research of swarm intelligence is important in both biology to gain newinsights about natural swarms and also in fields dealing with artificial swarms,such as swarm robotics, to obtain sophisticated design methodologies. The ideal1ool would allow to calculate fundamental features of swarm behavior, such asperformance, stability, and robustness, and would need only a few observedparameters in the case of natural swarm systems or a few designed parametersin the case of artificial swarms. We call this highly desired set of tools ‘swarmcalculus’ (calculus in its general sense). The underlying idea is to create a setof mathematical tools that are easily applied to a variety of settings in the fieldof swarm intelligence. In addition, these tools should be easily combined whichwould allow using them as mathematical building blocks for modeling. Thus,models will surely be an important part of swarm calculus.General properties and generally applicable models need to be found to ob-tain a general methodology of understanding and designing swarm systems.Today it seems that only few models exist that have the potential to be-come general swarm models. For example, swarm models in biology are par-ticularly distinguished by their variety (Okubo and Levin, 2001; Okubo, 1986;Vicsek and Zafeiris, 2012; Edelstein-Keshet, 2006; Camazine et al., 2001). Typ-ically, a specialized model is created for each biological challenge. It seems thatthe desire for models with wide applicability to a collection of natural swarmsis rather low in that community. In the field of artificial swarms, such as robotswarms, the desire for generality seems to be bigger, which is, for example, ex-pressed by several models in swarm robotics (Hamann, 2010; Berman et al.,2011; Prorok et al., 2011; Milutinovic and Lima, 2007; Lerman et al., 2005).The driving force for the creation of these models is to support the designof swarm robotic systems within a maximal range of applications. The focusof these models is on quantitative features of the swarm behavior, such as thedistribution of robots or required times for certain tasks. However, there is astruggle between the intended generality of the model and the creation of adirect mapping between the model and the actual description of the individualrobot’s behavior. A higher degree of generality is achievable if the demand fora detailed description of behavioral features is abandoned and focus is set onlyon high-level features such as overall performance or the macroscopic processof a collective decision. Such high-level models can be expressed by concisemathematical descriptions that, in turn, allow direct applications of standardmethods from statistics, linear algebra, and statistical mechanics. In this pa-per , we present two models of general properties of swarm systems concerningthe dependence of swarm performance on swarm density and the dependence ofcollective decisions on positive feedback. This paper is an extended version of Hamann (2012). The main extensions are the methodof deriving the probability of positive feedback based on observed decision revisions (Sec. 4.1),a discussion of additional methodology such as Markov chains, splitting probabilities, andmean first passage times (Sec. 4.3), and a comprehensive introduction of the Ehrenfest andthe Eigen urn models. Fundamentals of swarm performance and col-lective decisions
In this section we define the concepts of swarm performance and collectivedecision-making along with so-called urn models upon which the collective-decision model is based.
By ‘swarm performance’ we denote the efficiency of a swarm concerning a cer-tain task. For example, the swarm performance can be a success rate of howoften the task is accomplished on average, it can be the average speed of aswarm in collective motion etc. Here we are interested in the swarm per-formance as a function over ‘swarm density’, which is how many agents arefound on average within a certain area. For the following reason, the func-tion of swarm performance depending on swarm density cannot simply be alinear function. For a true swarm system, a very low density, which corre-sponds situations with only a few agents in the whole area, has to result inlow performance because there is neither a lot of cooperation between agentsbecause they seldom meet nor a significant speed-up. With increasing den-sity, the performance increases, on the one hand, because of a simple speed-up(e.g., two robots clean the floor faster than one) and, on the other hand, be-cause of increasing opportunities of cooperation (assuming that cooperation isan essential beneficial part of swarms). In natural swarms such increases in per-formance with increasing swarm size are revealed, for example, in productivitygains and also in the emergence of increased division of labor as an indicator forincreased cooperation (Jeanne and Nordheim, 1996; Karsai and Wenzel, 1998;Gautrais et al., 2002; Jeanson et al., 2007). Even superlinear performance in-creases are possible in this interval of swarm density and was reported for aswarm of robots (Mondada et al., 2005). For artificial swarms it was reportedthat at some critical/optimal density (Schneider-Font´an and Matari´c, 1996) theperformance curve will first level off and then decrease (Arkin et al., 1993) be-cause improvements in cooperation possibilities will be lower than the draw-back of high densities, namely interference (Lerman and Galstyan, 2002). Withfurther increase of the density, the performance continues to decrease as re-ported for multi-robot systems (Goldberg and Matari´c, 1997). Hence, swarmsgenerally face a tradeoff between cooperation, which is beneficial, and inter-ference, which is usually obstructive, however, has sometimes positive effectsboth within certain natural swarms and as a tool for designing swarm algo-rithms (Dussutour et al., 2004; Goldberg and Matari´c, 1997).In the following we report that many swarm systems not only show similarqualitative properties but show also similarities in the actual shapes of theirswarm performance over swarm size/density graphs (see function of swarm per-formance in Fig. 2(a) on page 9). Examples are the performance of foraging ina group of robots (Fig. 2(b) on page 9 and Fig. 10a in (Lerman and Galstyan,3002)), the activation dynamics and information capacity in an abstract cel-lular automaton model of ants (Figs. 1b and 1c in (Miramontes, 1995)), andeven in the sizes of social networks (Fig. 8b in (Strogatz, 2001)). A similarcurve is also presented as a hypothesis for per capita output in social wasps byJeanne and Nordheim (1996). The existence of this general shape was alreadyreported by Østergaard et al. (2001) in their expected performance model formulti-robot systems in constrained environments:
We know that by varying the implementation of a given task, we canmove the point of “maximum performance” and we can change theshapes of the curve on either side of it, but we cannot change thegeneral shape of the graph.
Traffic models of flow over density are related because traffic flow can be in-creased when cars cooperatively share streets but the flow decreases when thestreets are too crowded and cars interfere too much with each other. While the‘fundamental diagram’ of traffic flow (Lighthill and Whitham, 1955) is symmet-ric, more realistic models propose at least two asymmetric phases of free andsynchronized flow (e.g., Fig 3(b) in (Wong and Wong, 2002)). Actual measure-ments on highways show curves with shapes similar to Fig. 2(a) (e.g., see Fig. 6-4in (Mahmassani et al., 2009)). In these models, there exist two densities for agiven flow (except for maximum flow) similar to the situation here where wehave two swarm densities for each swarm performance (one smaller than theoptimal density and one bigger than the optimal density; the correspondingfunction that maps densities to performance values is surjective, not bijective).
In the context of swarm intelligence, collective decision-making is a process thatis distributed over a group of agents without global knowledge. Each agentdecides based on locally sampled data such as the current decision of its neigh-bors. There are many biological systems showing collective decision-making,for example, food source choice in honey bees (Seeley et al., 1991), nest siteselection in ants (Mallon et al., 2001), and escape route search in social spi-ders (Saffre et al., 1999). Collective decision-making systems are often modeledas positive-feedback systems that utilize initial fluctuations which are amplifiedand that way help to converge to a global decision (Deneubourg et al., 1990;Mallon et al., 2001; Nicolis et al., 2011). Interesting features of collective deci-sions are, for example, the speed-accuracy trade-off (Nicolis et al., 2011) or theinfluence of noise (Dussutour et al., 2009; Yates et al., 2009). Furthermore, itturns out that positive feedback is not always productive but can also generateirrational decisions (Nicolis et al., 2011).In the following, we limit our investigations to binary decision processesbecause they allow for a concise mathematical notation and, hence, allow amanageable application of mathematical methods. The investigated systemsare either inherently noisy (e.g., explicit stochastic processes within the agents’behaviors) or can validly be modeled as noisy processes (e.g., deterministic4haos with strong dependence on the initial conditions). We are not interestedin the quality of the final decision—that is the utility of choosing option A over option B or vice versa—and assume that there is no initial bias to one orthe other. In selecting an appropriate model for collective decisions, our mainconcern is simplicity while keeping focus also on the relation of how much isneeded as input to the model and how much is generated by it. We want tokeep the number of parameters small while achieving descriptions of qualitativeaspects of collective decisions as effects of the model. In order to obtain thesestandards, we choose a minimal macroscopic model which has only one statevariable describing the current status of the collective decision within the swarm(e.g., 80% for option A and consequently 20% for option B ).For simplicity, we view the asynchronous distributed process of collectivedecisions as a round-based game which allows only one agent at a time to eitherrevise its current decision or to convince a peer to revise its current decision. Therelation to natural systems is that decision events are serialized and intermediateperiods of time are ignored. The influence of this assumption to the steady statebehavior is considered to be low. For this purpose, we re-interpret well-knownurn models as models of collective decisions and extend them appropriately. Weuse simple models inspired by the urn model of Ehrenfest and Ehrenfest (1907)and by the urn model of Eigen and Winkler (1993). This urn model was originally introduced by Ehrenfest and Ehrenfest (1907) inthe context of dissipation, thermodynamics, statistical mechanics, and entropy.The dynamics of the model is defined as follows. Our urn is filled with N mar-bles. Say, initially all marbles are blue. Whenever we draw a blue one we replaceit with a red marble. If we draw a red marble we replace it with a blue one.Obviously, the two extreme states of either having an urn full of blue marblesor an urn full of red marbles are unstable. Similarly, this is true for all states ofunevenly distributed colors. To formalize this process, we keep track of how thenumber of blue marbles (without loss of generality) changes depending on howmany blue marbles were in the urn at that time. We can do this empirically orwe can actually compute the average expected ‘gain’ in terms of blue marbles.For example, say at time t we have B ( t ) = 16 blue marbles in the urn and atotal of N = 64 marbles. The probability of drawing a blue marble is therefore P B = = 0 .
25. The case of drawing a blue marble has to be weighted by − P R = = 0 .
75 which is weighted by +1. Hence, theexpected average change ∆ B of blue marbles per round depending on the currentnumber of blue marbles B = 16 is ∆ B ( B = 16) = 0 . −
1) + 0 . . B ( B ) = − · B
64 + 1 , (1)5 B ∆ B (a) Ehrenfest urn model -1 0 1 0 32 64 B ∆ B (b) Eigen urn model (note ∆ B (0) =∆ B (64) = 0) Figure 1: Expected average change in the urn models.which is plotted in Fig. 1(a). Hence, the average dynamics of this game is givenby B ( t + 1) = B ( t ) + ∆ B ( B ( t )).The recurrence B t = B t − − B t − + 1 can be solved by generating func-tions (Graham et al., 1998). For B = 0 we obtain the generating function G ( z ) = X t X k ≤ t (cid:18) (cid:19) k z t . (2)The t th coefficient [ z t ] of this power series is the closed form for B t . We get[ z t ] = B t = X k ≤ t (cid:18) (cid:19) k = 1 − (cid:0) (cid:1) t − . (3)Hence, for initializations B = 0 and the symmetrical case B = 64 thesystem converges in average rather fast to the equilibrium B = 32. The actualdynamics of this game is, of course, a stochastic process which can, for example,be modeled by B ( t + 1) = B ( t ) + ∆ B ( B ( t )) + ξ ( t ), for a noise term ξ .Several generalizations have been proposed for this urn model (Krafft and Schaefer,1993; Klein, 1956), however, these investigations mostly focus on mathemati-cally tractable variants. In the following we report generalizations, which focuson applications to feedback processes. The Ehrenfest model can be interpreted as an example of a negative feedbackprocess. Deviations from the fixed point B = 32 are corrected by negativefeedback (the predominant color will diminish on average). Eigen and Winkler(1993) reported a similar urn model to show the effect of positive feedback. Inthis model drawing a blue marble has the effect of replacing a red marble bya blue one and vice versa. The expected average change of blue marbles perround changes accordingly to 6 B ( B ) = ( B − , for B ∈ [1 , , else . (4)For a plot of ∆ B ( B ) see Fig. 1(b).While we still have an expected change of ∆ B = 0 for B = 32 as in theEhrenfest model, the fixed point B = 32 is unstable now as its surroundingdrives trajectories away from B = 32 and towards the stable fixed points B = 0and B = 64 respectively. In Sec. 4, we introduce a more general urn model thattakes the intensity of positive feedback as a parameter. This model can be usedto investigate collective decisions in swarms. Having identified the two main components (cooperation and interference) andthe typical shape of these graphs, we can define a simple model. The idea is tofit this model to empirical data for verification and predictions.
For a given bounded, constant area A the swarm density ρ is defined by theswarm size N according to ρ = N/A . Also a dynamic area A ( t ) could beassumed but throughout this paper we want to keep swarm density and swarmsize interchangeable based on the identity ρ = N/A . Although for a givenswarm density the swarm performance might be quantitatively and qualitativelydifferent for different areas, here we focus on describing such swarm–performancefunctions separately. We define the swarm performance Π depending on swarmsize N = ρA by Π( N ) = C ( N )( I ( N ) − d ) = a N b a exp( cN ) , (5)for parameters c < a , a > b >
0, and d ≥ d is subtracted to force a decreaseto zero (lim N →∞ I ( N ) − d = 0). The swarm performance depends on twocomponents, C and I . First, the swarm effort without negative feedback isdefined by the cooperation function (see also Fig. 2(a)) C ( N ) = a N b . (6)This function can be interpreted as the potential for cooperation in a swarmthat would exist without certain constraints, such as physical collisions or otherspatial restrictions. The same formula was used by Breder (1954) to modelthe cohesiveness of a fish school and by Bjerknes and Winfield (2010) to modelswarm velocity in emergent taxis. However, they used parameters of b < b >
1. In principle this is a major difference because b < b > I ( N ) = a exp( cN ) + d , (7)with d used for scaling (e.g., lim N →∞ I ( N ) = d ). The interference function canalso be interpreted as the swarm performance achievable without cooperation,that is, achievable swarm performance without positive feedback. The exponen-tial decrease of the interference function seems to be a reasonable choice because,for example, the Ringelmann effect according to Ingham et al. (1974) impliesalso a nonlinear decrease of individual performance with increasing group size;see also (Kennedy and Eberhart, 2001, p. 236). Nonlinear effects that decreaseefficiency with increasing swarm size are plausible due to negative feedback pro-cesses such as the collision avoidance behavior of one robot triggers the collisionavoidance behavior of several others in high density situations. Still, there aremany options of available nonlinear functions but best results were obtained withexponential functions. Also Fig. 10b in (Lerman and Galstyan, 2002) shows anexponentially decreasing efficiency per robot in a foraging task. To prove the wide applicability of this simple model we fit it to some swarmperformance plots that are available. We investigate four scenarios: foragingin a group of robots (Lerman and Galstyan, 2002), collective decision mak-ing (Hamann et al., 2012) based on BEECLUST (Schmickl and Hamann, 2011),aggregations in tree-like structures and reduction to shortest paths (Hamann,2006) similar to (Hamann and W¨orn, 2008), and the emergent taxis scenario(also sometimes called ‘alpha algorithm’, Nembrini et al. (2002); Bjerknes et al.(2007)).Given the data of the overall performance, the four parameters a , a , b , c of Π( N ) (Eq. 5) can be directly fitted. This was done for the three examplesshown in Figs. 2(b), 2(c), and 2(d). The equation can be well fitted to theempirical data (for details about the fitted functions see the appendix). In caseof the foraging scenario (Fig. 2(b)) we also have data of the efficiency per robot.We can use the model parameters a and c , that were obtained by fitting themodel to the overall performance, to predict the efficiency per robot, which is afunction that we suppose to be proportional to effects by interference. This isdone by scaling the interference function I ( N ) linearly and plotting it againstthe efficiency per robot. The result is shown in Fig. 2(b).We analyze the forth example, emergent taxis (Nembrini et al., 2002; Bjerknes et al.,2007), in more detail and, for this purpose, give a short description of the al-gorithm. The objective is to move a swarm of robots towards a light beacon.8 Π N swarm performance Πcooperationinterference (a) Model of swarm performance Π basedon cooperation and interference, examplesof swarm performance Π (Eq. 5), cooper-ation (Eq. 6), and interference (Eq. 7) de-pending on swarm size; a = 0 . a = 1, b = 2 . c = − . Π N efficiencyper robotgroup efficiency (b) Foraging in a group of robots, Eq. 5 fit-ted to group efficiency (upper solid line), pre-diction of interference (lower solid line, ef-ficiency per robot, linearly rescaled), datapoints extracted from Fig. 10 in Lermanet al. (Lerman and Galstyan, 2002); Π givesgroup/robot efficiency. Π N (c) Collective decision mak-ing (Hamann et al., 2012) based onBEECLUST (Schmickl and Hamann,2011). Π ρ (d) Aggregation in tree-like structures andreduction to shortest path (Hamann, 2006);Π gives the ratio of successful runs. Figure 2: Model of cooperation and interference and three scenarios with fittedperformance Π according to Eq. 5 9he robots are limited in their capabilities because they only have an omnidi-rectional beacon sensor that detects two states, beacon seen or not seen, butgives no bearing. If a robot does not see the beacon this is always because oneor several other robots are within the line of sight between this robot and thebeacon. Robots that see the beacon are called ‘illuminated’ and robots thatdo not see the beacon are called ‘shadowed’. The robots’ behavior is definedto differ depending on these two states. Shadowed robots have a shorter colli-sion avoidance radius than illuminated robots. Consequently if a shadowed andan illuminated robot approach each other, the illuminated robot will triggerits collision avoidance behavior while the shadowed robot will not be affecteduntil it gets even closer and triggers its own collision avoidance behavior aswell. This interplay between shadowed and illuminated robots generates a biastowards the beacon. In addition, the algorithm has a ‘coherence state’ whichaims at keeping the robots together within each others communication range.It is assumed that a robot is able to count the robots that are within range.Once this number drops below a threshold α , the robot will do a u-turn whichhopefully brings it back into the swarm. In the following investigations we setthe threshold at first to α = 0 which means we turn the coherence behavioroff. Later we follow (Bjerknes et al., 2007) and set the threshold to swarm size( α = N ) to enforce full coherence. Initially the robot swarm is distant from thebeacon and is randomly distributed while ensuring coherence. Interestingly, itis difficult to identify two separated behavior components of cooperation andinterference. The robots cooperate in generating coherence and in having shad-owed robots that approach illuminated robots which drives them towards thebeacon. However, collision avoidance behavior also has disadvantageous effects.In order to keep coherence the robots might be forced to aggregate too denselywhich might result in robots blocking each other.The following empirical data is based on a simple simulation of emergenttaxis. This simulation is noise-free and therefore robots move in straight linesexcept for u-turns according to the emergent taxis algorithm (a random turnafter regaining coherence was not implemented).First, we measure the performance that is achieved without cooperation.This is done by defining a random behavior that ignores any characteristic fea-ture of the actual emergent taxis algorithm. We set this parameter to α = 0in the simulation to obtain the cooperation-free behavior. Hence, no robot willever u-turn and they basically disperse in the arena. A simulation run is stoppedonce a robot touches a wall. The performance Π of the swarm is measured bythe total distance covered by the swarm’s barycenter multiplied by the swarmsize (i.e., an estimate of the sum of all distances that were effectively covered byeach robot). The performance obtained by this random behavior can be fittedusing the interference function of Eq. 7 (interpreting the interference function asa measurement of swarm performance without cooperation). The fitted inter-ference function and the empirically obtained data is shown in Fig. 3(a) labeled‘random’. The interference function does not drop to zero for large N ; thisbias towards the light source (positive covered distance) is due to the initialpositioning of the swarm closer to the wall that is farther away from the light10 Π N randomemergent taxis (a) Model fitted to data of random andemergent–taxis behavior; performance Πis the total distance covered by theswarm’s barycenter multiplied by theswarm size N . Π N (b) Fit of the cooperation function to theemergent–taxis performance based on thenarrow interval of N ∈ [20 ,
22] only, withfixed parameters of interference functionas shown in (a).
10 20 30 40-0.002 0 0.002 0.004 0.006 0 200 Nv f r e q u e n c y (c) Histogram of barycenter speeds v for differ-ent swarm sizes. Note the bimodality in theinterval N ∈ [16 , Figure 3: Performance of a random behavior and the actual self-organized emer-gent taxis behavior (also sometimes called ‘alpha algorithm’) (Nembrini et al.,2002; Bjerknes et al., 2007) with fitted model (Eq. 5); histogram showing twophases.source.In a second step, the full model of swarm performance Π (Eq. 5) is to befitted to the actual emergent taxis scenario. The data is obtained by setting thethreshold in our simple simulation of emergent taxis back to normal ( α = N ).The fitting is done by keeping the interference function fixed and we fit onlythe cooperation function (i.e., fitting a and b while keeping a , c and d fixed).The fitted swarm performance model is shown in Fig. 3(a) labeled ‘emergenttaxis’. This simple model is capable of predictions, if the interference functionhas been fitted and we fit the cooperation function only to a small interval of,for example, N ∈ [20 ,
22] (i.e., intervals close to the maximal performance).This is shown in Fig. 3(b). The implication is that if the interference functionis known as well as the optimal swarm size then the behavior within the otherintervals can be predicted.Note that the model operates on a single averaged value to describe the per-11ormance which does not fully catch the system’s behavior. At least in somescenarios, as here in the emergent taxis scenario, the performance does notjust continuously decrease due to continuously increasing interference. Instead,two coexisting phases of behaviors emerge: functioning swarms moving forwardand pinned swarms with extreme numbers of u-turns. In emergent taxis thisis shown, for example, by a histogram of barycenter speeds in Fig. 3(c). For
N <
15 the mean of a unimodal distribution increases with increasing N . Start-ing at about N = 15 a second phase emerges that shows slowly moving swarmsand generates a bimodal distribution. Hence, given the fully deterministic im-plementation of our simulation, there are two classes of initial states (robotpositions and orientations) that determine the two extremes of success or totalfailure. Consequently, swarm performance functions as shown in Fig. 3(a) needto be interpreted with care because they might indicate an average behaviorthat does actually not occur. Still, these swarm performance functions are use-ful if the values are interpreted relatively as success rates. That way Fig. 3(a)gives a good estimate of the frequency of the two phases. In other scenarios theinterference might truly increase continuously due to a qualitatively differentprocess, such as saturation of target areas with robots.The presented model of swarm performance has potential to be applicable tomany swarm systems. In the next section, a model is given for a subset of swarmsystems namely collective decision-making systems. The two models relate toeach other as in some cases they are both applicable to the same swarm system.A candidate for such a system could be a BEECLUST-controlled swarm. Theapplication of the swarm performance model to this system is given in Fig. 2(c).Data that supports, that an application of the following collective decision modelis likely, is given by Hamann et al. (2012). In such systems the effectivity of thecollective decision and the performance of the swarm are directly linked andconsequently the two reported models are, too. In the following, we investigate macroscopic models of collective decisions. Oneof the most general and at the same time simplest models of collective decisions isa model of only one state variable s ( t ), which gives the temporal evolution of theswarm fraction that is in favor of one of the options in a binary decision process.If we assume that there is no initial bias to either option (i.e., full symmetry),then we need a tie breaker for s = 0 .
5. A good choice for a tie breaker is noisebecause any real swarm will be noisy. The average change of s depending onitself per time (∆ s ( s ) / ∆ t ) is of interest. Given that the system should be ableto converge to either of the options at a time plus having the symmetric caseof s = 0 .
5, function ∆ s ( s ) / ∆ t needs to have at least three roots (∆ s ( s i ) / ∆ t = 0, s i ∈ S , | S | ≥
3) and consequently is at least a cubic function. Instead ofdeveloping a model, that defines such a function, we prefer a model that allowsthis function to emerge from a simple process. Once symmetry is overcome byfluctuations, swarm systems have a tendency to confirm and reinforce such a12reliminary decision due to positive feedback (say, for s = 0 . ǫ there is atendency towards s = 1). We define such a process depending on probabilitiesof positive feedback next. As discussed in the introduction, we define an urn model that was inspired bythe models of Ehrenfest and Ehrenfest (1907) and of Eigen and Winkler (1993).We use an urn model that has optionally positive or negative feedback dependingon the system’s current state and depending on a stochastic process. The urn isfilled with N marbles which are either red or blue. The game’s dynamics is turn-based. First, a marble is drawn with replacement followed by replacing a secondone determined by the color of the first marble. The probability of drawing ablue marble is implicitly determined by the current number of blue marbles inthe urn. The subsequent replacement of a second marble has either a positive ora negative feedback. Say, we draw a blue marble, we notice the color, and put itback into the urn. Then our model defines that with probability P a red marblewill be replaced by a blue one (i.e., a positive feedback event because drawinga blue one increased the number of blue marbles) and with probability 1 − P ablue one will be replaced by a red one (i.e., a negative feedback event becausenow drawing a blue one decreased the number of blue marbles). Hence, P givesthe probability of positive feedback.The analogy of this model to a collective decision making scenario is thefollowing. The initial drawing resembles the frequency of individual decisionsin the swarm over time within the turn-based model. This frequency is propor-tional to the current ratio s ( t ) of blue marbles in the urn. Consequently, weserialize the system dynamics and each system state s ( t ) has at most two pre-decessors ( s ( t −
1) = s ( t ) ± P is determined explicitly. We de-fine the determination whether positive feedback or negative feedback is ef-fective in a given system state s as a binary random experiment. The sam-ple space is Ω = { PFB , NFB } with PFB denoting positive feedback, NFBdenoting negative feedback, and we define a probability measure m , conse-quently m (PFB) + m (NFB) = 1 holds. Hence, the probability of positivefeedback is defined by m (PFB) = P ( s, ϕ ) (probability of negative feedbackis m (NFB) = 1 − P ( s, ϕ )) and for now we choose a sine function P ( s, ϕ ) = ϕ sin( πs ) . (8)Due to the symmetry P ( s, ϕ ) = P (1 − s, ϕ ) it is irrelevant whether s is set to theratio of blue marbles or the ratio of red marbles. Later we will find that similarbut different functions might be a better choice in certain situations (Sec. 4.4).The constant ϕ ∈ [0 ,
1] scales the amplitude of the sine function (see Fig. 4(a))and consequently defines the predominant ‘sign’ of the feedback and the prob-abilities of positive and negative feedback. The integral R P ( s, ϕ ) ds gives the13 s s s P ϕ = 0 ϕ = 0 . ϕ = 0 . ϕ = 0 . ϕ = 0 . (a) Examples of setting the probability ofpositive feedback for intensities of feed-back ϕ ∈ { , . , . , . , . } , squaresgive the redetermined values according toEq. 15. -1 0 1 0 0.5 1 s s s ∆ B a nd ∆ B ϕ = 0 ✁✁ ϕ = 0 . ϕ = 0 . ✁✁ ϕ = 0 . ϕ = 0 . (b) Average change of B ( t ) in terms of marbles,lines according to Eq. 9, squares give empiricaldata, number of samples is 8 × for eachpossible s , 64 marbles. Figure 4: Settings of the positive feedback probabilities and resulting averagechange in B ( t ) in the urn model over the ratio of marbles s .overall probability of having positive feedback independent of s . Negative feed-back is predominant for any s for ϕ < . s = 0 . ϕ > . P ( s, ϕ ) is plotted for different settings of ϕ in Fig. 4(a). There is maximumprobability for positive feedback for the fully symmetric case of s = 0 . s = 0 and s = 1 we have P ( s, ϕ ) = 0 because no positivefeedback is possible (either all marbles are already blue or all marbles are redand therefore no blue one can be drawn). For ϕ ≤ . ϕ ≤ . , ∀ s : P ( s, ϕ ) ≤ . s = 0 . B of bluemarbles B by summing over the four cases: drawing a blue or red marble,followed by positive or negative feedback, multiplied by the ‘payoff’ in terms ofblue marbles, respectively. Using the symmetry P ( s, ϕ ) = P (1 − s, ϕ ) we get∆ B ( s ) = sP ( s, ϕ )(+1) + s (1 − P ( s, ϕ ))( − − s ) P (1 − s, ϕ )( −
1) + (1 − s )(1 − P (1 − s, ϕ ))(+1)= 4( P ( s, ϕ ) − . s − . . (9)We defined P ( s, ϕ ) based on a trigonometric function but alternatively onecan choose, for example, a quadratic function P ( s, ϕ ) = ϕ (1 − s − . ) , (10)yielding a cubic function∆ B ( s ) = − s − .
5) + 4 ϕ ( s − . − ϕ ( s − . . (11)14n the following, however, we will use Eq. 9. Also note that by turning positivefeedback completely off ( P ( s, ϕ ) = 0) we obtain the Ehrenfest urn model again(cf. Eqs. 1 and 9) and by activating maximal positive feedback ( P ( s, ϕ ) = 1)we obtain the Eigen urn model (cf. Eq. 4).While in the above definition the positive feedback probability P ( s, ϕ ) isexplicitly given, it might be unknown in applications of this model. In theseapplications, the positive feedback itself might not be measurable but maybe thenumber of observed decision revisions of the agents. We introduce the absolutenumber of observed individual decision revisions from red to blue r b ( s ) over anygiven period and from blue to red r r ( s ). The ratio of red-to-blue revisions isdirectly related to the expected average change per round. Assuming payoffsof ±
1, the average change of blue marbles per round is obtained by the weightedsum r b ( s ) r b ( s ) + r r ( s ) · (+1) + r r ( s ) r b ( s ) + r r ( s ) · ( −
1) = ∆ B ( s ) , (12)which simplifies to r b ( s ) r b ( s ) + r r ( s ) = 0 . B ( s ) + 1) . (13)In addition, the ratio r b ( s ) r b ( s )+ r r ( s ) also directly relates to the positive feedbackprobability P ( s ) by considering the above mentioned four cases: drawing a blueor red marble, followed by positive or negative feedback. Red-to-blue revisionsare observed only in two of these four cases: drawing a blue marble followedby positive feedback and drawing a red marble followed by negative feedback.The summed probability of these two cases has to equal the ratio of red-to-bluerevisions, which is consequently interpreted as probability. The equation P ( s ) s + (1 − P (1 − s ))(1 − s ) = r b ( s ) r b ( s ) + r r ( s ) , (14)yields using the symmetry P (1 − s ) = P ( s ) P ( s ) = r b ( s ) r b ( s )+ r r ( s ) − s s − , for s = 0 . , r b ( s ) r b ( s ) + r r ( s ) ≤ max( s, − s ) , (15)an equation which is to be used to determine the positive feedback probabil-ity based on measurements of agents’ revisions. The pole at s = 0 . P ( s = 0 .
5) is reasonable. Any definition of P ( s = 0 . s = 0 .
5. We actually define P ( s = 0 .
5) in Eq. 8, whichis, however, only a simplification of notation. Also note that this mathematicaldifficulty has limited effect in applications because swarms are inherently dis-crete. For example, s = 0 . P ( s ) based on Eq. 15 and on measurements of r b ( s ) and r r ( s ). The values for P ( s = 0 .
5) cannot be given, as discussed, and15
0 0.5 1 0 32 64 0 0.05 0.1 ϕB Figure 5: Normalized histogram of blue marbles B over intensity of feedback ϕ after t = 2000 steps in the urn model, initialized to B (0) ∈ { , } , indicatinga pitchfork bifurcation at ϕ = 0 . s = 0 . B ac-cording to Eq. 9 to the empirically obtained average change ∆ B in terms ofnumber of marbles for the different settings of ϕ . The agreement between the-ory and empirical data is good (root mean squared errors of < × − ) asexpected. Two zeros s and s emerge additionally to s = 0 . ϕ > . s = π arcsin( ϕ ) and s = 1 − π arcsin( ϕ ). Positive values of ∆ B ( s ) for s < . s = 0 . s = 0 and vice versa for the other half( s > . ϕ . It shows a pitchfork bifurcation at ϕ = 0 .
5, whichis to be expected based on Fig. 4(b) and which is fuzzy because of the underlyingstochastic process. For ϕ > .
5, the curve defined by ∆ B ( s ) becomes cubic andgenerates two new stable fixed points while the former at s = 0 . Since positive and negative feedback is determined stochastically based on thecurrent global state in the urn model, it is straightforward to determine thenumber of replaced marbles stochastically, too. Instead of having a fixed ‘payoff’(number of replaced marbles) of +1 for positive feedback and − M of the event ‘having apayoff of +1’ (or − − M for a payoff of 0. Thusthe average payoff would be + M and − M for positive and negative feedbackrespectively. In addition, the payoff can be variant depending on the currentglobal state s , that is, we define a function M ( s ). It turns out that a definitionsimilar to the positive feedback probability P ( s, ϕ ) is useful here. We define the16ariant payoff by M ( s ) = c sin( πs ) + c , (16)for appropriately chosen constants c and c . M ( s ) defines the average overabsolute values of changes in s , similar to the diffusion coefficient in Fokker-Planck theory. Measurements of the diffusion coefficient in swarm systems werereported in Yates et al. (2009) and Hamann et al. (2010), which show low valuesfor s ≪ . s ≫ . s = 0 .
5. With M ( s ) defined by Eq. 16we have symmetry again ( M ( s ) = M (1 − s )) and hence the extension of Eq. 9is simply ∆ B ( s ) = 4 M ( s )( P ( s, ϕ ) − . s − . . (17) Note that the recurrence B t = B t − + ∆ B ( B t − ) is a trigonometric or cubicfunction based on Eq. 8 or Eq. 10, respectively. Thus, it is much more diffi-cult to be handled analytically and a concise result as for the Ehrenfest model(Eq. 3) cannot be obtained. When applying nonlinear equations for ∆ B we en-ter the domains of nonlinear dynamics and chaotic systems with all their knownmathematical intractabilities. Hence, in general we have to rely on numericalmethods.However, if we choose to investigate probability distributions ρ ( s, t ) insteadof single trajectories s ( t ), an interesting mathematical option is available. Thesteady state π ( s ) of the probability distribution ρ ( s, t ) over an ensemble ofrealizations of s ( t ) can be obtained analytically. We assume for simplicity thatthe current state of the collective decision system changes every step by exactlyone marble. The process defined by the urn model is memoryless, that is, it hasthe Markov property and we can define a Markov chain with N + 1 states. Wedefine the transition matrix T of the Markov chain by T , = 12 (cid:18) ∆ B (cid:18) N (cid:19) + 1 (cid:19) ,T i +1 ,i = 12 (cid:18) ∆ B (cid:18) i − N (cid:19) + 1 (cid:19) , for i ∈ [2 , N ] ,T i − ,i = 1 − (cid:18) ∆ B (cid:18) i − N (cid:19) + 1 (cid:19) , for i ∈ [2 , N ] ,T N,N +1 = 1 − (cid:18) ∆ B (cid:18) NN (cid:19) + 1 (cid:19) . (18)The steady state is then computed by determining the eigenvectors vTv = λ v . (19)Generally there are several eigenvectors but only one that has no changing signin its elements (all positive or all negative) which represents the equilibriumdistribution of the Markov chain π or − π respectively.17 σ a , b s ϕ = 0 . ϕ = . ϕ = 1 Figure 6: Splitting probabilities σ a,b between the two peaks in the steady stateprobability density π for varied positive feedback intensity ϕ ∈ { . , . , . . . , } .With the steady state at hand, several methods of statistical analysis areavailable. As indicated by the data obtained by simulation shown in Fig. 5,the equilibrium distributions are bimodal for ϕ > .
5. Hence, it is an optionto calculate splitting probabilities (Gardiner, 1985). The splitting probability σ a,b ( x ) gives the probability that the system initialized at s = x will reach thestate s = b before the state s = a . It is calculated by σ a,b ( x ) = Z xa π ( s ) − ds Z ba π ( s ) − ds ! − . (20)Note that this equation is based on a continuous distribution π ( s ) which can beobtained from a discrete distribution (e.g., the above equilibrium distributionbased on Markov theory) by fitting. Another option, which was used here, is toapply Fokker–Planck theory which allows to calculate a continuous equilibriumdistribution directly based on ∆ B ( s ) (Hamann et al., 2010). We set a and b tothe positions of the two peaks in the steady state probability density π that weobtain for ϕ = 1. Results for varied positive feedback intensity ϕ (while keeping a and b constant) are shown in Fig. 6. It is clearly seen that for ϕ < . a or b . This is becausethe steady state probability density π is unimodal for ϕ < .
5. Beginning with ϕ = 0 .
5, when the steady state probability density starts to be bimodal, andwith increasing ϕ the probability of switching from one peak to the other isdecreasing considerably which is an indicator for an effective collective decision.Another statistical property of Markov processes, especially those with bistablepotentials, is the mean first passage time. Of interest is the mean time to switchfrom the collective decision for option A ( s = s with ∆ B ( s ) = 0) to option B ( s = s with ∆ B ( s ) = 0) or vice versa due to symmetry. The switching timeis of particular interest in the context of swarm intelligence because it tellswhether a swarm is able to stay for a considerable time in a given state (e.g.,see Yates et al. (2009)). In case of collective motion, the mean switching time18 τ N Figure 7: Mean first passage times τ between the two collective decisions s and s with ∆ B ( s ) = ∆ B ( s ) = 0 over swarm size N , squares give values accordingto Markov theory, circles give measurements of the urn model, error bars givestandard deviations, lines are fitted functions τ ( N ) = a N b a exp( cN ).describes whether the swarm will be aligned long enough to cover a considerabledistance.Markov theory allows to calculate the mean first passage time using the tran-sition matrix T and specifying a target state i which is defined to be absorbing( T i,i = 1), that is, we convert the system into an absorbing Markov chain. Usingthe fundamental matrix M = ( I − Q ) − (with identity matrix I and Q is thetransition matrix of transitional states), a vector of mean first passage timesfor all transitional states is obtained by Mc with c is a column vector of all1’s (Grinstead and Snell, 1997). In addition, an estimate of the mean switchingtime can be determined numerically in the urn model by generating trajecto-ries s ( t ) and observing the switching behavior. This estimate is a lower bound(finite simulation time, finite number of samples, power-law distribution of firstpassage times). A comparison between theory and simulation along with fittedfunctions ( τ ( N ) = a N b a exp( cN )) is shown in Fig. 7. The switching timescales approximately exponentially with swarm size N . The fact that Eq. 5 isalso a good choice in fitting the mean first passage times can be interpreted asmore than mere coincidence. The mean first passage time may be seen as per-formance measure because longer times reflect a more stable swarm. However,for Eq. 5 we require c < c > Next we want to compare the data from our urn model (Fig. 4(b)) to data frommore complex models, such as the density classification scenario (Hamann and W¨orn,2007). In the following, we apply the more general definition of ∆ B (Eq. 17)but with an invariant payoff M ( s ) = c for a scaling constant c that scales the19 s ∆ s ( s ) / ∆ t (a) Density classification sce-nario (Hamann and W¨orn, 2007), change ofthe ratio of red robots for different timesduring simulation, squares give empiricaldata (from (Hamann et al., 2010)), lines arefitted according to Eq. 17. tϕ (b) Negative exponential function ϕ ( t ) =0 . − exp( − × − t ) fitted to feedbackintensities obtained from the density classi-fication scenario. Figure 8: Comparison of model and results from the density classification sce-nario (Hamann and W¨orn, 2007) and increase of positive feedback over time.average change for different payoffs.The density classification scenario (Hamann and W¨orn, 2007) is about aswarm of red and green agents moving around randomly. Their only interactionis constantly keeping track of those agents’ colors they bump into. Once anagent has seen five agents of either color it changes its own color to that it hasencountered most. Here, s gives the ratio of red agents. The name of this sce-nario is due to the idea that the swarm should converge to that color that wasinitially superior in numbers. It turns out that the averaged change ∆ s ( s ) / ∆ t (see Fig. 8(a)) starts with a curve similar to that of ϕ = 0 in Fig. 4(b) and thenconverges to a curve that is similar to that of ϕ = 0 .
75. That is, ∆ s ( s ) / ∆ t istime-variant and, for example, the above mentioned Markov model would haveto be extended to a time-inhomogeneous Markov model to achieve a better cor-relation with the measurements. Early in the simulation there is mostly negativefeedback forcing values close to s = 0 .
5. The negative feedback decreases withincreasing time, which results finally in positive feedback for s ∈ [0 . , . ϕ in Fig. 4(b), one can say that positive feedbackbuilds up slowly over time in the density classification scenario. By fittingEq. 17 to the data shown in Fig. 8(a) we get estimates for the feedback in-tensity ϕ . From the earliest and steepest line to the latest and only curvewith positive slope in s = 0 . ϕ ∈ [0 , , , . , . , . t ∈ [100 , , , , , ϕ according to our model. In Fig. 8(b), thedata points of feedback intensity ϕ obtained by fitting are shown along with a20 sP (a) Positive feedback probability as mea-sured in the simulation according toEq. 15 (squares) and fitted function ofEq. 21 (line). -0.001 0 0.001 0 0.5 1 s ∆ s ( s ) / ∆ t (b) Predicted average change based on themeasurements shown in (a) and Eq. 17 (line)compared to direct measurements in simula-tion (squares). s π ( s ) (c) Predicted steady state obtained by us-ing the measurements shown in (a) and byapplying Eqs. 17, 18, and 19 (line); directmeasurement in simulation (squares). Figure 9: Measured positive feedback probability, predicted and measured av-erage change ∆ s , and predicted and measured steady state for the density clas-sification scenario.function that was fitted to the data. This result supports the assumption of a‘negative exponential’ increase over time (1 − exp( − t )) of positive feedback inthis system as already stated in (Hamann et al., 2010).The positive feedback probability P ( s ) can also be measured via the ob-served decision revisions according to Eq. 15. This was done for the densityclassification scenario; the results are given in Fig. 9(a), which shows the mea-sured positive feedback probability P ( s ) for time t = 7900. The data was fittedusing the following axially symmetric function P ( s ) = c (cid:16) − c s (cid:17) , for s ≤ . c (cid:16) − c (1 − s ) (cid:17) , else , (21)for constants c and c , which turns out to be a better fit here than the sinefunction. Once the function is fitted, it can be used to calculate the expected21 s ∆ s ( s ) / ∆ t Figure 10: Model fitted to data from Fig. 3B of Yates et al. (2009) (local modelof swarm alignment in locusts) by ϕ = 1 (and c = 4 . × − ); data scaled to s ∈ [0 , s following Eq. 17 (only scaling constant c needs to be fitted).As seen in Fig. 9(b), we obtain a good fit (root mean squared error of < × − ). Indeed, our experience shows that the procedure of fitting P ( s ) insteadof ∆ s directly is much more accurate and needs fewer samples. Especially, itplaces the zeros ∆ s = 0 more accurate, which are important to predict thecorrect steady state of the system. Even with just 100 samples, good fits wereobtained, which indicates that this model could be applied to data from naturalswarm systems. The prediction of the steady state based on the measuredpositive feedback probability following Eqs. 17, 18, and 19 is shown in Fig. 9(c).In comparison to the measurements from simulation, this prediction shows areasonable agreement.Other examples showing similarities to the ϕ = 0 . ϕ = 1. In this paper, we have reported two abstract models of swarms with high gen-erality due to our long-term objective of creating a swarm calculus. The firstmodel focuses on the dependency of swarm performance on swarm density byseparating the system into two parts: cooperation and interference. It explainsthat an optimal or critical swarm density exists at which the peak performanceis reached. With the second model we describe the dynamics of collective deci-sion processes with focus on the existence and intensity of feedback. It gives anexplanation of how the typical cubic functions of decision revision emerge by anincrease of positive feedback over time.22he first model is simple and the existence of optimal swarm densities isa well-known fact. However, to the authors knowledge, no similar model com-bined with a validation by fitting the model to data from diverse swarm applica-tions was reported before except for the hypotheses stated by Østergaard et al.(2001). Despite its simplicity, the model has the capability to give predictionsof swarm performance, especially, if the available data, to which it is fitted, in-cludes an interval around the optimal density. That way this model might serveas a swarm calculus of swarm performance. In addition, we want to draw at-tention to the problem of masking special density-dependent properties by onlyinvestigating the mean performance. The example shown in Fig. 3(c) documentsthe existence of phases in swarm systems.The second model is abstract as well but has a higher complexity and ismore conclusive as it allows for mathematical derivations. Based on this urnmodel for positive feedback decision processes, the emerging cubic function ofdecision revisions can be derived (see Eq. 9). Hence, we generate the function ofdecision revision based on our urn model, which allows for an interpretation ofhow the function emerges while, for example, in (Yates et al., 2009) this functionis measured in a local model. Our model of collective decisions might qualify asa part of swarm calculus because those decision revision functions seem to be ageneral phenomenon in swarms.This model can also be used to predict probability density functions of steadystates in swarm systems. The workflow of measuring the positive feedback prob-ability P ( s ), fitting a function, using this function to calculate the expectedaverage change ∆ s , which can in turn be used to predict the expected proba-bility density function of the steady state (eigenvectors of transition matrix), isaccurate with comparatively small sample numbers.A model of notable similarity was published in the context of ‘sociophysics’(Galam, 2004). It is based on the assumption that subgroups form in collectivedecisions within which a majority rule determines the subgroup’s decision. Theaddition of contrarians, that is voters always voting against the current majoritydecision, generates dynamics that are similar to the reported observations innoisy swarms. Galam’s model is, however, focused on constant sizes of thesesubgroups and their combinatorics while in swarm intelligence these groups andtheir sizes are dynamic.A result of interest is also the particular function of the positive feedbackincrease over time (1 − exp( − t )) in the density classification task (see Fig. 8(b)).It is to be noted that this increase seems to be independent from respectivevalues of the current consensus s . Furthermore, extreme values, such as s ≈ s ≈
0, are not observed. An in-depth analysis of the underlying processes isbeyond this paper but we want to present two ideas. First, the final saturationphase (lim t →∞ ϕ = 0 .
8) is most likely caused by explicit noise in the simulation.The agent–agent recognition rate was set to 0.8 which keeps P ( s, ϕ ) small.Second, the initial fast increase of ϕ (after a transient which might also becaused by the simulation because agents revise their color only after a minimumof five agent–agent encounters) might be caused by locally emerging sub-groupsof homogeneous color within small areas that generate ‘islands’ of early positive23 P(s,φ)+
Figure 11: Time-variant feedback system; here for increasing probability ofpositive feedback.feedback. These properties might, however, be highly stochastic and difficultto measure. Time-variant positive feedback was also observed in BEECLUST-controlled swarms as reported before (Hamann et al., 2012). Hence, a feedbacksystem as given in Fig. 11 seems to be a rather common situation in swarmsystems. A and B are properties of a swarm system that form a feedbacksystem (e.g., amount of pheromone and number of recruited ants). C is a thirdswarm property, which is subject to a time-variant process and which influencesthe feedback of B on A . In terms of the above urn model we can mimic thissituation by saying A in Fig. 11 is the number of blue marbles (w.l.o.g.), B is theprobability of drawing a blue marble, P is the probability of positive feedback(i.e., this edge can also negatively influence A ), and C is an unspecified statevariable that increases positive feedback ( ϕ ) over time and is influenced byan additional, unknown process. This triggers the question of what C can beand how it influences the feedback process independent of the current swarmconsensus s .As seen in Fig. 10, we get maximally positive feedback ϕ = 1 for the dataof Yates et al. (2009) which has the effect that states of low alignment ( s ≈ . On of the main results reported in this paper is that generally applicable swarmmodels, that have simple preconditions, exist. The application of the swarmperformance model necessitates only a concept of swarm density and the appli-cation of the collective decisions model necessitates only a consensus variableof a binary decision. Although both models are simple, they have enough ex-planatory power to give insights into swarm processes such as the interplay ofcooperation and interference and the installation of positive feedback.The two presented models illustrate the methodology that can be applied24o find more models and to extend swarm calculus. The methodology is char-acterized by a combination of a heuristic approach and a simple mathematicalformalism. While the empirical part establishes a direct connection to applica-tions, the formalism allows the integration of superior mathematical methodssuch as Markov theory and linear algebra. That way there is reason for hopethat similar models might be found, for example, for swarm systems showingaggregation, flocking, synchronization, or self-assembly. The main benefit ofsuch models might be general insights in group behavior and swarms but alsodirect applications could be possible. It could be possible to implement variantsof these models on swarm robots if the global knowledge necessary for this kindof models can be substituted by local samplings. The models could also be usedas components and several such components could be combined to form special-ized, sophisticated models. The presented models could be combined to modela swarm showing collective decision–making, such as a BEECLUST-controlledswarm as pointed out above. Once a comprehensive set of such models has beencollected, research on swarms could also be guided by formalisms supplementalto empirical research. Hence, we contend that it is possible to generate a setof models and methods of general applicability for swarm science, that is, tocreate a swarm calculus.
A Details on curve fitting
All curve fitting was done with an implementation of the nonlinear least-squaresMarquardt-Levenberg algorithm (Levenberg, 1944; Marquardt, 1963) using gnu-plot 4.6 patchlevel 1 (2012-09-26) . A.1 Foraging in a group of robots
The data for the curve fitted in Fig. 2(b) is shown in Tab. 1.function P ( N ) = a N b a exp( cN )degrees of freedom 52root mean square of residuals 0.000146389parameter value asymptotic standard error a a b c -0.199589 +/- 0.002932 (1.469%)Table 1: Fitting data, foraging in a group of robots, Fig. 2(b). A.2 Collective decision making based on BEECLUST
The data for the curve fitted in Fig. 2(c) is shown in Tab. 2. see P ( N ) = a N b a exp( cN )degrees of freedom 22root mean square of residuals 0.0515291parameter value asymptotic standard error a a b c -0.0386915 +/- 0.002867 (7.409%)Table 2: Fitting data, collective decision making based on BEECLUST,Fig. 2(c). A.3 Aggregation in tree-like structures and reduction toshortest path
The data for the curve fitted in Fig. 2(d) is shown in Tab. 3. In this caseweighted fitting was used (values ρ = 0 . ρ = 0 . P < P ( ρ ) = a ρ b a exp( cρ )degrees of freedom 21root mean square of residuals 0.0924653parameter value asymptotic standard error a a b c -89.9857 +/- 6.985 (7.763%)Table 3: Fitting data, aggregation in tree-like structures and reduction to short-est path, Fig. 2(d). A.4 Emergent–taxis behavior
The data for the curve fitted in Fig. 3(a) is shown in Tab. 4. Fitting was donein two steps. First, the interference function I ( N ) was fitted. Second, theperformance function P ( N ) was fitted while keeping the parameters a and c fixed. A.5 Emergent–taxis behavior, narrow fit
The data for the curve fitted in Fig. 3(b) is shown in Tab. 5. The parameters a and c of the interference function I ( N ) as obtained in A.4 were reused. Theperformance function P ( N ) was fitted within the narrow interval of N ∈ [20 , a and c fixed.26unction I ( N ) = a exp( cN ) + d degrees of freedom 48root mean square of residuals 0.00479438parameter value asymptotic standard error a c -0.182333 +/- 0.007664 (4.203%) d P ( N ) = a N b a exp( cN )degrees of freedom 41root mean square of residuals 0.0403196parameter value asymptotic standard error a b P ( N ) = a N b a exp( cN )degrees of freedom 1root mean square of residuals 0.0180345parameter value asymptotic standard error a b A.6 Mean first passage times
The data for the curve fitted in Fig. 7 is shown in Tab. 6. Weighted fittingwas applied based on the measured standard deviation and weights scaled by √ τ theor respectively. A.7 Density classification
The data for the curve fitted in Fig. 8(a) is shown in Tab. 7. For times t ∈{ , , } , we set ϕ = 0 as otherwise the fitting would result in ϕ < A.8 Feedback intensities
The data for the curve fitted in Fig. 8(b) is shown in Tab. 8. Weighted fittingwas applied with zero-weight for data points of t < ϕ ( t ) = 0. Data points of t ≥ ≤ t < τ τ meas ( N ) = a meas1 N b meas a meas2 exp( c meas N )function theoretical τ τ theor ( N ) = a theor1 N b theor a theor2 exp( c theor N )degrees of freedom 7rms of residuals ( τ meas ) 0.0613924rms of residuals ( τ theor ) 1.35583parameter value asymptotic standard error a meas1 a meas2 b meas c meas a theor1 a theor2 b theor c theor A.9 Positive feedback probability
The data for the curve fitted in Fig. 9(a) is shown in Tab. 9.
A.10 Swarm alignment in locusts
The data for the curve fitted in Fig. 10 is shown in Tab. 10. We set ϕ = 1 asotherwise the fitting would result in ϕ >
1. Weighted fitting was applied valuesof s < .
17 and s > .
83 had double weight than values of 0 . ≤ s ≤ . References
Arkin, R. C., Balch, T., and Nitz, E. (1993). Communication of behavioral statein multi-agent retrieval tasks. In Book, W. and Luh, J., editors,
IEEE Con-ference on Robotics and Automation , volume 3, pages 588–594, Los Alamitos,CA. IEEE Press.Berman, S., Kumar, V., and Nagpal, R. (2011). Design of control policies forspatially inhomogeneous robot swarms with application to commercial polli-nation. In LaValle, S., Arai, H., Brock, O., Ding, H., Laugier, C., Okamura,A. M., Reveliotis, S. S., Sukhatme, G. S., and Yagi, Y., editors,
IEEE Inter-national Conference on Robotics and Automation (ICRA’11) , pages 378–385,Los Alamitos, CA. IEEE Press.Bjerknes, J. D. and Winfield, A. (2010). On fault-tolerance and scalability ofswarm robotic systems. In Martinoli, A., Mondada, F., Correll, N., Mer-moud, G., Egerstedt, M., Hsieh, M. A., Parker, L. E., and Støy, K., editors,
Proc. Distributed Autonomous Robotic Systems (DARS 2010) , pages 431–444,Berlin, Germany. Springer-Verlag. 28jerknes, J. D., Winfield, A., and Melhuish, C. (2007). An analysis of emergenttaxis in a wireless connected swarm of mobile robots. In Shi, Y. and Dorigo,M., editors,
IEEE Swarm Intelligence Symposium , pages 45–52, Los Alamitos,CA. IEEE Press.Breder, C. M. (1954). Equations descriptive of fish schools and other animalaggregations.
Ecology , 35(3):361–370.Camazine, S., Deneubourg, J.-L., Franks, N. R., Sneyd, J., Theraulaz, G., andBonabeau, E. (2001).
Self-Organizing Biological Systems . Princeton Univer-sity Press.Deneubourg, J.-L., Aron, S., Goss, S., and Pasteels, J. M. (1990). The self-organizing exploratory pattern of the argentine ant.
Journal of Insect Behav-ior , 3(2):159–168.Dussutour, A., Beekman, M., Nicolis, S. C., and Meyer, B. (2009). Noise im-proves collective decision-making by ants in dynamic environments.
Proceed-ings of the Royal Society London B , 276:4353–4361.Dussutour, A., Fourcassi´e, V., Helbing, D., and Deneubourg, J.-L. (2004). Opti-mal traffic organization in ants under crowded conditions.
Nature , 428:70–73.Edelstein-Keshet, L. (2006). Mathematical models of swarming and social ag-gregation.
Robotica , 24(3):315–324.Ehrenfest, P. and Ehrenfest, T. (1907). ¨Uber zwei bekannte Einw¨ande gegendas Boltzmannsche H-Theorem.
Physikalische Zeitschrift , 8:311–314.Eigen, M. and Winkler, R. (1993).
Laws of the game: how the principles ofnature govern chance . Princeton University Press.Galam, S. (2004). Contrarian deterministic effect on opinion dynamics: the“hung elections scenario”.
Physica A , 333(1):453–460.Gardiner, C. W. (1985).
Handbook of Stochastic Methods for Physics, Chemistryand the Natural Sciences . Springer-Verlag, Berlin, Germany.Gautrais, J., Theraulaz, G., Deneubourg, J.-L., and Anderson, C. (2002). Emer-gent polyethism as a consequence of increased colony size in insect societies.
Journal of Theoretical Biology , 215(3):363–373.Goldberg, D. and Matari´c, M. J. (1997). Interference as a tool for designingand evaluating multi-robot controllers. In Kuipers, B. J. and Webber, B.,editors,
Proc. of the Fourteenth National Conference on Artificial Intelligence(AAAI’97) , pages 637–642, Cambridge, MA. MIT Press.Graham, R., Knuth, D., and Patashnik, O. (1998).
Concrete Mathematics: AFoundation for Computer Science . Addison–Wesley, Reading, MA.29rinstead, C. M. and Snell, J. L. (1997).
Introduction to Probability . AmericanMathematical Society, Providence, RI.Hamann, H. (2006). Modeling and investigation of robot swarms. Master’sthesis, University of Stuttgart, Germany.Hamann, H. (2010).
Space-Time Continuous Models of Swarm Robotics Sys-tems: Supporting Global-to-Local Programming . Springer-Verlag, Berlin, Ger-many.Hamann, H. (2012). Towards swarm calculus: Universal properties of swarmperformance and collective decisions. In Dorigo, M., Birattari, M., Blum,C., Christensen, A. L., Engelbrecht, A. P., Groß, R., and St¨utzle, T., edi-tors,
Swarm Intelligence: 8th International Conference, ANTS 2012 , volume7461 of
Lecture Notes in Computer Science , pages 168–179, Berlin, Germany.Springer-Verlag.Hamann, H., Meyer, B., Schmickl, T., and Crailsheim, K. (2010). A model ofsymmetry breaking in collective decision-making. In Doncieux, S., Girard, B.,Guillot, A., Hallam, J., Meyer, J.-A., and Mouret, J.-B., editors,
From Ani-mals to Animats 11 , volume 6226 of
Lecture Notes in Artificial Intelligence ,pages 639–648, Berlin, Germany. Springer-Verlag.Hamann, H., Schmickl, T., W¨orn, H., and Crailsheim, K. (2012). Analysis ofemergent symmetry breaking in collective decision making.
Neural Computing& Applications , 21(2):207–218.Hamann, H. and W¨orn, H. (2007). Embodied computation.
Parallel ProcessingLetters , 17(3):287–298.Hamann, H. and W¨orn, H. (2008). Aggregating robots compute: An adaptiveheuristic for the Euclidean Steiner tree problem. In Asada, M., Hallam, J. C.,Meyer, J.-A., and Tani, J., editors,
The tenth International Conference onSimulation of Adaptive Behavior (SAB’08) , volume 5040 of
Lecture Notes inArtificial Intelligence , pages 447–456. Springer-Verlag.Ingham, A. G., Levinger, G., Graves, J., and Peckham, V. (1974). The Ringel-mann effect: Studies of group size and group performance.
Journal of Exper-imental Social Psychology , 10(4):371–384.Jeanne, R. L. and Nordheim, E. V. (1996). Productivity in a social wasp: percapita output increases with swarm size.
Behavioral Ecology , 7(1):43–48.Jeanson, R., Fewell, J. H., Gorelick, R., and Bertram, S. M. (2007). Emergenceof increased division of labor as a function of group size.
Behavioral Ecologyand Sociobiology , 62:289–298.Karsai, I. and Wenzel, J. W. (1998). Productivity, individual-level and colony-level flexibility, and organization of work as consequences of colony size.
Proc.Natl. Acad. Sci. USA , 95:8665–8669.30ennedy, J. and Eberhart, R. C. (2001).
Swarm Intelligence . Morgan Kauf-mann.Klein, M. J. (1956). Generalization of the Ehrenfest urn model.
Physical Review ,103(1):17–20.Krafft, O. and Schaefer, M. (1993). Mean passage times for triangular transi-tion matrices and a two parameter Ehrenfest urn model.
Journal of AppliedProbability , 30(4):964–970.Lerman, K. and Galstyan, A. (2002). Mathematical model of foraging in a groupof robots: Effect of interference.
Autonomous Robots , 13:127–141.Lerman, K., Martinoli, A., and Galstyan, A. (2005). A review of probabilisticmacroscopic models for swarm robotic systems. In S¸ahin, E. and Spears,W. M., editors,
Swarm Robotics - SAB 2004 International Workshop , volume3342 of
Lecture Notes in Computer Science , pages 143–152. Springer-Verlag,Berlin, Germany.Levenberg, K. (1944). A method for the solution of certain non-linear problemsin least squares.
Quarterly of Applied Mathematics , 2:164–168.Lighthill, M. J. and Whitham, G. B. (1955). On kinematic waves. II. A theoryof traffic flow on long crowded roads.
Proceedings of the Royal Society ofLondon , A229(1178):317–345.Mahmassani, H. S., Dong, J., Kim, J., Chen, R. B., and Park, B. (2009). Incor-porating weather impacts in traffic estimation and prediction systems. Tech-nical Report FHWA-JPO-09-065, U.S. Department of Transportation.Mallon, E. B., Pratt, S. C., and Franks, N. R. (2001). Individual and collectivedecision-making during nest site selection by the ant leptothorax albipennis . Behavioral Ecology and Sociobiology , 50:352–359.Marquardt, D. (1963). An algorithm for least-squares estimation of nonlinearparameters.
SIAM Journal on Applied Mathematics , 11(2):431–441.Milutinovic, D. and Lima, P. (2007).
Cells and Robots: Modeling and Controlof Large-Size Agent Populations . Springer-Verlag, Berlin, Germany.Miramontes, O. (1995). Order-disorder transitions in the behavior of ant soci-eties.
Complexity , 1(1):56–60.Mondada, F., Bonani, M., Guignard, A., Magnenat, S., Studer, C., and Flo-reano, D. (2005). Superlinear physical performances in a SWARM-BOT. InCapcarrere, M. S., editor,
Proc. of the 8th European Conference on Artifi-cial Life (ECAL) , volume 3630 of
Lecture Notes in Computer Science , pages282–291, Berlin, Germany. Springer-Verlag.31embrini, J., Winfield, A. F. T., and Melhuish, C. (2002). Minimalist coherentswarming of wireless networked autonomous mobile robots. In Hallam, B.,Floreano, D., Hallam, J., Hayes, G., and Meyer, J.-A., editors,
Proceedingsof the seventh international conference on simulation of adaptive behavior onFrom animals to animats , pages 373–382, Cambridge, MA, USA. MIT Press.Nicolis, S. C., Zabzina, N., Latty, T., and Sumpter, D. J. T. (2011). Collectiveirrationality and positive feedback.
PLoS ONE , 6:e18901.Okubo, A. (1986). Dynamical aspects of animal grouping: Swarms, schools,flocks, and herds.
Advances in Biophysics , 22:1–94.Okubo, A. and Levin, S. A. (2001).
Diffusion and Ecological Problems: ModernPerspectives . Springer-Verlag, Berlin, Germany.Østergaard, E. H., Sukhatme, G. S., and Matari´c, M. J. (2001). Emergentbucket brigading: a simple mechanisms for improving performance in multi-robot constrained-space foraging tasks. In Andr´e, E., Sen, S., Frasson, C.,and M¨uller, J. P., editors,
Proceedings of the fifth international conference onAutonomous agents (AGENTS’01) , pages 29–35, New York, NY, USA. ACM.Prorok, A., Correll, N., and Martinoli, A. (2011). Multi-level spatial modelsfor swarm-robotic systems.
The International Journal of Robotics Research ,30(5):574–589.Saffre, F., Furey, R., Krafft, B., and Deneubourg, J.-L. (1999). Collectivedecision-making in social spiders: Dragline-mediated amplification processacts as a recruitment mechanism.
J. theor. Biol. , 198:507–517.Schmickl, T. and Hamann, H. (2011). BEECLUST: A swarm algorithm derivedfrom honeybees. In Xiao, Y., editor,
Bio-inspired Computing and Communi-cation Networks . CRC Press.Schneider-Font´an, M. and Matari´c, M. J. (1996). A study of territoriality: Therole of critical mass in adaptive task division. In Maes, P., Wilson, S. W.,and Matari´c, M. J., editors,
From animals to animats IV , pages 553–561.MIT Press.Seeley, T. D., Camazine, S., and Sneyd, J. (1991). Collective decision-makingin honey bees: how colonies choose among nectar sources.
Behavioral Ecologyand Sociobiology , 28(4):277–290.Strogatz, S. H. (2001). Exploring complex networks.
Nature , 410(6825):268–276.Vicsek, T. and Zafeiris, A. (2012). Collective motion.
Physics Reports , 517(3-4):71–140.Wong, G. and Wong, S. (2002). A multi-class traffic flow model – an extensionof LWR model with heterogeneous drivers.
Transportation Research Part A:Policy and Practice , 36(9):827–841. 32ates, C. A., Erban, R., Escudero, C., Couzin, I. D., Buhl, J., Kevrekidis,I. G., Maini, P. K., and Sumpter, D. J. T. (2009). Inherent noise can fa-cilitate coherence in collective swarm motion.
Proc. Natl. Acad. Sci. USA ,106(14):5464–5469. 33unctions P ( s, ϕ ) = ϕ sin( πs ) M ( s ) = c ( c = 0)∆ B ( s ) = 4 M ( s )( P ( s, ϕ ) − . s − . c c c c ϕ c ϕ c ϕ c ϕ ϕ ( t ) = a − exp( bt )degrees of freedom 28root mean square of residuals 0.0173061parameter value asymptotic standard error a -0.000495857 +/- 2.064e-05 (4.161%) b -0.215755 +/- 0.01069 (4.956%)Table 8: Fitting data, feedback intensities, Fig. 8(b).34unction P ( s ) = c (cid:16) − c s (cid:17) , s ≤ . c (cid:16) − c (1 − s ) (cid:17) , elsedegrees of freedom 120root mean square of residuals 0.00562933parameter value asymptotic standard error c c P ( s, ϕ ) = ϕ sin( πs ) M ( s ) = c ( c = 0)∆ B ( s ) = 4 M ( s )( P ( s, ϕ ) − . s − . c2