Decentralised Coordination of Mobile Sensors Using the Max-Sum Algorithm

1 Introduction

In disaster response, and many other applications besides, the availability of timely and accurate information is of vital importance. Thus, the use of multiple mobile sensors for information gathering in crisis situations has generated considerable interest. 1 These mobile sensors could be autonomous ground robots or unmanned aerial vehicles. In either case, while patrolling through the disaster area, these sensors need to keep track of the continuously changing state of spatial phenomena, such as temperature or the concentration of potentially toxic chemicals. The key challenges in so doing are twofold. First, the sensors cannot cover the entire environment at all times, so the spatial and temporal dynamics of the monitored phenomena need to be identified in order to predict environmental conditions in parts of the environment that can not be sensed directly. Second, the sensors need to coordinate their movements to collect the most informative measurements needed to predict these environmental conditions as accurately as possible.

1 For example, one of the missions of both the Aladdin project (http://www.aladdinproject.org) and the Centre for Robot Assisted Search & Rescue (http://crasar.csee. usf.edu) is to use autonomous robots for information gathering in disaster response scenarios.

Recent work has addressed similar challenges by modelling the spatial and temporal dynamics of the phenomena using Gaussian processes (GPs) [Rasmussen and Williams, 2006]. GPs are a powerful Bayesian approach for inference about functions, and have been shown to be an effective tool for capturing the dynamics of spatial phenomena [Cressie, 1993]. This principled approach to modelling the environment has been used to compute informative deployments of fixed sensors [Guestrin et al., 2005], and informative paths for single [Meliou et al., 2007] and multiple mobile sensors [Singh et al., 2007].

However, the algorithms used to compute these informative deployments and paths are not suitable in our domain, since they are geared towards solving a one-shot optimisation problem in an off-line phase. Moreover, these algorithms are centralised. In hostile environments, this is undesirable, because it creates a single point of failure, thereby increasing the vulnerability of the information stream. Other work has employed on-line decentralised path planning using artificial potential fields to keep sensors in specific favourable formations [Fiorelli et al., 2006], or through multi-agent negotiation techniques to partition the environment and allocate the sensors to these partitions [Ahmadi and Stone, 2006]. However, in general, this work has used representations of the environment, that are less sophisticated than Gaussian processes, and are thus, less applicable for modelling complex spatial and temporal correlations.

To address this shortcoming, Low et al. [2008] combine these approaches and use Gaussian processes to represent the environment, and use Markov decision processes to compute non-myopic paths for multiple mobile sensors in an on-line fashion. Whilst such a non-myopic approach avoids the problem of local minima, it incurs significant computational cost (it is only empirically evaluated for systems containing just two sensors), and is again a centralised solution.

Thus, against this background, in this paper, we present a new on-line, decentralised coordination algorithm for teams of mobile sensors. This algorithm computes coordinated paths with an adjustable look-ahead, thus allowing the trade off between computation and solution quality, and uses GPs to represent generic temporal and spatial correlations of the phenomena. To this end, we represent each sensor as an autonomous agent. These agents are capable of taking measurements, coordinating their actions with their immediate neighbours, and predicting the state of the spatial phenomenon at unobserved locations. We then use the max-sum algorithm for decentralised coordination [Farinelli et al., 2008] to have the agents negotiate a joint plan by exchanging messages with their immediate neighbours. By applying max-sum in this manner, every agent controls its own movements using information it possesses locally, and the coordination mechanism is decentralised. We choose max-sum because it has been shown to generate good solutions to decentralised coordination problems, while limiting computation and communication. However, a standard application of max-sum is still too computationally costly within our particular domain. Thus, we introduce two novel and generic pruning techniques that speed up the max-sum algorithm, and hence, make decentralised coordination tractable in our application.

In more detail, in this paper we contribute to the state of the art in the following ways:

• We cast the multi-sensor monitoring problem as a decentralised constraint optimisation problem (DCOP), and present a new decentralised on-line coordination mechanism based on the max-sum algorithm to solve it.
• We present two novel, generic pruning techniques specifically geared towards reducing the number of function evaluations that is performed by max-sum. Thus alleviating a major bottleneck of this algorithm.
• We empirically show that a specific instantiation of our approach prunes 92% of joint moves for 5 sensors, and outperforms a greedy single step look-ahead algorithm by up to 50% in terms of root mean squared error.

The remainder of this paper is structured as follows. In section 2 we give a formal problem description. Section 3 describes how spatial phenomena are modelled. In section 4, we present our distributed algorithm, which we empirically evaluate in section 5, before concluding in section 6.

2 ProblemDescription

The problem formulation described in this section was inspired by [Meliou et al., 2007], and has been extended to deal with multiple sensors and limited local knowledge. Consider an environment in which M sensors monitor spatial phenomena that are modeled by a scalar field F : R3 → R, defined on one temporal and two spatial dimensions, at a finite set of locations V = {v1, v2 . . . } ⊂ R2 , and an indeterminate 2 number of of discrete time steps T = {t1, t2, ...}. To the measurement at location v ∈ V , and time t we associate a continuous random variable, Xv,t. The set of all random variables is denoted by X. The layout of the physical environment is given by a graph G = (V, E), where E encodes the possible movements between locations V . The locations accessible from v are denoted by adjG(v). Since it is generally not possible to visit all locations V during a single time step, each sensor selects an adjacent location at which to take a measurement at time step t + 1. Values at locations V are subsequently predicted with a statistical model using all measurements that the sensors have gathered so far. In order to do this, we model the scalar field F with a GP (see next section) that encodes both its spatial and temporal correlations.

Now, in order to select their movements, sensors need to be able to predict the informativeness of the samples that are collected along their paths with respect to the missing ones. Here, the informativeness of a set of samples O is quantified by a function f (O), that, depending on the context, can take on different forms [Meliou et al., 2007]. Our choice for this function f will be derived in the next section.

2 In uncertain and dynamic scenarios, the mission time is often not known beforehand.

Given this formalisation, we define the multi-sensor monitoring problem as follows. For every time step t, maximise f (Ot), where Ot = ∪Mi=1Oi t is the set of samples collected by all M sensors up to time step t at which the prediction is made. Moreover, while doing so, sensors can only communicate with and be aware of their immediate neighbours, such that no single point of control exists.

This problem is very challenging even for a single sensor. We therefore propose a distributed algorithm that computes paths with an adjustable look-ahead in Section 4, but first we discuss the way in which the spatial phenomena are modelled and derive function f .

3 Modeling the Spatial Phenomena

In order to predict measurements at unobserved locations, we model the scalar field F with a GP. Using a GP, F can be estimated at any location and at any point in time based on a set of samples collected by the sensors [Rasmussen and Williams, 2006]. In more detail, a single sample o of the scalar field F is a tuple 〈x, y〉, where x = (v, t) denotes the location and time at which the sample was taken, and y the measured value. Now, if we collect the training inputs x in a matrix X, and the outputs y in a vector y, the predictive distribution of the measurement at spatio-temporal coordinates x∗, conditioned on previously collected samples Ot = 〈X, y〉 is Gaussian with mean μ and variance σ2 given by:

where K(X,X′ ) denotes the matrix of covariances for all pairs of rows in X and X′ . These covariances are obtained by evaluating a function k(x, x′ ), called a covariance function, which encodes the spatial and temporal correlations of the pair (x, x′ ). Generally, covariance is a non-increasing function of the distance in space and time. For example, a prototypical choice is the squared exponential function where the covariance decreases exponentially with this distance:

where σf and l are called hyperparameters that model the signal variance and the length-scale of the phenomenon respectively. The latter determines how quickly the phenomenon varies over time and space 3 . If these hyperparameters are unknown before deployment of the sensors, they can be efficiently learnt on-line from collected samples using Bayesian Monte Carlo [Osborne et al., 2008].

One of the features of the GP is that the posterior variance in Equation 2 is independent of actual measurements y. This allows the sensors to determine the variance reduction that results from collecting samples along a certain path without the need of actually collecting them. Using this feature, we define the value f (O) to be the reduction in entropy that results at the coordinates of O after taking these samples 4 . This function exhibits the property of locality [Guestrin et al., 2005], that is exploited by our algorithm. This means that the correlation between two samples decreases rapidly (exponentially in the case of Equation 3) with increasing distance, such that samples that are far apart can be considered uncorrelated, and thus, mobile sensors that are far apart need not explicitly coordinate.

3 A slightly modified version of Equation 3 allows for different length-scales for the spatial and temporal dimensions of the process.

4 Decentralised Coordination

Since, in general, it is too expensive to perform on-line nonmyopic path planning, especially because the number of time steps in T (the mission time) is unknown beforehand, we present an algorithm that computes joint moves with a finite, adjustable look-ahead. An example of such a joint move for three sensors, consisting of a path of length 5 for each of them, is shown in Figure 1. The sensors run the algorithm every n time steps 5 to plan joint paths of length l ≥ n. The algorithm is flexible, in that it allows paths of various lengths, and additional constraints, to be considered. For instance, the sensors can coordinate over all possible paths of length l, or only those moves that keep them within their own partition of the environment.

Given any particular potential joint move, our algorithm must calculate the reduction in entropy, f (∪Mi=1Oi t ), that will result. In order to do so, we apply the chain-rule of entropies 6 to convert f (∪Mi=1Oi t ) to ∑Mi=1 f (Oi t ). This modified problem can then be cast as a decentralised constraint optimisation problem (DCOP), which lends itself to an agent-based solution paradigm. In this paradigm, sensors are modelled as agents that jointly plan their paths in order to maximise the value of the collected samples.

More formally, we denote the decision variable of agent i as pi, which takes values in the set Ai = {a1i , . . . , aqi i }, representing all qi moves that agent i is currently considering. The value of the samples that agent i collects along its path, given the movement of other sensors, is denoted by Ui(pi), where pi is the vector of decision variables on which agent i's utility depends (thus, pi ∈ pi will always hold). By the property of locality of f, this dependency relation is usually limited to a small set of nearby agents. Now, the agents will

4 Note that Guestrin et al. [2005] show that mutual information can also be used to define this value. However, the gain in doing so is small, and it is much more computationally expensive to evaluate.

5 For simplicity, and w.l.o.g. we assume that the speed of the sensors is 1 unit per time step.

collectively attempt to find joint move p∗ = [p∗1 , . . . , p∗M ], such that:

or, in other words, the joint move that maximises the total value obtained by the agents.

4.1 The Max-SumMessage-Passing Algorithm

The coordination problem encoded by Equation 4 is a DCOP, which can be solved by a wide range of algorithms. Unfortunately, many of these algorithms either compute the optimal solution at exponential cost, either in terms of the number or size of messages that are exchanged between agents (e.g. DPOP [Petcu and Faltings, 2005]), or require little local computation and communication, but produce approximate solutions (e.g. the Distributed Stochastic Algorithm [Fitzpatrick and Meertens, 2003]). However, there exists a class of algorithms usually referred to under the framework of the Generalised Distributive Law [Aji and McEliece, 2000], that can be used to obtain good approximate solutions. The maxsum message passing algorithm is one member of this class that is of particular interest here. This algorithm has been shown to compute better quality solutions than the approximate class with acceptable computation compared to representative complete algorithms [Farinelli et al., 2008].

In more detail, the max-sum algorithm operates on a factor graph: an undirected bipartite graph in which vertices represent variables pi and functions Uj. In such factor graphs, an edge exists between a variable pi and a function Uj iff pi ∈ pj, (i.e., pi is a parameter of Uj). Using the maxsum algorithm we exploit the fact that an agent's utility depends only on a subset of other agents' decision variables (locality), and that the global utility function is a sum of each agent's utility. Figure 2 shows an example factor graph that encodes Equation 4 for the coordination problem of Figure 1. In this example, the utility of agent 1 depends on its own action, and that of agent 2, so p1 = {p1, p2}. Similarly, p2 = {p1, p2, p3}, and p3 = {p2, p3}.

In yet more detail, using max-sum, each agent computes:

in a distributed way (i.e. based on local information and communication with direct neighbours). Agent i's optimal move p∗i is then obtained as follows:

Algorithm 1 Algorithm for computing pruning message from function Um to variable pn

• 1: compute Um(pn) ≤min p-n Um(pn, p-n)
• 2: compute Um(pn) ≥max p-n Um(pn, p-n)
• 3: send 〈Um(pn),Um(pn)〉 to pn

In order to do this, messages are passed between the functions Ui, and the variables in pi as described below 7 :

· From variable to function:

where αnm is a normalising constant that is chosen such that ∑pn Qn→m(pn) = 0, to prevent the messages from growing arbitrarily large.

· From function to variable:

Finally, ˜Ui(pi) is obtained by summing the most recent messages RUm′→pi (pi) received by pi. This computation is guaranteed to be exact when the factor graph is acyclic. However, since dependencies are usually mutual (i.e. agent i's action influences agent j's utility and vice versa), the factor graph will normally contain cycles. In this case, max-sum computes approximate solutions (i.e ˜Ui(pi) ≈ maxp-i ∑Mi=1 Ui(pi)). Despite this, however, there exists strong empirical evidence that max-sum produces good results even in cyclic factor graphs [Farinelli et al., 2008].

4.2 Speeding upMessage Computation

The straightforward application of max-sum to solve Equation 4 is not practical, because the computation of the messages from function to variable (Equation 8) is a major bottleneck. A naïve way of computing these messages for a given variable pn is to enumerate all joint moves (i.e. the domain of pm), and evaluate Um for each of these moves. Since the size of this joint action space grows exponentially with both the number of agents, and the number of available moves for each agent, the amount of computation quickly becomes prohibitive. This is especially true when evaluating Um is costly, as is the case in the mobile sensors domain 8 . Therefore, we introduce two pruning algorithms to reduce the size of the joint action space that needs to be considered. These algorithms are then applied to our mobile sensor domain, but they can just as easily be applied within other settings.

7 In what follows, we use adj(Um) to denote adjacent vertices of function Um in the factor graph, i.e., the set of variables in the domain pm of Um. Similarly, adj(pn) denotes the set of functions in which pn occurs in the domain.

8 Specifically, determining the value of a sample involves the inversion of a potentially very large matrix K(X, X) (Equation 2).

Algorithm 2 Algorithm for computing pruning messages from variable pn to all functions Um ∈ adj (pn)

• 1: if a new message has been received from all Um ∈ adj(pn) then 2: compute ⊥(pn) = PUm∈adj(pn) Um(pn) 3: compute ⊤(pn) = PUm∈adj(pn) Um(pn) 4: while ∃a ∈ An : ⊤(a) <max⊥(pn) do 5: An ←An \ {a} 6: end while 7: send updated domain An to each Um ∈ adj(pn) 8: end if

The Action Pruning Algorithm

The first algorithm attempts to reduce the number of moves each agent needs to consider before running the max-sum algorithm. This algorithm prunes the dominated states that can never maximise Equation 4, regardless of the actions of other agents. More formally, a move a′ ∈ An is dominated if there exists a move a∗ such that:

Just as with the max-sum algorithm itself, this algorithm is implemented by message passing, and operates directly on the variable and function nodes of the factor graph, making it fully decentralised:

• From function to variable: Function Um sends a message to pn, containing the minimum and maximum values of Um with respect to pn = an, for all an ∈ An. (see Algorithm 1).
• From variable to function: Variable pn sums the minimum and maximum values from each of its adjacent functions, and prunes dominated states. It then informs neighbouring functions of its updated domain (see Algorithm 2).

Using this distributed algorithm, functions continually refine the bounds on the utility for a given state of a variable, which potentially causes more states to be pruned. Therefore, it is possible that action pruning starts with a single move, and subsequently propagates through the entire factor graph.

Now, given the highly non-linear relations expressed in Equation 2, on which the agents' utility functions Um are based, it is very difficult to calculate these bounds exactly, without exhaustively searching the domain of xm for utility function Um. Needless to say, this would defeat the purpose of this pruning technique. Nonetheless, experimentation showed that by computing these bounds in a greedy fashion, a very good approximation is obtained. Thus, the lower bound Um(an) on a move an is obtained by selecting the neighbouring agents one at a time, and finding the move that reduces the utility of agent m's move the most. In a similar vein, the upper bound Um(an) is obtained selecting those moves of other sensors that reduce the utility the least.

The Joint Action Pruning Algorithm

Whereas the first algorithm runs as a preprocessing phase to max-sum, the second algorithm is geared towards speeding up the computation of the messages from function to variable (Equation 8), while max-sum is running. A naïve way of computing this message to a single variable pi is to determine the maximum utility for each of agent i's actions by exhaustively enumerating the joint domain of the variables in pm\{pi} (i.e. the Cartesian product of the domains of these variables), and evaluating the expression between brackets in Equation 8, which we denote by ˜RUm→pn (pm). This expression is the sum of the utility function Um and the sum of messages Q.

However, instead of just considering joint moves, we now allow some actions to be undetermined, and thus, consider partial joint moves, denoted by ˆa. By doing so, we can create a search tree on which we can employ branch and bound to significantly reduce the size of the domain that needs to be searched. In more detail, to compute RUm→pn (ain ) (a single element of the message from Um to variable pn) for a single state ain ∈ Ai in the domain of pn, we create a search tree T (ain ) as follows:

• The root of T (ain ) is a partial joint move ˆar = 〈∅, . . . , ∅, ain , ∅, . . . , ∅〉, which indicates that ain is assigned to pn, and the remaining variables are unassigned (denoted by ∅).
• The children of a vertex 〈a(1) 1 , . . . , a(k) k , ∅, . . .∅, ain ,∅, . . . ,∅〉 are obtained by setting the first unassigned variable pk+1 to each of its |Ak+1| moves.
• The leafs of the tree represent a (fully determined) joint move am (i.e. ∀i : pi ̸= ∅). In the tree, only leafs are assigned a value, which is equal to ˜RUm→pn (am).

The maximum value found in T (ain ) is the desired value. Now, in order to use branch and bound to find this value, we need to put bounds on the maximum value found in a subtree of T (ain ). These bounds depend on Um and the received messages Q. Now, in many cases we can put bounds on the maximum of the former, that is obtained by further completing a partial joint move ˆa′ in a subtree of T (ain ). The bounds on Um, combined with the minimum and maximum values of Q for ˆa′ (again, by further completing the partial joint move), gives us the desired bounds.

Figure 3 shows an example of a partially expanded search tree for computing a single element RUm→x3 (a13 ) of a message from function Um to variable p3. Given the lower and upper bounds on the maximum, subtree 〈a11 ,∅, a13 〉 can be pruned immediately after expanding the root. Similarly, subtree 〈a31 ,∅, a13 〉 is pruned after expanding leaf 〈a21 , a22 , a13 〉, which has the desired maximum value.

To compute these bounds on the maximum of Ui(ˆa) in the mobile sensor domain, note that partial joint move ˆa represents a situation in which only a subset of the agents have determined their move. Using this interpretation, we can obtain bounds as follows. The upper bound on this value is obtained by disregarding the agents that have not determined their action (i.e. agents i for which pi = ∅). Since the act of collecting a sample always reduces the value of other samples, disregarding the samples of these 'undecided' agents will give an upper bound on the maximum. To obtain a lower bound on the maximum, we use the locality property of f, which tells us that the interdependency between values of samples weakens as their distance increases. So, in order to calculate the lower bound, we move the undecided agents away from agent i's destination.

4.3 Ensuring Network Connectivity

In many situations, it is also important that the sensors maintain network connectivity in order to transmit their measurements to a base station. More importantly in the context of our algorithm, agents need to be able to communicate in order to negotiate over their actions. Not surprisingly, we can use the algorithm to accomplish this, by penalising disconnection from the network in the utility function Ui. To this end, we assume that every agent maintains a routing table that specifies which agents can be reached through each immediate neighbour. Thus, a move is only allowed if all agents will still be reachable through the remaining links. Otherwise, the agent risks disconnection from the network, in which case a large penalty is added to its utility function Ui.

5 Empirical Evaluation

To empirically evaluate our approach, we simulated five sensors on a lattice graph measuring 26 by 26 vertices. The data was generated by a GP with a squared exponential covariance function (see Equation 3) with a spatial length-scale of 10 and a temporal length-scale of 150. This means that the spatial phenomenon has a strong correlation along the temporal dimension, and therefore changes slowly over time. At every m time steps, the sensors plan their motion for the next l time steps (l ≥ m). In what follows, this strategy is referred to as MSm-l. Now, instead of considering all possible paths of length l from an agent's current position, which would result in a very high computational overhead, the action space is limited to the locations in G that can be reached in l time steps in 8 different directions, corresponding to the major directions on the compass rose. In the first experiment, we benchmarked MS1-1 and MS1-5 against four strategies:

• Random: Randomly moving sensors.
• Greedy: Sensors that greedily maximise the value of the sample collected in the next move without coordination.
• J(umping) Greedy: The same as Greedy, except that these sensors can instantaneously jump to any location.
• Fixed: Fixed sensors that are placed using the algorithm proposed in [Guestrin et al., 2005].

The averaged root mean squared error (RMSE) for 100 time steps is plotted in Figure 4(a). From this figure, it is clear that both MS strategies outperform the Greedy, Random, and Fixed strategies. Furthermore, the prediction accuracy of MS1-5 is comparable to that of JGreedy, whose movement is not restricted by graph G. Moreover, it shows that increasing the length of the considered paths from 1 to 5, reduces the RMSE by approximately 30%.

In the second set of experiments, we analysed the speed-up achieved by applying the two pruning techniques described in Section 4.2. Figure 4(b) shows the percentage of joint actions pruned plotted against the number of neighbouring agents. With 5 neighbours, the two pruning techniques combined prune around 92% of the joint moves. With such a number of neighbouring agents, the agents are strongly clustered, which occurs rarely in a large environment. However, should this happen, the utility function needs to be evaluated for only 8% of roughly 85 joint actions, thus greatly improving the algorithm's efficiency.

In the third experiment, we performed a cost/benefit analysis of various MSm-l strategies. More specifically, we examined the effect of varying m and l on both the number of utility function evaluations, and the resulting RMSE. Figure 4(c) shows the results. The results of MS1-1, MS2-2, MS4-4, MS5-5, and MS8-8 show an interesting pattern. Up to and including m = l = 4, both the number of function evaluations and the average RMSE decrease. This is due to the fact that planning longer paths is more expensive, but results in lower RMSE. However, for m, l > 4, the action space becomes too coarse (since only 8 directions are considered) to maintain a low RMSE. At the same time, the number of times the agents coordinate reduces significantly, resulting in a lower number of function evaluations. Finally, MS1-5 and MS4-8 provide a compromise; they compute longer paths, but coordinate more frequently. This leads to more computation compared to MS5-5 and MS8-8, but results in significantly lower RMSE, because agents are able to 'reconsider' their paths.

A video demonstrating the mobile sensors with the techniques from Section 4 can be found at www.youtube. com/ijcai09.

6 Conclusions

In this paper, we presented an on-line decentralised coordination algorithms for multiple sensors. We showed how the max-sum message passing algorithm can be applied to this domain in order to coordinate the motion paths of the sensors along which the most informative samples are gathered. We also presented two general pruning to speed up the max-sum algorithm. The first attempts to prune actions of the sensors that are not part of the optimal joint move. The second uses branch and bound to reduce the fraction of the joint action space that needs to be searched in order to compute the messages from functions to variables, which is the main bottleneck of the max-sum algorithm. We empirically showed that for 5 sensors, these techniques prune 92% of joint moves, thus significantly reducing the number of utility function evaluations, which are particularly expensive in the mobile sensor domain. Moreover, we showed that root mean squared error with which the spatial phenomenon is predicted by the MS15 strategy is approximately 50% lower compared to a greedy single step look-ahead algorithm. Our future work in this area is to extend the proposed approach by adopting techniques from sequential decision making to do non-myopic path planning, while keeping computational costs in check.

Acknowledgments This research was undertaken as part of the ALADDIN (Autonomous Learning Agents for Decentralised Data and Information Systems) Project and is jointly funded by a BAE Systems and EPSRC (Engineering and Physical Research Council) strategic partnership (EP/C 548051/1).