Abstract

Because sensory inputs typically have spike correlation, and local circuits have dense feedback connections, input spikes cause the propagation of spike correlation in lateral circuits; however, it is largely unknown how this secondary correlation generated by lateral circuits influences learning processes through STDP, or whether it is beneficial to achieve efficient spike-based learning from uncertain stimuli. To explore the answers to these questions, we construct models of feedforward networks with lateral inhibitory circuits and study how propagated correlation influences STDP learning, and what kind of learning algorithm such circuits achieve. We derive analytical conditions at which neurons detect minor signals with STDP, and show that depending on the origin of the noise, different correlation timescales are useful for learning. In particular, we show that non-precise spike correlation is beneficial for learning in the presence of crosstalk noise. We also show that by considering excitatory and inhibitory STDP at lateral connections, the circuit can acquire a lateral structure optimal for signal detection. In addition, we demonstrate that the model performs blind source separation in a manner similar to the sequential sampling approximation of the Bayesian independent component analysis algorithm. Our results provide a basic understanding of STDP learning in feedback circuits by integrating analyses from both dynamical systems and information theory.

Author Summary

For example, humans easily detect the mention of their names from across a noisy room, a phenomenon known as the cocktail party effect. Spike-timing-dependent plasticity (STDP) is a learning mechanism ubiquitously observed in the brain across various species and is considered to be the neural basis of such learning; however, it is still unclear how STDP enables efficient learning from uncertain stimuli and Whether spike-based learning offers benefits beyond those provided by standard machine learning methods for signal decomposition. To begin to answer these questions, we conducted analytical and simulation studies examining the propagation of spike correlation in feedback neural circuits. We show that non-precise spike correlation is useful for handling noise during the learning process. Our results also suggest that neural circuits make use of stochastic membrane dynamics to approximate computationally compleX Bayesian learning algorithms, progressing our understanding of the principles of stochastic computation by the brain.

Introduction

Despite the diversity and variability of input spike trains, neurons can learn and represent specific information during developmental processes and according to specific task requirements. Spike-timing-dependent plasticity (STDP) [1,2] is a candidate mechanism of neural learning. Extensive studies have revealed the type of information that a single neuron can learn through STDP [3—7]; however, the type of information that a population of neurons interacting with each other learns through STDP has not yet been determined. Understanding this extension from a single neuron to a population of neurons is crucial because a single neuron learns and represents only a limited amount of information that may be transmitted to it from thousands of inputs.

Previous theoretical results showed that neural circuits with lateral inhibition enhance signal detection [12,13] and improve learning performance [14—16]. Several simulation studies further revealed that neurons acquire receptive field [17—19] or spike patterns [20] through STDP by introducing lateral inhibition; yet, those studies were limited to simplified cases for which a large population of independent neurons was suggested to be sufficient [5,21,22]. Therefore, it remains unclear whether lateral inhibition plays a crucial role in STDP learning; in particular, the spike level effects of lateral inhibition remain elusive. Moreover, recent experimental results suggest that animals learn and discriminate mixed olfactory signals [23—25] or auditory signals masked by noise [26,27] , but it is still unknown how feedback interactions contribute to such learning.

We analyzed the propagation of spike correlations through inhibitory circuits, and revealed how such secondary correlations influence STDP learning at both feedforward and feedback connections. We discovered that the timescale of spike correlation preferable for learning depends on whether the noise is independent from any signal (random noise) or generated from the mixing of signals (crosstalk noise). We also found that excitatory and inhibitory STDP cooperatively shapes lateral circuit structure, making it suitable for signal detection. We further found a possible link between stochastic membrane dynamics and sampling process, which is necessary for neural approximation of learning algorithm of Bayesian independent component analysis (ICA). We applied our findings by demonstrating that STDP implements a spike-based solution in neural circuits for the cocktail party problem [26,28,29].

Results

Model

The external source layer represents the external environment or neural activity at sensory systems. The external layer also provides common inputs to the input layer and induces correlations in the neurons in the input layer. The input layer shows rate-modulated Poisson firing based on events at the external layer and external noise, which is approximated with the constant firing rate {no}. Subsequently, spike activity at the input layer projects to the output layer, which also receives inhibitory feedback from the lateral layer. Neurons in the lateral layers are excited by inputs from the output layer. We assumed that all neurons in the input layer and the output layer are excitatory, whereas lateral-layer neurons are assumed to be inhibitory. Although excitatory lateral interactions also exist in the sensory cortex, they are typically sparse [30] and weak [10] compared with inhibitory interactions; thus we concentrated on the latter. For the analytical treatment, the neurons in the output and lateral layers were modeled with a linear Poisson model. We first studied synaptic plasticity at the feed-forward connections (connections from the input layer to the output layer), while fixing lateral connections (i.e., connections from the output layer to the lateral layer and connections from the lateral layer to the output layer). For STDP, we used pairwise log-STDP (Fig 1B) [31], which replicates the experimentally observed long-tailed synaptic weight distribution [32,33]. We considered the case for information encoded in the correlated activity of input neurons [34,35], and fixed the average firing rate of all input neurons at the constant value UOX (See Table 1 and 2 for the list of variables and parameters). If the firing rate of input neuron i is of the neuron qiy, then common inputs from the external layer induce a temporal correlation proportional to Table 1. Definition of variables. yj(t) The spiking activity of output neuron j uk’(t) Membrane potential of inhibitory neuron k zk(t) The spiking activity of inhibitory neuron k wjix The synaptic weight of a feed-forward excitatory connection from ito j qifl Response probability of input neuron ito external source p 11X, 12X The correlation kernel functions used the gamma distribution With shape parameter kg = 3 in order to reproduce broad spike correlations typically observed in cortical neurons [36,37]. Synaptic weight dynamics by STDP is written as ods for details), the weight change of the feedforward connection WX can be approximated as Table 2. Parameter settings. Simulation time Neural population Neural subpopulation EPSP/IPSP time constants Synaptic weights Axonal delays Dendritic delays Synaptic (axonal) delays Correlation timescale Firing rates Learning rate Noise amplitude of plasticity STDP time windows Parameters for log-STDP Initial variance of synaptic weights LTD/LTP balance Where glx and gzx are scalar coefficients, C is the correlation matrix, and E is the identity matriX (see Eqs (25)—(30) for derivation). The first term describes the synaptic weight change directly caused by an input spike correlation and can be rewritten into the convolution of the temporal correlation and correlation kernel function XXI as where D = ZdXd+dY+dZ, and cy and 52 are EPSP/IPSP curves of output/ inhibitory neurons, respectively. This term primarily causes LTD as the sign flips through lateral inhibition

Previous simulation studies showed lateral inhibition has critical effects on excitatory STDP learning [17—19]; however, it has not yet been well studied how a secondary correlation generated through the lateral circuits influences STDP at feedforward connections, and it is still largely unknown how lateral inhibition functions with various stimuli in different neural circuits. For example, the correlation kernel of the feedback term eXhibits a delay as the signal propagates through the inhibitory circuit; yet, we do not know how much delay is permitted for effective learning or if realistic synaptic delays satisfy such a condition. Furthermore, it is also unknown what information a circuit can learn if there are several mixed signals with different amplitudes for which symmetry-breaking learning [5,39] is not valid. Therefore, using theoretical analysis and simulation, we first investigated the properties of the inhibitory kernel — Xé‘ (T; w) in STDP learning.

Lateral inhibition enhances minor source detection by STDP

On the other hand, if lateral inhibition is effective, different output neurons may acquire various components of the external structure. We first examined that point in a simple network model with two independent external sources (Fig 2A). In the model, each external source drives an independent subgroup of input neurons (we defined those input neurons as A-neurons and B-neurons), which project excitatory inputs to all of the output neurons. Here, we assume that source A drives input neurons with a higher probability than source

We refer to this as the minor source detection task below. Here, for lateral connections, we assumed that both excitatory-to-inhibitory (E-to-I) and inhibitory-to-excitatory (I-to-E) connections are well organized such that inhibition only works mutually between two output neuron groups (Fig 2A; blue lines are E-to-I and red lines are I-to-E connections. See also Eq (30) in Methods). The origin of these structured lateral connections is discussed later. When the network is excited by inputs from external sources, excitatory postsynaptic potential (EPSP) sizes of feedforward connections WX change according to STDP rules. Initially, in all output neurons, synaptic weights from A-neurons (blue triangles in Fig 2A) become larger because A-neurons are more strongly correlated with one another than B-neurons are. However, as learning proceeds, one of the output neuron groups becomes selective for the minor source B (Fig 2B). After 30 min, the network successfully learns both sources. If we focus on the peristimulus time histogram (PSTH) for the average membrane potential of output neurons aligned to external events, both neuron groups initially show weak responses to both correlation events, and yet the depolarization is relatively higher for source A than for source B (Fig 2C left). After 10 min of learning, both neuron groups show relatively stronger initial responses for source A, but group 1 shows a hyperpolarization soon after the initial response (Fig 2C middle). As a result, synaptic weights from A-neurons to group 1 become weaker, and group 1 neurons eventually become selective for the minor source B (Fig 2C right). The mean cross-correlation (see cross-correlation in Methods for details) between the external sources and the population activity of output neurons is maximized when the delay is approximately 10—15 ms (Fig 2E). If we fix the delay at 14 ms, then the cross-correlation gradually increases as the network learns both sources (Fig 2D). The same argument holds if mutual information is used for performance evaluation (green lines in Fig 2D and 2E). Interestingly, the network better detects the minor source when it is learned with a highly correlated source compared with when it is learned with another minor source (Fig 2F), because a highly correlated opponent source causes strong lateral inhibition on the output neurons, which enhances minor source learning. Similar results are also obtained for conductance-based leaky integrate-and-fire (LIF) neurons (81 Fig).

Lateral inhibition should be strong, fast, and sharp

Because both output excitatory neurons and lateral inhibitory neurons are bundled into groups, in the mean-field approximation, we can approximate M excitatory populations and N inhibitory populations into two representative output neurons and two inhibitory neurons. Similarly, input neurons can be bundled into three groups (A-neurons, B-neurons, and Background-neurons). In addition, we assumed that the synaptic connections from Background-neurons to output neurons are fixed because they showed little weight change in the simulation (orange lines in Fig 2B). In this approximation, by inserting Eq (32) into Eq (29), the mean synaptic weight changes of feedforward connections follow

dev Ma Ma d: g ZLanfv’ VEGA?

Z qquv’p — NaWZMa WYZL“ Wig", vi Z qquV/p v,=1 p v’=1 P

G IX and G2X are coefficients determined by synaptic delays, EPSP/IPSP (Inhibitory postsynaptic potential) shapes, and correlation structure, as shown in Eqs (3) and (4). By solving the self-consistency condition (Eq (34) in Methods), the firing rates of inhibitory neurons are approximated as

As a result, we found that when the mutual inhibition is weak (WI 2 10), the system has only one stable point at which w M is larger than M B (Fig 3A left). At this point, w2A is also larger than W23 (w2A = 9.64, W23 2 3.60; not shown in the figure), which means that both output neuron groups are specialized for the major source A (we call this state a winner-take-all state or T-state); however, if the inhibition is moderately strong (WI 2 21.5), two new stable fixed points and two unstable fixed points appear in the system (Fig 3A middle). In the stable point on the left, neuron group 1 picks up source B while neuron group 2 picks up source A (WZA = 12.52, W23 2 2.87). On the right-hand side, neuron group 1 selects source A while neuron group 2 selects source B (we denote those two states as winners-share-all states or S-states below). At the stable point in the middle, both groups detect source A (W1 A 2 WM 2 9.47, WI 3 2 W23 2 3.61). Note that because of the mutual inhibition, the synaptic weight from A-neuron is smaller when both groups learn A than it is when only group 1 learns A. For strong inhibition (W1 2 40.0), the stable point in the middle disappears, and the system is stable only when two neuron groups detect different sources (Fig 3A right). Simulation results confirm this analysis because strong inhibition indeed causes a winner-share-all state in which multiple neuron groups survive in competition [15] , whereas the network tends to show a winner-take-all learning when the inhibition is weak (Fig 3B). We measured the degree of winner-share-all/winner-take-all states by defining the specialization indeX wSI as

If two output groups are specialized for different sources, W31 becomes positive, whereas if two groups are specialized for the same source, W31 becomes negative. When the synaptic delay in the lateral connections is small, only S-states are stable, whereas at longer delays, both S-states and T-states are stable. In the simulation, the network typically grows toward the latter state in the bistable strategy (Fig 3C). Moreover, if we change the shape of the IPSP curve while keeping T213 2 5 TZA, for steep IPSP curves (i.e., both TZA and T23 are small), only the S-states are stable, whereas T-states also become stable for slower IPSPs (Fig 3D). Therefore, both analytical and simulation studies indicate that lateral inhibition should be strong, fast and sharp to detect higher correlation structure. Moreover, lateral inhibition does not need to be pathologically strong because the I/ E balance of Na wZ/waf 2 20% is sufficient to cause multistability.

Optimal correlation timescale changes depend on the noise source

Therefore, we next considered the effect of correlation timescales, especially on noise tolerance. In our current model, input neurons respond to external sources with input kernel qb(t) = l‘ge‘wt /20t3 (Fig 4A left), and so the correlation between input neuron i and l is given as

The correlation is precise When fit is small, Whereas it becomes broad at large values of 6t (Fig 4A right, Fig 4B). Because STDP causes homeostatic plasticity that does not depend on a correlation, as shown in the third term of Eq (5), in a more precise approximation, Eq (2) should be written as

Both g1X and g2X become smaller for a larger fit, but decreases in g2X are slower than those in gIX, and, as a result, K = g2X/g1X becomes larger for a longer correlation timescale (Fig 4C). This means that a longer temporal correlation is more suitable for the detection of multi-components. This is indeed confirmed in the simulation (Fig 4D). When 0t 2 0.5 and the minor component is slightly weaker than the major one (CA 2 0.36, C3 = 0.25), the minor component is no longer detectable. On the other hand, at 0t 2 2.0, the minor component is detectable even if the strength of the induced correlation is less than half

In the model above, the noise is provided through the spontaneous Poisson firing of input neurons as random noise (Fig 4E top, black dots are spikes caused by random noise). In reality, however, there would be crosstalk noise among input spike trains caused by the interference of external sources. We implemented this crosstalk noise by introducing non-diagonal components in the response probability matriX as where qs is the response probability to the preferred signal and qN is that to the non-preferred signal (Fig 4E bottom). We refer to this as the noisy source detection task below. To make a clear comparison, in the simulation of random noise, we kept qN = 0 and changed the spontaneous firing rate of the input neurons (no) to modify the noise intensity, whereas in simulation of crosstalk noise we removed random noise (i.e., no 2 0) and changed qN. For random noise, a smaller 6t enables better learning because a large g1X competes with the homeostatic force (Fig 4F). By contrast, for crosstalk noise, the performance is better at 6t 2 2.0 than at 6t 2 0.5 because strong lateral inhibition suppresses crosstalk noise (Fig 4G). Although for small noise regimens, the network performs better at 6t 2 0.5 than at 6t 2 2.0, but the difference is almost negligible. Therefore, to cope with crosstalk noise, the spike correlation needs to be broad, whereas a narrow spike correlation is better for random noise. We note that qualitatively the same arguments as above also hold for the eXponential kernel qbe(t) = e‘t/Qt /0t (83D and 83E Fig). However, the ratio of two coefficients (i.e., K6 2 gegx/ ge 1X) is typically smaller for this kernel than for the kernel we used throughout this study (83B and 83C Fig vs. Fig 4D) because lateral inhibition is less effective due to highly peaked spike correlation (83A Fig).

Excitatory and inhibitory STDP cooperatively shape structured lateral connections

This means that the network somehow knows a priori that the number of external sources is two; however, in reality, a randomly connected network should also learn such information. To test this idea, we introduced STDP-type synaptic plasticity in lateral excitatory connections and feedback inhibitory connections and investigated how different STDP rules cause different structures in the circuit.

For comparison, we also considered a model with random lateral connections in which all output neurons and inhibitory neurons are randomly connected with probability 0.5 (Fig 5A middle). When lateral connections are random, mean-field equations are modified as

In other cases, neurons are often organized into two groups with different population sizes. In such cases, for evaluating performance, we measured average weights from source A on the output neurons receiving stronger inputs from A-neurons than from B-neurons or Background-neurons. For randomly connected lateral inhibition, learning performance dropped significantly in noisy source detection (Fig 5B) and in minor source detection (Fig 5C); thus clustered connectivity is indeed advantageous for learning.

We first introduced Hebbian STDP for both E-to-I and I-to-E connections. With these learning rules, the lateral connections successfully learn a mutual inhibition structure (Fig 5D); however, this learning is achievable only when the learning of a hidden external structure is possible from the random lateral connections (magenta lines in Fig 5B and 5C; note that orange points are hidden by magenta points because they show similar behaviors in noisy cases), which means either when crosstalk noise is low or two sources have similar amplitudes. Nevertheless, once a structure is obtained in easy settings (qN = 0 or q A = qB), that network outperforms the network with random lateral connections in both noisy source detection (Fig 5E) and minor source detection (Fig 5F). In Fig 5E, we evaluated the performance of noisy source detection by first conducting STDP learning at qN = O, and then we terminated STDP and performed simulations at the various noise levels qN. Similarly, in the minor source detection task depicted in Fig 5F, we first performed STDP learning with q A = qB = 0.6, and then evaluated the performance for a smaller q B. STDP can also generate similar lateral connection structures When the total number of input sources is larger than two (82A and 82B Fig). Therefore, STDP at lateral connections helps signal detection by efficiently organizing the connection structure. We next studied the analytical conditions for learning of the clustered structure (see Analytic approach for STDP in lateral and inhibitory connections in Methods for details). The synaptic weight dynamics of lateral excitatory and inhibitory connections are approximately given as

From a linear analysis, we can expect that When gY1 is positive, E-to-I connections tend to be feature selective (see Eq (35) in Methods). Each inhibitory neuron receives stronger inputs from one of the output neuron groups and, as a result, shows a higher firing rate for the corresponding external signal. On the other hand, if g2 1 is positive, I-to-E connections are organized in reciprocal form, Where one of the reciprocal connections is enhanced and the other is suppressed (see Eq (36) in Methods). We can evaluate feature selectivity of inhibitory neurons by where QYA and QYB are the sets of excitatory neurons responding preferentially to sources A and B, respectively. Indeed, when the LTD time window is narrow, analytically calculated ng tends to take negative values (the green line in Fig 6A), and E-to-I connections organized in the simulation are not feature selective (the blue points in Fig 6A). By contrast, for a long LTD time window (i.e., when LTD is weakly spike-timing dependent), gY1 tends to take positive values, and E-to-I connections become clustered. In the simulation, WZ is also plastic, but as shown in Eq (10), the effect of WZ on the plasticity of WY is negligible in first-order approximations. Similarly, for I-to-E connections, we measure the degree of mutual inhibition (non-reciprocity) with

Note that the organized neuronal wiring patterns are not a pure product of the pre-post causality of STDP but the effect of spike correlations propagating through lateral inhibitory circuits. If the structural plasticity is merely caused by the pre-post causality, both ¢Y and (p2 can decrease with increases in the inhibitory population while maintaining the total synaptic weights because the causal effect becomes weaker as each synaptic weight becomes smaller [45]; however, in our simulations, the values of both quantities generally increased for larger inhibitory populations (82C Fig). Hebbian inhibitory STDP at lateral connections is not always beneficial for learning. For eXample, in minor source detection, if we use Hebbian inhibitory STDP, a slightly minor source is not detectable, whereas for anti-Hebbian STDP, a small number of neurons still detect the minor source because reciprocal connections from strong-source responsive inhibitory neurons to strong-source responsive output neurons inhibit synaptic weight development for the stronger source (Fig 6C).

Neural Bayesian ICA and blind source separation

To confirm this mechanism is indeed effective in realistic tasks, we applied the above method to blind source separation. We first examined the condition in which the network could capture external sources. We extended the previous network to include four independent sources mixed at the input layer (Fig 7A). In the present application, we used structured lateral connections because learning for clustered structures is difficult with noisy stimuli, as shown in the preceding section. The response probability matrix Q and correlation matrix C are given as namics follows WX % g‘IX WX C, we may expect that synaptic weight vectors converge to the ei-genvectors of the principal components; however, this was not the case in our simulations, even if we took into account the non-negativity of synaptic weights (see Fig 7B, where we renormalized the principal vectors to the region between 0 and 1). Instead, each weight vector converged to a column of the response probability matriX Q (Fig 7B, the left panel is the projection to the first two dimensions, and the right panel is the projection to the other two dimensions). This result implies that the network can extract independent sources, rather than principal components, from multiple intermixed inputs.

To this end, we compared the performance of the model with that of the Bayesian ICA algorithm, in which independence of external sources is treated as a prior [46,47]. In the algorithm, the learned mixing matrix may converge to its Bayesian optimal value estimated from a stream of inputs. Although we cannot directly argue the optimality of cross-correlations, if the mixing matrix is accurately estimated, external activity is also well inferred, and thus we can use the mean cross-correlation as a measure for the optimality of learning. In terms of discretized input activity X, the external source activity 8 and prior information I, we can express the conditional probability of the estimated response probability posterior P[Q|X, I] still depends on a prior given for S. If we introduce a prior that each external source follows an independent Bernoulli Process T/At L function is given as, Where

pf = 1 — <1 — rsAofi [1 — am: ask/szy] ,qsk = 2; 3 M + 1/2>At12exp[—<k + mam/at].

We approximated this Bayesian ICA algorithm by a sequential sampling source activity instead of calculating the integral over all possible combinations in the estimation of the log-pos-terior of the response probability matriX Q. In this approximation, the learning rule of the estimated response probability matriX Q obeys where Y is the sampled sequence, and pik(Y1‘k'1) is the sample based approximation of pik in the previous equation. This rule has spike-timing and weight dependence similar to those seen in STDP (Fig 7D). Although the performance of STDP is much worse than the ideal case (when the true Q is given), this performance is similar to that for the sample-based learning algorithm discussed above (Fig 7C). Therefore, the network detects independent sources if crosstalk noise is not large. We further studied the response of the models for the same inputs and 1 found that the logarithm of the average membrane potential ufi = W Z well approxi-,Ll jEQP‘ mates the log-posterior estimated in Bayesian ICA, even in the absence of a stimulus (Fig 7E). This result suggests that in the STDP model, eXpected external states are naturally sampled through membrane dynamics that are generated through the interplay of feedforward and feedback inputs.

We created “sensory” inputs by miXing four artificially created auditory sequences (Fig 8A and SI Auditory File). In the auditory corteX, various frequency components of a sound, particularly high-frequency components, are represented by specific neurons typically organized in a tonotopic map structure [48], whereas low-frequency components are eXpected to be perceived as a change in sound pressure. Furthermore, populations of neurons in the primary auditory corteX are known to synchronize the relative timing of their spikes during auditory stimuli and provide correlated spike inputs for higher cortical areas in which the auditory scene is fully analyzed and perceived [49,50]. We modeled these features by assuming that input neurons have a preferred frequency {fi} defined as and auditory inputs are provided as time-dependent response probabilities, Which follow

In this representation, each sound source is represented by correlated spikes of neural populations (right panel of Fig 8C). Even if signals have overlapping frequency components {651110)}? blind separation is possible as long as {aql(t)}q are independent and have sharp rising profiles sufficient to cause spike correlations. After learning, four output neuron groups successfully detected changes in the sound pressure of the four original auditory signals (colored lines in @ 8B) by correctly identifying the input neurons that encoded the signals. Therefore, STDP rules implemented in a feedforward neural network With lateral inhibition serve as a spike-based solution to the blind source separation or cocktail party effect problem.

Discussion

We showed that a population of neurons could learn multiple signals with different strengths or mixed levels. In addition, we found that to perform learning from signals corrupted with random noise, the timescale of the input correlations needed to be in the range of milliseconds, whereas the timescale was broader for crosstalk noise, which may eXplain why the spike correlation of cortical neurons often exhibits a large jitter (approximately 10 ms) [36,37]. We also investigated the functional roles of STDP at lateral excitatory and inhibitory connections to demonstrate that

Our results also suggested that anti-Hebbian plasticity was helpful for learning from minor sources and implied that different STDP rules at lateral connections induced different algorithms at feedforward connections. Furthermore, we derived an STDP-like online learning rule by considering an approximation of Bayesian ICA with sequence sampling. This result suggested that lateral inhibition adjusted the membrane potentials of postsynaptic neurons so that their spiking processes accurately performed sequence sampling. We also demonstrated that this mechanism was applicable to blind source separation of auditory signals.

Noise characteristics and correlation timescales

Neurons in the lateral geniculate nucleus show strong spike correlations [42,51], while correlations in V1 [36,52] or higher Visual areas [37] are less precise. Our results indicate the interesting possibility that these differences may reflect the different characteristics of the noise with which the various cortical areas need to contend. At an early stage of sensory processing, the major noise component may be environmentally produced background noise from various sources; thus precise spike correlation is beneficial at this stage for noise reduction during signal detection and learning (Fig 4G). By contrast, in higher sensory cortices, crosstalk noise accumulated through signal propagation in circuits may form the primary noise source, so less precise spike correlation is preferable (Fig 4H). It would be intriguing to examine whether lower and higher cortical areas similarly change the strength of spike correlations for other sensory modalities.

STDP in E-to-I and I-to-E connections

We showed that in a feedback circuit, Hebbian inhibitory STDP preferred winner-take-all while anti-Hebbian inhibitory STDP tended to cause winner-share-all (see Fukai and Tanaka 1997 for winner-share-all) at eXcitatory neurons (Fig 6D). This result indicates that different inhibitory STDP imposes different functions for excitatory STDP, which suggests that a neural circuit may select optimal inhibitory STDP for a specific purpose or strategy of learning, and this may differ across regions and be modified by neuromodulators. A recent study showed that inhibitory plasticity even directly influences the plasticity at excitatory synapses of the postsynaptic neuron [62]. In such cases, algorithm selection would play a more important role than it did for the standard STDP implemented in our model.

In addition, the authors suggested the candidate mechanism for this change in activity is STDP at E-to-I synapses. Our results examining E-to-I STDP confirmed this configuration of inhibitory cells modulated by plasticity at feedforward excitatory connections (Fig 5D, 82A and 82B Fig). In our model, although inhibitory neurons are not directly projected from input sources, as excitatory neurons learn a specific input source (Fig 5D, left panel), inhibitory neurons acquire feature selectivity through Hebbian STDP at synaptic connections from those excitatory neurons (Fig 5D, middle panel). Furthermore, our results indicate an important function of these feature-selective inhibitory neurons. Once an adequate circuit structure is learned and inhibitory connections are organized into a feature-selective pattern, even if the input to the network becomes noisy or faint, the network can still robustly detect signals (Fig 5E and SF). This robustness would be useful for spatial learning, as contextual information is often uncertain.

STDP and Bayesian ICA

First, output neurons were able to detect hidden external sources, without capturing principal components (Fig 7B). Previous results suggest that for a single output neuron, an additional homeostatic competition mechanism is necessary to detect an independent component [7,22]. In addition, when information is coded by firing rate, homeostatic plasticity is critically important, because STDP itself does not mimic Bienenstock-Cooper-Munro learning [18]. However in our model, information was encoded by correlation, and mutual inhibition naturally induced intercellular competition so that intracellular competition through homeostatic plasticity was unnecessary. Moreover, our analytical results suggested the reason that independent sources are detected. To perform a principal components analysis using neural units, the synaptic weight change needs to follow where LT[] means lower triangle matrix [@675]. This LT transformation protects principal components caused by the lateral modification from higher order components; however in our model, because all output neurons receive the same number of inhibitory inputs Eq (2), all neurons are decorrelated with one another and develop into independent components.

In the model used by those authors, the synaptic weight matrix is treated as a hyper parameter and estimated by considering the maximum likelihood estimation of input spike trains. By contrast, in the Bayesian ICA framework, the mixing matrix (corresponding to synaptic weight matrix) is treated as a probabilistic variable. Using this framework, we needed to calculate an integral over all possible source activities in the past to derive stochastic gradient descendent; however, as shown in Fig 7C, the stochastic learning was well performed by employing an approximation with sequential sampling. Moreover, we naturally derived an adequate LTP time window from the response kernel of input neurons to external events (Fig 7D). We also found that STDP self-organized a lateral circuit structure that performed better than a random global inhibition (Fig 5E and SF). Mathematically, to perform sampling from a probabilistic distribution, we first needed to calculate the occurrence probability of each state; however, in a neural model, membrane potentials of output neurons approximately represent the occurrence probability through membrane dynamics. In machine learning methods, integration over possible source activities is often approximated using Markov chain Monte Carlo (MCMC) sampling [68]. Interestingly, a recent study showed that a recurrent network performed MCMC sampling [69,70] , suggesting that our network may perform a more accurate sampling in the presence of recurrent excitatory connections.

Suboptimality of STDP

Previous theoretical results suggest that STDP can modulate synaptic weights in a way that optimizes information transmission between pre and postsynaptic neurons [71,72]. In the Bayesian ICA framework, blind source separation can be formulized as an optimization problem, but, in this case, the problem itself is ill-defined because optimality does not guarantee the true solution. In addition, local minima are often unavoidable for online learning rules. Nevertheless, the problems faced by the brain are often ill-defined, and suboptimality is inevitable [73]. Because we performed both nonlinear dynamics-based and machine learning-based analyses, we can offer some insights regarding the origins of local minima in stochastic gradient descendent learning. In the initial state, synaptic weights are typically homogeneously distributed, and this state is often locally stable. As a result, the homogeneous stable point is more likely to be selected in learning (Fig 2C and 2D) than the non-homogenous, more desirable, points; however, introducing additional noise may change this situation. Indeed, in Fig 4B and Fig 7C, the performance of the model was improved by adding a small amount of noise to input activities, although the improvement was not significant; however, because a large amount of noise is harmful for computations and stable learning, the benefit of noise addition is highly limited, and the brain may recruit other mechanisms for near optimal learning.

Neural mechanism of blind source separation

The mechanism underlying the cocktail party effect remains elusive [26,28,29] , although several solutions have been proposed [74,75]. An effective solution for this problem is ICA [76—78], and the neural implementation of the algorithm has been studied by several authors [14,18,79,80]. Our study extended these results through a rigorous analytical treatment on biologically plausible STDP learning of spiking neurons, and our analyses enabled us to discover interesting functions of correlation coding. Moreover, by explicitly modeling inhibitory neurons, we found that STDP at E-to-I and I-to-E connections cooperatively organized a lateral structure suitable for blind source separation. In addition, we successfully extended a previous model for the formation of static visual receptive fields [18,19] to a more dynamic model in an auditory blind source separation task. In realistic auditory scene analysis, the frequency spectrum of acoustic signals is first analyzed in the cochlea, where each frequency component is the mixture of sound components from independent sources. Components belonging to the same source may be separated and integrated by downstream auditory neurons for the perception of the original signal. These frequency components can be considered a mixed signal in the ICA problem [81]; thus even if signals are mixed in frequency space, if the amplitudes of the signals are temporally independent, blind separation is still achievable. In the neural implementation of the problem, if two frequencies are commonly activated in the same signal, neurons representing those frequencies show spike correlation under the presence of the signal; thus the learning process is naturally achieved by STDP learning. These results indicate an active role of spike correlation and STDP in efficient biological learning.

Methods

Model

Based on the previous study [7] , we constructed a network model with one external layer and three layers of neurons (Fig 1A). The first layer is the external layer that corresponds to external stimuli or the sensory system’s response to these stimuli. For simplicity, we approximated the actiVity of external sources using a Poisson process with the constant rate 1180. If we define the Poisson process with rate r as 6(1’), the actiVity of the external source [,4 at time tis written as sfl(t) = 6 (see Table 1 for the list of variables). Neurons in the input layer fire spikes in response to actiVity in the external layer. By assuming a rate-modulated Poisson process, the spiking activity of the input neuron 1' follows p where no is the instantaneous firing rate defined with rf = vff — E qmvj, q,” is the response [421 probability for the hidden external source {4, and qb(t) = t2e‘t/9t /20t3 is the response kernel for each external event. In most theoretical studies, cross-correlations give an exponential decay or a delta function [5,38], but here we used a response kernel that produces broader correlations (Fig 4A right), because the actual correlations observed in the cortex are usually not sharply peaked [36,37]. For instance, for the exponential kernel qbe(t) = e‘t/Qt /0t, correlations show a peaked distribution even if the timescale parameter fit is several milliseconds (83A Fig). Because of the common inputs from the external layer, input neurons show highly correlated activity, which enables population coding of the hidden structure. Although here we explicitly assumed the presence of the external layer, these analytical results can also be applied for arbitrary realization of a spatiotemporal correlation. Output neurons are modeled with the Poisson neuron model [5,38,45] in which the membrane potential of neuron j at time tis described as where wjiX and wjkz are the EPSPs/IPSPs of input currents from input neuron x,- and lateral tory connections are clin and dij. For feedforward excitatory connections, the synaptic delay clle is given by the sum of the axonal delay dija and dendritic delay dijd, whereas for inhibitory connections, we assume for simplicity that the delay is purely axonal. The response of the output neuron follows yj(t) = 6 Similarly, inhibitory neurons in the lateral layer show tic delay of the lateral connection dij. The synaptic delay of the excitatory lateral connection is also approximated as the axonal delay. The spiking activity of the inhibitory neurons is given with zk(t) = 6 For analytical tractability, we use a linear response curve gE(u) = u and g1<u> = u.

For most of this study, we focused on synaptic plasticity in the feedfor-ward connection WX, with fixed lateral synaptic weights WY and W2. When the timing of the spikes at the cell bodies of pre and postsynaptic neurons is tpre and tpost, spike timings at the and ddji. For every pair of time and tspost, synaptic weight change is given With For the synaptic weight dependence of STDP, we considered a pairwise log-STDP [31] in Which LTP/LTD follows where E is a Gaussian random variable. The log-weight dependence well replicates eXperimentally observed synaptic weight distributions [32,33] and is suggested to have an important function in memory modulation [82]. Analytical treatment below is applicable to other types of synaptic weight dependence, yet in the additive STDP (i.e. fi,(w) 2 CP and fd(w) 2 Cd), the mean-field equation typically does not have any stable fixed point. In addition, under the multiplicative STDP in which LTD has a linear rather than a logarithmic dependence on synaptic weight, strong correlation is often necessary to induce salient LTP [31]. The coefficients CP 2 1 and C d = Cptjf/tjf are chosen so that total LTP and LTD are balanced around the referential synaptic weight. The STDP at E-to-I connections and I-to-E connections is similarly defined. For simplicity, we assume that synaptic delays are solely axonal (i.e., dz]. = 61,3“,di J. = 61,373.“), and the change in synaptic weight does not depend on the synaptic weight. To maintain the balance between LTP

2 LC? 2 ’yZCi‘L'i/‘L'i ,nZ = 0.317wf/ We also modify constant (initial) synaptic weights to way 2 50.0 and woZ = 25.0, and bounded synaptic weights with wYmax = 100.0 and wzmax = 50.0. In this normalization, the total lateral inhibition takes the same value as that in the non-plastic model at the initial state. Time windows are defined as TPY = Id = 11,2 = 1,5 = 20.0 ms. In Fig 6C, anti-Hebbian STDP was calculated by for Q = Y or Z. Similarly, the correlation detector type of STDP in 82 Fig was defined as The anti-correlation detector was calculated by changing the sign of above equations.

In the main text, we performed all simulations With a linear Poisson model for analytical purposes, although we also confirmed those results With a conductance-based LIF model (81 Fig). In the LIF model, the membrane potentials of excitatory neurons follow Where ngE and ngI are excitatory and inhibitory conductances, respectively, and 1‘5 and tks are the spike timings of input neuron i and lateral neuron k. Similarly, for inhibitory neurons in the lateral layer,

Both cross-correlation and mutual information behave as they do in the Poisson model, but the performance is slightly better, possibly because the dynamics are deterministic (Fig ID and IE, SIB and SIC Fig); however, membrane potentials show different responses for correlation events (SID Fig) because output neurons are constantly in high-conductance states, so that correlation events immediately cause spikes. As a result, membrane potentials drop to the Vref, and the average potential goes down. Interestingly, after neuron groups detect different signals, a preferred signal initially causes hyperpolarization due to firing, but, subsequently, a non-preferred signal causes hyper-polarization due to lateral inhibition (Fig ID right). The PSTH of firing shows that the behavior of the membrane potential in the Poisson model is similar (Fig IC and SIE Fig). This is natural, because in the linear Poisson model, the firing rate has linear relationship with the membrane potential, whereas in LIF model relationship between the average membrane potential and firing rate is highly nonlinear. Bayesian ICA. If discretized with At, the time series of the external source activity is writ-formation I, the joint probability of sources S and the estimated response probability matriX Therefore, by considering marginal probability, By considering maximum likelihood estimation for a given prior P[S|I] , Q can be optimally estimated [46,47]. In our problem setting, by assuming that external signals are independent, and input neurons respond to signals With a Bernoulli process, Where

pf = 1 — <1 — rsAofi [1 — am: ask/w] ,qsk = 2; 3 M + 1/2>At12exp[—<k + mam/at].

Therefore, log-likelihood becomes By taking gradient descendent,

H615.

Therefore, we need to calculate the integral over all possible combinations of sources in the past to obtain stochastic gradient descendent; however, such a calculation is computationally difficult and incompatible With neural computation. Instead, we used sequential sampling of Where Note in the above equations, x is given as a fixed value and not a random variable. Under this sample-based approximation, the stochastic gradient descendent follows For Fig 7C, we discretized the activity of hidden sources and input neurons With 5 ms bins, and performed learning With a learning rate 1] using the sample sequence Y. For the ideal case, we performed sequential sampling from the true response probability Q. at the connection ql-y caused by an output spike yk'k # = 1 for Xik = 1 is written as In the absence of the input spike (Xik = 0), an output spike yk'k’y = 1 causes LTD in total dependence similar to those in STDP. In Fig 7D, we plotted Aqifl (gm: 0.1, 0.3, 0.5). Blind source separation. In the blind source separation task, we created the original and AqiTD for different 511.” source by calculating high-frequency and low-frequency components separately. First, the spectrum of the signal q at a high frequency was defined as where fl,” is a characteristic frequency of signal q, and kflm are the harmonics of that frequency. The standard deviation was defined as 0W 2 kazf for 0'sz = 20Hz. Low-frequency components were directly given as an exponential oscillation as below. flu is a characteristic frequency, and 6qu is the delay. By combining these two components, the amplitude of a mixed sound is given as Summation over frequency f is performed using 400 representative values that correspond to the tuned frequency of each input neuron:

Simulations were calculated using the Runge-Kutta method, with a 0.05 ms time step. Initial synaptic weights were randomly chosen with W13?) 2 w3(1 —|— min

Analytical consideration of synaptic weight dynamics

: 166t3 Average synaptic weight velocity. The synaptic weight dynamics defined above can be rewritten as ij 7 7 of time and also using a stochastic Poisson process, synaptic weight change follows aptic weight dynamics can be analytically estimated. Because the spike probability linearly depends on the membrane potential in our model, cross-correlation follows Since we define cross-correlation among input neurons as the first term is written as This result is consistent With that in previous studies [5,38,45]. The analysis can be extended to the cross-correlation between an input neuron and a lateral inhibitory neuron as Theoretically, expansion over a lateral connection should be performed infinite times to obtain the exact solution, but at each expansion, the delay caused by synaptic delay dZ+dY and EPSP/IPSP rise times is accumulated so that the effect on correlation rapidly becomes small, especially when the original input cross-correlation C(t) is narrow; however, even if C(t) is broad, the effect for learning is bounded by the STDP time window. Therefore, higher order terms practically influence weight dynamics only through firing rates, so that by applying the approximation the last term can be obtained. In general, an is not analytically calculable, but by considering the balanced condition, it can be estimated. Therefore, the second term is given as Therefore, if we denote

average synaptic weight dynamics satisfy

The first two terms are Hebbian terms that depend on correlation by FX1 and FXZ, Whereas the remainders are homeostatic terms. In all terms, synaptic weight dependence is primarily caused by joi and not by other synapses. By inserting the explicit representation of correlation into the equation above, PM and FX2 can be rewritten as

P P ji) : vSGf(W;§)Zqiuqlu7 : VSG§(WJ§)ZquqZfl7 “=1 M=1 < 1-,.)

2[ dsF<w;§,s>[ drax<r>[ dt'¢<t'>¢<t' — (r — s + 2%)), —oo 0 max(0,r—s+2dXd) < 1,) 2[ dsF<w;§,s>[0 dram 0 dqay<q>[0 dr'am

Where t” = r+q+r’-s+2dxd+dz+dy. Note that G IX and G2X do not depend on any indexes of the neurons, except for synaptic weight dependency, and so the two values are considered basic constants that decide how correlation shapes learning. The correlation kernel 11X was derived from

If the correlation structure C(s) is simply organized, further analytical consideration is possible. In the two-source model shown in Fig 2A, lateral connections are structured nonreciprocally, and EPSP/IPSP sizes are constants. The synaptic weight matrices are written as Therefore, the original L x M differential equations can be reduced into 2 x 2 equations of representative neurons as The firing rates of inhibitory neurons can be approximated as Therefore, by solving the simultaneous equations for VIZ and v22, This analytical approach is applicable only when the synaptic weight change is sufficiently slow relative to the neural dynamics. Also, because we ignored the variance in the synaptic weights, numerically the accuracy is limited. Analytic approach for ST DP in lateral and inhibitory connections. Using a similar calculation as above, synaptic weight development of the lateral connections is given as Where 00 where = [ dr5X(r). The meaning of these equations is made clear by summarizing the 0 correlation propagation in the diagrams (82D i—iii Fig). In the diagram, blue wavy lines represent intrinsic correlation, and arrows are synaptic connections. To estimate how a blue correlation influences STDP at a red arrow, we need to determine all the major trajectories in which the correlation reaches pre and postsynaptic neurons. In the linear Poisson framework, for a given trajectory, the propagation of a correlation is calculated by simply using integrals as above. From this diagram, we can safely assume that ng and gY3 are negligibly smaller than gYI, because trajectories (ii) and (iii) are secondary correlations and also contain synaptic delays. In this approximation, we additionally assume that w and the eigenvalue is es lY — Because the eigenvector develops by eXp[ch1Y — Wfift], When ng is positive, the E-to-I connections are more likely to be structured in a way that the inhibitory neurons become feature selective. On the other hand, if that value is negative, such structure may not be obtained. Note that (1, -1, -1, 1) is not the principal Similarly, for inhibitory connections We approximated with only two terms because the third term is negligible (82D iv—vi Fig).

If we assume WY 2 < > , and gZZ = 0, then the synaptic weight change follows

Awf1 — Awf2 = Aw;2 — Aw§1 = csglz — wf)2(w§ — This means that if glz is positive, reciprocal connections are enhanced (or inhibitory connections to the neurons coding a similar feature are enhanced), whereas for negative gII, inhibitory connections develop non-reciprocal-ly (i.e., lateral connections function as mutual inhibition between output excitatory neuron groups).

Although it is difficult to study all combinations of STDPs, we still provide analytical insights by investigating the behaviors of g1Y and gIZ. 82E Fig shows the behaviors of four different types of STDP. This indicates that the anti-correla-tion detector type of E-to-I STDP [53] tends to cause non-feature-selective lateral connections. In addition, under the anti-coincidence detector type of I-to-E STDP [55] , mutual inhibition structures would be preferred; however, the implication of our analytical method is limited, and further study will be necessary to fully understand the functions of the various types of STDP.

Evaluation of the performance

We evaluated the performance by measuring the mean cross-correla-tion between the external sources and population activity of the output neurons. For time bin At kAt is a set of output neurons coding a source 1/. For these, cross-correlation is defined as

We used TC 2 10 ms for the analysis. Correspondence between k=1 sources and output groups are arbitrary, and so the learned correlation should be given between sources and output groups. For example, When For the models With randomly connected lateral inhibition and (e+i) STDP, we defined output neuron j as belonging to QHY if for octh = 1.5, and the cross-correlation was calculated based on QHY. Mutual information. Based on the discretized hidden external source/ output neuron activity syk, yv k, we defined the binary variables statement X. Therefore, mutual information can be defined as

Supporting Information

(A) Synaptic weight developments at the feedforward connection. (B) Cross-correlation and mutual information calculated for various delays. Both values were calculated by averaging five independent simulation results. (C) Development of two values for the simulation shown in (A). (D) PSTH of the membrane potential calculated for gray areas in (A). (E) Peristimulus time histogram (PSTH) of the firing probability for the same simulation. (EPS)

(A, B) Synaptic weight development When the number of external inputs is three (A) and four (B). Thick lines represent averages over all synapses, and thin lines represent individual synaptic weights. Colors represent detected sources for output neurons (left) and inhibitory neurons (middle right). (C) Relationship between the number of inhibitory neurons and the lateral structure. (D) Propagation of structure. i to iii correspond to lateral excitatory connections, and iv to vi correspond to feedback inhibitory connections. (E) Analytic results for various types of STDP.

(A) Cross-correlations among input neurons responding to the same source calculated from simulated

independent signals.

Acknowledgments

We thank Matthieu Gilson and Florence Kleberg for helpful discussions and comments on the manuscript.

Author Contributions

Performed the experiments: NH. Analyzed the data: NH. Wrote the paper: NH TF.