Behind any complex system in nature or engineering, there is an intricate network of interconnections that is often unknown. Using a control-theoretical approach, we study the problem of network reconstruction (NR): inferring both the network structure and the coupling weights based on measurements of each node’s activity. We derive two new methods for NR, a low-complexity reduced-order estimator (which projects each node’s dynamics to a one-dimensional space) and a full-order estimator for cases where a reduced-order estimator is not applicable. We prove their convergence to the correct network structure using Lyapunov-like theorems and persistency of excitation. Importantly, these estimators apply to systems with partial state measurements, a broad class of node dynamics and internode coupling functions, and in the case of the reduced-order estimator, node dynamics and internode coupling functions that are not fully known. The effectiveness of the estimators is illustrated using both numerical and experimental results on networks of chaotic oscillators.

Understanding the interplay between the network topology of interconnected units (or components) and their emergent behavior is a fundamental problem in network science.^{1–3} For instance, the topology of interactions between different institutions in financial networks is central to predicting behaviors that arise during a crisis.^{4} Network topology affects coordination in human ensembles^{5} and the stability of power networks.^{6} Moreover, recent developments on controllability,^{7} observability,^{8} and control^{9} of complex networks require the topology to be known. In many applications, however, including mobile sensor networks, biological networks, and social networks, the network topology may be unknown or time-varying.^{10–12} Uncovering the topology of such networks is fundamental to understanding and controlling their complex behavior.^{11,12} Using adaptive control theory, we study the problem of inferring the network topology of complex systems based on measurements of each node’s activity. In particular, we construct two different estimators, a reduced-order estimator and a full-order estimator, that can infer the network topology in real time.

## I. INTRODUCTION

Discovering the underlying network topology of complex systems is a fundamental problem in fields ranging from biology to engineering.^{11,12} Consider, for example, a gene regulatory network. Even if the individual dynamics of gene expression are well understood, determining the tangled network of interactions is still a significant challenge.^{13,14} As another example, the network topology of a power grid might change over time. For instance, renewable sources depending on weather conditions might be connected/disconnected or transmission failures might occur.^{15} A pressing challenge is to design appropriate methods to infer (or reconstruct) the topology of a network in real time.

We consider a general network model of $N>1$ nodes, each described by a set of nonlinear ODEs of the form

where $xi\u2208Rn$ is the state variable of the $i$th node and $N={1,\u2026,N}$ is the finite set of indexes. The state could represent the abundance of species in ecological networks^{16} or voltages and currents in electrical networks.^{17} The measurable “output” of a node is $yi$, a linear combination of the states specified by the output matrix $Ci\u2208Rmi\xd7n$. The function $fi:Rn\u21a6Rn$ is continuous and represents the intrinsic dynamics of each node. The constants $aij$ are the entries of a weighted adjacency matrix $A=[aij]$, where $aij\u22600$ indicates that there is an interconnection link from node $j$ to node $i$ through a continuous coupling function $hij:R2n\u21a6Rn$. We assume that the diagonal elements $aii$ are zero, as such terms can be absorbed into $fi$.

We are interested in the network reconstruction (NR) problem, i.e., inferring the coupling weights $aij$ in real time using measurements of each node’s activity $yi$.

Several methods have been proposed for NR, and they can be classified in two major categories: “model free” and “model based.” In the model-free approach, information theory can be exploited for NR without requiring full knowledge of the network model (1).^{12} In particular, interconnection links are inferred using metrics based on conditional,^{18} transfer,^{19} or causal entropies.^{20} These techniques are not able to precisely reconstruct the coupling weights associated with each link, however, and they typically require a large amount of data and computational resources,^{19} making their real-time implementation cumbersome. Similarly, the network reconstruction problem has also been addressed within the context of graphical models where the network topology is assumed to have locally treelike graphs^{21} or the node states are normally distributed.^{22}

In contrast, model-based techniques make use of model (1). For instance, in Ref. 13, the network is linearized around a structurally stable equilibrium point, and the network is reconstructed using classic regression techniques. Extensions of this strategy have been proposed for cases where the network solutions are periodic^{23} or can be driven to an equilibrium point by adding linear feedback controllers at every node of the network.^{24} (For a more comprehensive list of NR techniques, see Ref. 11 and references therein.)

More recently, Ref. 25 proposed different metrics inspired by the Kuramoto model to infer coupling links, while in Ref. 26, the NR problem is presented and solved as an underdetermined regression model based on the fact that all $x\u02d9i$ are known or can be computed from measurements $xi$. A related approach is described in Ref. 27. However, none of these approaches is able to exactly infer the interconnection links and their associated coupling weights. In fact, they lack a proper convergence analysis, and they often rely on strong assumptions like special steady-state solutions, local approximations (linearization), full state measurements ($yi=xi$), and perturbing or controlling all the nodes of the unknown network. The network reconstruction problem is also known as network tomography^{28,29} where particular diffusion models are considered such as the Combine-Then-Adapt (CTA) diffusion mechanism.^{29} A similar diffusion process is considered in Ref. 30 where the goal is to reconstruct the topology associated with the graph Laplacian.

Copy-synchronization^{31–34} is a model-based approach to the NR problem where exact convergence can be guaranteed. This method consists of an estimated “copy” (or slave) network with adaptive coupling weights, and the goal is to drive the error between the copy network and the unknown (or master) network to zero. Often, but not always, copy-synchronization strategies require full knowledge of the network dynamics ($fi$ and $hij$), full state measurements ($yi=xi$), and particular coupling functions (e.g., linear diffusive coupling).^{31,32,34,35} In addition, the convergence of these strategies has not been fully understood within the context of complex network systems (see, for instance, the supplementary material in Ref. 36 and page 18 in Ref. 11).

To overcome these limitations, in this paper, we interpret network copy-synchronization as a control-theoretic adaptive observer technique, where the copy network is an adaptive observer that reconstructs the master network. This formulation allows us to derive novel NR methods guaranteeing exact convergence while making fewer assumptions about the network and its measurements, as described below. We wish to emphasize that the proposed adaptive estimators are similar to gradient-descent-based algorithms that are known to track slowly varying coupling weights.^{37–39}

### A. Related work and contributions of this paper

Using a copy network for NR was originally proposed in Ref. 31. As pointed out in Ref. 32, however, convergence to the real network structure fails when all the trajectories of the unknown network converge to each other, i.e., synchronization. The authors state that a sufficient condition to guarantee exact convergence to the real network topology is that all the node trajectories should be linearly independent (LI) on the zero-error manifold, that is, when the error between the copy network and the real one is zero. More recently, extensions of the copy-synchronization approach satisfying the LI condition were given in the context of outer-synchronization^{34} (where the copy network can be an oscillator of lower dimension, but full state measurements are required), linearly coupled networks with partial state measurements,^{35} and systems with delays and uncertainties.^{40–42} The requirement of LI on the zero-error manifold is a significant limitation of these approaches, however, since it is only a sufficient condition that might be very difficult, if not impossible, to test for networks with a large number of components.^{34} Reference 33 instead uses the persistent excitation (PE) condition, which is simpler to test. Only networks with diffusive couplings and full state measurements were addressed, however.

Despite the fact that there is a direct relation between network copy-synchronization and adaptive observers, until now the relationship has not been made explicit, even though master-slave synchronization has been shown to be closely related to observers in the context of single chaotic systems.^{43,44}

In this paper, we exploit adaptive observer theory to study the NR problem for a more general class of networks, including networks with nonlinear heterogeneous nodes and coupling functions, networks where full state sensing is not available, and networks where the node dynamics and coupling functions are not fully known. Compared to the existing literature, the contribution of the paper is threefold:

Section III: We introduce the “reduced-order copy-synchronization estimator,” where the dynamics of each node are projected to a single state variable. It is always possible to design such an estimator if all state variables of the unknown network are measurable. Due to the smaller number of estimator state variables, the reduced-order estimator has lower computational complexity than existing copy-synchronization estimators, making the approach more suitable for large networks. In addition, reduced-order estimators can be designed for systems where only a subset of states can be measured, and it is not necessary to have full knowledge of the node dynamics or coupling functions. Using the PE condition, we derive necessary and sufficient conditions guaranteeing convergence of the copy network to the unknown network.

Section IV: In cases where a reduced-order estimator cannot be constructed, we extend the “full-order copy-synchronization estimator” approach to (i) a more general class of networks than previously considered (e.g., networks with nonidentical coupling functions and non-Lipschitz node dynamics) and (ii) cases where only partial state measurements are available. Convergence to the unknown network structure is guaranteed using PE rather than LI. Preliminary results for networks with diffusive couplings were presented in Ref. 45.

Sections V and VI: We validate the new real-time NR estimators on simulated networks as well as experimental networks of chaotic electrical oscillators (Chua’s circuits). We find that, in practice, there is a trade-off between the accuracy of the estimation and the convergence rate that depends on the estimator parameters.

## II. PRELIMINARIES AND PROBLEM STATEMENT

### A. Notation

The $N\xd7N$ identity matrix is written $IN$ and $0$ is a matrix of zeros of appropriate dimension. The notation $\Vert \u22c5\Vert $ denotes the Euclidean norm for vectors and the spectral $L2$-norm for matrices. A block diagonal matrix $D$ with diagonal blocks $d1,\u2026,dN$ is indicated by $D={diag}{d1,\u2026,dN}$ and $\lambda min(A)$ denotes the minimum eigenvalue of a square matrix $A$. $A\u227b0$ indicates that the matrix $A$ is positive definite. $A\u2ab0B$ indicates that the matrix $A\u2212B$ is positive semidefinite.

A signal $\omega (t)$ satisfies the “PE condition” if it excites all the modes of the system.^{38} This condition has been widely used in the context of adaptive control and adaptive observers.^{39,44,46} (The reader is referred to Chap. 6 of Ref. 38 for more details about the PE condition.)

### B. Problem formulation

Let $x:=[x1\u22a4,\u2026,xN\u22a4]\u22a4$ be the stack vector of all network states. Then, (1) can be written as

where $ai\u2208RN\u22121$ is the stack vector of unknown coupling weights and $Hi(x)\u2208Rn\xd7(N\u22121)$ is the regressor matrix, i.e.,

for all $i,k,j\u2208N$ such that $k,j\u2260i$ and $j>k$. In this paper, we address the problem of designing adaptive strategies for updating an estimate $a^i$ of the unknown network topology $ai$ based on measurements $yi$. The network (3) can be seen as a dynamical system with linearly parametrized unknowns,^{37,38} so we can consider estimators of the form

where $a^i$ are the estimated coupling weights and $a^0,i$, $\zeta 0,i$ are the initial conditions. Functions $\chi i(y,\zeta i)$ and $qi(\zeta i,y,a^i)$ are the “estimator functions,” where $\chi i$ represents the adaptation law for reconstructing the coupling weights. In what follows, we will show that NR is guaranteed by properly choosing these functions.

*for all*$i\u2208N$.

### C. The original copy-synchronization approach

The original copy-synchronization approach,^{31} which this paper extends, consists of estimators (6) and (7) with functions

where $\mu i$ and $\kappa i$ are positive constant parameters representing the adaptation gain and control strength, respectively. Note that this copy-synchronization approach requires (i) the dimension of each estimator state $\zeta i$ to be the same as the node dynamics of the unknown network ($s=n$), (ii) full state measurements, i.e., $Ci=In$, and (iii) full knowledge of the intrinsic node dynamics $fi$ and coupling functions $hij$.

In Sec. III, we show that the dimension of the estimator state $\zeta i$ can often be reduced to one by using a one-dimensional projection. Partial rather than full state measurements are considered, and most importantly, full knowledge of the intrinsic node dynamics and coupling functions is not required.

## III. REDUCED-ORDER COPY-SYNCHRONIZATION

To reconstruct the network structure (the coupling matrix $A$), it is often not necessary to estimate all the nodes’ states or to have full knowledge of the node dynamics and coupling functions. We can project each node’s dynamics to a one-dimensional space, meaning we only need to estimate this scalar value to reconstruct the coupling weights. The ability to do this relies on the existence of a particular kind of projection operator, which is discussed below.

### A. One-dimensional projection

### B. Reduced-order estimator

Note that the projected dynamics of the $i$th node can be seen as a first-order system with $N\u22121$ linearly parametrized unknowns. We can directly make use of classic observer-based estimators^{37,38} to obtain the following result.

$(b)$ A stronger sufficient but not necessary condition for asymptotic NR is that all $\varphi i(yi)$ have continuous second partial derivatives$,$ all $hij(xi,xj)$ have continuous Jacobians$,$ and all $wi$ are PE $($see Definition II.1).

From Assumption II.1, we have that all $xi$ are bounded. Then, by continuity of the functions $\u2207\varphi i$ and $hij$, we have that all $wi$ are bounded as well. Using Theorem A.1 in the Appendix with $A=\u2212\kappa i$ and the fact that all $wi$ satisfy (16), we can conclude that $limt\u2192\u221e\Vert ei(t)\Vert =\Vert a~i(t)\Vert =0$ for all $i\u2208N$.

In contrast to the original full-order copy-synchronization approach, the reduced-order estimators (14) and (15) use partial state measurements and the dimension of the estimator state $zi$ is one. Indeed, each estimator $i$ only needs $N$ state variables to reconstruct the entire network, rather than $n+(N\u22121)$ as in (10). In addition, both the node dynamics $fi$ and coupling functions $hij$ are not required to be fully known, as illustrated in the following example.

Theorem III.1 provides two conditions guaranteeing asymptotic NR using (16) and (2), respectively. Condition (16) is necessary to guarantee convergence to the real network topology. Similar to the PE condition$,$ (16) measures the richness of the signals $wi$. In practice, however$,$ the PE condition (2) is easier to verify than (16),^{38} but extra conditions on the projection function $\varphi i(yi)$ and coupling functions $hij$ are required. Indeed$,$ in our experimental results, we show that condition (2) can be easily tested $[$e.g.$,$ see Fig. 8(c) $]$.

In the case of full state measurements ($Ci=In$ for all $i\u2208N$), as usually considered in the NR literature,^{31–34,41,42} a one-dimensional projection that is output-compatible $($see Definition III.1) always exists. In some cases of partial state measurements$,$ however$,$ such a projection might not exist. In such cases$,$ a projection onto a $pi$-dimensional space may be possible instead$,$ i.e.$,$ a function $\varphi i:Rmi\u21a6Rpi$ with $1<pi<n$. This also results in a reduced-order estimator.

## IV. FULL-ORDER COPY-SYNCHRONIZATION

In Sec. III, we showed that the state of each estimator can be reduced to one by using the reduced-order observer (14) and (15) under Assumption III.1. In some cases, however, there may not exist a valid projection of the dynamics to a lower-dimensional space (see also Remark III.2), and full-order estimation is required. In this section, we derive a full-order copy-synchronization estimator that is convergent for a broader class of networks than previously studied in the copy-synchronization literature.

We assume that the intrinsic node dynamics can be written as the sum of linear and nonlinear terms, i.e.,

where $Ai\u2208Rn\xd7n$ and $\xi i:\Omega \u2286Rn\u21a6Rn$. Many nonlinear systems can be written in this form, e.g., Lorenz oscillators, Chua’s circuit, and others.^{47} Moreover, we rewrite the coupling functions as $hij(xi,xj)=\Gamma ih~ij(xi,xj)$, with $\Gamma i\u2208Rn\xd7n$ being a square matrix. For instance, in the example discussed in Sec. V B, the matrix $\Gamma i$ represents the fact that nodes are coupled though the first state variable only. This matrix is often called the “inner coupling matrix.”^{2} We will show that $\Gamma i$ must satisfy “the matching condition” [see Eq. (23)] to guarantee stability of the full-order copy-synchronization estimator.

With these definitions, classic adaptive estimators^{37,38} can be used for NR,

where $H~i(\zeta ):=[h~ik(\zeta i,\zeta k),\u2026,h~ij(\zeta i,\zeta j)]$ for all $i,k,j\u2208N$ such that $k,j\u2260i$ and $j>k$. The $n\xd7mi$ matrices $Li$ and $Mi$ should be properly designed according to the following result.

Consider the unknown network (1) where $fi(xi)$ can be written as in (20)$,$ Assumption II.1 holds$,$ and the elements of the network adjacency matrix $A$ satisfy $|aij|\u2264amax$. Then, estimators (21) and (22) asymptotically reconstruct the network topology if the following conditions are satisfied$:$

$(C5)$ all matrices $Hi\u22a4(\zeta )=H~i\u22a4(\zeta )\Gamma i\u22a4$ are PE $($see Definition II.1).

For the sake of clarity, we split the proof in three steps as follows:

Condition (C4) of Theorem IV.1 guarantees $V\u02d9\u22640$, implying $ei$ and $a~i$ are bounded for all $i\u2208N$, and $limt\u2192\u221eV(t)=V\u2217<\u221e$.

Note that the linear part of the vector field is crucial for stabilizing the estimator$,$ as with the classical Luenberger and adaptive observers with Lipschitzian nonlinearities.^{48–50} Finding appropriate matrices satisfying condition $(C1),$ however, might not be an easy task. This is a common problem in the literature of adaptive observers.^{50,51} The existence of matrices $Li$ can be studied using the Kalman-Yakubovich-Popov lemma on strictly positive real functions$,$^{52} and its solutions can be found using standard software for Linear Matrix Inequalities (LMI) $($see the example in Sec. VI C 2).

## V. EXAMPLES

### A. Networks of Lorenz oscillators

We consider the network reconstruction problem in (1) where the intrinsic node dynamics are Lorenz systems given by^{2}

where $\sigma =10$, $\omega =\u22128/3$, and $\rho $ is an unknown constant. We consider the nonlinear coupling functions in (1) to be given by $hij(xi,xj)=[tanh\u2061((vi\u2212vj)/\epsilon i),0,0]\u22a4$, with $\epsilon i$ being a positive constant for all $i\u2208N$. In addition, we assume that only the states $vi$ and $qi$ are accessible for measurements for each node $i\u2208N$; that is, $yi=Cxi$, where $C\u2208R2\xd73$ is given in (19).

To implement the reduced-order estimators (14) and (15), we need to find functions $\alpha i$ and $\beta ij$ satisfying Assumption III.1. We start by choosing the one-dimensional function $\varphi i(yi)=(1/2)vi2\u221210qi$, whose gradient is given by $\u2207\varphi i\u22a4(yi)=[vi,\u221210]$. Then, using Eqs. (11) and (12), we find that

Note that the functions $\alpha i$ and $\beta ij$ depend only on the output states, satisfying Assumption III.1. Also, similar to Example III.1, there is no need to know the unknown constant $\rho $.

To test the reduced-order estimator, we consider the unknown network topology depicted in Fig. 1 where $N=8$ ($N={1,2,3,4,5,6,7,8}$). We choose $\rho =28$ and $\epsilon i=1.5i/(i+1)$ for all $i\u2208N$ for the node dynamics and coupling functions. In addition, we select $\kappa i=\kappa =100$ and $\mu i=\mu =10$ for all $i\u2208N$.

To verify the PE condition of Theorem III.1, we define $Qi=\u222btt+Twi(yi(\tau ))wi(yi(\tau ))\u22a4d\tau $, and we calculate the minimum eigenvalue $\lambda min(Qi)$ setting $T=10s$ for all $i\u2208N$. The time evolution of the estimated network topology and $\lambda min(Qi)$ are shown in Figs. 2(a) and 2(b), respectively. Note that $\lambda min(Qi)$ is always positive, thus fulfilling the PE condition of Theorem III.1 and guaranteeing convergence to the unknown coupling weights.

### B. Networks of Chua’s circuits

As a second example, we consider networks of chaotic circuits interconnected through identical resistors.

#### 1. Node dynamics

The dynamics of each node are given by the electrical oscillator depicted in Fig. 3(a). The element $D$ represents a modified nonlinear Chua’s diode whose voltage-current relation is a linear function as shown in Fig. 3(b). The circuit realization of $D$ is depicted in Fig. 3(c). Using Kirchhoff’s laws, we find the circuit dynamics are given by^{54}

where $vC1$ and $vC2$ are the voltages across capacitors $C1$ and $C2$, respectively, $IL$ is the current through the inductor $L$, and $INet$ is an incoming current from interconnections with neighboring circuits. In addition, $ID(vC1)$ is the diode current given by

where

All parameter values are reported in Table I.

Symbol . | Parameter . | Value . | Unit . |
---|---|---|---|

C_{1} | Capacitance | 10.05 | nF |

C_{2} | Capacitance | 100.65 | nF |

L | Inductance | 18.5 | mH |

R | Resistance | 1796 | Ω |

R_{1} | Resistance | 219.79 | Ω |

R_{2} | Resistance | 219.83 | Ω |

R_{3} | Resistance | 2200 | Ω |

R_{4} | Resistance | 22 | kΩ |

R_{5} | Resistance | 22 | kΩ |

R_{6} | Resistance | 3296.25 | Ω |

V_{sat} | Voltage | 14 | V |

Symbol . | Parameter . | Value . | Unit . |
---|---|---|---|

C_{1} | Capacitance | 10.05 | nF |

C_{2} | Capacitance | 100.65 | nF |

L | Inductance | 18.5 | mH |

R | Resistance | 1796 | Ω |

R_{1} | Resistance | 219.79 | Ω |

R_{2} | Resistance | 219.83 | Ω |

R_{3} | Resistance | 2200 | Ω |

R_{4} | Resistance | 22 | kΩ |

R_{5} | Resistance | 22 | kΩ |

R_{6} | Resistance | 3296.25 | Ω |

V_{sat} | Voltage | 14 | V |

#### 2. Network model

We consider $N>1$ identical interconnected Chua’s circuits. The variables $vC1,i$ and $vC2,i$ represent the voltages across $C1$ and $C2$, and $IL,i$ represents the current through the inductor $L$ for the $i$th circuit. The interconnections between circuits are through identical resistors $Rlink$. The overall dynamics of the electrical network are given by

where $Sij=1$ if there is an interconnection between circuits $i$ and $j$, while $Sij=0$ otherwise. Letting $xi=[vC1,i,vC2,i,IL,i]\u22a4$, the network of Chua’s circuits can be rewritten as in (1) with $fi(xi)=Axi+\xi i(xi)$, where

The unknown coupling weights are given by $aij=Sij/Rlink$ and the coupling functions $hij(xi)=(vC1,j\u2212vC1,i)/C1$. We assume that both $vC1,i$ and $vC2,i$ for all $i\u2208N$ are accessible for measurements, i.e., $yi=Cxi$, where $C$ is given in (19).

### C. Estimator design and numerical simulations

#### 1. Case a: Reduced-order estimator

We first design the reduced-order estimators (14) and (15) by choosing the function $\varphi i(yi)=vC1,i$ for all $i\u2208N$. Using (49) and (11) and (12), we find that $\alpha (yi)=(1/(RC1))(vC2,i\u2212vC1,i)$ and $\beta ij=(\u22121/C1)ID(vC1,i)$, and Assumption III.1 is satisfied. We assume that the coupling resistor $Rlink$ is unknown but belongs to the interval $Rlink\u2208[20,70]k\Omega $, implying $|aij|\u22641/(20k\Omega )=50\mu S$. Note that all estimates of $aij=Sij/Rlink$ take very small values (on the order of $10\u22125$) if there is an interconnection link between nodes $i$ and $j$. Although this is not an issue for numerical simulations, in practice, there is not exact convergence due to unavoidable model mismatches and noise (see Sec. VI). This induces oscillations around the true value and would make difficult to distinguish between links that are zero or those which are in the order of $10\u22125$. To overcome this issue, we replace the state $a^i$ in the left hand of (15) by $a^i/ar$, where $ar$ is a positive constant representing the rescaling of the estimate. The estimate in (14) will asymptotically converge to $a^i\u2192arai=arSij/Rlink$. By setting $ar=100k\Omega $, we have that the rescaled estimates would converge to either zero or values between $(100/70,100/20)=(1.4286,5)$ for nonexisting ($Sij=0$) or existing ($Sij=1$) links, respectively. This makes it easier to distinguish between links even if there are oscillations.

For the numerical simulations, we consider the case where the unknown network topology is a line [see Fig. 7(a)] and $Rlink=50k\Omega $. The time evolution of the reduced-order estimator using the rescaled variable with $ar=100k\Omega $ is shown in Fig. 4(a) for $\kappa i=\kappa =50000$ and $\mu i=\mu =20000$. Note that network reconstruction is attained provided the PE condition of Theorem III.1 is satisfied. Indeed, the minimum eigenvalue of $Qi=\u222btt+Twi(yi(\tau ))wi(yi(\tau ))\u22a4d\tau $ for $T=3\mu s$ depicted in Fig. 4(b) is always positive. We have chosen these values of $\kappa $ and $\mu $ in order to compare the simulations with the full-order copy-synchronization approach and with experimental results. We show in Sec. VI that there is a trade-off between the rate of convergence and accuracy of the estimation.

#### 2. Case b: Full-order estimator

To design the estimators (21) and (22), we use Theorem IV.1. We first start by noticing that the coupling functions can be rewritten as $hij(xi)=(vC1,j\u2212vC1,i)/C1=\Gamma h~ij$ with $\Gamma ={diag}{1,0,0}$ and $h~ij=(1/C1)(xj\u2212xi)$ for all $i,j\u2208N$. Because the coupling functions $h~ij$ are linear, the Lipschitz condition (C3) of Theorem IV.1 is easily verified with $Lij=Lmax=1/C1$.

Next, assuming all nodes share identical parameters, we need to find matrices $Pi=P$, $Ri=R$, and $Mi=M$ $i\u2208N$ such that the first two conditions of Theorem IV.1 are satisfied. Using the software cvx, we solve the LMI of condition (C1),

and

Next, we show that $\xi (\u22c5)$ satisfies the QUAD condition (25). Using the fact that $\xi (\u22c5)$ is a piecewise smooth linear function, we have that $(z1\u2212z2)\u22a4P(\xi (z1)\u2212\xi (z2))\u2264\rho (z1\u2212z2)\u22a4(z1\u2212z2)$, for all $z1,z2\u2208R3$, with $\rho =\u2212m0/C1=7.5415\xd7104$. From the earlier calculation, we have that $Lmax=1/C1$ and $|aij|\u2264amax=50\xd710\u22126$. Plugging $Mmax=1$, $L=Lmax$, and $N=4$ into the right-hand side of the inequality in condition (C4) of Theorem IV.1, we get $\gamma >\rho +amaxMmaxL(1+N)=1.003\xd7105$. It follows that $\gamma =1.3672\xd7105$ and condition (C4) is fulfilled. As also done in the previous case, we rescale the estimation variable $a^i$ in (22) by substituting $a^i/100$ k$\Omega $. The time evolution of the estimated network topology is shown in Fig. 5(a) for $\mu i=\mu =20000$. Similar to the previous examples, we calculate $\lambda min(Qi(t))$, where $Qi=\u222btt+THi(\tau )\u22a4Hi(\tau )d\tau $. We find that setting $T=3\mu s$, all $\lambda min(Qi(t))$ are lower-bounded by a positive constant as shown in Fig. 5(b), satisfying the PE condition of Theorem IV.1.

## VI. EXPERIMENTAL VALIDATION

To carry out the experimental validation of the network reconstruction strategies, we use the experimental setup of interconnected Chua’s circuits shown in Fig. 6. The setup was developed at the University of Naples Federico II, Naples, Italy, and is fully described in Ref. 55. It consists of (1) an electrical module of Chua’s circuits, (2) a data acquisition system, (3) software developed in MATLAB for real-time monitoring of the network states, and (4) an external storage unit to record all the network output signals. The experimental parameters of the real circuits are those reported in Table I.

For the sake of simplicity, the experiments were conducted on networks of four nodes with three different topologies, as depicted in Fig. 7. The unknown coupling resistance is $Rlink=49.027k\Omega $. We obtained measurements of $vC1,i$ and $vC2,i$ for each node with a sampling rate of $13.3\mu s$. We use those signals to drive the estimators designed in Sec. V.

### A. Case a: Reduced-order estimator

We first consider the case where the unknown network topology is a line [see Fig. 7(a)]. We use the reduced-order estimator designed in Sec. V C where $\varphi i(yi)=vC1,i$ and the rescaled variable $a^i/100k\Omega $.

The time trajectories of the estimated network topology and their corresponding $\lambda min(Qi)$ are shown in Figs. 8(a) and 8(b), respectively. Note that $\lambda min(Qi)$ is always positive, satisfying the PE condition of Theorem III.1 and suggesting that the signals are rich enough to reconstruct the network. Indeed, the estimates $a^i$ converge with a low residual error [0.2472, see Fig. 8(b)]. This error is due to unavoidable mismatches between the mathematical model and the real experimental model and also from small delays and quantization errors, which are added during the signal acquisition.

Unlike the numerical simulations reported in Sec. V, we use a larger value of $\kappa $ and a smaller value of the adaptation gain $\mu $. In Fig. 4(a), convergence is attained faster (in 0.005 s) by setting $\mu =20000$. In practice, however, the measurements $yi$ have small quantization errors, and high values of $\mu $ amplify them, adversely affecting the network estimation. To illustrate this trade-off between accuracy and rate of convergence, we plot the estimation error $\Vert a~\Vert $, varying the adaptation gain $\mu $ (see Fig. 9) when the unknown network topology is a ring. Note that for lower values of $\mu $, convergence is slower and the estimation error is small, while the error increases for larger values of $\mu $.

A similar conclusion can be made for all three network topologies (see Fig. 7) as shown in Fig. 10(a). Note that the steady-state error increases for larger values of $\mu $. Moreover, we also observe that the error can be decreased by choosing larger values of $\kappa $ as shown in Fig. 10(b).

### B. Case b: Full-order estimator

## VII. CONCLUSIONS

Using a control theoretical approach, we have studied the inference problem of an arbitrary network structure based on measurements of each node’s activity. By recasting the problem in the context of adaptive observers, full-order and reduced-order estimators were proposed. Moving beyond existing methods in the literature, these estimators are effective for a broader class of node dynamics and coupling functions, arbitrary coupling weights, and partial state measurements. We derived necessary and sufficient conditions guaranteeing asymptotic network reconstruction, and we showed that under appropriate projections, the number of estimator variables can be substantially reduced, making the estimation suitable for large networks. Numerical and experimental results validate the effectiveness of the proposed strategies.

In the examples in this paper, the chaotic dynamics of the network nodes ensure that the PE condition is satisfied. For some complex network dynamics, the PE condition (see Definition II.1) may not be satisfied, as in the case that all node trajectories converge to each other (synchronization). The node trajectories are linearly dependent and the PE condition is not satisfied, making it impossible to reconstruct the network.^{32} Recent studies indicate that network symmetries represent a fundamental obstacle to network reconstruction.^{36} This issue can be sidestepped if the transient dynamics before synchronization are rich enough^{56} or if we break the network symmetry by injecting perturbation signals on a fraction of nodes. Finding the minimal number and location of such inputs is not a trivial task, however, since the only available information is the node state trajectories. We are currently investigating this issue.

## ACKNOWLEDGMENTS

This work was supported by the Office of Naval Research (ONR) under Grant No. N00014-13-1-0331 and by the Army Research Lab. In addition, we wish to thank the Department of Electrical Engineering and Information Technology, University of Naples Federico II, Italy, and, in particular, to Professor Mario di Bernardo for providing hosting in his lab, and Professor Massimiliano de Magistris and Professor Carlo Petrarca for providing access to their experimental setup. Additionally, we acknowledge the Young Investigator Training Program (YITP)—ISCAS2018 for providing the scholarship that funded the stay of Daniel Alberto Burbano Lombana in Italy. The dataset and codes used in the examples will be provided upon request.

### APPENDIX: STABILITY OF ADAPTIVE SYSTEMS

#### (Theorem 2.17 and Corollary 2.3 in Ref. 38)

#### (Ref. 46)

## REFERENCES

*et al.*, “

*2018 IEEE International Symposium on Information Theory (ISIT)*(IEEE, 2018), pp. 1839–1843.

*American Control Conference (ACC)*(IEEE, 2018), pp. 3398–3403.

*Chua’s Circuit Implementations: Yesterday, Today and Tomorrow*, World Scientific Series on Nonlinear Science Series A Vol. 65 (World Scientific, 2009).

*International Conference on Nonlinear Dynamics of Electronic Systems*(Springer, Cham, 2014), pp. 203–210.