Leader Selection in Stochastically Forced Consensus Networks

noise_corrupted_k40 

Fu Lin, Makan Fardad, and Mihailo R. Jovanovic

June 2013

Matlab Files
Download all Matlab files

Presentation
2013 Allerton Talk

Paper
2014 IEEE TAC Paper

Purpose

This website provides a Matlab implementation of algorithms that quantify suboptimality of solutions to two combinatorial network optimization problems. For undirected consensus networks, these problems arise in assignment of a pre-specified number of nodes, as leaders, in order to minimize the mean-square deviation from consensus. For distributed localization in sensor networks, one of our formulations is equivalent to the problem of choosing a pre-specified number of absolute position measurements among available sensors in order to minimize the variance of the estimation error.

Two combinatorial problems involving graph Laplacian

Diagonally strengthened graph Laplacian

Let L = L^T in {bf R}^{n times n} be the graph Laplacian of an undirected connected network. A diagonally strengthened graph Laplacian is obtained by adding a vector kappa in {bf R}^n with positive elements to the main diagonal of L,

 	bar{L} ;=; L ;+; {rm diag} , (kappa).


For connected graphs, the diagonally strengthened Laplacian bar{L} is a positive definite matrix.

We are interested in selecting a pre-specified number N_l of diagonal elements of L such that the trace of the inverse of the strengthened graph Laplacian is minimized. This problem can be expressed as

         begin{array}{lrcl}         {rm minimize}         &         J(x)         & = &         {rm trace}         left(         (L ,+, {rm diag} , (kappa) , {rm diag} , (x) )^{-1}         right)         [0.25cm]         {rm subject~to}         &         x_i         & in &         {0,1},         ~~~~~         i ; = ; 1,ldots,n         [0.15cm]         &                 displaystyle{sum_{i , = , 1}^n} , x_i         & = &          N_l         end{array}         ~~~~~~~         ({rm LS1})

where the graph Laplacian L and the vector kappa with positive elements are problem data, and the Boolean-valued vector x with cardinality N_l is the optimization variable.

This combinatorial optimization problem is of interest in several applications. For example, in leader-follower networks subject to stochastic disturbances, the problem ({rm LS1}) amounts to assigning N_l nodes as leaders (that have access to their own states) such that the mean-square deviation from consensus is minimized. We refer to ({rm LS1}) as the noise-corrupted leader selection problem. Also, it can be shown that the problem of choosing N_l absolute position measurements among n sensors in order to minimize the variance of the estimation error in sensor networks is equivalent to the problem ({rm LS1}).

Principal submatrices of graph Laplacian

Let L_f be a principal submatrix of the graph Laplacian L obtained by deleting a pre-specified number N_l of rows and columns. For undirected connected graphs, L_f is a positive definite matrix.

The problem of finding a principal submatrix L_f of L such that {rm trace} , (L_f^{-1}) is minimized can be expressed as

         begin{array}{lrcl}         {rm minimize}         &         J_f(x)         & = &         {rm trace}         ,         (L_f^{-1})         [0.25cm]         {rm subject~to}         &         x_i         & in &         {0,1},         ~~~~~         i ; = ; 1,ldots,n         [0.15cm]         &                 displaystyle{sum_{i , = , 1}^n} , x_i         & = &          N_l         end{array}         ~~~~~~~         ({rm LS2})

where the Boolean-valued vector x with cardinality N_l is the optimization variable. The index set of nonzero elements of x identifies the rows and columns of the graph Laplacian L that need to be deleted in order to obtain L_f.

Similar to {(rm LS1)}, the combinatorial optimization problem ({rm LS2}) arises in several applications. For example, in leader-follower networks, the problem ({rm LS2}) amounts to assigning N_l nodes as leaders (that perfectly follow their desired trajectories) such that the mean-square deviation from consensus is minimized. We refer to ({rm LS2}) as the noise-free leader selection problem.

Performance bounds on the global optimal value

For large graphs, it is challenging to determine global solutions to the combinatorial optimization problems ({rm LS1}) and ({rm LS2}). Instead, we develop efficient algorithms that quantify performance bounds on the global optimal values. Specifically, we obtain lower bounds by solving convex relaxations of ({rm LS1}) and ({rm LS2}) and we obtain upper bounds using a simple but efficient greedy algorithm. We note that the greedy approach also provides feasible points of ({rm LS1}) and ({rm LS2}).

Lower bounds from convex relaxations

Since the objective function J in ({rm LS1}) is the composition of a convex function {rm trace} , (bar{L}^{-1}) of a positive definite matrix bar{L} with an affine function bar{L}(x) ,=, L ,+, {rm diag} , (kappa) , {rm diag} , (x), it follows that J is a convex function of x. By enlarging the Boolean constraint set x_{i} in {0,1} to its convex hull x_{i} in  [0,1], we obtain the following convex relaxation of ({rm LS1}):

         begin{array}{lrcl}         {rm minimize}         &         J(x)         & = &         {rm trace}         left(         (L ,+, {rm diag} , (kappa) , {rm diag} , (x) )^{-1}         right)         [0.25cm]         {rm subject~to}         &         x_i         & in &         [0,1],         ~~~~~         i ; = ; 1,ldots,n         [0.15cm]         &                 displaystyle{sum_{i , = , 1}^n} , x_i         & = &          N_l         end{array}         ~~~~~~~         ({rm CR1})

Schur complement can be used to cast ({rm CR1}) as a semidefinite program (SDP). Thus, ({rm CR1}) can be solved using general-purpose SDP solvers with computational complexity of order n^4. In contrast, computational complexity of a customized interior point method that we develop is of order n^3.

In contrast to ({rm LS1}), the objective function J_f in ({rm LS2}) is a nonconvex function of x. By a change of optimization variables, we identify the source of nonconvexity in the form of a rank constraint. Removal of the rank constraint and relaxation of the Boolean constraints can be used to obtain the following convex relaxation of ({rm LS2}):

         begin{array}{lrcl}         {rm minimize}         &         J_f(Y,y)         & = &         {rm trace}         left(         ( L circ Y          ,+,          {rm diag} left( {bf 1} - y right) )^{-1}         right)         ,-,          N_l         [0.25cm]         {rm subject~to}         &         y_i         & in &         [0,1],         ~~~~~         i ; = ; 1,ldots,n         [0.15cm]         &         Y_{ij}         & in &         [0,1],         ~~~~~         i,j ; = ; 1,ldots,n         [0.15cm]         &         displaystyle{sum_{i,=,1}^n}          ~ y_i & = & N_f,         ~~~~~         displaystyle{sum_{i,j,=,1}^n}          ~Y_{ij}          ~ = ~          N_f^2         [0.25cm]         &         Y & succeq &         0         end{array}         ~~~~~~~         ({rm CR2})

Here, the symmetric matrix Y in {bf R}^{n times n} and the vector y in {bf R}^n are the optimization variables, circ is the elementwise multiplication of two matrices, {bf 1} in {bf R}^n is the vector of all ones, and the integer N_f is given by N_f = n - N_l.

Similar to ({rm CR1}), the convex relaxation ({rm CR2}) can be cast as an SDP and solved using general-purpose SDP solvers. However, since this approach requires the computational complexity of order n^6, we develop an alternative approach that is well-suited for large problems. Our approach relies on the use of the alternating direction method of multipliers (ADMM) which allows us to decompose ({rm CR2}) into a sequence of subproblems that can be solved with O(n^3) operations.

Upper bounds using a greedy algorithm

With lower bounds on the global optimal value resulting from the convex relaxations ({rm CR1}) and ({rm CR2}), we use a greedy algorithm to compute upper bounds. This algorithm selects one leader at a time by assigning the node that provides the largest performance improvement as the leader. Once this is done, an attempt to improve a selection of N_l leaders is made by checking possible swaps between the leaders and the followers. In both steps, substantial improvement in algorithmic complexity is achieved by exploiting structure of the low-rank modifications to Laplacian matrices.

Detailed description of the aforementioned algorithms (i.e., the greedy algorithm, the customized interior point method for ({rm CR1}), and the ADMM-based algorithm for ({rm CR2})) are provided in our paper

Description of the main function

Matlab Syntax
>> [Jlow,Jup,x] = leaders(L,Nl,kappa,flag);

Our Matlab function leaders.m takes the graph Laplacian L in {bf R}^{n  times n}, the positive integer N_l < n, the vector kappa in {bf R}^n with positive elements, and the indicator variable

  • {rm flag} ; = ; 1  –  for the noise-corrupted leader selection formulation ({rm LS1});

  • {rm flag} ; = ; 2  –  for the noise-free leader selection formulation ({rm LS2}).

leaders.m returns the lower and upper bounds on the global optimal value of the leader selection problem ({rm LS1}) or ({rm LS2}), along with the Boolean-valued vector x which represents a feasible point for the selection of leaders.

Acknowledgements

This project is supported in part by the National Science Foundation under

  • CAREER Award CMMI-0644793

  • Award CMMI-0927509

  • Award CMMI-0927720.