Exponential distance distribution of connected neurons in simulations of two-dimensional in vitro neural network development

Abstract

The distribution of the geometric distances of connected neurons is a practical factor underlying neural networks in the brain. It can affect the brain’s dynamic properties at the ground level. Karbowski derived a power-law decay distribution that has not yet been verified by experiment. In this work, we check its validity using simulations with a phenomenological model. Based on the in vitro two-dimensional development of neural networks in culture vessels by Ito, we match the synapse number saturation time to obtain suitable parameters for the development process, then determine the distribution of distances between connected neurons under such conditions. Our simulations obtain a clear exponential distribution instead of a power-law one, which indicates that Karbowski’s conclusion is invalid, at least for the case of in vitro neural network development in two-dimensional culture vessels.

Graphical abstract

Keywords

distance distribution / connected neurons / development / exponential / power-law / neural networks / complex systems

1 Introduction

Complicated functions of the brain are largely attributed to the complex topological connections of a neural network (NN) that consists of about10¹¹ neurons and10¹⁴ synaptic links. The geometric distance distribution of connected neurons in an NN is pertinent to its dynamic properties because it is related to real properties such as message transmission [1, 2], energy cost [3], and response time [4, 5]. The measurement of such connections can be made at three levels: between brain regions [6, 7], between voxels, which are artificially divided blocks of brain tissue, and between neurons, which are tiny hence difficult to discern. Experimental results have exhibited the small-world topology of structural connections at the first two levels. At the microscopic level, nodes in the network are neuronal somas, while the links are synapses between a tip of the axon of one neuron and a terminal of a dendrite of another neuron. Large-scale in vivo measurements of the topological properties of an NN and the geometric distances between neurons are still challenges for researchers, although small-scale precise observation has been successfully performed on C. elegance [8] and mammalian visual cortexes [9].

In contrast to in vivo experiments, significant achievement has been made on the in vitro development of NNs. For instance, Ito [10] found that synchronized bursts of an NN on a multi-electrode array in a cultured vessel required an initial plating density of at least 250 cells/mm² for Neuron Culture Medium and 500 cells/mm² for DMEM/serum. They also found that the final densities of surviving neurons at one month decreased greatly compared with the initial plating densities and stabilized in denser cultures. In another experiment, the densities of antibody-labeled synaptic terminals in both types of cultures increased gradually from 7 to 21–28 days in vitro (DIV). They did not increase further at 35 DIV and tended to become saturated [11], which serves as a metric for numerical simulations to match.

From theoretical aspect, Karbowski [12] presented an optimal wiring principle and plateaus regarding the degree of separation for cortical neurons. He derived this result based on scaling theory from statistical physics using two neuron linking rules: minimal total axonal length and minimal total energy consumption during information transport, and let distance distribution p(r) reside at the balance point between these two factors. By assuming that neurons in the cortex can be reached in just a few steps (i.e., the small-world property, which is now widely accepted), he obtained neural distance distribution p(r)\sim {(r/\sigma )}^z for \frac{r}{\sigma }\gg 1, i.e., the power-law decay for large distances with z=-2\frac{1-\alpha }{1-2\beta } and p(r)\to 1 for r\to 0. The quantity σ denotes a microscopic length characterizing a neuron’s size or the extent of local intra-cortical connections, where and , so that . Karbowski’s power-law decay for distance distribution p(r) between neurons has not been verified up to now, as direct in vivo measurement in human brains is not available. The purpose of the present work is therefore to check its validity in an indirect way.

Based on a simple model, we simulate the in vitro development of rat NNs in culture vessels. Considering the death of isolated neurons and synapses formation process between connected neurons, we fit the model to the saturation time of the synapses number with the results of experiments done by Ito [10, 11] by tuning the simulation parameters that describe various conditions. We then determine the distance distributions in the final networks that fit the development and saturation periods in Ito’s experimental results well, so as to compare the Euclidean distance distributions obtained from the simulations with the power-law decay in Karabowski’s prediction. We found that distance distribution is exponential rather than power law and thus Karbowski’s prediction is invalid.

2 Two-dimensional model for in vitro NN development

To mimic in vitro NN development, we set up a simple model with the following rules. (i) Initial setting. We assume a cultured square vessel of size L with N neurons randomly scattered on its bottom plate. All coordinates of such nodes are recorded. Based on Ito’s experiments, the initial densities D=N/L^2 of the neurons are set to 500/mm² and 1000/mm², respectively, to ensure the formation of a final NN. The total culture time is 35 DIV. (ii) Growth of axons and dendrites. We choose the soma radius of rat neurons to be 3\times {10}^3 nm [13]. Each axon of a neuron grows at a speed that is uniformly distributed in the range (0,V_{\max }) in an arbitrary direction, where V_{\max } is an adaptable parameter to fit the experimental results and is assumed to be proportional to the vessel size. The axon grows 0–5 branches after two weeks, and each branch is shorter than 1/5th the length of the axon, which corresponds with information in neuroscience textbooks. Each dendrite grows hierarchically with 2–4 branches every three days after an initial growth of the first few days (3–5 DIV). We set its growth velocity to 130 nm/day, which is much smaller than V_{\max } used for the axons [15]. In addition, the terminal of a branch can generate 2–4 new branches every three days after it has been generated, and this process is constrained by the maximal length, 5000 nm, of such branches. The center angle of the dendrites is opposite to that of the axon, and they span an angle of {240}^{\circ }, which is tunable. (iii) Formation of synapses: When a tip of an axon is less than 100 nm from a terminal of a dendrite, they form a synapse, because any synapse is actually a junction with a gap [14]. (iv) The number of surviving isolated neurons M(t) follows the rule [3]:

\begin{array}{}& M(t)=M(0)\exp (-\lambda t),\end{array}

where M(0) is the number of remaining rat neurons isolated outside the NN and λ is a decay constant in units of 1/DIV [3]. In this phenomenological model, we have adjustable parameters V_{\max }[15] and λ to fit the saturation time of the number of synapses in rat NN development in vitro [11].

3 Results of numerical simulations

3.1 Fitting synapse number saturation time in the development process

To ensure that the initially scattered neurons finally form a connected NN and the number of synapses saturates within the range of DIV reported in experimental results, we need to carefully adjust a parameter set that consists of decay coefficient λ and maximal growth velocity V_{\max }. The simulation results are not always good for arbitrary parameter values. For example, if we chose \lambda =0.1 and V_{\max }=0.04L, L=4.0 mm, with D=500/mm², no saturation could be found in the simulations. Saturation time in the in vitro development experiment was fitted well by tuning parameters λ and V_{\max } in the simulations. For a growth velocity in the range (0,V_{\max }), where V_{\max }=0.08L, the number of saturated synapses is close to the empirical values using \lambda =0.2 for culture vessels of sizes L=2.0,4.0, and 8.0 mm. Typical results are shown in Figs. 1(a) and (b). With an extended range of parameters, say \lambda =0.16 or 0.18/DIV, similar results were also obtained. Typical results are shown in Figs. 1(c) and (d). One can see that the saturation times S_a are 30,29,29, and 28 DIV for the parameters in Figs. 1(a-d), respectively, which all agree with experimental data.

3.2 Fitting the stabilization time of the number of surviving neurons in the development process

The number of surviving neurons decreases with time and arrives at a stable value when the number of synapses is saturated, which is obtained using the same modeling parameter set in our simulations. Therefore, stabilization time S_n for surviving neurons is estimated to be around S_a, which is supported by our simulations. Typical examples are shown in Figs. 2(a) and (b) using the parameters D=1000/mm², and V_{\max }=0.08L for \lambda =0.2/DIV [Fig. 2(a)] and \lambda =0.16 and 0.18/DIV [Fig. 2(b)], respectively. Stabilization time S_n for the surviving neurons was found to be S_n=27, 29, and 28 DIV for D=1000 cells/mm², \lambda =0.2/DIV, and V_{\max }=0.08L/DIV for L=2.0,4.0, and 8.0 mm, respectively. Furthermore, S_n=29 and 31 DIV for the same values of D and V_{\max }, L=4.0 mm, and \lambda =0.16/DIV, and 0.18/DIV, respectively. One can see that the simulated stabilization time values are all close to those obtained in experiments, although stabilization time S_n is not always the same as synaptic number saturation time S_a.

3.3 Calculation of distance distribution between physically connected neurons in final NNs

Given the simulation results for which both synaptic saturation time S_a and neuron number stabilization time S_n fit well with experimental values, we can count the number of all possible distance values r between all pairs of connected neurons in a NN after its saturation (stabilization) and calculate the distribution function of such distances. Simulation results exhibit quite a definite exponential distance distribution: p(r)\sim {\mathrm{e}}^{-\alpha r} for sufficiently large distance r, where α is a constant that relies on parameters of the culture vessels. Typical distributions of p(r) are shown in Figs. 3(a)–(c). With \lambda =0.2 and V_{\max }=0.08L for sizes L=2.0,4.0 and 8.0 mm, we see an exponential decay tendency versus distance r, where r\gg \sigma, for a fixed initial density of neurons D=1000 cells/mm². Moreover, we see that larger numbers of neurons in the simulations lead to a clearer exponential distribution, which reveals the effect of statistical laws for finite systems. Similar results appear for the cases \lambda =0.16 [Fig. 3(b)] and \lambda =0.18/DIV [Fig. 3(c)], each with L=4.0 mm and L=8.0 mm, respectively. All these results convince us that the distance distributions of NNs growing in vitro are exponential instead of the power-law of Karbowski for the condition r\gg \sigma. In fact, such an exponential trend is obscured when growth velocity increases beyond a certain level (V_{\max }>0.12L), which limits the parameter range of the present model. Actually, the sizes of rat somas can change during development [13], but we use a simplified case with a steady size at the present stage. Moreover, when we checked our simulation range (0,V_{\max }) with V_{\max }=0.0808L for L=4.0 mm and L=8.0 mm, we found that velocities in these ranges correspond well with previously reported measured results [15]. In addition, limiting simulation parameter λ in the range (0.16,0.20) was inspired by experiments [3]. These practical parameters ensure successful S_a and S_n fitting and extended results p(r) of the proposed phenomenological model.

All results simulated using parameters that mimic in vitro conditions in the culture vessels are listed in Table 1. The sizes of the culture vessels and the time needed to stabilize neuron numbers are not available from the reference. In this work, we increased sizes L while keeping the density D of neurons invariant. In this way, we repress fluctuations because of the finite size effect.

4 Conclusions and discussion

To summarize, we obtained the necessary parameter values for both the range of axon growth velocity and decay rate of surviving neurons by fitting simulation results to the saturation (stabilization) time obtained from in vitro development in two-dimensional culture vessels. We then obtained an exponential distance decay law for connected neurons in NNs under such growth circumstances. The main point of this model is that surviving neurons grow their axons and dendrites ends to form synapses while neurons remaining isolated outside the network die out at a constant rate. Our simulation results do not support Karbowski’s power-law decay function for physically connected neurons, although we cannot disprove it, as we do not know practical range of the parameters in the brain, and this function was derived against a quasi-three dimensional background for the whole brain. However, because it is parameter-independent, which is quite different from what we have seen in experiments, the present simulation results indicate that it is invalid, at least for the development of NNs in two-dimensional culture vessels.

Actually, the exponential decay distribution of distances between individual neurons in vivo has been revealed by both empirical investigation [8] and theoretical models [16-18]. The present work simply contributes a phenomenological model for in vitro development in culture vessels. Indeed, nothing other than exponential decay can be obtained this way, as it contains a definite length scale (L). However, L-based axonal maximal growth velocity resides within the range of experimental results [15]. The exponential decay of the distance distribution actually implies a small-world network, in which the wiring cost, defined as the sum of the interneuronal distances, plays an important role in the actual neuronal network [8]. It is characterized by dense local clustering (compared with that in scale-free networks) or the cliquishness of connections between neighboring nodes, yet a short path length between any (distant) pair of nodes is caused by the existence of a relatively few long-range connections. This is a model for the organization of anatomical and functional brain networks because a small-world topology can support both segregated/specialized and distributed/integrated information processing. Moreover, small-world networks are economical, tending to minimize wiring costs while supporting high dynamical complexity [19]. However, how the movable neurons in the brains of infants possibly rearrange [18] their positions to reach permanent sites, or if neurons in vitro really undergo such a rearrangement, remains a mystery for us to explore.

References

[1] M. J. Chacron, L. Maler, and J. Bastian, Electroreceptor neuron dynamics shape information transmission, Nat. Neurosci. 8(5), 673 (2005)

[2] L. Agnati, L. Santarossa, S. Genedani, E. Canela, G. Leo, R. Franco, A. Woods, C. Lluis, S. Ferré, and K. Fuxe, On the nested hierarchical organization of CNS: Basic characteristics of neuronal molecular networks, Comput. Neurosci. 3146, 24 (2004)

[3] E. Bullmore and O. Sporns, Complex brain networks: Raph theoretical analysis of structural and functional systems, Nat. Rev. Neurosci. 10(3), 186 (2009)

[4] C. L. Leveroni, M. Seidenberg, and A. R. Mayer, Neural systems underlying recognition of familiar and newly learned daces, J. Neurosci. 20(2), 878 (2000)

[5] G. Shahaf and S. Marom, Learning in networks of cortical neurons, J. Neurosci. 21(22), 8782 (2001)

[6] X. Liang, J. H. Wang, and Y. He, Human connectome: Structural and functional brain networks, Chin. Sci. Bull. 55(16), 1565 (2010)

[7] P. Hagmann, M. Kurant, X. Gigandet, P. Thiran, V. J. Wedeen, R. Meuli, and J.P. Thiran, Mapping human whole-brain structural networks with diffusion MRI, PLoS One 2(7), e597 (2007)

[8] Y. Y. Ahn, H. Jeong, and B. J. Kim, Wiring cost in the organization of a biological neuronal network, Physica A 367, 531 (2006)

[9] R. Kötter, Online retrieval, processing, and visualization of primate connectivity data from the CoCoMac database, Neuroinformatics 2(2), 127 (2004)

[10] D. Ito, H. Tamate, M. Nagayama, T. Uchida, S. N. Kudoh, and K. Gohara, Minimum neuron density for synchronized bursts in a rat cortical culture on multi-electrode arrays, Neuroscience 171(1), 50 (2010)

[11] D. Ito, T. Komatsu, and K. Gohara, Measurement of saturation processes in glutamatergic and GABAergic synapse densities during long-term development of cultured rat cortical networks, Brain Res. 1534, 22 (2013)

[12] J. Karbowski, Optimal wiring principle and plateaus in the degree of separation for cortical neurons, Phys. Rev. Lett. 86(16), 3674 (2001)

[13] M. Miller and A. Peters, Maturation of rat visual cortex (II): A combined Golgi-electron microscope study of pyramidal neurons, J. Comparative Neurology 203(4), 555 (1981)

[14] B. Hayes and A. Roberts, Synaptic junction development in the spinal cord of an Amphibian Embryo: An electron microscope study, Z. Zellforsch. 137, 251 (1973)

[15] C. G. Dotti, C. A. Sullivan, and G. A. Banker, The establishment of polarity by hippocampal neurons in culture, J. Neurosci. 8(4), 1454 (1988)

[16] M. Kaiser, C. C. Hilgetag, and A. van Ooyen, A simple rule for axon outgrowth and synaptic competition generates realistic connection lengths and filling fractions, Cereb. Cortex 19(12), 3001 (2009)

[17] M. Ercsey-Ravasz, N. T. Markov, C. Lamy, D. C. Van Essen, K. Knoblauch, Z. Toroczkai, and H. Kennedy, A predictive network model of cerebral cortical connectivity based on a distance rule, Neuron 80(1), 184 (2013)

[18] M. Kaiser and C. C. Hilgetag, Nonoptimal component placement, but short processing paths, due to long-distance projections in neural systems, PLOS Comput. Biol. 2(7), e95 (2006)

[19] D. S. Bassett and E. Bullmore, Small-world brain networks, Neuroscientist 12(6), 512 (2006)

Rights & permissions

Higher Education Press and Springer-Verlag Berlin Heidelberg 2017