Eigenvalue spectra of asymmetric random matrices for multi-component neural networks
This paper focuses on large neural networks whose synaptic connectivity matrices are randomly chosen from certain random matrix ensembles. The dynamics of these networks can be characterized by the eigenvalue spectra of their connectivity matrices. In reality, neurons in a network do not necessarily behave in a similar way, but may belong to several different categories. The first study of the spectra of two-component neural networks was carried out by Rajan and Abbott. In their model, neurons are either ’excitatory’ or ’inhibitory’, and strengths of synapses from different types of neurons have Gaussian distributions with different means and variances. A surprising finding by Rajan and Abbott is that the eigenvalue spectra of these types of random synaptic matrices do not depend on the mean values of their elements. In this paper we prove that this is true even for a much more general type of random neural network, where there is a finite number of types of neurons, and their synaptic strengths have correlated distributions. Furthermore, using the diagrammatic techniques, we calculate the explicit formula for the spectra of synaptic matrices of multi-component neural networks.
In neuroscience, interconnections of neurons are often represented by synaptic matrices whose elements are drawn from a certain random matrix ensemble [1, 2]. Knowing the distribution of eigenvalues of these random matrices is very important in studying spontaneous activities and evoked responses of the network. To calculate the eigenvalue distribution of these matrices, it is often necessary to work with asymmetric (non-hermitean) random matrix theory, which has been successfully applied to many fields of physics and interdisciplinary sciences, e.g. the phase diagram of QCD [3, 4], nuclear decay and resonances in multichannel chaotic scattering  and neural networks[1, 2, 6].
A prominent result of asymmetric random matrix theory is Girko’s circle law . In its variation with partial symmetry, the circle becomes an ellipse . These classic results, however, can not be directly applied to realistic neural network models where neurons do not behave in the same way [1, 8, 9]. Assume there are number of neurons and let be the synaptic matrix. In the model of Rajan and Abbott , there are number of neurons which are ’excitatory’, and all others are ’inhibitory’. To model this neural network, elements in columns of are sampled from a Gaussian distribution with mean and variance and elements in the remaining columns of are Gaussian variables with mean and variance . Therefore, the synaptic matrix has the structure , where is drawn from the real Ginibre ensemble  such that , , and , where and are identity matrices of dimension and , respectively. is a constant matrix whose elements are the mean strength of the synapses. Since there are two types of synapses, every row of is identical and in each row, the first elements are equal to and the remaining elements equal to . In particular, is chosen to be in a ’balanced’ situation such that [13, 14]. To confine eigenvalues inside a unit circle, a second constraint  is introduced which requires that the strengths of the synapses attached to each neuron independently sum to zero. It is found in  that, in the limit , modifying the mean strengths of excitatory and inhibitory synapses has no effect on the eigenvalue spectra of the synaptic matrices. Therefore, the spectrum of is identical to that of .
It is natural to wonder why the ’mean’ strength matrix has no effect on the spectra. Moreover, in real biological neural systems, several different types of neurons may connect each other to form a multi-component network [8, 9]. Distributions of synaptic strengths of different types of neurons are distinct  and could be non-Gaussian [10, 11]. Dynamics of this network therefore depend on properties of each type of neuron. It is interesting to find whether or not the eigenvalue density of this type of networks depends on the mean value of each individual type of synapse. These questions are addressed in section 2. One of the main results of this paper shows that even without the second constraint in  and when synaptic strengths have certain non-Gaussian distributions, the spectrum of the network still does not depend on the mean synaptic strengths.
Aside from biological motivations, the eigenvalue problem of random plus fixed matrices has been a research topic in both random matrix theory and condensed matter physics for a long time [15, 16, 17, 18]. A different point of view of the problem in this paper is: how is the density function of large random matrix perturbed by the rank-1 constant matrix . Note that random matrix is not of Wigner type, nor are its elements independently and identically distributed (iid). Therefore this paper provides new results to similar problems studied in [17, 18].
Furthermore, finding the eigenvalue density of random matrices of the form , where is drawn from a random matrix ensemble and is a fixed matrix, has been an interesting topic in random matrix theory and mesoscopic physics . When is drawn from the circular unitary ensemble (CUE), an exact result is given in , and the large-N limit is calculated in . In section 3, we calculate the density function of where belongs to the real Ginibre ensemble using the method introduced in [22, 23]. Discussion and remarks are made in the last section.
2 Synaptic strength of non-Gaussian distributions
Let the dimensional real matrix be the synaptic matrix of an -neuron network. Assume there are types of neurons and the i-th type of neuron has a population of , . Define a constant diagonal matrix
where is the -dimensional identity matrix. Let be an -dimensional row vector with the following form,
where is the mean strength of the synapses from neurons of the i-th type. Define the dimensional ’mean’ matrix , whose rows are all equal to . The synaptic matrix in our model takes the form
where is an dimensional real random matrix drawn from the ensemble
where is an arbitrary function and is the normalization constant. The case corresponds to Ginibre ensemble where elements of are statistically independent Gaussian variables. By symmetry, the mean of vanishes, . The variance of is determined by , which is normalized to be , i.e. .
Synaptic matrix defined in Eq.(2.3) has column blocks. Each block corresponds to one type of neuron. By construction, elements of have the following statistical properties,
i.e. the i-th type of synaptic strength has variance and mean . The matrix has a similar column block structure as the one defined in . Following , we also choose to put the synapses at the ’balanced’ situation, i.e.
We want to know the eigenvalue density of in the limit , with fixed ’s. Let , the density is related to the Green’s function as
In the above formula, means averaging over the ensemble Eq.(2.4) of matrix . We write the Green’s function as to emphasize it is not analytic on a 2-dimensional region of the -plane, more details can be found in . This region is called the support of the density function since on which we have . Since we will be dealing with both and dimensional matrices, to remove ambiguity, we use as the trace operator for matrices. We will work on asymmetric random matrices with the methods introduced in [22, 23]. For consistency, we adopt the notation convention of  in the remaining of this paper. Define a dimensional matrix
Introducing the self-energy , we have
we have the relation
It is impossible to calculate explicitly for arbitrary . But for our purpose it is sufficient to know its basic structure. Without loss of generality, assume , so that we can expand as
We expand the higher order terms in Eq.(2.14) and use the quadratic term to calculate the ensemble averages, denoted by . Let . Then because of the presence of matrix , we have
The self-energy can be written in terms of the cumulants of , i.e. , . And it is well known that in the limit , to leading order in , each of these cumulants is the sum of all connected plannar diagrams with external ’s and ’s. All diagrams which contribute to are also planar diagrams, as shown in Fig.1. By Eq.(2.15), the self-energy has the following structure
where scalars and are functions of and and are determined by and is defined in Eq.(2.1). In appendix A, we prove that . This fact, together with Eq.(2.7), completes the proof that the eigenvalue spectrum of is identical to that of the random matrix , as discovered in  when belongs to Ginibre ensemble.
In Fig.2, we compare the eigenvalue spectra of and . In both cases is drawn from a non-Gaussian ensemble. All spectra are generated by Monte-Carlo simulations. Since , where , it is sufficient to show the dependence of eigenvalue density function on radius . We find these functions match quite well. In comparison, we replace with a constant matrix which does not have the column structure, and find the density function is rather different.
Next, we test our result on the two-component Gaussian network of . Denote the weight matrix by and follow the brief discussion in the introduction. We fix and choose the variance of the Gaussian distributions for excitatory and inhibitory synaptic strengths to be and , respectively. The only constraint on the mean strength matrix is the ’balance’ condition. As in the previous example and already pointed out in , the bulk of the spectra of with different matrices are almost identical. From numerical simulations, it appears that when has larger elements, there are more eigenvalues outside the support of the spectrum, see Eq.(3.22) for the definition. To show this is a finite-size effect, we calculate the percentage of ’outliers’, for different matrix size . It shows in Fig.3 that, for all values of , as , percentages of ’outliers’ for approach to that of .
3 Synaptic strength of Gaussian distribution
In this section, we calculate the eigenvalue density of multi-component Gaussian network. Assume there are types of neurons in the network and synaptic strengths have different Gaussian distributions. From the previous section we know the density functions of and are identical when has the column structure and satisfies the ’balance’ condition, even without the additional constraint of . Therefore, spectrum of this network is the density function of the following random matrix
where is drawn from the real Ginibre ensemble, i.e. in Eq.(2.4), and is defined in Eq.(2.1). The case is solved in  with the method in . For , we find the technique developed in [22, 23] is more convenient. Define an operator which, when acts on an matrix , gives . By Eq.(2.15), the equation for the one particle irreducible (1PI) self-energy is
Eq.(3.20) has multi-number of solutions. The correct one for our problem is the one satisfying the boundary condition
The boundary of spectrum is determined by the transition point , which corresponds to the circle with radius , such that
The disk region defined by is the support of the spectrum. Off the support, we always have . In the case , the formula in Eq.(3.22) gives the same result for spectrum boundary obtained in  by solving a saddle-point equation. From (2.11), the Greens’s function is given by the following formula
From Eq.(3.20), when on the support of spectrum, we get
Finally, by Eq.(2.7), we get the eigenvalue density of
Introduce the notation , where is a constant. From Eq.(3.25), we find the eigenvalue density at the centre and boundary of the spectrum
When , from Eqs.(3.20) and (3.25) we easily recover the well known result for Ginibre ensemble . For , choosing the solution for quadratic equation Eq.(3.20) satisfying condition Eq.(3.21), then from Eq.(3.25), we successfully recover the results obtained in . For large , it is hardly possible to have an analytic solution for Eq.(3.20). But it is very simple to find the numerical solution for this algebraic equation. It turns out that there is always only one solution on , which is just what we need according to the boundary conditions.
In Fig.5, we compare the density function in Eq.(3.25) with numerical simulations for synaptic strengths with different variances ’s but the same ’s. In Fig.6, we let ’s change but keep ’s fixed. In both cases, we observe very good match between numeric data and analytical results. The only significant deviation happens near . In fact, this deviation already appears when and is shown due to finite-size effect .
In the first part of this paper, we show that modifying the mean strengths of synapses of a neural network does not change the density function of synaptic matrix even when there are several types of neurons and the strengths of their synaptic connections have correlated distributions.
In Eq.(2.4), the ensemble of random matrix is chosen to be invariant so that all elements of random matrix have the same distribution. Differences between different types of neurons are introduced only by and . In fact, we can draw the synaptic matrix from more general ensembles. As long as Eq.(2.15) holds and has the column block structure, eigenvalue spectra of will not depend on .
We therefore prove that the density functions of large random matrices described by Eq.(2.3)-(2.4) are not changed by perturbations of the rank-1 matrix . This type of random matrices are not of Wigner type or have iid elements as in [17, 18].
It is its structure that makes irrelevant to the eigenvalue density function. In reality, we may need to choose the mean value of synaptic connections to be of the same order of their fluctuations, i.e. . But this is not necessary in our proof. If the ’balance’ condition is not imposed, the eigenvalue spectra will be identical to the ’balanced’ case except the eigenvalues at zero will be shifted [1, 17].
In the second part of this paper we calculate the density function of random matrices of the form , where belongs to Ginibre ensemble. These matrices describe random networks with multiple independent components. We find closed formulas for the eigenvalue density at both the centre and the boundary of the spectrum in terms of variances of synaptic strengths.
When is drawn from the ensemble in Eq.(2.4) and , we know by the Single-Ring Theorem [24, 25, 26] that the support of the eigenvalue spectrum is either a disk or an annulus. It will be interesting to find out whether or not the Single-Ring Theorem still holds when is diagonal but not proportional to . Eq.(3.22) shows when has Gaussian distribution the support of the spectrum is always a disk of radius , but never an annulus. This indeed agrees with the Single-Ring Theorem. Clearly, to prove the Single-Ring Theorem for , where belongs to the general ensemble defined in Eq.(2.4), we need to take different approaches. Work on this topic is currently in process.
I am grateful to Prof. R.T. Scalettar for helpful communication on numerical methods and giving me his Monte-Carlo code. I also thank Dr. K. Rajan for drawing my interest to this research. This work was supported by the Swartz Foundation.
In this section we show that in the large limit, due to the structure of , the Green’s function defined in Eq.(2.11) equals to defined similarly for , i.e.
Here we used the fact that the (12)-element of is proportional to identity matrix and the following formula from linear algebra
Step 2. Let be a -dimensional matrix with all elements equal to 1 and let for . Then has the following -block structure
where , , , and , for . This claim is proved by induction. First, notice that has the same -block structure as in (A.30) except all its parameters are of order 1. Using the fact
and (A.29), we find that for , indeed has the properties described by (A.30). Then assume the claim is true for . By straightforward calculation using (A.29) and (A.31), we find the claim is also true for .
- K Rajan and L.F. Abbott, Phys. Rev. Lett. 97, 188104 (2006).
- H. Sompolinsky, A. Crisanti and H.J. Sommers, Phys. Rev. Lett, 61 259 (1988).
- M.A. Stephanov, Phys. Rev. Lett. 76 4472 (1996).
- J.J.M. Verbaarschot and T. Wettig, Ann. Rev. Nucl. Part. Sci. 50 343 (2000).
- J.J.M. Verbaarschot, H.A. Weidenmüller and M.R. Zirnbauer, Phys. Rep. 129 367 (1985).
- H.J. Sommers, A. Crisanti, H. Sompolinsky and Y. Stein, Phys. Rev. Lett. 60 1895 (1988).
- V.L. Girko, Theory Probab. Its Appl. 29 694 (1984).
- G.M. Shepherd et al, The Synaptic Organization of the Brain (Oxford University Press, 2004).
- C. Holmgren, T. Harkany, B. Svennenfors and Y. Zilberter, J. Physiol. 551.1 139 (2003).
- B. Barbour, N. Brunel, V. Hakim and J.P. Nadal, Trends in Neurosciences, 30 No.12 622 (2007).
- S. Song, P.J. Sjöström, M. Reigl, S.B. Nelson and D.B. Chklovskii, PLoS Biol. 3 507 (2005).
- J. Ginibre, J. Math. Phys. 6 440 (1965).
- M.N. Shadlen and W.T. Newsome, Curr. Opin. Neurobiol. 4 569 (1994).
- T.W. Troyer and K.D. Miller, Neural. Comput. 9 971 (1997).
- E. Brézin, S. Hikami, and A. Zee, Phys. Rev. E 51 5442 (1995).
- A. Zee, Nucl. Phys. B 474 726 (1996).
- T. Tao, arXiv:1012.4818v4 [math.PR], Probab. Theory and Related Fields (to appear).
- A. Pizzo, D. Renfrew and A. Soshnikov, arXiv:1103.3731v4 [math.PR], Ann. Inst. Henri Poincaré, B (to appear).
- K. Zyczkowski and H.-J. Sommers, J. Phys. A 33 2045 (2000).
- Y. Wei and Y.V. Fyodorov, J. Phys. A 41 502001 (2008).
- E. Bogomolny, J. Phys. A 43 335102 (2010).
- R.A. Janik, M.A. Nowak, G. Papp and I. Zahed, Nucl. Phys. B 501 603 (1997).
- J. Feinberg and A. Zee, Nucl.Phys. B 504 579 (1997).
- J. Feinberg and A. Zee, Nucl. Phys. B 501 643 (1997).
- J. Feinberg, R. Scalettar and A. Zee, J. Math. Phys. 42 5718 (2001).
- J. Feinberg, J. Phys. A 39 10029 (2006).