Dynamical and Coupling Structure of Pulse-Coupled Networks in Maximum Entropy Analysis.

Maximum entropy principle (MEP) analysis with few non-zero effective interactions successfully characterizes the distribution of dynamical states of pulse-coupled networks in many fields, e.g., in neuroscience. To better understand the underlying mechanism, we found a relation between the dynamical structure, i.e., effective interactions in MEP analysis, and the anatomical coupling structure of pulse-coupled networks and it helps to understand how a sparse coupling structure could lead to a sparse coding by effective interactions. This relation quantitatively displays how the dynamical structure is closely related to the anatomical coupling structure.

1. Introduction

Binary-state networks—each node in one sampling time bin is binary-state—arise from many research fields, e.g., gene regulatory modeling and neural dynamics [1,2,3]. Statistical distributions of network states are essential to encode information [4,5,6,7,8]. For example, with statistical distributions of network states, experimental studies show that rats can perform awake replays of remote experiences in the hippocampus [9]. The number of network states of n binary-state nodes, , exponentially grows as the node number, n, increases, which creates a challenge of characterizing the probability distribution of network states. Many works effectively characterize the distribution of network states in various systems, e.g., a network of ∼100 neurons [10], with a low-order maximum entropy principle (MEP) analysis [10,11,12,13,14,15,16,17,18]—a method with few (far less than ) non-zero effective interactions (see a precise definition in Equation (1)) constrained by low-order statistics. In the MEP analysis, the dynamical structure of the network, i.e., how nodes interact with each other in the recording of dynamical data, is characterized by effective interactions. This dynamical structure has been used to study the functional connectivity of networks [16,19]. For example, experimental studies show that the second-order effective interaction map of the retina is sparse and dominated by local overlapping effective interaction modules [19]. In this sense, those effective interactions can be regarded as a sparse coding of the information that encoded in the state distribution. Here, the sparseness is defined as the ratio of the number of non-zero effective interactions to (the number of total effective interactions). Several studies show that high-order effective interactions could be important for characterizing observed distributions of network states [6,10,20]. Although high-order effective interactions are required, the number of required non-zero high-order effective interactions are very small in these experimental studies [10]. What leads to the sparsity of effective interactions remains to be clarified.

We address how a sparse anatomical coupling structure (In the following, for simplicity, we use “coupling structure” instead of “anatomical coupling structure”.) could lead to a sparse coding by effective interactions. Network dynamical structure often closely relates to the underlying coupling structure [21]. For example, when the input of each node is independent from others, (i) high-order (≥2) effective interactions are zero in a network of no connections, and (ii) high-order effective interactions are large in a dense and strong connected excitatory network. To efficiently encode information, a realistic system often incorporates a coupling structure with certain features [22,23], e.g., sparsity, small-world, or scale-free. However, it is still unclear how the coupling structure affects the dynamical structure of effective interactions.

In this work, we consider a general class of pulse-coupled networks. The state of each node is binary-state, i.e., active when the node sends pulses to its child nodes, otherwise, silent. We first establish the connection between the coupling structure and the number of non-zero effective interactions in the full-order MEP analysis (constrained by all moments) through an observed Fact—which is independent from node dynamics. We then examine the observed Fact by numerical simulations. Through our analysis, we can estimate the number of non-zero effective interactions for a given coupling structure when the external input of each node is independent from one another. Our results show that a sparse network could lead to a lot of vanishing high-order effective interactions. For illustration, we verify our results by estimating the number of non-zero effective interactions for each order in a network with Erdos–Renyi connection structure, in which our estimated number is much smaller than , the number of all possible kth-order effective interactions. Our results establish a connection between the dynamical structure and the network coupling structure. This connection provides an insight into how a sparse coupling structure can lead to a sparse coding scheme. In this work, we would mainly use neural networks as examples for illustration, while our results apply to general binary-state networks.

2. Results

In the following analysis, we use binary vector to represent the state of n nodes within the sampling time bin labeled by l. To obtain correlations up to the mth-order requires to evaluate all , where , , and is defined by for any function and is the total number of sampling time bins in the recording. The mth-order MEP analysis is to find the desired probability distribution for n nodes by maximizing the entropy subject to correlations up to the mth-order (). To solve this optimization problem [24], one introduces Lagrangian multipliers for each constraint, that is, for the constraint of and for the constraint of = 1. The optimization problem is where and is the state of the th node in the network state V. Since the entropy S is a convex function of , the unique distribution can be solved by taking a derivative of the above objective function with respect to each as follows: where, following the terminology of statistical physics, we call a kth-order effective interaction (), the partition function Z is the normalization factor. Equation (1) is referred to as the mth-order MEP distribution. When , the constants can be determined by all constraints through a commonly used iteration method [13]. When , all order moments of the experimentally observed distribution are used as constraints; thus, the above optimization problem has only one feasible solution, that is, the experimentally observed distribution [25,26].

First, we discuss the relationship between effective interactions and the statistical distribution of network states. By taking the logarithm of both sides of Equation (1) for , we can get a set of linear equations of all-order effective interactions for all states V. Since is the same as the experimentally observed distribution [25], we can obtain the effective interactions in in terms of the experimentally observed distribution [27]. For example, , we can obtain and , where represents the probability of the network state . By applying , we have and . Our earlier study has shown a recursive structure among effective interactions, that is, the st-order effective interaction can be obtained as follows [27]: first, we switch the state of the st node in from silent to active to obtain a new term , e.g., from to ; then, we subtract from the new term to obtain , i.e.,

Without loss of generality, we randomly select two nodes labeled by 1 and 2. By the recursive relation, any kth-order effective interaction that includes node 1 and 2 can be expressed as the summation of terms with the following basic form:

Note that we use superscript b to denote “basic”. For any i and j (), the states of node l () in is silent. is a function of states of nodes ’s () by replacing the silent state of every node l () in with variable . For example, and . We can observe that if nodes 1 and 2 are independent conditioned on all other nodes, i.e., , any effective interaction containing these two nodes is zero.

Next, we would show what kind of coupling structure could entail the conditional independence of two nodes. Here, we define some notations. In any sampling time bin with state , , we denote as node i’s input from the outside of the network, denote as the input from the node i to node j, denote as the set of all child nodes of node i, denote , denote as the probability of event e, denote .

For n pulse-coupled nodes with binary-state dynamics on a network with a coupling structure where H is a possible state of nodes in

We justify our two assumptions as follows: to avoid the influence of correlation in external inputs when we are studying the relation between the dynamical structure and the coupling structure, we assume that the external input of each node is independent from one another, i.e., assumption (a). The second assumption implicates a Markov-like property; that is, for a connected pair of pulse-coupled nodes in an equilibrium state, the pulse from the parent node to the child node only depends on the state of the parent node but is independent from inputs imposed on the parent node. For example, in neural networks, a neuron sends out spikes only when this neuron is active, regardless of what inputs are imposed on the neuron.

The argument for the conclusion in Equation (4) is as follows. By assumption (a), node i and node j can be dependent only through the coupling structure . When we are considering how node i and node j affect each other by changing their states through the coupling structure , we can consider a simplified coupling structure, , which ignores those connections that are independent from states of node i and node j, i.e., and . , i.e., any other node k, its state is fixed when we are considering the conditional probability in Equation (4). By assumption (b), for node k’s any child node l, the input from node k to node l is independent from and . Thus, the connections started from those nodes in are fixed for different states of and . Therefore, is a simplified coupling structure that only keeps those connections that originated from node i and node j in . In , any connection only exists in either sub-network or sub-network . Under the condition , i.e., they neither are connected nor share any common child node, sub-network and sub-network are two isolated sub-networks. and cannot affect each other by changing their states through the coupling structure , that is, node i and j are independent conditioned on the states of all other nodes.

Figure 1 displays an example to illustrate our observed Fact. The coupling structure is shown in Figure 1a. We focus on node 1 and node 2, where they neither are connected nor share any child node. When the state of other nodes (black) are fixed, all outputs from black nodes can be ignored in the simplified coupling structure , as shown in Figure 1b. Node 1 and node 2 respectively belong to two separate sub-networks. Therefore, nodes 1 and node 2 are independent conditioned on the state of all other nodes.

Figure 1

Structure vs. simplified structure.

Based on the recursive structure of effective interactions and the observed Fact, we reach the following conclusion: with the two assumptions in the observed Fact, for a group of nodes , if there exists at least one pair of nodes that neither are connected nor share any child node, effective interaction is zero.

The system we would use to examine our conclusion is an integrate-and-fire (I&F) network, a general pulse-coupled network, with both excitatory and inhibitory nodes [21]. For the ith node, the dynamics of its state variable with time scales are governed by where and are the reversal values of excitation (ex) and inhibition (in), respectively. is the background input with magnitude f and time scale , is a Poisson process with rate , is the Heaviside function, is the excitatory input from other jth excitatory nodes, and is the inhibitory input from other jth inhibitory nodes. The jth excitatory (inhibitory) node evolves continuously according to Equation (5) until it reaches a firing threshold . That moment in time is referred to as a firing event (say, the kth spike) and denoted by ( Then, is reset to the reset value () and held for an absolute refractory period of Each spike emerging from the jth excitatory (inhibitory) node causes an instantaneous increase () in (), where and are the excitatory and inhibitory coupling strengths, respectively. The model (5) describes a general class of physical networks [1,3,21,28,29]. To be intuitive, the I&F model described in Equation (5) can be understood through a resistor–capacitor circuit. Each neuron is a leaky circuit, which consists of a capacitor with dimensionless capacitance as 1 in parallel with two resistors. is the voltage. on the left-hand side is the current which passes through the capacitor. On the right-hand side, is the leaky current. The first resistor has a reversal potential with conductance . Therefore, the second term on the right-hand side obtained by the Ohm’s law is the current that passes through the first resistor. Note that the conductance is affected by inputs. It is similar for the third term on the right-hand side, i.e., the current that passes through the second resistor.

In the first example, two excitatory and two inhibitory I&F nodes form a ring coupling structure (Figure 2a). For any pair of nodes, say, node i and j, we compute , where H is a state vector of other two nodes. By our observed Fact, the conditional independent pairs are and , and other pairs are categorized as dependent pairs. In Figure 2b, the strengths of of independent pairs (green) are almost two orders of magnitude smaller than those of dependent pairs (red). We then shuffle spike trains of each node. We similarly compute for 10 different pieces of shuffled data. Blue dots and cyan dots in Figure 2b are results of all shuffled data of dependent pairs and independent pairs, respectively. The strength of of independent pairs (green)—computed from the observed data—are within the statistical error of shuffled data. We then solve effective interactions in the full-order MEP analysis for this ring network. As shown in Figure 2c, the effective interaction strengths of independent pairs ( and are within the statistical error of shuffled results (red). Since every high-order (≥3) effective interaction includes at least one independent pair of nodes, as predicted, the strengths of all high-order effective interactions are within the statistical error of shuffled results as shown in Figure 2d.

Figure 2

Anatomical structure vs. effective interactions of integrate-and-fire networks. Each row shows a numerical case. In the first column, black arrows and red arrows represent excitatory and inhibitory connections, respectively. In the second column, red and green dots are the strengths of of dependent and independent pairs, respectively. Blue dots and cyan dots are the strengths of of dependent and independent pairs from ten shuffled spike trains, respectively. Each dot is for one . The third and fourth columns display absolute effective interaction strengths (blue bars). The corresponding node indexes for each effective interaction are shown in the abscissa. The mean and standard deviation of absolute strengths of each effective interaction of ten shuffled spike trains are also displayed by garnet bars. The simulation time for each network is . The time bin size for analysis is [12,13]. Independent Poisson inputs for each network are and . The firing rate of each node is about . Parameters are chosen [28] as , , , , , , , and , .

In the second example in the second row in Figure 2, results are similar in that dependent pairs and independent pairs can be identified through our observed Fact, and the strength of any effective interaction that includes the independent pair of nodes (node 1 and node 3) is within the statistical error of shuffled data. In this example, as shown in Figure 2h, is very small, i.e., within the statistical error of shuffled results. However, does not include the independent pair of nodes (node 1 and node 3); thus, the theoretical estimation of zero-strength effective interactions misses . This example indicates that our theoretical estimation result may give rise to an upper bound of the number of non-zero effective interactions. In contrast, for a network of all excitatory nodes with the same coupling structure as the one in Figure 1e, is significantly larger than zero (not shown). Since the estimation of the strength of high-order effective interactions involves the estimation of high dimensional joint probability distribution, a very long recording due to the curse of dimensionality constrains us from examining for a large network.

Based on the relation between the coupling structure and effective interactions, the number of non-zero high-order effective interactions can be small in a sparsely connected network compared with , which is the number of all possible kth-order interactions. For example, we estimate the number of each-order non-zero effective interactions in a network with an Erdos–Renyi connection structure. We generate 1000 Erdos-Renyi random networks of 100 nodes (the same connection probability but different random samples). In each network, for any pair of node i and node j, we assign the connection from node i to node j by the following rule: we generate a random number from the uniform distribution on ; if the number is smaller than , we assign a connection from node i to node j. As shown in Figure 3, the number of non-zero kth-order () effective interactions is much smaller than (too large to be shown). The number of high-order effective interactions (order higher than 11th) almost vanishes (order higher than 20th not shown). In this example, the sparseness, i.e., the ratio of the number of non-zero effective interactions to , is less than .

Figure 3

Non-zero effective interactions in the Erdos-Renyi random networks. We generate 1000 Erdos-Renyi random networks of 100 nodes (the same connection probability but different random samples). The connection probability between two nodes is . The number of non-zero effective interaction is plotted against effective interaction order. The mean and standard deviation are respectively shown by the black line and shaded area.

3. Conclusions and Discussion

In summary, we have established a relation between effective interactions in MEP analysis and the anatomical coupling structure (“coupling structure” for simplicity in the following) of pulse-coupled networks to understand how a sparse coupling structure could lead to a sparse coding by effective interactions. For example, the sparseness of the case in Figure 3 is less than . Since effective interactions characterize how nodes interact with each other in the recording of dynamical data, i.e., network dynamical structure, our study quantitatively displays how the dynamical structure closely relates to the coupling structure.

Even though high-order effective interactions are often much smaller compared with low-order ones [27], it is still unclear why small high-order effective interactions do not accumulate to have a significant effect in a large network [10,30]. For example, it has been shown that MEP distribution with sparse low-order effective interactions—non-zero effective interactions are sparse and vanish when the order is higher than the eighth-order—can well capture the state distribution of 99 ganglion cells in the salamander retina responding to a natural movie clip or natural pixel [10]. In this study, we show that a large amount of effective interactions vanish in a sparse coupling structure, thus rationalizing the absence of the accumulation of high-order interactions for a large network.

Finally, we point out that some important issues remain to be elucidated in the future. First, we have ignored correlations in external inputs when estimating the number of non-zero effective interactions. Correlated inputs can induce non-zero high-order effective interactions [31]. It still needs to be considered how the statistics of inputs affect the sparsity of effective interactions. Second, current algorithms for estimating non-zero effective interactions (not limited to the second-order) for a large network (e.g., ∼100 nodes) are very slow, e.g., Monte Carlo based methods [30,32]. Our undergoing work is developing a fast algorithm that exploits the sparsity of effective interactions. We have seen an indication that the algorithm can work well for an I&F network with sparse coupling structure; however, that work has yet to be fully verified to be conclusive.