Efficient Network Reconstruction from Dynamical Cascades Identifies Small-World Topology of Neuronal Avalanches

Cascading activity is commonly found in complex systems with directed interactions such as metabolic networks, neuronal networks, or disease spreading in social networks. Substantial insight into a system's organization can be obtained by reconstructing the underlying functional network architecture from the observed activity cascades. Here we focus on Bayesian approaches and reduce their computational demands by introducing the Iterative Bayesian (IB) and Posterior Weighted Averaging (PWA) methods. We introduce a special case of PWA, cast in nonparametric form, which we call the normalized count (NC) algorithm. NC efficiently reconstructs random and small-world functional network topologies and architectures from subcritical, critical, and supercritical cascading dynamics and yields significant improvements over commonly used correlation methods. With experimental data, NC identified a functional and structural small-world topology and its corresponding traffic in cortical networks with neuronal avalanche dynamics.

1. Determine the rule by which the n a source nodes will be chosen for every target node.
In our current implementation in reconstructing networks from neuronal avalanches, we took the approach developed in [29,30] in which the events were binned with a time window ∆t chosen as a mean of the time-interval distribution from successively activated nodes in the network. In this particular case, only nodes in the preceding time bin are treated as the potential source nodes. Thus, the raster becomes a N t × N sparse matrix, where N t = T exp /∆t is the total number of time bins.
2. Determine the number of propagation steps, N p , i.e. number of the incidences in which each of the two successive time-bins contains at least one active node. This number should be greater than the number of potential links in the network (N (N −1)) ( Figure  4A) 3. Evaluate the NC ij factor for all links, according to Equation 18. If highly supercritical dynamics is suspected, the more general Eqs. 14 and 15 should be used. ij , obtain the reconstructed network topology and weights according to Eqs. 20 and 21.
Additional steps might include validating the reconstruction by varying the length of the raster records used, to ensure that the reconstruction is not dependent on N p . The implementation of this algorithm was done in C and in Matlab. The pre-compiled version of this code will be made available upon publication at http://mscl.cit.nih.gov/spaj (link to PWAnetrec). Step; a simple Bayesian estimation step consisting of a single target node at a particular time instance and a subset of potential source nodes.

Mathematical Notation
O -generic symbol designating an observation, i.e. a pattern of activations N -generic symbol designating a network N c -instance of network topology (adjacency matrix). For STES it is a particular configuration by which source nodes connect to the target node ( Figure 1E) n a -number of active source nodes considered in STES �n a -average n a over all STES in a given experiment n c -number of active source nodes that connect to the target node in a particular N c being considered p l -probability that a given link l exists in IB; also a link prior p t -threshold for p l in IB; link exists if p l ≥ p t p b -uniform prior probability for link existence used in PWA p(N c ) -prior probability of a given network configuration N c p(O|N c ) -likelihood term (here called dynamics term) in Bayesian reconstruction p c D -dynamics term for a particular configuration N c Π -generic designation for posterior probabilites p Π -generic term for posterior probability for links p Π l -posterior value of p l obtained for link l within IB iteration using Equation 6 Π c (n c , n a ) -the posterior probability of a particular network configuration N c , having n c existing links out of n a possible links P(n c , n a ) -compounded the dynamics and the prior terms, p D (n c )p r (n c , n a ) N R -network obtained through reconstruction from the dynamics N NC R -network obtained by the normalized-count (NC) approach N IB R -network obtained by the Iterative Bayesian (IB) approach n max a -cut-off value for n a in IB reconstruction above which the corresponding STES is skipped N p -number of propagation steps in a time-binned raster, i.e. number of times that both of the successive time-bins contain at least one active site N STES -total number of STES used in reconstruction t -index into individual time bins, or STES t i -time of an event "i" occurred, or for a particular STES, an event that occurred on the i th node A i -amplitude of an event "i", or the one occurring on the i th node p ij -activation probability for the link i → j; a fixed probability that an active source node i will activate its target node j in the avalanche dynamics p F (...) -function that represents a-priori knowledge of the continuous time branching process dynamics, and describes how the timing of the events (in most cases only the timing differences, t j − t i ) and amplitudes of the events A affect the original activation probabilities, p ij . p F ij -predicted probability for the i → j activation, obtained using a-priori knowledge contained in p F (...) p d -mean-field approximation to our branching process dynamics, p d = �p ij , or just a uniform activation probability p c d -critical value of the probability p d w ij -weight of a directed link i → j a ij -binary indicator of the existence of a link i → j m (m BA , m OHO ) -number of edges added in growing networks, OHO and BA m 0 -initial number of nodes in growing networks K -number of the nearest neighbors (2K) that a node connects to in WN network w Π -weighting factor in PWA σ d -the branching parameter, which determines the dynamical regime Π norm -the normalization term for the posterior probability {s l } -source node index that connects to the target through the l th link W norm -normalization factor for PWA weighted measures V norm -normalization factor for a correlation measure ∆t -the width of a time window at which events were binned z,σ A -parameters determining the heterogeneity of the branching process initiation N -number of nodes in a network k d -degree of a particular node C -clustering coefficient of a network �k d -average node degree for the whole network �d -average node-to-node distance C ER -C of the randomized network using ER randomization C DSPR -C of the randomized network using DSPR randomization p S -edge density, or sparsity of network, defined as p S = #links/(N (N − 1)) p ER -probability that a link exist in Erdös-Rényi network p WN -probability that a random link exist in Watts-Newman network p ext -in our simulations, the probability that the node activation is caused by noise FC ij -frequency of successive activations between the nodes i and j WC ij -weighted measure of successive activations between the nodes i and j NC ij -scalar measure assigned to each link according to Equation 18  -reconstruction error relative to all possible links, N (N − 1) E U -the percent difference between two topologies, quantified as the total number of differences divided by the number of links that exist in either of the two (range: 0-100%) ρ I , ρ U -correlation between two architectures using the links that exist in both networks (intersect, ρ I ), or in either network (union, ρ U ) E ER D , ρ ER I , ρ ER U -measures E U , ρ I , and ρ U obtained by comparing the ER randomized networks