Sensor Selection Cost Optimization for Tracking Structurally Cyclic Systems: a P-Order Solution

TL;DR

This paper presents a polynomial-time graph-theoretic method for optimizing sensor measurement selection to minimize costs while ensuring observability in large-scale structurally cyclic systems.

Contribution

It introduces a novel polynomial-time solution for sensor selection cost optimization in structurally cyclic systems using graph theory.

Findings

The proposed method guarantees minimal sensing cost under observability constraints.

The solution operates with a computational complexity of O(m^3).

Applicable to large-scale systems due to polynomial-time algorithm.

Abstract

Measurements and sensing implementations impose certain cost in sensor networks. The sensor selection cost optimization is the problem of minimizing the sensing cost of monitoring a physical (or cyber- physical) system. Consider a given set of sensors tracking states of a dynamical system for estimation purposes. For each sensor assume different costs to measure different (realizable) states. The idea is to assign sensors to measure states such that the global cost is minimized. The number and selection of sensor measurements need to ensure the observability to track the dynamic state of the system with bounded estimation error. The main question we address is how to select the state measurements to minimize the cost while satisfying the observability conditions. Relaxing the observability condition for structurally cyclic systems, the main contribution is to propose a graph theoretic approach to solve the problem in polynomial time. Note that, polynomial time algorithms are suitable for large-scale systems as their running time is upper-bounded by a polynomial expression in the size of input for the algorithm. We frame the problem as a linear sum assignment with solution complexity of O(m3).

Figures

Tables

Table 1. Table 1: Table of Notation.

$A$	system matrix
$𝒜$	structured system matrix
$H$	measurement matrix
$\mathscr{H}$	structured measurement matrix
$\nu$	system noise
$\eta$	measurement noise
$k$	discrete time index
$x$	system state
$y$	measurement
$c$	state-sensor cost matrix
$𝒞$	SCC-sensor cost matrix
$𝒵$	assignment matrix
$\tilde{c}$	pseudo-cost
$n$	number of states
$m$	number of sensors/measurements
$𝒥$	Jacobian matrix
$𝒳$	set of state nodes
$𝒴$	set of sensor nodes

Equations

$$\dot{\underline{x}}=A\underline{x}+\underline{\nu }$$

$$\underline{x}(k+1)=A\underline{x}(k)+\underline{\nu }(k)$$

$$y_i=H_i\underline{x}+{\eta }_i$$

$$y_i(k)=H_i\underline{x}(k)+{\eta }_i(k)$$

$$\underset{\mathcal{H}}{\min }$$

$$(A,H)-observability,$$

$${\mathcal{H}}_{ij}\in \{0,1\}$$

$$\underset{\mathcal{H}}{\min }$$

$$(\mathcal{A},\mathcal{H})-observability,$$

$${\mathcal{H}}_{ij}\in \{0,1\}$$

$$\underset{\mathcal{H}}{\min }$$

$$(\mathcal{A},\mathcal{H})-observability,$$

$$\sum _{i=1}^m{\mathcal{H}}_{ij}\le 1$$

$$\sum _{j=1}^n{\mathcal{H}}_{ij}=1$$

$${\mathcal{H}}_{ij}\in \{0,1\}$$

$${\mathcal{C}}_{ij}=\min \{c_{il}\},x_l\in SC{C}_j,i,j\in \{1,\dots ,m\}$$

$$\underset{\mathcal{Z}}{\min }$$

$$\sum _{j=1}^m{\mathcal{Z}}_{ij}=1$$

$$\sum _{i=1}^m{\mathcal{Z}}_{ij}=1$$

$${\mathcal{Z}}_{ij}\in \{0,1\}$$

Full text

1 Introduction

Sensors and sensing devices are widespread in everyday use and are involved in many aspect of human life. Nowadays, sensors are advanced beyond the physical world and even are introduced in online social networks. The emerging notion of IoT and the so-called Trillion Sensors roadmap further motivates sensor and actuator implementation in many physical systems and cyber networks A few examples are: in ecosystems and environmental monitoring (May, 1972), security and vulnerability of social networks (Pequito et al., 2014; Doostmohammadian & Khan, 2014a), eHealth and epidemic monitoring (Nowzari et al., 2016), Dynamic Line Rating (DLR) in smart power grids (Ilic et al., 2010; Hug et al., 2015), etc. In these large-scale applications the cost of sensing is a challenge. The cost may represent energy consumption, the economic cost of sensors, and even the additive disturbance due to, for example, long distance communication in wireless sensor networks. The rapidly growing size of IoT and sensor networks motivates minimal cost sensor-placement solution for practical applications.

There exist different approaches toward sensor selection optimization. In (Joshi & Boyd, 2009), authors study sensor selection for noise reduction. This work introduces combinatorial problem of selecting $k$ out of $m$ sensors to optimize the volume of probabilistic confidence ellipsoid containing measurement error by adopting a convex relaxation. Authors in (Pequito et al., 2014), consider minimum sensor coverage for dynamic social inference. Their idea is to find minimum sensor collection to ensure generic social observability. The authors show the relaxation lies in set covering category and is generally NP-hard111Note that, the NP-hard problems are believed to have no solution in time complexity upper-bounded by a polynomial function of the input parameters. to solve. In (Zejnilovic et al., 2013) the source localization problem under observability constraints is addressed. The authors aim to find the minimal possible observers to exactly locate the source of infection/diffusion in a network. They state that this problem is NP-hard and propose approximations to solve the problem. In another line of research, distributed optimization is discussed in (Wang & Elia, 2010), where dynamic feedback algorithms robust to disturbance is proposed to minimize certain cost function over a sensor network. Optimal sensor coverage with application to facility allocation is studied in (Sayyaadi & Moarref, 2011). The authors propose a distributed deployment protocol as a local optimal solution in order to assign resources to group of mobile autonomous sensors under certain duty to capability constraints. Minimal actuator placement ensuring controllability is discussed in (Tzoumas et al., 2015). Authors provide P-order approximations to a generally NP-hard problem by considering control energy constraints. In (Doostmohammadian & Khan, 2014a), the minimal number of observers for distributed inference in social networks is discussed. Similarly, minimizing the number of actuators for structural controllability following specific rank constraints is addressed in (Commault & Dion, 2015).

This paper studies minimum sensor placement cost for tracking structurally cyclic dynamical systems (see Fig.1). In general, state measurements are costly and these costs may change for different sensor selection. This is due to various factors, e.g. sensor range and calibration, measurement accuracy, embedding/installation cost, and even environmental conditions. Therefore, any collection of state-sensor pairs may impose specific sensing costs. The main constraint, however, is that not all collection of sensor measurements provide an observable estimation. Observability determines if the system outputs convey sufficient information over time to infer the internal states of the dynamic system. With no observability, no stable estimation can be achieved and the tracking error covariance grows unbounded. In conventional sense, observability requires algebraic tests, e.g. the observability Gramian formulation (Bay, 1999) or the Popov-Belevitch-Hautus test (Hautus, 1969). These methods are computationally inefficient especially in large-scale systems. In contrast, this paper adopts a structural approach towards observability. The methodology is irrespective of numerical values of system parameters, while the structure is fixed but the system parameters vary in time (Dion et al., 2003). Indeed, this approach makes our solution practically feasible for Linear Structure-Invariant (LSI) systems. This may arise in linearization of nonlinear dynamics, where the linearized Jacobian is LSI while the values depend on the linearization point. As it is known the observability and controllability of Jacobian linearization is sufficient for observability and controllability of the original nonlinear dynamics (Liu et al., 2011; Slotine & Li, 1991).

Another structural property of the system is the system rank. In particular, for full rank systems the adjacency graph of the system (system digraph) is structurally cyclic, i.e. there exist a disjoint union of cycles covering all state nodes in system digraph. An example is dynamic system digraphs including self loops, which implies that for every system state the $1$st-derivative is a function of the same state (referred to as self-damped systems in (Chapman & Mesbahi, 2013)). An example arises in systems representing ecological interactions (May, 1972, 1973) among species, where intrinsic self-dampening ensures the eco-stability around equilibrium state. In (Zhao et al., 2015), authors extended the case to network of coupled $d$th-order intrinsic dynamics. In such case, the self-dynamic impose a cyclic subgraph in the large-scale structure of the system digraph. In addition, in distributed estimation and sensor network literature typically it is assumed that the system matrix is invertible (Battistelli et al., 2012) and therefore system is full-rank. In (Liu et al., 2013; Arcak & Sontag, 2006) authors study the stability and observability of cyclic interconnected systems resulted from biochemical reactions, while in (Nowzari et al., 2016) authors analyze the self-damped epidemic equations integrated in social and human networks. These examples motivates the study of structurally cyclic systems in this paper.

Contributions: Towards sensor selection cost optimization, this paper, first, considers LSI dynamic systems. In such systems the parameters vary in time while the system structure is unchanged as in Jacobian linearization of nonlinear dynamic systems. This is also the case, for example, in social systems with invariant social interactions, and in power systems with fixed system structure but time-varying parameters due to dynamic loading. Second, the observability constraint for cost optimization is framed as a selection problem from a necessary set of states. We relax the observability constraint by defining the equivalent set of states necessary for observability. This is a novel approach towards cost optimization for estimation purposes. Third, the optimization is characterized as a Linear Sum Assignment Problem (LSAP), where the solution is of polynomial complexity. In this direction, the NP-hard observability optimization problem, as reviewed in the beginning of this section, is relaxed to a P-order problem for the case of structurally cyclic systems. Note that, for general systems, this problem is NP-hard to solve (see (Pequito et al., 2014) for example). The relaxation in this paper is by introducing the concept of cost for states measured by sensors and the concept of structural observability. Further, the definition of the cost is introduced as a general mathematical concept with possible interpretation for variety of applications. Particularly, note that this P-order formulation is not ideal but practical, and this work finds application in monitoring large-scale systems such as social systems (Pequito et al., 2014; Doostmohammadian & Khan, 2014a), eco-systems (May, 1972), and even epidemic monitoring (Nowzari et al., 2016). To the best of our knowledge, no general P-order solution is proposed in the literature for this problem.

Assumptions: The following assumptions hold in this paper:

(i) The system is globally observable to the group of sensors.

(ii) Number of sensors is at least equal to the number of crucial states necessary for observability, and at least one state is accessible/measurable by each sensor.

Without the first assumption no estimation scheme works, and there is no solution for optimal observability problem. In the second assumption, any sensor with no access to a necessary state observation is not a player in the optimization game. Other assumptions are discussed in the body of the paper.

The outline of the paper is as follows. In Section 2, relevant structural system properties and algorithms are reviewed. Section 3 states the graph theoretic approach towards observability. Section 4 provides the novel formulation for the cost optimization problem. Section 5 reviews the so-called assignment problem as the solution. Section 6 states some remarks on the results, motivation, and application of this optimal sensor selection scheme. Section 7 illustrates the results by two academic examples. Finally, Section 8 concludes the paper.

Notation: We provided a table of notations in Table 1 to explain the terminologies and symbols in the paper.

2 Structured System Theory

Consider the state of linear system, $\underline{x}$, evolving as222The underline notation represents a vector variable.:

[TABLE]

and in discrete time as:

[TABLE]

where $\underline{x}\in\mathbb{R}^{n}$ is the vector of system states, and $\underline{\nu}\sim\mathcal{N}(0,V)$ is independent identically distributed (iid) system noise. Consider a group of sensors, indexed by $y_{i},~{}i=1,\ldots,m$ each taking a noise-corrupted state measurement as:

[TABLE]

or in discrete time as:

[TABLE]

where $H_{i}$ is a row vector, $y_{i}\in\mathbb{R}$ is the sensor measurement, and $\eta_{i}\sim\mathcal{N}(0,Q_{i})$ is the zero-mean measurement noise at sensor $i$.

Let $\mathcal{A}\sim\{0,1\}^{n\times n}$ represent the structured matrix, i.e. the zero-nonzero pattern of the matrix $A$. A nonzero element implies a system parameter that may change by time, and the zeros are the fixed zeros of the system. Similarly, $\mathcal{H}\sim\{0,1\}^{m\times n}$ represents the structure of measurement matrix $H$. A nonzero entry in each row of $\mathcal{H}$ represents the index of the measured state by the corresponding sensor. This zero-nonzero structure can be represented as a directed graph $\mathcal{G}_{sys}\sim(\mathcal{X}\cup\mathcal{Y},\mathcal{E})$ (known as system digraph). Here, $\mathcal{X}$ is the set of state nodes $\{{x}_{1},\ldots{x}_{n}\}$ each representing a state, and $\mathcal{Y}$ is the output set $\{{y}_{1},\ldots{y}_{m}\}$ representing the set of sensor measurements. The nonzero entry $\mathcal{A}_{ij}$ is modeled by an edge ${x}_{j}\rightarrow{x}_{i}$. The set $\mathcal{E}=(\mathcal{X}\times\mathcal{X})\cup(\mathcal{X}\times\mathcal{Y})$ is the edge set. Edges in $\mathcal{E}_{xx}=\mathcal{X}\times\mathcal{X}$ represent the dynamic interactions of states in $\mathcal{G}_{A}\sim(\mathcal{X},\mathcal{E}_{xx})$, and edges $\mathcal{E}_{xy}=(\mathcal{X}\times\mathcal{Y})$ in $\mathcal{G}_{xy}=(\mathcal{X}\cup\mathcal{Y},\mathcal{E}_{xy})$ represent the flow of state measurement information into sensors. It is clear that $\mathcal{G}_{sys}=\mathcal{G}_{A}\cup\mathcal{G}_{xy}$. Define a path as a chain of non-repeated edge-connected nodes and denote $\xrightarrow{path}\mathcal{Y}$ as a path ending in a sensor node in $\mathcal{Y}$. Define a cycle as a path starting and ending at the same node.

Similar to the structured matrices and the associated digraphs for linear systems, one can define a digraph for the Jacobian linearization of the nonlinear systems, also referred to as inference diagrams (Liu et al., 2013). In the nonlinear case the structure of the system digraph is related to the the zero-nonzero structure of the Jacobian matrix $\mathcal{J}$. For system of equations $\dot{\underline{x}}=f(\underline{x},\underline{\nu})$ if $\mathcal{J}_{ij}=\frac{\partial f_{i}}{\partial x_{j}}$ is not a fixed zero, draw a link ${x}_{j}\rightarrow{x}_{i}$ in the system digraph. This implies that ${x}_{i}$ is a function of ${x}_{j}$ and state ${x}_{j}$ can be inferred by measuring ${x}_{i}$ over time. Following this scenario for all pairs of states and connecting the inference links the system digraph $\mathcal{G}_{A}$ is constructed. It should be mentioned that we assume the nonlinear function $f$ is globally Lipschitz, and therefore the system of equation has a unique solution and the Jacobian matrix is defined at all operating points.

The properties of the system digraph and its zero-nonzero structure are closely tied with the generic system properties. Such properties are almost independent of values of the physical system parameters. It is known that if these specific properties of the system hold for a choice of numerical values of free parameters, they hold for almost all choices of system parameters, where these system parameters are enclosed in nonzero entries of the system matrix. Therefore, the zero-nonzero structure of the system and the associated system digraph ensures sufficient information on such generic properties. In general, efficient structural algorithms are known to check these properties while the numerical approach might be NP-hard to solve (Dion et al., 2003). An example of generic properties are structural rank ($\mathcal{S}$-rank) and structural observability and controllability (Dion et al., 2003; Doostmohammadian & Khan, 2014a).

2.1 Structurally Cyclic Systems

The following definition defines structurally cyclic systems:

Definition 1**.**

A system is structurally cyclic if and only if its associated system digraph includes disjoint family of cycles spanning all nodes (Van der Woude, 1999).

There exist many real-world systems which are structurally cyclic. As mentioned in Section 1, any complex network/system governed by coupled $d$th-order differential equations and randomly weighted system parameters is structurally cyclic (see (Zhao et al., 2015) for more information). Such structures may arise in biochemical reaction networks (Arcak & Sontag, 2006; Liu et al., 2013), epidemic spread in networks (Nowzari et al., 2016), ecosystems (May, 1972, 1973) and even in social networks where each agent has intrinsic self-dampening dynamics represented as self-loop in social digraph (see the example in (Pequito et al., 2014)).

There are efficient methods to check if the graph is cyclic and includes a disjoint cycle family, namely matching algorithms. A matching of size $m$, denoted by $\mathcal{M}_{m}$, is a subset of nonadjacent edges in $\mathcal{E}_{xx}$ spanning $m$ nodes in $\mathcal{X}$. Define nonadjacent directed edges as two edges not sharing an end node. Define a maximum matching as the matching of maximum size in $\mathcal{G}_{A}$, where the size of the matching, $m$, is defined by the number of nodes covered. A perfect matching is a matching covering all nodes in the graph, i.e. $\mathcal{M}_{n}$ where $n=|\mathcal{X}|$.

Lemma 1**.**

A system of $n$ state nodes is structurally cyclic if and only if its digraph includes a perfect matching $\mathcal{M}_{n}$.

Proof.

The detailed proof is given in (Murota, 2000). ∎

The size of the maximum matching is known to be related to the structural rank of the system matrix defined as follows:

Definition 2**.**

For a structured matrix $\mathcal{A}$, define its structural rank ($\mathcal{S}$-rank) as the maximum rank of the matrix $A$ for all values of non-zero parameters.

Note that, the $\mathcal{S}$-rank of $\mathcal{A}$ equals to the maximum size of disjoint cycle family, where the size represents the number of nodes covered by the cycle family (Harary, 1962). For structurally cyclic systems the maximum size of the cycle family equals to $n$, number of system states. This implies that for structurally cyclic system $\mathcal{S}\mbox{-rank}(\mathcal{A})=n$. This result can be extended to nonlinear systems as $\mathcal{S}\mbox{-rank}(\mathcal{J})=n$ at almost all system operating points.

As an example, consider a graph having a random-weighted self-cycle at every state node. In this example, every self-cycle is a matching edge, graph contains a perfect matching and therefore is cyclic. On the other hand, these self-loops imply that every diagonal element in the associated structured matrix, $\mathcal{A}$, is nonzero. Having random values at diagonal entries and other nonzero parameters the determinant is (almost) always nonzero and system is structurally full rank 333This can be checked simply by MATLAB, considering random entries as nonzero parameters of the matrix. The probability of having zero determinant is zero.. This is generally true for network of intrinsic $d$th-order self dynamics (instead of $1$st order self-loops) as addressed in (Zhao et al., 2015). In the coupled dynamic equations a $d$th-order dynamic represents a cyclic component in system digraph. The intra-connection of these individual dynamics construct the digraph of large-scale physical system. Assuming time-varying system parameters the system remains structurally full rank.

3 Graph Theoretic Observability

Observability plays a key role in estimation and filtering. Given a set of system measurements, observability quantifies the information inferred from these measurements to estimate the global state of the system. This is irrespective of the type of filtering and holds for any estimation process by a group of sensors/estimators. Despite the algebraic nature of this concept, this paper adopts a graph theoretic approach towards observability. This approach is referred to as structural observability and deals with system digraphs rather than the algebraic Gramian-based method. The main theorem on structural observability is recalled here.

Theorem 1**.**

A system digraph is structurally observable if the following two conditions are satisfied:

• •

Every state, is connected to a sensor via a directed path of states, i.e. $x_{i}\xrightarrow{path}\mathcal{Y}$, $i\in\{1,\ldots,n\}$.

• •

There is a sub-graph of disjoint cycles and output-connected paths that spans all state nodes.

Proof.

The original proof of the theorem is available in the work by Lin (Lin, 1974) for the dual problem of structural controllability, and more detailed proof is available in (Reinschke, 1988). The proof for structural observability is given in (Liu et al., 2013). ∎

The conditions in Theorem 1 are closely related to certain properties in digraphs. The second condition holds for structurally cyclic systems, since all states are included in a disjoint family of cycles. The first condition can be checked by finding Strongly Connected Components (SCCs) in the system digraph (Doostmohammadian & Khan, 2011). Recall that a SCC includes all states mutually reachable via a directed path. Therefore, the output-connectivity of any state in SCC implies the output-connectivity of all states in that SCC, and consequently, this satisfies the first condition in Theorem 1. By measuring one state in every SCC, all states in that SCC are reachable (and observable), i.e. $x_{i}\in SCC_{l}$ and $x_{i}\xrightarrow{path}\mathcal{Y}$ implies $SCC_{l}\xrightarrow{path}\mathcal{Y}$. This further inspires the concept of equivalent measurement sets for observability stated in the following lemma.

Lemma 2**.**

States sharing an SCC are equivalent in terms of observability.

Proof.

This is directly follows from the definition. Since for every two states $x_{i}$ and $x_{j}$ we have $x_{j}\xrightarrow{path}x_{i}$, then following Theorem 1 having $x_{i}\xrightarrow{path}\mathcal{Y}$ implies that $x_{j}\xrightarrow{path}\mathcal{Y}$. See more details in the previous work by the first author (Doostmohammadian & Khan, 2011). ∎

Observationally equivalent states provide a set of options for monitoring and estimation. This is of significant importance in reliability analysis of sensor networks. These equivalent options are practical in recovering the loss of observability in case of sensor/observer failure (Doostmohammadian & Khan, 2014b). In order to explore states necessary for observability, we partition all SCCs in terms of their reachability by states in other SCCs.

Definition 3**.**

• •

Parent SCC: is a SCC with no outgoing edge to states in other SCCs, i.e. for all $x_{i}\in SCC_{l}$ there is no $x_{j}\notin SCC_{l}$ such that $x_{i}\rightarrow x_{j}$.

• •

Child SCC: is a non-parent SCC (a SCC having outgoing edges to other SCCs), i.e. there exist $x_{i}\in SCC_{l}$ and $x_{j}\notin SCC_{l}$ such that $x_{i}\rightarrow x_{j}$.

Lemma 3**.**

Parent SCCs do not share any state node.

Proof.

The above lemma is generally true for all SCCs. The proof is clear; if two components share a state node they in fact make a larger component. ∎

Following the first condition in Theorem 1, the given definitions inspire the notion of necessary set of equivalent states for observability.

Lemma 4**.**

At least one measurement/sensing from every parent SCC is necessary for observability.

Proof.

This is because the child SCCs are connected to parent SCCs via a direct edge or a directed path. Therefore, $\mathcal{Y}$-connectivity of parent SCCs implies $\mathcal{Y}$-connectivity of child SCCs. In other words, $x_{i}\in SCC_{l}$, $x_{i}\rightarrow x_{j}$ and $x_{j}\in SCC_{k}\xrightarrow{path}\mathcal{Y}$ implies $SCC_{l}\xrightarrow{path}\mathcal{Y}$. See detailed proof in the previous work by the first author (Doostmohammadian & Khan, 2011). ∎

This further implies that * the number of necessary sensors for observability equals to the number of parent SCCs in structurally cyclic systems*. In such scenario, it is required to assign a sensor for each parent SCC in order to satisfy the observability condition. For more details on SCC classification and equivalent set for observability refer to (Doostmohammadian & Khan, 2011, 2014a).

4 Cost Optimization Formulation

In sensor-based applications every state measurement imposes certain cost. The cost may be due to, for example, maintenance and embedding expenses for sensor placement, energy consumption by sensors, sensor range and calibration, and even environmental condition such as humidity and temperature. In this section, we provide a novel formulation of the minimal cost sensor selection problem accounting for different sensing costs to measure different states. Contrary to (Pequito et al., 2014), the final formulation in this section has a polynomial order solution as it is discussed in Section 5.

Problem Formulation 1**.**

Assume a group of sensors and a cost $c_{ij}$ for every sensor $y_{i},i\in\{1,\ldots,m\}$ measuring state $x_{j},j\in\{1,\ldots,n\}$. Given the cost matrix $c$, the sensor selection cost optimization problem is to minimize sensing cost for tracking the global state of the dynamical system (1) (or discrete time system (2)). Monitoring the global state requires observability conditions, leading to the following formulation:

[TABLE]

where $A$ and $H$ are, respectively, system and measurement matrix, and $\mathcal{H}$ represents the $0-1$ structure of $H$.

In this problem formulation, $\mathcal{H}$ is the $0-1$ pattern of $H$, i.e. a nonzero element $\mathcal{H}_{ij}$ represents the measurement of state $x_{j}$ by sensor $y_{i}$. First, following the discussions in Section 2, the observability condition is relaxed to structural observability.

[TABLE]

Notice that in this formulation $(\mathcal{A},\mathcal{H})$-observability implies the structural observability of the pair $(A,H)$. Primarily assume that the number of sensors equals to the number of necessary measurements for structural observability. In control and estimation literature (Doostmohammadian & Khan, 2014a; Commault & Dion, 2015), this is addressed to find the minimal number of sensors/actuators . This consideration is in order to minimize the cost. Notice that extra sensors impose extra sensing cost, or they take no measurements and play no role in estimation. Therefore, following the assumptions in Section 1, the number of sensors, at first, is considered to be equal to the number of necessary measurements for observability (i.e. number of parent SCCs). This gives the following reformulation of the original problem.

Problem Formulation 2**.**

Considering minimum number of sensors for observability, the sensor selection cost optimization problem is in the following form:

[TABLE]

The added conditions do not change the problem. The constraint $\sum_{i=1}^{m}\mathcal{H}_{ij}\leq 1$ implies that all states are measured by at most one sensor, and $\sum_{j=1}^{n}\mathcal{H}_{ij}=1$ implies that all sensors are responsible to take a state measurement. Notice that, in case of having, say $N$, sensors more than the $m$ necessary sensors for observability, this condition changes to $\sum_{j=1}^{N}\mathcal{H}_{ij}\leq 1$ to consider the fact that some sensors are not assigned, i.e. they take no (necessary) measurement.

Next, we relax the observability condition following the results of Section 3 for structurally cyclic systems. Revisiting the fact that parent SCCs are separate components from Lemma 3, the problem can be stated as assigning a group of sensors to a group of parent SCCs. For this formulation, a new cost matrix $\mathcal{C}_{m\times m}$ is developed. Denote by $\mathcal{C}_{ij}$, the cost of assigning a parent set, $SCC_{j}$, to sensor $y_{i}$. Define this cost as the minimum sensing cost of states in parent $SCC_{j}$:

[TABLE]

This formulation transforms matrix $c_{m\times n}$ to matrix $\mathcal{C}_{m\times m}$. This transfers the sensor-state cost matrix to a lower dimension cost matrix of sensors and parent SCCs. Further, introduce new variable $\mathcal{Z}\sim\{0,1\}^{m\times m}$ as a structured matrix capturing the assignment of sensors to parent SCCs. Entry $\mathcal{Z}_{ij}$ implies sensor indexed $i$ having a state measurement of SCC indexed $j$, and consequently $SCC_{j}\xrightarrow{path}y_{i}$. Recalling that sensing all parent SCCs guarantee observability (see Lemma 4), the problem formulation can be modified accordingly in a new setup as follows.

Problem Formulation 3**.**

For structurally cyclic systems, having a set of $m$ sensors to be assigned to $m$ parent SCCs, the sensor selection cost optimization is given by:

[TABLE]

In this formulation, the new constraint $\sum_{j=1}^{m}\mathcal{Z}_{ij}=1$ is set to satisfy sensing of all parent SCCs as necessary condition for observability. The formulation in (9) is well-known in combinatorial programming and optimization. It is referred to as Linear Sum Assignment Problem (LSAP) (Pentico, 2007). It is noteworthy that the three statements in this section represent the same problem and the differences stem from mathematical relaxations and observability consideration.

The above formulation is one-to-one assignment of sensors and parent SCCs. By changing the first constraint to $\sum_{j=1}^{m}\mathcal{Z}_{ij}\geq 1$ we allow more than one parent SCC to be assigned to each sensor. This is the generalization to the primary assumption of assigning only one state to each sensor. For the second constraint in (9), considering $\sum_{i=1}^{m}\mathcal{Z}_{ij}\geq 1$ implies that more than one sensor may be assigned to a parent SCC. This adds redundancy in sensor selection and consequently increases the cost, and thus should be avoided. On the other hand, $\sum_{i=1}^{m}\mathcal{Z}_{ij}\leq 1$ violates the necessary condition for observability as some of the SCCs may not be assigned and tracked by sensors.

5 Linear Sum Assignment Problem (LSAP)

The novel formulation of sensor selection problem proposed in the Problem Formulation 3 is known to be a classical optimization problem referred to as the Assignment problem. Assignment problem is widely studied as many problems, e.g. in network flow theory literature, are reduced to it. The problem deals with matching two sets of elements in order to optimize an objective function. Linear Sum Assignment Problem (LSAP) is the classical problem of assigning $m$ tasks to $m$ agents (or matching $m$ grooms with $m$ brides, $m$ machines/companies to $m$ jobs, etc.) such that the matching cost is optimized (Pentico, 2007). The LSAP is mathematically similar to the weighted matching problem in bipartite graphs. This problem is also called one-to-one assignment as compared to one-to-many assignment problem in which one agent is potentially assigned to more than one task. There have been many solutions to this problem. From the original non-polynomial solution to later polynomial-time primal-dual solutions including the well-known Hungarian method. The Hungarian Algorithm, proposed by Kuhn (Kuhn, 1955) and later improved by Munkres, is of complexity order of $\mathcal{O}(m^{4})$ with $m$ as the number of tasks/agents. The algorithm was later improved by (Edmonds & Karp, 1972) to the complexity order of $\mathcal{O}(m^{3})$. The algorithm is given in Algorithm 1.

Other than these original solutions, recently new linear programming methods to solve the classical one-to-one LSAP and variations of this original setting are discussed. To name a few, distributed assignment problem based on a game-theoretic approach is proposed in (Zavlanos et al., 2008). Sensors/Agents are assigned to tasks relying only on local information of the cost matrix. The complexity of the algorithm is $\mathcal{O}(m^{3})$ in the worst case scenario. In (Bertsekas, 1981) a new algorithm is proposed whose average complexity matches Edmonds Hungarian method in large-scale. All these solutions can be applied to solve the Problem Formulation 3 and its variant, for example even when the sensing costs are changing. However, in terms of performance the Edmonds’ Hungarian algorithm (Edmonds & Karp, 1972) is more practical and used in programming softwares like MATLAB. The algorithm by (Zavlanos et al., 2008) is practical in distributed setting while the algorithm by (Bertsekas, 1981) is as practical as (Edmonds & Karp, 1972) only in large-scale applications. Note that the focus of this paper is on the polynomial complexity of such algorithms to be practical in large-scale application. Therefore, although other non-polynomial solutions to LSAP may exist, they are not of interest in large-scale sensor selection optimization.

Note that, in the LSAP the cost matrix has to be a complete $m$ by $m$ matrix. However, in practical application some states may not be measured by some sensors (not realizable by some sensors). This may, for example, caused by mismatch in range/calibration of the sensor and what is required for the state measurement. In the sensor-state cost matrix, $c$, this simply implies that some entries are not defined. For this unmeasurable states, the cost is infinite, in application a large enough cost (pseudo-cost $\tilde{c}_{ij}$) can be given. By introducing $\tilde{c}_{ij}$ and having a complete cost matrix, the LSAP problem can be solved using anyone of the polynomial methods mentioned in this section. Notice that if the optimal cost from LSAP in (9) is greater than the pseudo cost $\tilde{c}_{ij}$, the sensor selection has no feasible solution. A possible explanation is that at least one parent SCC is not realizable by any sensor, implying the assignment of a pseudo cost by LSAP. In case the feasible solution exists, no non-realizable state is assigned and the LSAP gives the optimal feasible solution in polynomial time.

6 Remarks

This section provides some remarks to further illustrate the results, motivation, and application of the polynomial order sensor selection solution proposed in this work.

Remark 1**.**

The main motivation on this paper is to find a polynomial order solution to optimize sensor selection problem for cyclic systems. Notice that in general, as mentioned in the introduction and literature review, the problem is NP-hard to solve, see for example (Pequito et al., 2014) and references therein. However, we showed that if system is cyclic there exist a polynomial order solution for sensor selection optimization. Note that this is significant in large-scale system monitoring as polynomial order algorithms are practical in large-scale applications because their running time is upper-bounded by a polynomial expression in the size of input for the algorithm. Examples of such large-scale cyclic systems are given in the introduction.

Remark 2**.**

It should be mentioned that LSI dynamics are practically used in state estimation and complex network literature, see (Liu et al., 2011, 2013; Doostmohammadian & Khan, 2014a; Dion et al., 2003; Davison & Wang, 1973, 1974; Ji & Egerstedt, 2007) and references therein. The motivation behind structured system theory is that this approach holds for systems with time-varying parameters while the system structure is fixed. This is significant in system theory as in many applications the system non-zero parameters change in time while the zero-nonzero pattern of the system matrix is time-invariant. Indeed, this structural analysis deals with system properties (including observability and system rank) that do not depend on the numerical values of the parameters but only on the underlying structure (zeros and non-zeros) of the system (Dion et al., 2003; Davison & Wang, 1973, 1974). It is known that if a structural property holds for one admissible choice of non-zero elements/parameters it is true for almost all choices of non-zero elements/parameters and, therefore, is called generic property (Dion et al., 2003). Another motivation is in linearization of nonlinear systems where the nonlinear model is linearized over a continuum of operating points, see (Doostmohammadian & Khan, 2014a; Liu et al., 2011). In this case the structure of the Jacobian matrix is fixed while the matrix elements change based on the linearization point, therefore implying the LSI system model. In this direction, the observability/controllability of LSI model implies the observability/controllability of the nonlinear model (Doostmohammadian & Khan, 2014a; Liu et al., 2011) and therefore the results of the LSI approach leads to conclusions on the nonlinear model.

Remark 3**.**

Based on the mentioned features of LSI model in Remark 2, the structural observability almost always implies algebraic observability, therefore LSI relaxation in Problem Formulation 2 almost always holds. Further, for structurally cyclic systems the problem can be exactly framed as a LSAP, and therefore the relaxation in Problem Formulation 3 and polynomial order solution is exact for cyclic systems.

Remark 4**.**

Note that the SCC decomposition is unique (Cormen et al., 2001) and therefore the cost matrix $\mathcal{C}$ and the formulation in (9) are uniquely defined.

Remark 5**.**

While this work focuses on sensor selection and observability, the results can be easily extended to the dual problem of controllability and particularly input/actuator selection. In this case the problem is to choose among the possible inputs to direct/control the dynamical system to reach the desired state with optimal cost. Note that the only mathematical difference is that the constraint in Problem Formulation 1 and 2 changes to $(\mathcal{A},\mathcal{H})$-controllability and the same graph theoretic relaxation holds. Because of duality the problem changes to assigning sensors optimally to Child SCCs resulting the same formulation as in Problem Formulation 3 where the solution is known via Hungarian algorithm. In fact, in the context of control of networked systems, this problem is also known as the so-called leader selection. In this problem using LSI model the idea is to determine the control leaders in structured multi-agent system. In this case the cost, for example, may represent energy consumption by agents. For more information on this subject we refer interested readers to (Fitch & Leonard, 2013; Lin et al., 2011, 2014; Clark et al., 2014).

7 Illustrative Examples

This section provides academic examples to illustrate the results of the previous sections.

Example 1: Consider a dynamical system with the associated digraph given in Fig.2. Every node represents a state of the system and every edge represent the dynamic interaction of two states. For example an edge from $x_{3}$ to $x_{1}$ and the self-loop on $x_{1}$ implies $\dot{x_{1}}=a_{13}x_{3}+a_{11}x_{1}$. Assuming a nonlinear dynamic system, the same link represents a possible nonlinear interaction function $\dot{x_{1}}=f_{1}(x_{3},x_{1})$ where the Jacobian linearization is in the form $\dot{x_{1}}=\frac{\partial f_{1}}{\partial x_{3}}x_{3}+\frac{\partial f_{1}}{\partial x_{1}}x_{1}$.444Notice that having self-cycle at every node implies that the diagonal entries of the Jacobian matrix, $\frac{\partial f_{i}}{\partial x_{i}}$, are non-zero and the Jacobian is structurally full rank. Such terminology holds for all state nodes and edges in the system digraph and relates the system digraph to the differential equation governing the dynamic phenomena.

A group of sensors are needed to track the $15$ states of this dynamic system. If a state is measurable by a sensor, the measurement is associated with a cost $c_{ij}$. However, not every state is measurable by every sensor. In this example, we assume the pseudo-cost as $\tilde{c}_{ij}=max(c_{ij})*m*n$. This prevents the assignment algorithm to assign these non-realizable sensor-state pairs. Among the realizable states, some states are not necessary to be measured. This is defined based on the SCC classification discussed in Section 3. In the particular example of Fig. 2, the inner component and the outer components have no outgoing edges to other SCCs; therefore, $\{x_{1},x_{2},x_{3}\}$, $\{x_{9},x_{10}\}$, $\{x_{11},x_{12},x_{13}\}$, and $\{x_{14},x_{15}\}$ are parent SCCs. The other components, $\{x_{4},x_{5},x_{6}\}$ and $\{x_{7},x_{8}\}$, are child SCCs. For observability, each parent SCC is needed to be tracked by at least one sensor. The selection of which state to be measured in each parent SCC is cost-based; in every parent SCC, each sensor measures the state with minimum cost. This takes the problem in the form given in Fig.3 and in the form of (9). Then, Hungarian method in Section 5 is applied to solve this LSAP.

For numerical simulation, in this system graph example we consider uniformly random costs $c_{ij}$ in the range $(0,10)$. Number of sensors equals to $m=4$ that is the number of parent SCCs. The non-realizable states are defined randomly with probability $50\%$, i.e. almost half of the states are not measurable by sensors 555We should mention that this is only for the sake of simulation to check the algorithm. In real applications if the measurable states are not observable no sensor selection optimization algorithm can provide an observable estimation of the system, therefore in real applications it is usual to assume that at least one observable solution exist for the problem otherwise no sensor selection and estimation scheme works.. In the assignment algorithm, the pseudo-cost of a non-realizable state is considered $\tilde{c}_{ij}=max\{c_{ij}\}*m*n$, which is certainly more than $\sum_{i=1}^{m}\sum_{j=1}^{n}c_{ij}$. Fig. 4 shows the cost of all realizable observable and non-observable assignments. In this figure, for the sake of clarification the indexes are sorted in ascending and descending cost order respectively for observable and non-observable assignments. Among the realizable state-sensor pairs, if the selected sensors do not measure one state in each parent SCC, this sensor selection is not observable (violating the Assumption (i) in Section 1). Among the observable selections, the optimal sensor selection has the minimum cost of $10.88$, which matches the output of the LSAP using the Hungarian method. Note that, the naive solution in Fig. 4 has complexity $\mathcal{O}(m!)$, and is only provided for clarification and checking the results of the LSAP solution.

Example 2: In order to show the advantage of using the proposed approach over the existing methods we provided another example. This example, shown in Fig.5, is similar to the example given in (Pequito et al., 2014), where the graph represents a social dynamic system. In such social digraphs each node represents an individual and the links represent social interaction and opinion dynamics among the individuals (Doostmohammadian & Khan, 2014a; Friedkin, 2006). According to (Doostmohammadian & Khan, 2014a; Friedkin, 2006), the social system is generally modeled as LSI system where the social interactions (as the structure) are fixed while the social influence of individuals on each-other change in time. This is a good example stating the motivation behind considering LSI model in this work. We intentionally presented this example to compare our results with (Pequito et al., 2014). As (Pequito et al., 2014) claims for such example there is no polynomial order solution to solve the problem while here we present a sensor selection algorithm with polynomial complexity of $\mathcal{O}(m^{3})$. Similar to the Example 1, a group of sensors (referred to as information gatherers in social system (Pequito et al., 2014)) are required to monitor $20$ states of the social system. Notice that having a cycle family covering all states the system is Structurally full-rank. The measurement of each state by each sensor is associated with a cost $c_{ij}$ and if not measurable the cost is assigned with $c^{\prime}_{ij}=max\{c_{ij}\}*m*n$. Applying the DFS algorithm one can find the SCCs and Parent/Child classification in $\mathcal{O}(m^{2})$ as shown in Fig.5. Then, the assignment problem in Fig.6 and in the form of (9) can be solved using the Hungarian method in $\mathcal{O}(m^{3})$.

Again we consider uniformly random costs $c_{ij}$ in the range $(0,10)$ for numerical simulation. Since there are $5$ parent SCCs in this social digraph we need $m=5$ social sensor. Among the sensor-state pairs, the non-realizable states are defined randomly with probability $30\%$ with pseudo-cost $\tilde{c}_{ij}=max\{c_{ij}\}*m*n$. In Fig. 7 the costs of all realizable observable and non-observable assignments are shown, where among the observable cases the optimal sensor selection has the minimum cost of $7.048$. As expected, this value matches the output of the proposed solution, i.e. the Hungarian algorithm for the LSAP in equation (9).

8 Conclusion

There exists many efficient algorithms to check the matching properties and structural rank of the system digraph, namely Hopcraft-Karp algorithm (Hopcroft & Karp, 1973) or the Dulmage-Mendelsohn decomposition (Dulmage & Mendelsohn, 1958) of $\mathcal{O}(n^{2.5})$ complexity. Moreover, the SCC decomposition and the partial order of SCCs is efficiently done in running time $\mathcal{O}(n^{2})$ by using the DFS algorithm (Cormen et al., 2001), or the Kosaraju-Sharir algorithm (Hopcroft, 1983). As mentioned earlier in Section 5, the LSAP solution is of complexity of $\mathcal{O}(m^{3})$. This gives the total complexity of $\mathcal{O}(n^{2.5}+m^{3})$ to solve the sensor coverage cost optimization problem. In dense graphs typically the nodes outnumber parent SCCs; assuming $m\ll n$, the complexity of the algorithm is reduced to $\mathcal{O}(n^{2.5})$. In case of knowing that the system is structurally cyclic, e.g. for self-damped systems, the running time of the solution is $\mathcal{O}(n^{2})$.

In practical application using MATLAB, the sprank function checks the structural rank of the system (i.e. the size of maximum matching in the system digraph). System is structurally cyclic if sprank(A) equals to $n$, the size of the system matrix. To find the partial order of SCCs the straightforward way (but not as efficient) is to use dmperm function. This function takes the system matrix $A$ and returns the permutation vectors to transfer it to upper block triangular form and the boundary vectors for SCC classification. The function assignDetectionsToTracks solves the assignment problem using Munkres’s variant of the Hungarian algorithm. This function takes the cost matrix and the cost of unassigned states/sensors as input, and returns the indexes of assigned and unassigned states/sensors as output.

As the final comment, recall that we consider sensor cost optimization only for structurally cyclic systems. For systems which are not structurally cyclic, other than parent SCCs, another type of observationally equivalent set emerges, known as contraction (Doostmohammadian & Khan, 2014a). Number of contractions equals to the number of unmatched nodes in the system digraph which in turn equals to system rank deficiency. Contractions and Parent SCCs determine the number of necessary states for system observability. The key point is that, unlike Parent SCCs which are separate sets, contractions may share state nodes with each other and with parent SCCs(Doostmohammadian & Khan, 2014a). This implies that the problem cannot be generally reformulated as LSAP and Problem Formulation 3 is only valid and exact for structurally cyclic systems. In general systems, particularly in structurally rank-deficient systems, a combination of assignment problem and greedy algorithms may need to be applied, which is the direction of future research.

Acknowledgement

The first author would like to thank Professor Usman Khan from Tufts University for his helpful suggestions and feedback on this paper.