Nicolas Fontbonne

Nicolas Maudet

Nicolas Bredeche

Cooperative Co-Evolution and Adaptive Team Composition for a Multi-Rover Resource Allocation Problem

Keywords: ad hoc autonomous agent teams, multi-agent systems, marginal contribution, team composition, multi-robots, cooperative co-evolutionary evolutionary algorithms (CCEA), evolutionary computation, evolutionary robotics

In this paper, we are interested in ad hoc autonomous agent team composition using cooperative co-evolutionary algorithms (CCEA). In order to accurately capture the individual contribution of team agents, we propose to limit the number of agents which are updated in-between team evaluations. However, this raises two important problems with respect to (1) the cost of accurately estimating the marginal contribution of agents with respect to the team learning speed and (2) completing tasks where improving team performance requires multiple agents to update their policies in a synchronized manner. We introduce a CCEA algorithm that is capable of learning how to update just the right amount of agents' policies for the task at hand. We use a variation of the El Farol Bar problem, formulated as a multi-robot resource selection problem, to provide an experimental validation of the algorithms proposed.

Introduction

When multiple individuals get together to solve a task, it is sometimes difficult to identify who is actually contributing, and who is not. This is especially problematic when the benefits are equally shared among individuals, including with free-riders who invest a minimal amount of effort. Nature abounds from such examples and various strategies evolved to mitigate the detrimental cost of free-riding, such as partner choice or reputation [START_REF] Noë | Biological markets: supply and demand determine the effect of partner choice in cooperation, mutualism and mating[END_REF][START_REF] Stuart A West | Social semantics: altruism, cooperation, mutualism, strong reciprocity and group selection[END_REF].

The problem of identifying the marginal contribution of individuals has also been studied extensively in cooperative game theory [START_REF] Shoham | Multiagent systems: Algorithmic, gametheoretic, and logical foundations[END_REF]. However, exact methods such as computing the Shapley value [START_REF] Shapley | A value for n-person games. Contributions to the Theory of Games II[END_REF] require strong assumptions (e.g. ability to replay coalitions) and unrealistic computation time, which have led to a flourishing research into the design of approximate methods [START_REF] Tumer | Learning sequences of actions in collectives of autonomous agents[END_REF][START_REF] Wolpert | Optimal payoff functions for members of collectives[END_REF][START_REF] Wolpert | An introduction to collective intelligence[END_REF]. The basic idea of such methods is to identify the individual's contribution by computing the difference between the group performance with and without this very individual (e.g. by removing it or by replacing it with an individual with a default strategy).

In this research, we are interested in ad hoc autonomous agent teams where agents must act together without pre-coordination [START_REF] Stone | Ad hoc autonomous agent teams: collaboration without pre-coordination[END_REF], which implies that the environment is non-stationary as all agents learn in parallel. This requires using methods that can only alter individuals' strategy, with neither a default strategy being known nor the possibility to remove temporarily one particular individual.

Such problem settings have been explored in evolutionary computation for multiagent systems, and notably with Cooperative Co-Evolutionary Algorithms (CCEA) which were first introduced in [START_REF] Mitchell | A cooperative coevolutionary approach to function optimization[END_REF][START_REF] Mitchell | Cooperative coevolution: An architecture for evolving coadapted subcomponents[END_REF] and largely explored since then (see [START_REF] Ma | A survey on cooperative co-evolutionary algorithms[END_REF] for a review of variants and applications). In particular, CCEA have been explored in setups involving multiple robotic agents in tasks such as exploration and foraging [START_REF] Gomes | Novelty-driven cooperative coevolution[END_REF] and environment monitoring [START_REF] Rahmattalabi | D++: Structural credit assignment in tightly coupled multiagent domains[END_REF][START_REF] Rockefeller | Multi-level fitness critics for cooperative coevolution[END_REF][START_REF] Zerbel | The power of suggestion[END_REF].

In this paper, we address the problem of isolating team members so as to identify their contribution within the collective. On the one hand, one could allow only a single agent to learn at a given time, making it possible to measure accurately its contribution as other agents' strategies would remain stationary. On the other hand, several (or all) of the agents' policies could be iterated at the same time, possibly speeding up learning thanks to parallelization.

Balancing between providing accurate estimation of an individual contribution and parallelization of learning actually depends on the context at hand. When rewards are sparse and depend on a single individual's behavioural innovation, it is preferable to bootstrap learning with large-scale exploration. However, whenever team performance increases it is more efficient to turn towards a more conservative search so as to retain previous improvements. Finally, one less obvious situation arises when the synergy between individuals is required to improve performance, for example when two robots are required to open a door, none of which would gain any benefit by acting alone.

The rest of the paper is organized as follows. Section 2 presents the problem of team composition in CCEA. Section 3 presents two CCEA algorithms that enable to tune the number of learning agents within one learning step. The two algorithms differ with respect to how the balance between contribution estimation accuracy and learning parallelization is set: fixed, or self-adaptive. Section 4 presents the problem used for evaluation, which is a modified version of the famous Bar El Farol problem [START_REF] Bredeche | Inductive reasoning and bounded rationality[END_REF] formulated as a multi-robot resource allocation problem where coordination is required (i.e. several resources must be harvested and harvesting is extremely beneficial it the optimal fixed number of robots is met). Results are presented in Section 5 for both the ad hoc version of and the self-adaptive versions of the algorithm.

Cooperative Coevolutionary and Team Composition

In its most simple setup, cooperative co-evolutionary algorithms (CCEA) rely on a collection of independent evolutionary algorithms, with each dedicated to optimizing the policy of one particular agent of the team [START_REF] Jong | Evolutionary computation: a unified approach[END_REF][START_REF] Agoston | Introduction to Evolutionary Computing[END_REF]. Each independent algorithm works to improve the performance of one agent's control parameters using an assessment of the team performance.

Each algorithm maintains a population of parameter sets, which define candidate policies for the agent this algorithm is in charge of. At each generation, performance assessment is computed for various teams. Then, each instance of the CCEA tries to improve its agent's performance by using classic evolutionary operators of selection, mutation and recombination.

The problem faced by CCEA is thus a black-box optimization problem, with the additional twist that evaluation is for the whole team, and optimization is performed at the individual level, thus implying a weak link between team evaluation and the actual individual behaviour. Defining θ as the vector containing the parameters provided by each algorithm of the CCEA, F as the fitness function used to assess team performance and f the fitness value, we have:

F : team parameters θ -→ fitness value f
Figure 1 illustrates the learning loop of a simple multi-agent black-box optimization procedure for cooperative co-evolution. A given algorithm in charge of a particular agent i will provide policy parameters θ i for this agent to be evaluated in a team. These parameters will be evaluated alongside parameters provided by the other agents. The team is then evaluated and a fitness value is returned.

Environment Learning algorithm:

Fig. 1: The learning loop of interaction. All agents submit their own parameters independently for evaluation. They are evaluated at the same time on the environment or task. This return the fitness f of the whole team. Finally this feedback is used by all agent to update their parameters submission for next iteration.

After an evaluation, each agent has to evaluate if the new set of parameters proposed θ i has contributed to the team fitness in a negative, positive or neutral manner. It is necessary to extract from the fitness f = F (θ), the part which depends only on the parameters of the agent i, F i (θ i ).

From an agent's viewpoint, an increase in team performance may be due to others, and may even shadow a decrease in the very contribution the agent performs. In order to mitigate the intrinsically noisy fitness evaluation due to team heterogeneity (team composition changes over time), multiple evaluations of the same set of policy parameters will be performed for a given agent, so that different versions of θ i can be ranked and further selected to create new candidate policies for the next generation. However, obtaining an exact assessment of one individual's contribution to the team remains elusive unless all other individuals follow static policies. Considering teams formed of n team agents, with each agent's evolutionary algorithm maintaining a population of p candidate policies, complexity would be O(p n ) at each generation.

In order to provide results in a reasonable time, CCEA generally relies on a partial evaluation of agents' policy contribution, by evaluating a subset of all possible team compositions at each iteration. Though such implementation breaks down complexity, CCEA algorithms have been shown to have a tendency to prematurely converge to stable states because of a deceptive fitness landscape created by the choice of collaborators for evaluation [START_REF] Funes | Intransitivity revisited coevolutionary dynamics of numbers games[END_REF][START_REF] Panait | Theoretical convergence guarantees for cooperative coevolutionary algorithms[END_REF]. Several methods have been proposed to address these issues, including novelty-based rewards to escape local minima [START_REF] Gomes | Novelty-driven cooperative coevolution[END_REF] or automatically merging populations when agents' behaviours are similar to address scalability issues [START_REF] Gomes | Dynamic team heterogeneity in cooperative coevolutionary algorithms[END_REF][START_REF] Gomes | Challenges in cooperative coevolution of physically heterogeneous robot teams[END_REF]. However, the contribution of one specific agent remains approximated rather than precisely measured.

Cooperative Co-Evolutionary Algorithms with Limited Team Composition Update

In the simplest case, the global fitness F (θ) is the sum of each agent's individual fitness F i (θ i ) for the current team:

F (θ) = i∈agents F i (θ i )
With θ i the policy parameters for agent i, and θ = (θ 0 , . . . , θ i , . . . , θ N -1 ), i.e. the team policy parameters composed of the policies of N agents.

Whenever a single agent updates its policy parameters, the variation in overall fitness F (θ t+1 ) -F (θ t ) will be equal to the variation in the fitness due to the change in behaviour of the agent concerned δF i . This can be written as:

F (θ t+1 ) -F (θ t ) = F -i (θ -i t+1 ) + F i (θ i t+1 ) -F -i (θ -i t ) -F i (θ i t ) = F i (θ i t+1 ) -F i (θ i t ) = δF i (1)
With F -i (θ -i ) the performance of all individuals minus the agent i, assuming θ -i is stationary between t and t + 1. Though an exact value for the contribution of agent i remains unavailable, δF i gives a proxy which provides sufficient information to measure both the direction and amplitude of the change due to agent i's new policy.

Assuming agents are independent, the above equation holds true and can be used as long as only one agent's policy is changed at a time. However, this assumption incurs two important limitations:

-The computational cost of iterating over a single agent's strategy at a time is high (see previous Section), and there is a trade-off between the quality of one agent's contribution estimation and the expected gain at the level of the team (e.g. whenever a single robot is needed to significantly improve team performance, trying with all robots is relevant); -The task may require coupling between the agents' behaviour, and any team fitness improvements may require that several agents change their policy parameters simultaneously (e.g. moving a heavy object can only be done with two robots). If not, a CCEA can get stuck on a local minimum if we assume independence between agents and change only one agent at a time.

In order to address these issues and still retain the benefit of accurate estimation of the agents' contribution, we propose a CCEA algorithm where it is possible to modulate the number of agents that are updated in-between team evaluations. We use a collection of (1+1)-GA algorithm where each (1+1)-GA algorithm i provides the policy parameters θ i for its corresponding agent i, and the whole team is evaluated using all agents with their policy parameters, i.e. θ = (θ 0 , . . . , θ N -1 ).

Algorithm 1: CC-(1+1)-GA k fixed

Introducing k mutants per iteration 1 k ← number of team members to be updated; 2 N ← total number of agents; 3 θparent ← parameters initialisation; 4 fparent ← F (θparent); 5 for gen < nb max generation do 6 ID ← randomly sample k agents;

7 θ ID child ← mutate(θ ID parent ); 8 f child ← F (θ ID child , θ -ID parent ); 9 if f child ≥ fparent then 10 θ ID parent ← θ ID child ; 11 fparent ← f child ; 12 end 13 end
Algorithm 1 details the complete CCEA, which runs multiple instance of (1+1)-GA in parallel, which we will refer to as the CC-(1+1)-GA k fixed algorithm from now on. Each (1+1)-GA algorithm maintains a population of two individuals [START_REF] Beyer | Evolution strategies -a comprehensive introduction[END_REF], a parent θ i parent and a child θ i child . Both are candidate policy parameters for agent i. The parent is replaced when the child fares better, and a new child is created by applying mutation on the new parent. Whenever a child fails to outperform its parent, it is replaced by a new child mutated from the current parent. The mutation operator depends on the problem (e.g. Gaussian mutation, bit-flip, uniform draw).

At each new iteration, k agents are drawn and randomly changed in the team, with 0 < k ≤ N . The k parameter is used to tune the amount of renewal k/N for the team composition in-between iterations of the CCEA algorithm. The k new team members are kept only if they provide an increase in the team fitness. Therefore, the challenge is to find the most efficient size k of the number of team agents to be modified at each CCEA steps. So far, k is fixed beforehand by the user, and may benefits from prior knowledge on the task regarding possible required coupling between agents, in terms of number of agents to change simultaneously to reach the global optimum in terms of team efficiency.

However, such prior knowledge may not be available and a relevant value of k not only depends on the problem (e.g. some problems may require coupling between agents, others may not), but also on the current state of the optimization (e.g. broad initial search steps vs. refined tuning near the optimal solution). In order to address this, we propose the CC-(1+1)-GA k adaptive , where the k parameter is learned during the course of evolution (see Algorithm 2). We propose to choose the number of team members to be updated by using the adversarial bandit learning algorithm EXP3 (Exponential-weight algorithm for Exploration and Exploitation [START_REF] Auer | The nonstochastic multiarmed bandit problem[END_REF]) that tracks the success rate of various possible values for k so far, which means both exploiting the current best value and exploring alternate values. The goal of the adaptation mechanism is to converge to the best possible value for k for the context at hand, i.e. the value that leads to the largest increase of fitness, whether through rare but important increases or through small but frequent increases.

As described in Algorithm 2, we define a set of J possible values for k (k 0 , ..., k J-1 ), each associated with a weight W (k j ) monitoring the success rate of a particular k j . Lines 10-13 of the algorithm detail which k j is selected for a particular iteration. We compute the probability distribution of each k j which depends on the weight W (k j ) and the parameter γ of the algorithm. γ → 1 favours exploration (i.e. the choice of k j will be almost uniform). On the contrary, γ → 0 favours exploitation, taken into account the importance of the weights W (k j ). The fitness gain is normalized between [0, 1] (Line 20) and then used to update the weight W (k j ) (Line 21).

The Multi-Rover Resource Selection Problem

We define a problem that is a variation of the well-known El Farol bar problem [START_REF] Bredeche | Inductive reasoning and bounded rationality[END_REF][START_REF] Wolpert | An introduction to collective intelligence[END_REF] where each individual must choose a day to go to the bar among M possible choices with the criterion of not being too numerous each days. In our setup, we consider the problem where N independent robots must spread over M resources, and where each resource has an optimal capacity c in terms of number

f child ← F (θ ID child , θ -ID parent ); 17 if f child ≥ fparent then 18 θ ID parent ← θ ID child ; 19 fparent ← f child ; 20 R ← tanh ( f child f parent ); 21 W [j] ← W [j] exp( γR P [j]J ); 22 end 23 end
of robots necessary to optimally harvest it. This is illustrated in Figure 2, which provides an example where each robot chooses a resource.

Team performance f is obtained by adding each resource's satisfaction ϕ c . For each resource, its satisfaction ϕ c depends on the number of robots r who chose it. This satisfaction is described by the following equation:

ϕ c (r) = M r exp( -r c ) if r = c r exp( -r c ) else.
(

) 2
where r represents the amount of robots on the resource, M the total number of resources, and c controls the optimal number of robots required for the resource. The satisfaction function diverges from the original formulation of the El Farol Bar problem as the best team performance always implies that the number of robots per resource must be optimal (exactly c robots per resource), even if it implies some resources are left aside or only partially filled. The satisfaction boost for the r = c case ensures that filling a maximum number of resources with the c robots is the optimal strategy. An example of such function with c = 10 is plotted on figure 3. The fitness of a team is then the sum of the satisfaction for all resources:

f = m∈[0,M -1] ϕ c (r m )
Where r m is the number of robots at resource m. The robots must coordinate to optimally fill a maximum number of resources.

To increase the value of this function, it is then necessary to move individuals from crowded resources to resources with fewer robots, up to the extent that resources with the exact number of robots are favoured.

The number of robots N , the number of resources M , and the optimal number of robots per resource c can be modified to change the structure of the problem. In the next section, different instances of this problem are used to study various properties of the Algorithms we proposed in the previous Section. In particular, it is possible to set up the problem so that either single or multiple changes in the team composition may always yield too few or too many permutations in the team distribution over resources for team performance to increase.

Results

Experimental setting

We use the Multi-Rover Resource Selection problem, with different number of resources M , number of agents N , and optimal number of robots per resource c. The three setups used are:

-Setup 1 with N = 120, M = 300, and c = 30. There are many resources, each requiring a large number of robots. The maximum performance could be reached by a team of exactly M × c = 300 × 30 = 9000 agents. Given the limited number of agents, they must spread over a few of the resources (N/c = 120/30 = 4 resources) so that the team reaches optimal performance; -Setup 2 with N = 900, M = 300, and c = 3. The number of robots involved makes it possible to reach the optimal team performance value for this setup (N = M × c) if all agents are uniformly spread over the resources; -Setup 3 with N = 60, M = 7, and c = 10. There exists several configurations of pairing agents and resources which are local optima and cannot be escaped by updating only one agent in the team. Figure 4 gives an example of a sub-optimal configuration for which using k = 1 is detrimental.

When the algorithm gets into such a configuration, all possible updates of a unique agent will decrease the team fitness. Escaping such a local optimum requires either exploring new configurations at the cost of a (hopefully temporary) decrease in team performance (see [START_REF] Gomes | Novelty-driven cooperative coevolution[END_REF] for example using novelty search in CCEA, which is out-of-scope of the current paper) or modifying several agents at the same time (which is possible with > 1).

In the following, we use both the CC-(1+1)-GA k fixed algorithm with either k = 1, 10 or 30, and the CC-(1+1)-GA k adaptive algorithm (using EXP3) for learning the value of k in {1, 10, 30}. All experiments are replicated 32 times. Mean and standard deviation for all algorithm variants are traced. Evaluations is used on the x-axis to provide a fair comparison in terms of computational effort.

Fixed vs. Adaptive Methods for Team Composition Update

Starting with the three variants of the CC-(1+1)-GA k fixed Algorithm, we can observe different learning dynamics depending both on the value of k and the setup at hand.

In the first setup (Fig. 5(a)), we observe a clear benefit for using k = 10 and k = 30 during the first iterations. But this initial gain in performance does not allow it to converge faster when using k = 1. In particular, a value of k = 30 is extremely deleterious for the convergence as it fails to reach the optimum value within the allocated evaluation budget. This tendency is even more visible in the Fig. 4: The resource selection problem has local minimums that can't be escaped by mutating only one agent. In this example, 60 agents must spread on 7 resources by being 10 on 6 of them. In the state where 4 resources are selected by 10 agents, 2 are selected by 7, and 1 by 6, modifying the selection of one agent can only decrease the fitness of the system. The outcome is rather different in the third setup (Fig. 5(c)) as using k > 1 allows to reach better performances and prevent the algorithm to get stuck on a local optimum. Indeed, the algorithm becomes stalled when using k = 1, the structure of the problem making it impossible to improve team performance without considering coupled synergies when updating team members.

Figures 5(a) and (b) show that the CC-(1+1)-GA k adaptive Algorithm follows the curves of the best performing algorithms using a fixed value of k. Figure 5(c) also shows that the adaptive algorithm is able to adapt to a situation where the CC-(1+1)-GA k fixed algorithm would fail because of its fixed k value (here, using k = 1 withholds convergence to an optimal team composition). Overall, dynamically modulating the number of policies updated in the team composition always results in performance curves closely matching the best out of the algorithmic variants using a fixed value of k. This remains true even when the best variant with a fixed k value is outperformed by another variant with a different value of k, which confirms the relevance of the adaptive algorithm to act as an anytime learning algorithm. In other words, the CC-(1+1)-GA k adaptive Algorithm presents the best choice when the problem and the evaluation budget are not known.

Dynamics of Adapting the Number of Team Agents to Update

We analyze how the CC-(1+1)-GA k adaptive Algorithm is changing the value of k throughout evolution for the three setups at hand. Figure 6 represents the median value of k over the 32 repetitions for each of the three experimental setups. We observe that the algorithm switches between the different values for k, and follows different dynamics depending on the setup.

In the first setup, the method slightly favours k = and k = 10 after a few iterations of exploration. This bias is consistent with the performances observed for k fixed, where the k = 30 version is less efficient (see Fig. 5). In the second setup, the value of k decreases during the learning process to stabilize at k = 1, allowing for some fine-tuning of team composition. The third setup displays somewhat different dynamics for the value of k, quickly switching from one value to the other. The difference in performance between the different group sizes is not large enough to make a radical choice, and the method chooses k uniformly at each step without impacting the performances.

Sensitivity of meta-parameters

As described in Algorithm 2, CC-(1+1)-GA k adaptive uses two meta-parameters, which are:

-K = (k 0 , ..., k J-1 ), the set of possible values for k; -the egalitarianism factor γ that determines at each step whether k should be chosen at random (uniform sampling), or selected with respect to the weights of the k values, obtained from the cumulative fitness gain for each particular value of k. The value of γ balances between exploitation and exploration, and the EXP3 algorithm for multi-armed bandit problems has been extensively studied elsewhere [START_REF] Auer | The nonstochastic multiarmed bandit problem[END_REF][START_REF] Richard | Reinforcement learning, Second Edition: An introduction[END_REF];

In the previous section, these meta-parameters were fixed as follow: the set of possible k was limited to {1, 10, 30}, the egalitarianism factor γ was set to 0.1. Figure 7 shows the sensitivity of the algorithm with respect to the different meta-parameters. From top to bottom, the sensitivity to γ and the set of possible k. The general conclusion from this study is that the algorithm is robust and remains a relevant choice for anytime learning. Learning curves remain close to what has been shown previously, with some exceptions for extreme values. In particular, we can observe that:

-γ does not have a significant impact on the algorithm, provided that it is not too small nor too high to efficiently modulate the exploration and exploitation of k's -the algorithm is somewhat also sensitive to the cardinal of the set of possible values for k. When there are too many possibilities to explore, the evaluation of each choice takes more time and is, therefore, less accurate if the context changes too fast. The effect of this exploration can especially be observed for the second setup where the algorithm is less accurate when the value for k can be chosen among 30 possible values (k ∈ [START_REF] Bredeche | Inductive reasoning and bounded rationality[END_REF]30]).

Conclusion

In this paper, we present a cooperative co-evolutionary algorithm (CCEA) that implements a collection of (1+1)-GA algorithms, each endowed with the task to optimize the policy parameters of a specific agent while performance is assessed at the level of the team. Our algorithm acts on team composition by continuously updating a limited number of team agents, depending on the task at hand and the level of completion. Therefore, the main contribution of this paper is to describe an algorithm with a self-adapting team composition update mechanism used throughout learning. We showed that modulating through time the number of new policies added to the current team makes it possible to provide efficient anytime learning, without requiring a priori knowledge on the problem to be solved. Moreover, we show that the algorithm can deal with problems where coupling between agents during learning is mandatory to improve team performance.

Experimental validation was conducted using a variant of the El Farol Bar problem, a famous problem in collective decision making, which was modified to capture a multi-agent resource selection problem. Our algorithm is indeed also relevant for multi-robotic setups, which have been recently studied using various CCEA algorithms [7-9, 15, 16, 25], and future works are currently being conducted in this direction.

Fig. 2 :

 2 Fig.2:The N = 60 robots are represented here as little rovers that each must choose between M = 7 resources. Here the selected agent chooses resource 5. When all robots have made their choices, the satisfaction for each resources are computed and summed to obtain the fitness f of the team.

Fig. 3 :

 3 Fig. 3: Satisfaction function with c = 10 of a given resource, depending on the number of robots that picked it

10 0 10 1 10 2 10 3 10 4 setup 3 Fig. 5 :

 435 Fig. 5: Performance of the CC-(1+1)-GA k fixed algorithm with either k = 1, k = 10 and k = 30, as well as CC-(1+1)-GA k adaptive with k ∈ {1, 10, 30}. Performance f is plotted as mean (solid lines) and standard deviation for the three setups considered with respect to the number of fitness evaluations. Curves are plotted on a x-log scale. There are 32 replications for each algorithm and for each experiment.

10 0 10 1 10 2 10 3 10 4

 4

Fig. 6 :

 6 Fig. 6: Median value of k over the 32 repetitions for the first (left), second (middle) and third setups (right). Curves are plotted on a x-log scale.

10 0 10 1 10 2 10 3 10 4 Fig. 7 :

 47 Fig.7: Sensitivity of the algorithm to meta-parameters. Each column represents one of the different experimental setups. On the rows, one of the meta-parameter is fixed and the other one is varying. On the top, the set of k are fixed but γ varies. At the bottom, γ is fixed but the set of k varies.

  Algorithm 2: CC-(1+1)-GA k adaptive Replacing a varying number of team agents per iteration 1 K ← table of possible number of team members to update simultaneously ; 2 W ← table of weights for each k; 3 P ← table of probability for each k; 4 k ← number of team members to be updated ; 5 N ← total number of agents; 6 θparent ← parameters initialization; 7 γ ← real ∈ [0, 1], parameter for the EXP3 algorithm ; 8 fparent ← F (θparent); 9 for gen < nb max generation do

	10	for j = 1, . . . , J do
	11	P [j] ← (1 -γ) W [j] J i=1 W [i] + γ J
	12	end
	13	kj← random draw in K[] with probabilities P [];
	15	θ ID child ← mutate(θ ID parent );
	16

14

ID ← randomly sample k agents;

Acknowledgements

This work is funded by ANR grant ANR-18-CE33-0006.


Leo Cazenille

Nicolas Bredeche

José Halloy

Automated optimization of multilevel models of collective behaviour: application to mixed society of animals and robots

Keywords: collective behaviour, decision-making, multi-level modelling, evolutionary algorithms, biohybrid systems

Animal societies exhibit complex dynamics that require multi-level descriptions. They are difficult to model, as they encompass information at different levels of description, such as individual physiology, individual behaviour, group behaviour and features of the environment. The collective behaviour of a group of animals can be modelled as a dynamical system. Typically, models of behaviour are either macroscopic (differential equations of population dynamics) or microscopic (such as Markov chains, explicitly specifying the spatio-temporal state of each individual). These two kind of models offer distinct and complementary descriptions of the observed behaviour. Macroscopic models offer mean field description of the collective dynamics, where collective choices are considered as the stable steady states of a nonlinear system governed by control parameters leading to bifurcation diagrams. Microscopic models can be used to perform computer simulations or as building blocks for robot controllers, at the individual level, of the observed spatial behaviour of animals. Here, we present a methodology to translate a macroscopic model into different microscopic models. We automatically calibrate the microscopic models so that the resulting simulated collective dynamics fit the solutions of the reference macroscopic model for a set of parameter values corresponding to a bifurcation diagram leading to multiple steady states. We apply evolutionary algorithms to simultaneously optimize the parameters of the models at different levels of description. This methodology is applied, in simulation, to an experimentally validated shelter-selection problem solved by gregarious insects and robots. Our framework can be used for multi-level modelling of collective behaviour in animals and robots.

Introduction

A "mixed society" is defined as a group of robots and animals that are able to integrate and cooperate: each robot is influenced by the animals, but can, in turn, influence the behaviour of the animals and of the other robots. Individuals, natural or artificial, are perceived as equivalent, and the collective decision process results from the interactions between natural and artificial agents [START_REF] De Schutter | Animal-robots collective intelligence[END_REF][START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF][START_REF] Gribovskiy | Towards Mixed Societies of Chickens and Robots[END_REF][START_REF] Halloy | Towards bio-hybrid systems made of social animals and robots[END_REF]. Robots are useful for a number of reasons [START_REF] Floreano | Bio-inspired artificial intelligence: theories, methods, and technologies[END_REF][START_REF] Mitri | Using robots to understand social behaviour[END_REF]: validating models in silico [START_REF] Garnier | Self-organized aggregation triggers collective decision making in a group of cockroach-like robots[END_REF], inducing stimuli to observe animal feedback [START_REF] Sempo | Integration of an autonomous artificial agent in an insect society: experimental validation[END_REF][START_REF] Gribovskiy | Towards Mixed Societies of Chickens and Robots[END_REF][START_REF] Landgraf | A biomimetic honeybee robot for the analysis of the honeybee dance communication system[END_REF][START_REF] Cazenille | How mimetic should a robotic fish be to socially integrate into zebrafish groups?[END_REF], modulating animal collective behaviour [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF], etc. Earlier research works have used robots in mixed societies to study individual and collective animal behaviours: robots have been mixed with a large and growing number of animal species; including cockroaches in [START_REF] Sempo | Integration of an autonomous artificial agent in an insect society: experimental validation[END_REF][START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF], chicks in [START_REF] Gribovskiy | Towards Mixed Societies of Chickens and Robots[END_REF][START_REF] Gribovskiy | Automated analysis of behavioural variability and filial imprinting of chicks (g. gallus), using autonomous robots[END_REF], honeybees in [START_REF] Landgraf | A biomimetic honeybee robot for the analysis of the honeybee dance communication system[END_REF], fruit flies in [START_REF] Zabala | A simple strategy for detecting moving objects during locomotion revealed by animal-robot interactions[END_REF], fish in [START_REF] Polverino | Zebrafish (danio rerio) behavioural response to bioinspired robotic fish and mosquitofish (gambusia affinis)[END_REF][START_REF] Landgraf | Blending in with the shoal: Robotic fish swarms for investigating strategies of group formation in guppies[END_REF][START_REF] Cazenille | How mimetic should a robotic fish be to socially integrate into zebrafish groups?[END_REF][START_REF] Butail | Influence of robotic shoal size, configuration, and activity on zebrafish behavior in a free-swimming environment[END_REF][START_REF] Butail | Fish-robot interactions in a freeswimming environment: Effects of speed and configuration of robots on live fish[END_REF][START_REF] Bierbach | Insights into the social behavior of surface and cave-dwelling fish (poecilia mexicana) in light and darkness through the use of a biomimetic robot[END_REF][START_REF] De Lellis | Model-based feedback control of live zebrafish behavior via interaction with a robotic replica[END_REF][START_REF] Lei | Computational and robotic modeling reveal parsimonious combinations of interactions between individuals in schooling fish[END_REF][START_REF] Romano | Robot-fish interaction helps to trigger social buffering in neon tetras: The potential role of social robotics in treating anxiety[END_REF] and multi-species societies in [START_REF] Bonnet | Robots mediating interactions between animals for interspecies collective behaviors[END_REF]. See also [START_REF] Romano | A review on animal-robot interaction: from bio-hybrid organisms to mixed societies[END_REF] for a recent review on mixed societies.

Groups of animals are able to reach collective consensus when presented with mutually exclusive alternatives. Over the years, scientists have compiled a large collection of dynamics observed in collective decision-making systems based on experimental observations. These systems can be complex and it can be challenging to build models that appropriately describe the observed behaviours.

Animal societies are systems with a very large parameter space. They can be modelled in numerous ways, using information about individual physiology, individual behaviour, group behaviour and features of the environment [START_REF] Mermoud | Design, Modeling and Optimization of Stochastic Reactive Distributed Robotic Systems[END_REF][START_REF] Mondada | A general methodology for the control of mixed natural-artificial societies[END_REF][START_REF] Schmickl | Biomimetic and bioinspired design of collective systems[END_REF]. The collective behaviour of a group of animals can be viewed as a dynamical system, that exhibits dynamics at several levels of organization (hierarchical organization). One of the difficulties in the modelling process is to find the appropriate levels of description.

In such systems, complementary approaches to modelling (macroscopic vs. microscopic; analytical vs. simulation) can be used: different models deliver the data necessary for the robot design process, provide explicit and analytic descriptions of observed collective behaviour, yield predictions that may be used for the modulation of the collective behaviour of the society, and ease the development of robot controllers.

Models describing dynamical systems can typically be categorized into two groups, describing two different levels of abstraction: macroscopic and microscopic (Fig. 1). There are many studies, mainly in physics, examining methods and applications for both groups of models and the relations between them. Macroscopic models describe the system at the population level [START_REF] Camazine | Self-organization in biological systems[END_REF]. They formalize the dynamics of the system mathematically, but they generally cannot describe the state of individual agents. They cannot be used directly to drive the behaviour of agents in simulation, or to drive the behaviour of robots in experiments. Microscopic models explicitly describe the state of each individual agent (e.g. agent-based models of flocking, like the Vicsek model [START_REF] Vicsek | Novel type of phase transition in a system of self-driven particles[END_REF]). They can capture the individual behaviours and their relations with the environment and, moreover are easier to implement into robotic controllers.

These two kind of models offer complementary descriptions of the system. In this context, collective choices can be described by the stable steady states of a nonlinear system and are governed by control parameters leading to bifurcation diagrams. These diagrams give the mean field asymptotic solutions of the system. Microscopic models can be used to simulate the observed spatial behaviour of animals [START_REF] Jm Amé | Collegial decision making based on social amplification leads to optimal group formation[END_REF][START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF][START_REF] Cazenille | Multi-objective optimization of multi-level models for controlling animal collective behavior with robots[END_REF].

However, working with several models at different levels of abstraction can be difficult and requires appropriate modelling frameworks and methodologies [START_REF] Mondada | A general methodology for the control of mixed natural-artificial societies[END_REF][START_REF] Schmickl | Biomimetic and bioinspired design of collective systems[END_REF].

In particular, these models must be designed and calibrated to all exhibit the same individual and collective dynamics. This poses the problem of navigating between models of different levels of abstraction [START_REF] Mondada | A general methodology for the control of mixed natural-artificial societies[END_REF][START_REF] Mermoud | Design, Modeling and Optimization of Stochastic Reactive Distributed Robotic Systems[END_REF].

Objectives

In this article, we present a methodology that automates the calibration of microscopic target models on the basis of a reference macroscopic model, so that the dynamics of the microscopic model can be described with the same bifurcation diagram as the dynamics of the macroscopic model. The scientific question that we address in this paper is thus the following: how to automatically calibrate models at different level of description to exhibit the same collective dynamics, at all levels of description, for specific sets of parameter values corresponding to a bifurcation diagram leading to multiple steady states.

These tasks are usually performed empirically or manually by different authors, sometimes for the same experimental system. Here, our methodology aims to formalise and automate these tasks and at the same time to optimise the parameters of the models according to the experimental data and the different collective behaviours observed.

The proposed method is general and does not depend on the specific details of the models presented here. The selected models serve as a case study because for this experimental system all three types of models exist in the literature.

We propose to automate the calibration of microscopic models, using information both at the macroscopic level (a pre-existing macroscopic model) and the microscopic level (pre-established knowledge of the animal individual behaviour). The objective of this approach is to automatically calibrate this new microscopic model (i.e. optimize the parameters) to exhibit collective dynamics that fit the predictions of the macroscopic model, with the added ability of accurately simulating the microscopic interaction between individuals.

Methodologies to calibrate microscopic models directly from experimental data, by using optimization algorithms, were already presented in [START_REF] Mermoud | Design, Modeling and Optimization of Stochastic Reactive Distributed Robotic Systems[END_REF][START_REF] Correll | System identification of self-organizing robotic swarms[END_REF]. Our presented methodology improves upon these works by optimizing a microscopic model to correspond to the dynamics exhibited by an entire bifurcation diagram.

This means that in this type of system, several states can coexist for the same parameter values. Moreover, the calibration takes into account the variation of one or more bifurcation parameters which make the system's states appear or disappear, i.e. different types of coexisting collective behaviour.

Our methodology calibrates models with generalization capabilities, enabling them to exhibit different dynamics for different experimental parameters.

In the following, we apply our methodology to the collective decision-making problem described in [START_REF] Jm Amé | Collegial decision making based on social amplification leads to optimal group formation[END_REF][START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF], where a group of cockroaches must reach a consensus on a preferred resting site (a shelter ). These papers introduced an experimentally validated ordinary differential equations (ODE) model of cockroach shelter-selection dynamics.

Here, starting from a mean-field ODE model (macroscopic), we use experimental data on individual cockroach behaviour from [START_REF] Jeanson | A model of animal movements in a bounded space[END_REF][START_REF] Garnier | Self-organized aggregation triggers collective decision making in a group of cockroach-like robots[END_REF] as a-priori microscopic information. We show how our method can be used to calibrate a target model using these two sources of information. We consider two target models: a Markov-Chain (MC) agent-based microscopic model [START_REF] Cazenille | Multi-objective optimization of multi-level models for controlling animal collective behavior with robots[END_REF] and an agent-based Hybrid model, combining macroscopic and microscopic information, that was already used with manually defined parameters in [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF]. These models drive the behaviour of virtual agents in simulation. We use evolutionary algorithms to automatically calibrate the parameters of the MC and Hybrid models. They are validated by comparing their shelter-selection dynamics to those exhibited by the MF model.

In [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF], Halloy et al integrate robots into a group of cockroaches to modulate their collective behaviour. Here, we consider this problem in simulation, with cockroaches agents and a small number of robotic agents. The cockroach agents are driven by the MC and Hybrid models optimized previously in animal-only simulations. The robotic agents are driven by MC and Hybrid models with human-calibrated parameters. We show that it is possible to program the robots to modulate the collective behaviour of the whole society. As such, we show that our methodology could help the design of robotic controllers to modulate the collective behaviour of societies in biohybrid systems (societies of animals and robots).

Macroscopic models can convincingly describe collective dynamics, but cannot be implemented directly in robotic controllers. Robot controllers are intricately microscopic, as they describe the behaviour of individual agents. One major challenge must be overcome to design appropriate robotic behaviour in mixed-societies of animals and robots: how to go from the collective decision dynamics observed in animals to an algorithmic implementation in robots. In previous studies on mixed societies, this process has been carried out empirically. For example, in [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF] collective decision-making in cockroaches is modulated using robots. The authors used observation both to build a macroscopic model and to program the robot behaviour by tinkering. Although the results are promising, designing the robot behaviour proved very challenging, suggesting that automation would be highly beneficial.

The relevance of our approach is not limited to the field of animal collective behaviour and biohybrid systems, as it tackles the problem of automatically moving between models at different levels of description (from macroscopic to microscopic), a key problem in the modelling of nonlinear dynamical systems. It is especially relevant to the design of mixed-societies robotic controllers. Moreover, in a large portion of the literature, the calibration of model parameters is only done for specific solutions -typically only one state of the system. Here, with our methodology we can automatically calibrate all models for a set of states of the system corresponding to a bifurcation diagram. Macroscopic and microscopic models complement each other, describing dynamics at different levels. Microscopic models can include spatial information about the agents, which enables them to be used to simulate the modelled behaviour.

Microscopic Level

Methods

We simulate the experimental set-up from [START_REF] Jm Amé | Collegial decision making based on social amplification leads to optimal group formation[END_REF][START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF] (Fig. 2): this set-up is composed of a circular arena with two identical shelters (resting sites). Each shelter is sufficiently large to host the entire insect group. Two species of cockroaches are considered: P. americana and B. germanica. The cockroaches choose collectively to rest under one of these shelters [START_REF] Jeanson | Selforganized aggregation in cockroaches[END_REF][START_REF] Jm Amé | Collegial decision making based on social amplification leads to optimal group formation[END_REF]. Individuals have no a priori information about the shelters occupation and spatial position, and decide only between staying under a shelter and leaving it to search for another. The cockroaches tend to aggregate under the shelters. This set-up is well adapted to the study of collective decision-making because it allows to quantitatively analyse the interplay of social and environmental mechanisms leading to collective choices. Group-living animals have to choose between alternative resource sites. In this context, a central question includes determining which individuals induce the decision, when and how [START_REF] Jm Amé | Collegial decision making based on social amplification leads to optimal group formation[END_REF][START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF].

We consider three models of cockroach collective behaviour in a shelter selection problem (Table 1): a macroscopic mean-field model (MF model), a microscopic markov chain model (MC model), and a hybrid model combining macroscopic and microscopic levels of abstraction. All three models can handle time-discrete data. The MF model does not include extended spatial information on individuals, while the MC and Hybrid models include explicit spatial information. Because the latter also include a microscopic component, they can be implemented as robotic controllers. A classification of models according to their level of abstraction can be found in [START_REF] Mermoud | Design, Modeling and Optimization of Stochastic Reactive Distributed Robotic Systems[END_REF].

Name

Experimentally Validated Wall-Following Behaviour Constant Speed Mean Field (MF) yes [START_REF] Jm Amé | Collegial decision making based on social amplification leads to optimal group formation[END_REF][START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF] no yes Markov Chain (MC) yes [START_REF] Jeanson | A model of animal movements in a bounded space[END_REF][START_REF] Garnier | Self-organized aggregation triggers collective decision making in a group of cockroach-like robots[END_REF] yes no Hybrid partially [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF] no yes

Table 1: Comparison of the models studied. The MF model is a global description of the problem. The MC model is an agent-based model using a Markov chain representation. The Hybrid model combines macroscopic information (nonlinear propensities drawn from the MF model) and spatial information (with an approach similar to that of the MC model).

Mean field description: Ordinary differential equation model

Halloy et al. [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF] describe a mathematical model of the collective dynamics of mixed groups of cockroaches and robots in a shelter-selection problem (from [START_REF] Jm Amé | Collegial decision making based on social amplification leads to optimal group formation[END_REF]). This model was designed to take in account the following experimental facts: (i) individuals explore their environment by moving randomly, and randomly reach the shelters; (ii) they rest in shelters according to their quality (in this set-up, it is mainly determined by darkness); (iii) the presence of conspecifics influences their behaviour, through social amplification of their resting time; and (iv) no long-range (across shelters) interactions occur among individuals. This model describes mixed groups of animals and robots, in set-ups with two shelters; animals and robots exercise an equivalent influence on the collective decision-making process, and exhibit similar and homogeneous behaviour.

The model is used as a quantitative explanation as well as an overall (macroscopic)

guidance for the design of the robot and its controller. The model takes in account the fact that robots and insects do not occupy the same surface.

The following set of ODEs represents the evolution of the number of individuals in each shelter (and outside), in a set-up with two shelters:

dx i dt = x e µ i 1 - x i + ωr i S i
Probability of animals to join site i

-x i θ i 1 + ρ x i +βr i S i n
Probability of animals to leave site i

(1)

dr i dt = r e µ ri 1 - x i + ωr i S i Probability of robots to join site i -r i θ ri 1 + ρ r γx i +δr i S i
nr Probability of robots to leave site i

(2)

C = x e + x 1 + x 2 , R = r e + r 1 + r 2 , M = R + C (3)
We use a mean field description of the system, instead of a exact representation, to take the fluctuations of the system into account. Table 2 describes the parameters of this ODE model. Variables x i and r i represent respectively the numbers of cockroaches and robots present in shelter i; and x e and r e the numbers outside the shelters. Parameters C and R are respectively the total numbers of cockroaches and robots (state variables). Parameter M is the total number of agents (cockroaches and robots). Parameter S i is the carrying capacity of shelter i. The parameter ω corresponds to the surface of one robot expressed as a multiple of the surface of one cockroach. Equations 1 and 2 take into account the probabilities of animals and robots to join or leave (corresponding to 1/mean resting time) a site. The parameter µ i is the maximal kinetic constant of entering the shelters for insects; µ ri is the equivalent parameter for robots. The parameter θ i is the maximal probability of leaving a shelter for insects (θ ri for robots). The parameters ρ and n characterise the influence of the insect conspecifics (ρ r and n r for robots). When both shelters are identical (as it is the case in this study), the parameters describing them are equal: S 1 = S 2 ; µ 1 = µ 2 ; µ r1 = µ r2 ; θ 1 = θ 2 ; θ r1 = θ r2 . Parameters γ, β and δ are respectively the influence of cockroaches on robots, of robots on cockroaches, and of robots on robots. The greater they are, the greater the mutual influences. The influence of animals on animals is equal to 1.0, and is not considered in [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF]: the assumption is made that this parameter is imposed by biology, and cannot be changed in experiments. However, parameters γ, β and δ can be modulated by changing the design of the robots, either in term of hardware or control (behaviour). In [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF], the robots are coated with a pheromone, as the interaction dynamics of cockroaches societies is mainly chemotactile. A higher concentration of pheromone corresponds to a higher value of β.

Because of crowding effects, the probability of an individual joining a shelter decreases with its level of occupancy. We define a measure σ = S/C, corresponding to the sites' carrying capacity as a multiple of the total population. When no robots are present (R = 0) and only animals are considered, two different dynamics are observed. The bifurcation point is close to σ = 0.8 for P. americana, and σ = 1.0 for B. germanica. Before the bifurcation point (0.4 ≤ σ < 0.8 for P. americana, 0.4 ≤ σ < 1.0 for B. germanica), only one configuration exists, corresponding of an equipartition of the individuals (x 1 /C = x 2 /C = 1/2, x e = 0). After the bifurcation point (σ > 0.8 for P. americana, σ > 1.0 for P. americana), two stable configurations exist, corresponding to all individuals in one of the shelters (either [START_REF] Jm Amé | Collegial decision making based on social amplification leads to optimal group formation[END_REF]. Only results with a population of 50 cockroaches are represented in Fig. 3, but similar dynamics are observed with different population sizes. [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF] and [START_REF] Jm Amé | Collegial decision making based on social amplification leads to optimal group formation[END_REF]. We only consider the case where M = 50. In set-ups with two shelters, the MF, MC, and Hybrid models have 18, 46 and 20 parameters respectively. The parameter values used for the MC and Hybrid models were obtained through the calibration process described in Sec. 3 and Supplementary S1. The animal influence on animals is equal to 1.0, and is kept constant in [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF]: we assume that this parameter is imposed by biology and cannot be changed in experiments. All parameters of the MC and Hybrid models are optimized, using the method described in Supplementary S1, to exhibit the collective dynamics described in the reference MF model. Shelters are in the central zone. When an agent enters a shelter, it has a probability of stopping for a random duration before exiting the shelter. The probability of stopping under a shelter depends on the number of neighbours present in the shelter, and can differ for each shelter. Only 10 neighbours are considered in our experiments. In set-ups with two shelters, this model has 46 parameters per population.

x 1 ≈ 0, x 2 ≈ 1, x e ≈ 0 or x 1 ≈ 1, x 2 ≈ 0, x e ≈ 0)

Markov chain model

We use the Markov chain described in [START_REF] Cazenille | Multi-objective optimization of multi-level models for controlling animal collective behavior with robots[END_REF] as an agent-based model of cockroach and robot behaviour. This model is inspired by the agent-based aggregation models in [START_REF] Jeanson | A model of animal movements in a bounded space[END_REF][START_REF] Garnier | Self-organized aggregation triggers collective decision making in a group of cockroach-like robots[END_REF] that describe the collective behaviour of cockroaches in a similar set-up.

Cockroaches tend to follow the walls of the arena when they are already close to them. The model defines two zones in the arena. The ring area that borders the walls of the arena is called the peripheral zone, while the rest of the arena is labelled as the central zone. In the peripheral zone, agents follow a wall-following behaviour for a random number of time steps (the mean time is denoted τ c,exit ). In the central zone, agents follow a random-walk behaviour, with trajectories composed of a recurring alternation of straight lines (of randomly chosen length, with a mean length of l c ) and rotations (of randomly chosen angles, with a geometric mean of a c ). The shelters are all in the central zone. We do not model the actual trajectories of cockroaches.

When agents enter a shelter, they have a probability of stopping (parameter s c,i,n ) for a random duration (parameter τ c,i,n ) before moving away from the shelter. Similarly to [START_REF] Garnier | Self-organized aggregation triggers collective decision making in a group of cockroach-like robots[END_REF], this probability depends on the number of agents present under the shelter, as cockroaches are gregarious during their resting period. However, in this model (as opposed to [START_REF] Jeanson | A model of animal movements in a bounded space[END_REF][START_REF] Garnier | Self-organized aggregation triggers collective decision making in a group of cockroach-like robots[END_REF]), the probability of stopping when under a shelter differs between shelters: this model is more general, and can be used to describe more complex behaviours with asymmetric decision-making dynamics [START_REF] Cazenille | Multi-objective optimization of multi-level models for controlling animal collective behavior with robots[END_REF]. When the cockroaches are not under a shelter, their movements are not influenced by the presence or absence of neighbours.

Figure 4 represents the Markov chain used in this model. The relevant model parameters are found in Table . 2.

We use two different parameter sets of the MC model to describe either cockroaches or robot behaviour. However, all cockroaches are considered to exhibit an homogeneous behaviour (all cockroach agents share the same parameters). Similarly, all robots agents share the same parameters. We make this choice for simplicity reasons: the main message of this paper can be explored without the need to take individual variability into account.

Hybrid model

Here, we introduce a model of the collective behaviour of cockroaches using information at both macroscopic and microscopic levels of abstraction. We call this multi-level model the 'Hybrid' model. This model was already used with manually defined parameters in [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF], but was not formally described previously. This hybrid model was done to facilitate the development of the behavioural architecture of the robots [START_REF] Asadpour | Robot-animal interaction: Perception and behavior of insbot[END_REF][START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF]. Compared to the MC model, which presents a biomimetic description of the insects trajectories, the Hybrid model is a compromise between biomimetism and ease of implementation as robotic controller. The robot control architecture is a behaviour-based controller [START_REF] Rc Arkin | Behavior-based robotics[END_REF] composed of a multi-level collection of behaviours. Each behavioural building block can take inputs from the robot sensors and/or from other behavioural building blocks, and send outputs to the robot actuators and/or to other behaviours. The behaviours are arranged in a hierarchy in which the behaviours on the higher levels integrate or arbitrate the ones on the lower levels. At the higher level the Hybrid model is used as a building block that takes into account the speciality of the agents, and thus allows to build the robot controller. The Hybrid model is a crossover between the macroscopic MF model, which easily describes collective behaviour and site occupation, and the microscopic MC model, which details the spatio-temporal behaviour of single agents. As such, it is a multi-agents model (like the MC model), but it also takes into account macroscopic information (like the MF model). For ease of implementation, the hybrid model does not include a wall-following behaviour, and only considers simple arenas with no distinctions between central zone and peripheral zone. The agents can have two states: moving, or resting under a shelter. In contrast to the MC model, the agents move with a constant speed v.

The Hybrid model builds on the MF model, introducing several parameters from the MF model (Table 2). The Hybrid model has a smaller dimensionality than the MC model: in set-ups with two shelters, the MC models has 46 parameters (43 without counting wall-following dynamics), while the Hybrid model has 20 parameters. It allows the Hybrid model to be easier to calibrate than the MC model.

Figure 6 describes the Hybrid model using a Markov chain representation of the behaviour of a single agent. When the agents are not under a shelter, they follow a random-walk behaviour (microscopic behaviour). As in the MC model, this randomwalk behaviour involves trajectories composed of a recurring alternation of straight lines of randomly chosen length, with a mean length of l, and rotations of randomly chosen angles, with a geometric mean of a. When agents enter a shelter, they stop, and have a probability of leaving the shelter at each subsequent time-step. This probability, taken from Eq. 1 and 2, is computed using macroscopic information. For the cockroaches agents, this probability of leaving the shelter is defined as:

θ i 1 + ρ x i +βr i S i n (4)
For the robotic agents, this probability of leaving the shelter is defined as:

θ ri 1 + ρ r γx i +δr i S i nr (5)
This behaviour can be described as macroscopic, as it requires information about the density of agents under the shelter. This combination of microscopic and macroscopic components makes it a multi-level (or hybrid) model.

As in the MC model, we use two different parameter sets of the hybrid model to describe either cockroaches or robot behaviour. However, all cockroaches are considered to exhibit an homogeneous behaviour (all cockroach agents share the same parameters). Similarly, all robots agents share the same parameters.

Figure 5 presents examples of the trajectories of single cockroaches in a simulation with a population of 50 cockroaches with the Hybrid model. 2 gives the parameters of the model. The model describes two kinds of behaviour: when the agents are not under a shelter, they will exhibit a random-walk behaviour, following a recurring alternation of straight lines and rotations. This behaviour can be described as microscopic because agents use only local information to determine their course of action. When agents encounter a shelter, they stop. At each subsequent timestep, the stopped agent has a probability of

θ i 1+ρ x i +βr i S i
n (for cockroaches) or

θ ri 1+ρr γx i +δr i S i
nr (for robots) of leaving the shelter and returning to random-walk behaviour. This behaviour can be described as macroscopic, as it requires information about the density of agents under the shelter. This combination of microscopic and macroscopic components makes it a multi-level (or hybrid) model.

Models calibration

To use the MC and Hybrid models describing animal behaviour in simulation, we must calibrate them to exhibit the same decision-making dynamics as the MF model. As the MF model is parametrized using experimental data, it allows the MC and Hybrid models to accurately describe the (macroscopic) site-selection dynamics of the cockroaches. The calibration process is described in Fig. 7. We optimize the parameters for the individual cockroaches in the MC and Hybrid models. The optimization algorithm will identify interesting solutions (i.e. sets of model parameters). The performance of each proposed solution will then be tested ("evaluated"), by first computing a collection of simulations where agents are driven by the model with the parameters sets associated with the solution, and second quantifying how the agents behaviors in those simulations match the dynamics of the reference (MF) model. The latter is achieved as follows. We consider simulations with differing values of the number of robots (N ) and of the sites' carrying capacities (S). This allows us to compute the bifurcation diagrams of each tested solution (w.r.t. σ = S/N ). This bifurcation diagram will finally be compared to the bifurcation diagram of the reference (MF) model. The performance score ("Fitness") of each individual corresponds to the difference between the solution diagram and the reference diagram measured by the Hellinger distance. The optimization algorithm will iteratively search for highperforming solutions, based on this score.

Table 2 lists the parameters of these two optimized models. Instances of the MC and Hybrid models using these parameters are simulated for different values of σ. This yields bifurcation diagrams for each optimized individual, similar to those in Fig. 3.

As there is little a priori information about the parameter space, and as it is relatively high-dimensional, we use the state-of-the-art CMA-ES evolutionary optimization method [START_REF] Auger | A restart CMA evolution strategy with increasing population size[END_REF] to optimize the parameters of the MC and Hybrid models. To evaluate the difference between two parameter sets, we use a distance metric between the two resulting bifurcation diagrams. This method is described in the Supplementary S1.

Results

Models Calibration

Our goal is to find parameter sets of the MC and Hybrid models describing animal and robot behaviour so that the resulting collective dynamics, observed in simulations, fit the solutions of the MF model describing animal behaviour. We use use the methodology presented in Sec. 2.4. We consider two types of simulations, for both P. americana and B. germanica cockroach species. The first type describes a purely biological system, with only 50 cockroaches (either P. americana or B. germanica) and no robots. It is used as the biological reference case. The second type is devoted to biohybrid groups made up of 45 cockroaches (either P. americana or B. germanica) and 5 robots. The number of robots is kept small to reflect the settings used in a mixed-society experiment [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF], where a minority of robots can control the whole mixed group behaviour. The parameters sets of models describing robot behaviour are chosen empirically.

We consider populations of 50 individuals. Similar results are observed with populations of 16 and 100 (results not shown).

Figure 8 shows to the distribution of agents in the two shelters, using parameters from the best-performing optimized individuals after 100 optimization runs. Panels A and C show results from simulations with 50 cockroaches and no robots. Panels B and D show results from simulations with 45 cockroaches and 5 robots. Only results from the bifurcation diagram at selected values of σ are shown. More generally, results before the bifurcation point (σ < 0.8) are similar to results at σ = 0.4, and results after the bifurcation point (σ ≥ 0.8) are similar to results at σ = 1.2.

Both the MC and Hybrid models can be optimized to approximate correctly the decision-making dynamics described by the MF model, as shown in Fig. 8. Our methodology can generate many different parameter sets for the MC and Hybrid models. Optimized parameters that produce the collective dynamics described by the MF model can be associated to highly variable agent behaviour. In the MC model, the parameter d, the diameter of the central zone of the arena, is optimized: when this parameter is very close to the diameter of the peripheral zone, the resulting agents do not exhibit any wall-following behaviour. In the MC and Hybrid models, the parameters that influence stopping behaviour (s c,i,n , τ c,i,n , θ i ) vary less than the other parameters, with only a few islands of relevant values in the explored ranges.

D hellinger (P, Q) = 1 √ 2 d i=1 (
√ P i - √ Q i ) 2
is the Hellinger distance [START_REF] Deza | Dictionary of distances[END_REF]. This approach is described in detail in the Supplementary S1

We show that simulations performed with 45 cockroaches and 5 robots exhibit the same dynamics as the simulations of groups with 50 cockroaches and no robots (Fig. 8). In this case the robots are governed by the same behavioural models as the insects, but do not have the same parameter sets as those used to describe the natural behaviour of the cockroaches. The detailed microscopic behaviours of the robots, e.g. trajectories and movement patterns, can be very different from the microscopic behaviours of the animals. Nevertheless, we show that our methodology can be used to optimize the parameters of robot behavioural models in biohybrid systems to mimic correctly the decision-making dynamics of the animals.

Modulation of collective behaviours

We also test if our methodology could be applied to calibrate robot models so that they can modulate the collective dynamics of a group of cockroaches and robots. We consider simulations, for both P. americana and B. germanica cockroach species. In these simulations, cockroach behaviour is described by the MC or Hybrid models calibrated to have the same collective dynamics as the reference MF model. Robot behaviour is also described by the MC or Hybrid models, but with different parameter sets. We find relevant parameters of models describing robot behaviour empirically, for different proportion of robots, and for different values of σ. These examples of results of the modulation of the mixed-society of cockroaches and robots are found in Fig. 9. This shows that our methodology could be also applied to optimize parameter sets of models describing robot behaviour for a modulation task. Additionally, Fig. 9 also presents examples of modulation of a mixed-society described by the MF model. In this case, we find empirically interesting robot behaviour related parameters of the MF model. These preliminary results suggest that a larger number of robots than 5 may be needed to modulate a population of 50 agents. We will investigate how robots models can be calibrated automatically to modulate the collective dynamics of a mixed-society in a subsequent study.

Discussion

We tackle the problem of moving between models of different levels of abstraction in the context of animal collective decision-making. Animal collective behaviour can be described macroscopically (analytical description of the behaviour of the population) or microscopically (explicit description of the behaviour and states of individuals and their interactions with the environment). The two types of models are complementary. Our methodology enables translation from one to the other: we automatically optimize the parameters of microscopic target models on the basis of a reference macroscopic model from the literature. We apply this methodology to the cockroach shelter-selection problem described in [START_REF] Jm Amé | Collegial decision making based on social amplification leads to optimal group formation[END_REF][START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF]. The Mean Field macroscopic model used as a reference is described in [START_REF] Halloy | Social integration of robots into groups of cockroaches to control self-organized choices[END_REF]. We consider two target models, both agent-based. The MC model [START_REF] Cazenille | Multi-objective optimization of multi-level models for controlling animal collective behavior with robots[END_REF] is a microscopic model inspired by the literature on individual cockroach behaviour [START_REF] Jeanson | Selforganized aggregation in cockroaches[END_REF][START_REF] Garnier | Self-organized aggregation triggers collective decision making in a group of cockroach-like robots[END_REF]. The Hybrid model uses both macroscopic and microscopic information. Both the MC and Hybrid models can be used both to replay the behaviour of animals in simulation and for implementation as robot controllers. We automatically generate the parameters of the MC and Hybrid models for cockroach agents, calibrating them to display the collective behaviour and site-selection dynamics described in the Mean Field model.

The MC and Hybrid models presented in this study can directly be implemented into a robot. However, here our approach does not explain how to translate them automatically into robotic controllers. This question was tackled in [START_REF] Lopes | Supervisory control theory applied to swarm robotics[END_REF] by using formal methods and supervisory control theory to automatically generate robot controller code and to validate it so that it translates into robot behaviours matching a given formal specification.

Here, we use a user-defined metric (the bifurcation diagrams) to compare the results in simulation of two different pairs of models and parameter sets, during the optimisation process. As a result, the obtained models and parameter sets only approximate the collective dynamics of individuals, and not the behaviours (trajectories) of each individual. This could be further improved by considering the calibration of these models as a multi-objective process where the first objective would cater to the collective dynamics, and additional objectives would cater to the individual behaviour of the individuals. Alternatively, one could design a methodology similar as the like in [START_REF] Li | Turing learning: a metric-free approach to inferring behavior and its application to swarms[END_REF], where no metric is specified, by co-evolving simultaneously models of robot behaviours and classifiers of the resulting behaviour in simulation. These classifiers would be trained to identify whether or not the resulting behaviours of the optimised models are distinct from the behaviours from reference experiments.

More generally, complex systems exhibit multi-level dynamics (hierarchical organisation), with both global and local behavioural patterns. Recent studies have investigated the micro-macro link : the relationship between macroscopic and microscopic descriptions of multi-level behavioural dynamics [START_REF] Hamann | A framework of space-time continuous models for algorithm design in swarm robotics[END_REF][START_REF] Reina | A quantitative micro-macro link for collective decisions: the shortest path discovery/selection example[END_REF][START_REF] Martinoli | Modeling swarm robotic systems: A case study in collaborative distributed manipulation[END_REF][START_REF] Yamins | Towards a theory of local to global in distributed multi-agent systems (i)[END_REF]. This problem also applies to the design of swarm group robotic controllers [START_REF] Hamann | A framework of space-time continuous models for algorithm design in swarm robotics[END_REF][START_REF] Martinoli | Modeling swarm robotic systems: A case study in collaborative distributed manipulation[END_REF][START_REF] Lerman | A review of probabilistic macroscopic models for swarm robotic systems[END_REF][START_REF] Vigelius | Multiscale modelling and analysis of collective decision making in swarm robotics[END_REF][START_REF] Mermoud | Real-time automated modeling and control of self-assembling systems[END_REF]. Our methodology is a first step toward the automatic generation of controllers for robots in a mixed society of animals and robots. Mixing animals and robots can be useful for the study of animal behaviour, and even to modulate their individual or collective behaviour.

Few works in the literature on animals and robotics attempt to tackle the problem of transitioning from models of one level of abstraction (reference model) to another level of abstraction (target model). Moreover, these studies have generally considered the transition from microscopic to macroscopic models [START_REF] Martinoli | Modeling swarm robotic systems: A case study in collaborative distributed manipulation[END_REF][START_REF] Degond | Continuum limit of self-driven particles with orientation interaction[END_REF][START_REF] Degond | A hierarchy of heuristicbased models of crowd dynamics[END_REF]. The transition methodology adopted by these studies is incremental, and relies on the creation of intermediate models, dealing with both macroscopic and microscopic information, and that share some parameters with both the reference and the target models. In [START_REF] Martinoli | Modeling swarm robotic systems: A case study in collaborative distributed manipulation[END_REF], this methodology is applied to go from a microscopic model to a macroscopic model of the behaviour of a swarm of autonomous robots in a collaborative task. The resulting model outperforms human-calibrated macroscopic models. In [START_REF] Degond | Continuum limit of self-driven particles with orientation interaction[END_REF], a time-continuous kinetic mean field version of the Couzin-Vicsek model is obtained from its discrete microscopic version. In [START_REF] Degond | A hierarchy of heuristicbased models of crowd dynamics[END_REF], continuous macroscopic models of pedestrian behaviour are obtained from discrete microscopic agent-based models. Little work in the literature has investigated how to automate and generalize the transition between models at different levels of abstraction. The transition process can also be more challenging if the reference and target models have no (or few) common parameters, or if their formulation is too different.

Our approach can be applied to other modelling problems involving a transition between different levels of description: either from a macroscopic one to a microscopic one, or the other way around, provided that the considered models at several levels of abstraction already represent the same kind of dynamics. This approach is illustrated here with three models that we designed to represent the same studied system: a macroscopic model that deals with the dynamics of site occupancy, and two models with microscopic components that mirror this problem by describing how a focal agent can move from a resting state under a site to a random walk in the middle or at the edge of the arena. To generalise to other systems, it is necessary to first design several models that would represent a particular system, but at different levels of abstraction.

Then our approach can be used to calibrate the parameters of these models so that they exhibit the exact same dynamics, with respect to a given bifurcation diagram from a given reference model (that can be at any level of abstraction).

In particular, it may be possible to apply our approach to model, calibrate, and modulate the collective behaviour of other species. For instance, it could be used to calibrate models of fish behaviour to match what would be observed in collective decision making experiments. The work of Couzin et al. [START_REF] Iain D Couzin | Uninformed individuals promote democratic consensus in animal groups[END_REF] is a good example of a fish collective behavioural study where models at different levels of abstraction are introduced to describe various features of the observed experiments. In the latter, the authors investigate consensus in fish groups, and asses whether a strongly opinionated minority can exert its influence on group movement decisions. The authors propose two models of the observed behaviours: an adaptive-network model (macroscopic and analytical) and a spatially explicit model (microscopic and in simulations). The authors present the phase diagram of the studied system based on the results of the adaptivenetwork model. Our approach could have been applied here to calibrate other models (for instance the spatially explicit model proposed by the authors) to directly match the dynamics of this phase diagram. Another example would be the work of Calovi et al. [START_REF] Daniel S Calovi | Swarming, schooling, milling: phase diagram of a data-driven fish school model[END_REF], that presents an in-depth analysis of three prominent collective decision-making dynamics of fish: swarming, schooling and milling. That study relies on a fish collective dynamics model from the literature [START_REF] Gautrais | Deciphering interactions in moving animal groups[END_REF] and the authors compute its phase diagram. In this setting, our approach could be used to calibrate another model (microscopic and/or macroscopic) to exhibit the dynamics of a reference phase diagram. These two studies are similar to the present article in the sense that they focus on dynamics with a co-existence of three states (two stables and one unstable states) that correspond to collective choices. All three works rely on bifurcation or phase diagrams to present the transition between these states. These diagrams could in turn be used as reference of our methodology to calibrate any models describing these dynamics.

In our simulation, both the MC and hybrid models exhibit relatively similar performances -however their differences in design imply different challenges to implement them as robotic controllers in experiments ("sim-to-real" task). Both models can easily be implemented and calibrated on robots, as the only major difference between the two models is the type of parameters used. Our algorithm does not use complex features in its implementation on robots: it is only required that each robot knows the number of robots in its neighbourhood. In this case, the implementation of the Hybrid model and the MC model have approximately the same level of complexity. In a robot implementation, the parameters of the behavioural models often need to be adjusted to obtain the same results as in simulation. In the case of the MC model, all parameters will be related to the local interactions of the robot, which may make calibration easier than for the hybrid model (which deals with both macroscopic and microscopic information); however, the Hybrid model has fewer parameters than the MC model, which may make the Hybrid model easier to calibrate in robotic experiments. Overall difficulties will depend on the problem and personal preferences.

Another study could include an application of this methodology to more complex set-ups, with more than two shelters and more than two population types. Our methodology could also be extended by generating microscopic Markov Chain models from scratch, without a priori structural knowledge (i.e. the type and number of states).

Figure 2 :

 2 Figure 2: Experimental set-up used in [2]. It includes two identical shelters and both cockroaches (P. americana, approximate size: ∼ 4cm, surface: 600mm 2 , or B. germanica, size: ∼ 0.25cm, surface: 3mm 2 ) and robots (surface: 1230mm 2 in P. americana set-ups, 6.15mm 2 in B. germanica set-ups) in a circular arena (diameter: 1m for P. americana, 14cm for B. germanica).

Figure 3 :

 3 Figure 3: Bifurcation diagrams and distributions of 50 (C) B. germanica (Panels A and B) or P. americana (Panels C and D) cockroaches in the first shelter as a function of σ [28]. The bifurcation diagrams are represented as bi-dimensional histograms of the results of using 1000 resolutions for each parameter set. In the bifurcation diagrams (Panels A and C), the (greyscale) colour intensity of each bin of the histograms corresponds to the frequency of observed experiments. The diagrams are symmetric for all tested values of σ, so only one shelter is represented. The bifurcation point is close to σ = 0.8 for P. americana, and σ = 1.0 for B. germanica. Before the bifurcation point (0.4 ≤ σ < 0.8 for P. americana, 0.4 ≤ σ < 1.0 for B. germanica), only one configuration exists, an equipartition of the individuals between the two shelters (x 1 /C = x 2 /C = 1/2, x e = 0). After the bifurcation point (σ > 0.8 for P. americana, σ > 1.0 for P. americana), two stable configurations exist, with all individuals concentrated in one of the two shelters (either x 1 ≈ 0, x 2 ≈ 1, x e ≈ 0 or x 1 ≈ 1, x 2 ≈ 0, x e ≈ 0).

Figure 4 :

 4 Figure 4: Markov chain model of the cockroach individual behaviours. The arena contains two zones: the peripheral zone (where agents display a wall-following behaviour), and the central zone (where agents display a random-walk behaviour). Shelters are in the central zone. When an agent enters a shelter, it has a probability of stopping for a random duration before exiting the shelter. The probability of stopping under a shelter depends on the number of neighbours present in the shelter, and can differ for each shelter. Only 10 neighbours are considered in our experiments. In set-ups with two shelters, this model has 46 parameters per population.

Figure 5 :Figure 6 :

 56 Figure 5: Two example trajectories of a one simulated cockroach, using the Hybrid model, in a population of 50 cockroaches. The arena contains two shelters. Each grey line represents the (random-walk) trajectory of one agent. These trajectories are not meant to fit the natural trajectories of actual cockroaches: we designed our models to reproduce qualitatively the observed random exploration. The opacity of the line reflects simulation time. The timeframe of all simulations is 8 hours.

Figure 7 :

 7 Figure 7: Workflow of the automated calibration of models by optimization. The optimized bifurcation diagram and the reference bifurcation diagram are both converted to one-dimensional histograms, by normalizing the sum of all bin values to 1.0. The optimizer will maximizes the fitness, which is computed by the formula: F itness = 1.0 -D hellinger (B optimized /N u , B ref erence /N u ) where N u is the number of columns in the bifurcation diagrams (10) and B optimized and B ref erence are one-dimensional histogram versions of the respective bifurcation diagrams. The term N u is a normalization term.

Figure 8 :

 8 Figure 8: Proportion of 50 agents in the first shelter for chosen values of σ, using three different models: MF, MC and Hybrid. Panels A and C are obtained from simulations using 50 cockroaches (Panel A: P. americana ; Panel C: B. germanica). Panels B and D are obtained from simulations using 45 cockroaches and 5 robots (Panel B: P. americana ; Panel D: B. germanica). Results for MF, MC and Hybrid models results are shown respectively in red, green, and blue. The bifurcation point is close to σ = 0.8 for P. americana, and σ = 1.0 for B. germanica. The other σ parameter values chosen are before the bifurcation point (σ = 0.4 for P. americana , σ = 0.6 for B. germanica), and just after the bifurcation point (σ = 1.2 for P. americana , σ = 1.8 and σ = 2.0 for B. germanica). The best sets of optimized model parameters are used, after 100 runs of optimization. The diagram is symmetric for all tested values of σ, so only one shelter is represented. Calibrated versions of the MC and Hybrid models behave similarly to the MF model: (1) Before the bifurcation point (0.4 ≤ σ < 0.8 for P. americana, 0.4 ≤ σ < 1.0 for B. germanica), only one configuration exists, an equipartition of the individuals between the two shelters (x 1 /N = x 2 /N = 1/2, x e = 0); (2) After the bifurcation point (σ > 0.8 for P. americana, σ > 1.0 for P. americana), two stable configurations exist, with all individuals in only one of the shelters (either x 1 ≈ 0, x 2 ≈ 1, x e ≈ 0 or x 1 ≈ 1, x 2 ≈ 0, x e ≈ 0).

Figure 9 :

 9 Figure 9: Examples of modulation of the biohybrid group behaviour when robots are optimized to change the behaviour of cockroaches (A: P. americana, B: B. germanica). (red: MF model, green: MC model, blue: Hybrid model). Results in red, green and blue are the final states corresponding to the change of steady states induced by the robots. Results in dark red, dark green, and dark blue correspond to the reference results, from experiments with only insects and no robots (from Fig. 8). Values of σ are chosen around the bifurcation point (P. americana: σ = 0.8, B. germanica: σ = 1.0), and just after the bifurcation point (P. americana: σ = 1.2, B. germanica: σ = 1.8). Results before the bifurcation points are not shown..

Table 2 :

 2 Parameters in the MF, MC, and Hybrid models. The parameter values for the MF model are from

Acknowledgement

The experiments presented in this paper were carried out using the Grid'5000 testbed, supported by a scientific interest group hosted by Inria and including CNRS, RENATER and several universities as well as other organizations (see https://www.grid5000.fr).