A multi-agent model for growing spiking neural networks

09/21/2020 ∙ by Javier Lopez Randulfe, et al. ∙ 0

Artificial Intelligence has looked into biological systems as a source of inspiration. Although there are many aspects of the brain yet to be discovered, neuroscience has found evidence that the connections between neurons continuously grow and reshape as a part of the learning process. This differs from the design of Artificial Neural Networks, that achieve learning by evolving the weights in the synapses between them and their topology stays unaltered through time. This project has explored rules for growing the connections between the neurons in Spiking Neural Networks as a learning mechanism. These rules have been implemented on a multi-agent system for creating simple logic functions, that establish a base for building up more complex systems and architectures. Results in a simulation environment showed that for a given set of parameters it is possible to reach topologies that reproduce the tested functions. This project also opens the door to the usage of techniques like genetic algorithms for obtaining the best suited values for the model parameters, and hence creating neural networks that can adapt to different functions.



There are no comments yet.


This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.


Artificial Intelligence has looked into biological systems as a source of inspiration. Although there are many aspects of the brain yet to be discovered, neuroscience has found evidence that the connections between neurons continuously grow and reshape as a part of the learning process. This differs from the design of Artificial Neural Networks, that achieve learning by evolving the weights in the synapses between them and their topology stays unaltered through time.

This project has explored rules for growing the connections between the neurons in Spiking Neural Networks as a learning mechanism. These rules have been implemented on a multi-agent system for creating simple logic functions, that establish a base for building up more complex systems and architectures. Results in a simulation environment showed that for a given set of parameters it is possible to reach topologies that reproduce the tested functions.

This project also opens the door to the usage of techniques like genetic algorithms for obtaining the best suited values for the model parameters, and hence creating neural networks that can adapt to different functions.


I would like to express my most sincere gratitude to my supervisor Leon Bonde Larsen, who has always been a source of support, ideas, and resources. His dedication motivates the people around and unquestionably increased the quality of my work.

I would also like to thank the other students and researchers at the MMMI institute, whose feedback during this project is much appreciated.

Thanks to my classmates Katrine, Job, Anne, and Sergi for joining in the work sessions and making the duty of the thesis more pleasing. Special thanks to Job and Katrine for listening and proof reading my work repeatedly, and offering very precious tips. Last but not least, I would also like to thank Leo, as her company during those days turned invaluable.

Finally, I would like to thank my friends and family back home.


1.1 Neuromorphic engineering

In order to understand the aim and motivation of this project, it is paramount to understand what neuromorphic engineering is, as the project is highly influenced and motivated by the ideas that build up this field.

Neuromorphic engineering is a field that has been continuously growing in the last few decades. Its main goal is the replication of biological brains with technology, under the assumption that those systems have been evolving and adapting to the real world during thousands or millions of years and therefore they are highly reliable and optimised (Liu and Wang, 2001).

It is true that already existing digital systems and computer science have overtaken the biological brain in several aspects, specially those concerning speed for solving arithmetic and logical tasks. However, artificial computing systems are still very far from biological brains regarding energy efficiency, and their performance is very poor when dealing with the vast whole extent of the real world (Indiveri and Horiuchi, 2011). It is actually nowadays a big struggle for engineering to develop systems that can solve problems which require a big amount of computation without spending large amounts of energy (A comparison is offered in Fig. 1.2). Also, the task turns non-feasible if it requires interacting with most of the information available in a real environment.

Figure 1.2: Comparison of the power density and clock frequency of different processors developed during the last 40 years, as well as the human brain. Image obtained from (Merolla et al., 2014).

This motivation has led to the exploration of the idea of replicating some of the concepts that rule biological systems. Thus, neuromorphic engineering has put focus in the development of new hardware systems that are closer to the physical structure of biological systems, as well as software solutions that implement algorithms which resemble the behaviours of the biological brain (Soman et al., 2016).

Regarding the introduction of new computing paradigms, perhaps the most noticeable breakthrough has been the appearance of Spiking Neural Networks (SNNs). Whereas traditional Artificial Neural Networks (ANNs) calculate the value of the system outputs based on the value of the inputs in the same discrete time, SNNs assess the time evolution of the inputs i.e. The value of an output not only depends on the current value of the inputs, but also in their activity in the past. This new paradigm is closer to biological postulates, and the mathematical models describing it are based on the theories that explain the behaviour of biological neurons.

SNNs allow, on the one hand, to drastically reduce the amount of computational units (neurons), due to the fact that one neuron can theoretically encode a much bigger amount of information, as the value at each time step is relevant for the outcome of the system. On the other hand, the processing times can be considerably accelerated, as fewer neurons need to be assessed in order to determine the value of the outputs. This feature is enhanced too by the asynchronous nature of these networks, as they do not need to wait for a synchronization clock. This means that neurons do not have to be evaluated synchronously at the same discrete time. Actually, all of them can work as independent units that react based on incoming pulses, no matter when these happen.

Regarding the development of hardware solutions closer to biological structures, some initiatives have appeared recently, such as the TrueNorth circuit developed by IBM (Merolla et al., 2014), or the EU Human Brain Project. The latter involves, among others, the development of the SpiNNaker board (Furber et al., 2014). All of them have in common the implementation of systems that are composed by small computing units that work asynchronously and independent of each other, introducing a big degree of parallelism between different cores, and where the concept of memory is represented by the current state of the internal neurons. In Fig. 1.3 is depicted a schematic of the structure followed by these systems.

Figure 1.3: Conceptual structure of the TrueNorth circuit, which ressembles to the SpiNNaker board. They are formed by multiple processing systems that work asynchronously and in parallel. Image obtained from (Merolla et al., 2014).

Their structure is opposed to the traditional and predominant Von Neumann architecture for electronic and digital systems, where the different computing entities share a common system memory that they need to access synchronously. The main benefit of biological structures is a highly optimization of the resources needed for making calculations or operations, as they will consume only the minimum amount of computational units and reduce the read/write operations to the minimum extent. The hope is to drastically reduce the energy consumption from the current order of megawatts for super computing systems to the few tens of watts used by the human brain.

1.2 Brief introduction to the biological neuron

Neurons are specialized cells which carry out a central role in the nervous systems of animals. They are electrically excitable, feature that they use for receiving, processing, and transmitting information, mostly to other neurons. They are grouped in large populations, forming clusters that are able to originate what is known as intelligence and, in the case of mammals, most of the neurons are allocated in the brain, which is the centre of the nervous system.

In this section, a brief description of the biological neurons is offered, as well as some of their most relevant details. A more detailed description of the working of biological neurons and neural circuits is offered in section 2.1. In any case, a deep study on the biological neuron is out of the scope of this project, and since the study of biological neurons has been a very relevant topic in science for more than a century, the reader may look for further information in more specialized literature if interested (Gerstner et al., 2014; Dayan and Abbott, 2001).

Neuron structure:

Physically, the structure of a neuron is typically decomposed into three main components, which play different roles (See Fig. 1.4):

  • Soma: It is the main body of the cell, and it contains the main organelles. It is covered by a membrane that can be charged with electric potential. When this potential is high enough the neuron gets excited, and generates what is known as a spike.

  • Dendrites: These are filaments organized in structures known as dendritic trees. They are sensitive to incoming signals. Thus, their role is to receive information and propagate it to the soma of the neuron

  • Axon: A neuron typically has only one axon, and its function is to propagate the electric pulses generated in the neuron. The stem of the axon is denominated axon shaft, and at its tip is located the growth cone, which leads the elongation and movement of the axon through the neural circuit. It can travel as far as 1 meter in the case of humans

Figure 1.4: Schematic of the biological neuron, with its main components. Image obtained from Wikipedia (https://en.wikipedia.org/wiki/Neuron).
Neuron synapse and synaptic plasticity:

The synapse is the connection between two neurons, where one of them is transmitting information and the other one is receiving it. The synapse is driven by chemical and electrical forces, as the information is transmitted by the change of the chemical composition in the junction between both neurons, which leads to a modification in the electrical potential of the receptor neuron.

The chemical processes involved in the synapse are complex. Activity in a neuron leads to a modification of its state, which leads to some of its organelles to start releasing chemical components that travel through the axon shaft to the synapse and eventually reach the membrane of the receptor neuron. These components are denominated neurotransmitters, and the process by which they flow within the axon to the synapse and activate the connection is called neurotransmission.

Traditionally, it was believed that the wiring between neurons and the strength of their synapses were formed during early stages of animal life (i.e. during childhood) and stayed mostly unchanged afterwards. This idea actually reinforced the conception that learning is a process that happens mostly during the first years of animals’ life, being drastically reduced during their adulthood. However, research found evidences against this idea, raising the concept of synaptic plasticity. Synaptic plasticity is a property of neurons which makes their synapses dynamic. This implies that the synaptic weight is evolving through time, as well as the dendritic structure of the neurons (Feldman, 2009). Hence, it is nowadays assumed that learning is produced throughout the whole life of the animals (although less intense from adulthood) due to the aforementioned synaptic plasticity.

Hebbian learning:

The previous paragraph introduces the concept of synaptic plasticity, and how the synaptic weights evolve through time, provoking as well the creation and destruction of new synapses between neurons. In order to create a neuronal model replicating this behaviour, it is crucial to have an equation or set of equations that provides quantitative values, and Hebb’s rule follows this purpose.

Hebbian learning is a concept introduced half a century ago which states that the synaptic strength between 2 neurons is modified based on the correlation between the firing times of the 2 neurons involved. It offers a simple approach for explaining how synapses between neurons get stronger, based on the time differences between their spikes. A simple way of summarizing the rule would be as follows:

If neuron A and neuron B are connected by a synapse, and neuron A spikes repeatedly before neuron B spikes, the synapse from neuron A towards neuron B becomes stronger. Conversely, if neuron A spikes repeatedly after neuron B spikes, the synapse from neuron A towards neuron B becomes weaker.

Moreover, neuroscience has found evidences showing that the short-term fast growth of neuron spines is lead by a sort of Hebbian learning process (Feldman, 2009). This means that the spines of a neuron tend to grow towards other neurons that tend to spike before the actual neuron in a correlated way.


Neurotaxis is the travelling of neurons’ dendritic spines and axons through the neural circuits in order to reorganize the neural structure, being learning one of its main purposes (the other main purpose is the restoration of damaged areas of the nervous system after injuries). Biology explains neurotaxis mostly by the appearance of chemical potentials which provoke the attraction or repulsion of the neuron filaments. It has been proven that the growth cone at the tip of the axon shaft is sensitive to the proteins present in the neuron environment, reacting in different ways depending on the type of protein. Thus, the concentration of these proteins creates the moving behaviour leading to the desired structure (Dickson, 2002).

However, neurons are sensitive as well to electric potentials, and the movement of filaments towards their destination could be explained as well by the appearance of electric forces acting on them (galvanotaxis). Although literature supports that chemotaxis is the main process behind neurotaxis, there are plenty of experiments showing that dendrites are sensitive to electric potentials, and neurons emit as well electric pulses when they get excited (Patel and Poo, 1982). Thus, despite neurotaxis can not be fully explained if chemical mechanisms are not taken into account, considering only electric forces for explaining the movement of dendrites offer a simple representation of the movement of the spines within the neural circuits.

1.3 Project motivations

One of the main motivations of this project was the idea that ANNs have drawbacks that can be hardly solved without approaching their development with a wholly different paradigm. This was the main inspiration leading to the inception and development of SNNs, lead by the knowledge learned by neuroscience.

Traditional neural networks suffer from intrinsic limitations, mainly for processing large amount of data or for fast adaptation to a changing environment. Several characteristics […] are strongly restrictive compared with biological processing in natural neural networks. (Paugam-Moisy and Bohte, 2012)

Neuroscience has found evidence showing that neurons follow a Hebbian-like learning pattern. Moreover, Hebb’s rule also explains how dendrites grow and create new connections.

This phenomenon is believed to be paramount in the short-term learning in animal brains (Feldman, 2009). In fact, fast spine growth does a fine tuning of dendritic trees and thus makes the layout of the circuits to adapt to external events. This happens after the main structure of a neural network has been already created

With this in mind, the aim of this project was to establish a set of mathematical rules for modeling SNNs whose neuron connections are dynamic.

As mentioned in section 1.1, neuromorphic engineering follows the idea of creating asynchronous independent computing units (i.e. neurons). Therefore, the implementation and testing of the neuron model has been done by using multi-agent systems, where a whole network can be distributed in several units which follow the same rules with a big degree of independence to each other.

All this lead to the formulation of the main premise for the current project:

There is a set of rules that defines the evolution of the structure of SNNs and makes them produce intelligent systems.

There are many ways of defining intelligence, depending on the perspective, scope, or field in which the definition is formulated. In this project, intelligence has been considered to be “the capacity of a system to solve a logical problem”. It is a simple definition, and lacks many details that have to be taken into consideration for a complete and broad application of the concept. In any case, it is enough for formulating the set of statements which give shape to the hypothesis that have been tried to validate or reject in this project.

The hypothesis proposed to be validated or rejected in this project were the following:

  1. A set of rules exists for making neurons grow in an SNN.

    1. The established growth rules are determined by a set of parameters.

    2. Growth can be obtained by applying Hebbian learning.

    3. Growth depends on the spatial distribution of the neurons.

    4. Neurons can also show behaviours inhibiting growth.

  2. There are logical problems that can be solved by an SNN with a certain topology.

    1. It is possible to propose a logical problem that has at least one solution.

    2. There is a quantitative value for determining how successful a solution is. It can range from a true-false boolean to a percentage score.

  3. The growth rules can be optimized by using a Machine Learning algorithm.

  4. For a given function , an SNN exists with a minimum number of neurons, , that is able to grow to a topology able to solve the function.

  5. The topology of an SNN can change from one structure to another based on its environment and the problem it has to solve.

This project has mainly focused on the two first points in the previous list (Although points 1.c, 1.d, 2.b have only been partially tested). Point 3 was initially included in the scope of this project, but the goals were narrowed afterwards, and research on this area has been left for future works. The last two points are very relevant, as their veracity could allow this networks to be reused for different purposes without the need of going through the design stage i.e. it would not be required a designer to tune the network until the correct setup was reached. In any case, work on this direction is still far, and more efforts are needed in the previous points first.

1.3.1 Main goals

Based on the aforementioned motivation outline, the main goals that had driven the efforts made during the project are the following:

In-depth study of neurotaxis, and review of the current knowledge:

As the final task is to implement SNNs that can grow dynamically, and it is known that biological brains present this feature, it was paramount to investigate and understand the current knowledge in neuroscience about neurotaxis.

Establishment of a mathematical model for the growth of neuron spines:

Once understood how neurotaxis works in biological brains, a set of rules that simulate this behaviour was proposed. These rules were designed with the trade off between creating a faithful model of the biological neurons and shortcoming the huge complexity that they present.

Testing of the previous model in SNNs by adapting it to existing spiking models:

The created set of rules that implements Hebb’s rule for simulating neurotaxis can not be put into practice without a mathematical model representing the firing of the neuron. This is due to the fact that Hebb’s rule establishes the growth of neuron connections based on the timings of the firing of the neurons involved. Therefore, it can only exist if neurons fire along time.

Once this and the previous points were implemented, it was obtained the design of a system capable of being tested autonomously without the need of other algorithms.

Design of a multi-agent model for implementing the aforementioned rules:

Due to the nature of neural networks, the designed model was implemented by creating a multi-agent system. This way, the system could be split into independent computational units, which would make it easier for handling networks with a big amount of neurons.

A multi-agent system was prepared, and the simulation environment RANA was chosen. In chapter 3 is explained the design, and the main features of the aforementioned environment.

Establish the grounds for the application of genetic algorithms for tuning the growth model:

The project involves the design of a model for replicating the behaviour of biological neural networks, and it was tested for implementing simple topologies and functions.

As it will be explained later in this report, the designed model contains several parameters that can drastically alter the performance of the implemented neural networks. Setting the parameters to different values can lead to a wide variety of topologies and logic functions.

Therefore, this is presented as a potential system to be optimized by the usage of a complementary learning mechanism that can be responsible of the tuning of the mentioned parameters. Namely, it is believed that evolutionary algorithms can be used for obtaining an adequate set of parameters that lead to the implementation of a desired function.

1.4 Initial considerations

As mentioned in the previous sections, the goal of this project was to establish growth rules for the dendrite spines of the neurons in Spiking Neural Networks (SNNs). Moreover, this rules were based on the knowledge provided by neuroscience and, more specifically, inspired by Hebbian learning.

Starting from this motivation and the initial hypothesis, a set of considerations regarding the scope and implications of the project were made, and they are summarized within this section.

1.4.1 Problem statement

The problem that this project dealt with was finding a set of rules that created intelligence in SNNs by dynamically modifying their topologies.

Let us formalize the problem by setting an SNN formed by k neurons, so the set of all the neurons in the system is . If neuron i has j input connections , whose values vary during time following the set of functions , there is a function dependent of that defines the value of the neuron’s output . Moreover, its value depends on how the signals have evolved through time:


The previously mentioned system has to comply with the following specifications:

  • The system is formed by 1 to n neurons, where .

  • The connections of the neurons may change during their lifetime i.e. The amount of input connections does not need to stay constant.

  • Every neuron follows the same rules, although they may be tuned with different parameters values.

  • The neurons work independently to each other, in parallel, and asynchronously.

  • The state of the neurons depend on the evolution through time of their different inputs (eq. 1.1).

1.4.2 The benefits of using Multi-Agent systems

This project deals with the design of a set of rules able to reproduce some of the behaviours observed in biological neural circuits. The final goal is to produce intelligence by creating SNNs that follow the aforementioned rules.

In order to implement these rules, a multi-agent system (MAS) was designed. A MAS is a system where the computation is distributed among a certain amount of units called agents. These agents possess some sort of intelligence, and they can operate both autonomously or they may require some degree of cooperation between themselves.

According to literature (Ferber and Weiss, 1999), the applications of MAS cover different purposes, among which this project puts focus on two of them:

  • Multi-Agent Simulation: This application of MAS creates simulations of the observed behaviours (mostly in natural systems), with the goal of reproducing and validating the theories and inferences made about those systems. In the case of this project, it was simulated systems reproducing neural circuits by implementing some of the rules established in the field of neuroscience. With the application of MAS is possible to modify and simplify deductions made in neuroscience, and therefore test if the reduced set of rules is able to produce intelligence and reproduce the behaviours observed in biological neural circuits.

  • Problem solving: This other application of MAS focuses on splitting and distributing into several agents the computational effort that has to be done for doing a task. This way, it is possible to obtain a better organisation of the tasks, and to optimize its solution.

1.4.3 The problem of using rules with time dimension for implementing discrete time logic functions

Traditional non-spiking Artificial Neural Networks (ANNs) tipically implement logic functions of different complexities, depending on the size of the network. This means simple ANNs can implement basic logic functions such as OR or AND gates.

Following a reductionist approach, the initial idea for developing and testing the neuron growth in SNNs was to develop rules that, once applied to a given network, were able to create a network which could implement simple logic functions.

Furthermore, it was also among the initial goals for the developed project that the obtained rules could be applied in a general way to any kind of network. Then, after the learning process, the network would end up with a certain topology that would implement the desired function based on the value of the input and output signals during the learning period. This means that the same growth rules applied to the same initial network conditions would end up generating different topologies based on the signal values during the learning process.

However, the neuron model implementing growth rules designed so far can not meet this last requirement for one of the most basic networks i.e. a 3-neurons network (see Fig. 1.5). According to the previous statement, the developed rules should be able to generate a network implementing both an OR or an AND function, depending on how the signals behaved during the learning process.

Figure 1.5: Spiking neural network formed by 2 input neurons connected to a third neuron. The blue dots represent the neurons’ somas, whereas the red ones represent the dendrites’ segments. The arrows indicate the information direction flow.

First of all, for a given static set of parameters, the previous network can only implement one specific function once it gets wired as shown, due to the absence of synaptic weights. For example, for a certain set of parameters the network provokes the output neuron to spike when at least 2 synapses were detected within a short period of time, as shown in Fig. LABEL:fig:leaky_time_plot. Furthermore, if the rate of the train of input spikes is lower, it would be necessary to get more input pulses in order to reach the threshold voltage. In any case, the plotted function would correspond to an activity detector which, without considering the time dimension, is equivalent to an AND gate i.e. the output neuron would only spike when the 2 input synapses trigger simultaneously.

On the other hand, if the parameter is modified so it is big enough to make the neuron spike after one single input synaptic pulse, then the network would become an OR gate if the time dimension is ignored.

The absence of specific weights for each synapse, plus the fact that one single synapse firing several consecutive times has the same effect as several synapses firing once at similar times, makes the neurons agnostic to which input is sending pulses. Actually, these neurons behave as activity detectors, as their outputs are functions of how many synapse pulses have been received in the short term, independent of which were the inputs providing these synaptic pulses.

Perhaps it is not valid to assume that SNNs can work as logic units at single time instants. In fact time dimension is crucial to understand their behaviour, and as spiking neurons implement functions that depend on the activity through time of their inputs, trying to implement classic logic functions that ignore time evolution of signals is somewhat preposterous.

Several spines connecting the same 2 neurons:

The designed model allows spines from a neuron to grow towards any other soma in its neighbourhood based on the Hebbian growth rules. Moreover, once a spine joins 2 neurons, new spines have the restriction of not being able to grow towards the same neuron again.

However, the drawback of not having synaptic weights may be solved by allowing several spines to connect to the same destination neuron. This way, when the input neuron spikes the membrane potential of the destination neuron would be increased by , where is the number of input spines from the same neuron.

Therefore, the rule limiting spines to grow towards the same neuron can be modified so the ANNs concept of synaptic weight could be applied i.e. more spines connecting the neurons increase the weight between both neurons.

Synaptic weight based on Hebbian-based synaptic strength:

Another alternative that can eliminate this problem is the usage of the Hebbian learning for establishing the strength of the synapse once the neurons are connected. If this alternative were to be implemented, the effect of electric pulse in the receptor would be diminished if the pulses of both neurons follow an anti-Hebbian pattern. On the other hand, neurons following a Hebbian pattern would imply a stronger synapse and therefore a higher increment of the membrane potential after incoming pulses in the synapse.

1.4.4 Neuron scope versus network scope

One main question that emerged during the early stages of the project was which should be the main direction that the development of the research should focus on. Whereas the design and implementation of a model of the neuron growth was the main target of the project, 2 directions in order to test the model and validate were depicted.

On the one hand, this project could be focused on testing the robustness of a single neuron in different situations and environment conditions. Since the outcome of the Hebbian-based growth is neuron spines growing towards other neurons with correlated firing sequences, this axiom can be validated by setting scenarios with pairs of neurons with different correlation levels in their spikes. Moreover, different levels in the environment noise can be tried out, as well as strong disturbances coming from specific spatial points.

On the other hand, the neuron model can be implemented in large networks by preparing tests with several neurons. Therefore, the performance of the model can be evaluated in different network layouts, and its behaviour can be assessed. Some benchmark networks with known outcomes can be prepared, so it can be measured the success rate of the model on the different scenarios.

In the end, a trade-off between both approaches was followed, though closer to the first option. Hence, it was mostly assessed the robustness of the Hebbian-based growth by setting single pairs of neurons and small pre-defined layouts. Moreover, its performance was evaluated within networks of bigger complexity as well.

The implementation of a neuron reservoir for testing the performance of the introduced rules js a very interesting test-bench. If successful, it can be a first step towards a tool that could be used in engineering for implementing intelligent systems able to solve logic tasks. Despite this may be an interesting later stage of the project, the current project is focused on the creation of a neuron model and testing the mathematical properties of the model, as well as its robustness against different situations. Moreover, this project also tries to get insight on how such model can benefit from the application of a MAS for deploying the growth rules.

1.5 Related work

1.5.1 Self-organizing maps

A self-organizing map (SOM), also named Kohonen map after their inventor, is a variant of unsupervised learning in ANNs. It consists in the mapping of an input vector into a N-dimensional output, typically a 2-D map. Through an iterative process, random input vectors are chosen, and the weights whose distance is closer to the input, as well as their neighbours, are modified

(Kohonen, 1990). In Fig. 1.6 is shown a graphical explanation of how this algorithm work. The weights between the input and the output layers evolve based on similarity between different elements in the input. At the end, inputs with similar values are clustered together and mapped to close output neurons.

Figure 1.6: Graphical example of a SOM. The n input variables are fed into the system, and these are mapped into the output layer. Image obtained from

Therefore, this learning creates ANNs where the outputs form clusters with their closest inputs. The distance between the weights and the inputs is normally calculated with the Euclidean distance.

The weights are not a means of deciding how to activate output neurons, but a feature defining the proper output neuron i.e. they are the value identifying the output neuron. This is, they define their relative distance to the inputs in an Euclidean space.

This family of ANNs is used for classification, as it allows to cluster inputs without former knowledge about them.

There are popular modifications to this algorithms that typically include a traditional ANN setup before the output layer for optimizing the values fed to the SOM (Thurau et al., 2003).

The main similarity between SOMs and the proposed SNNs is that SOMs follow to some degree the postulate of Hebbian learning: Although the concept of spiking does not apply to SOMs, the connection between their neurons evolve based on the similarity between their values, which is represented by the weights. In this specific ANN, the weights are a way of representing the characteristic feature of the output neurons, opposed to the traditional role of representing the connection strength.

However, there are several and relevant differences between SOMs and the proposed networks:

  • Despite SOMs get some inspiration from Hebbian learning, neurons in an SOM do not spike, and therefore the original definition of Hebbian learning can not be applied to them.

  • SOMs do a mapping between input and output neurons by evolving the weights between these two layers, whereas the proposed design modifies the structure of the network itself, creating and destroying connections between the neurons.

  • A SOM is mainly focused in the task of classification, whereas the proposed design serves as a learning mechanism for general purpose SNNs

  • SOMs are time-driven systems, and they do not take into account the time dimension, opposed to the networks implemented in this project, which are event-driven. In any case, this is one of the main differences of SNNs and ANNs, and it is already covered in the introduction of this document.

1.5.2 Recurrent neural networks

Traditional ANNs obtain output values based on the values of the input vector at a specific discrete time. This approach turns disadvantageous when the solution to a problem depends on the recent states of the system. For example, when analysing the meaning of a sentence in speech recognition, doing an assessment of the previous words makes the problem easier than evaluating the meaning of each word separately.

This shortcomming in ANNs led to the design and usage of Recurrent Neural Networks (RNNs). They are neural networks very similar to the vanilla ANNs, being their most relevant feature the existence of neurons that feed their output to their own input vector as well (Graves et al., 2008). Therefore, when a neuron in an RNN evaluates its new state, it takes into account its former state (See Fig. 1.7

). This is a way of giving memory to the neurons, and making their outputs dependent not only in the current input vector, but also in the states they were in previous steps. This architecture is normally called long short-term memory (LSTM).

Figure 1.7: Representation of a neuron feeding back its state in an RNN. It can be observed that on each step, a neuron takes as inputs both the input vector and its former value. Image taken from https://machinelearning-blog.com.

Although this approach may look similar to SNNs, they are very different approaches for solving a problem. It is important to realize that an RNN feeds back the former states of some of its neurons, and use that information as an additional input to the system. On the contrary, a neuron in an SNN does not feed back its output, but calculates its state based on the recent activity on all its inputs i.e. An SNN analyses the history of its inputs for calculating its output, whereas an RNN only uses the current value of the inputs, as well as the history of its previous outputs.

1.6 Outline

The rest of this document is divided in 4 chapters covering the different aspects relevant to the project.

In chapter 2 is offered a review of the background of the project. It has been divided in 3 main fields: On the one hand, it is summarized some of the knowledge acquired in the last century in the field of neuroscience that is related to the growth and learning mechanisms of the neurons. The discoveries regarding spine mechanics and neurotaxis are relatively recent, whereas neural learning is a topic that has been generating relevant literature since the early 20th century. On the other hand, there is a survey of SNNs, which cover the main trends in the last decades. Finally, the last section offers a brief overview of state-of-the-art multi-agent systems.

In chapter 3 is offered a thorough explanation of the proposed set of rules for neuron growth in SNNs, as well as an adapted spiking model that is compatible with those rules. The mathematical model is designed on the grounds of Hebbian learning, which is further developed by translating its results into actual spine movement. Fundamental mechanics are used for creating the kinematics of the spine, provoked by the forces applied by the environment. This set of rules depends mostly on the firing times of the neurons, reason why a spiking model was needed. It was chosen the Leaky integrate and fire model, which was modified so it could be adapted to the designed growth rules. Finally, it has been included an overview of the application of multi-agent systems for the project, explaining the different types of agents implemented and how they interact with each other.

Chapter 4 introduces the different experiments that have been performed with the designed model, which basically consist in different network layouts that allow to assess different aspects of the model performance. Afterwards, it is offered an analysis of the different parameters of the model, showing what is their impact on the system performance, and how they could be altered for obtaining different results.

Finally, chapter 5 offers a final conclusion of the project, which further work can be done in the future, and what can be expected from the designed model.

2.1 Neuroscience theory

2.1.1 The learning process and the role of connections growth

Hebbian learning was originally an hypothesis that explained that the synaptic connections between neurons get stronger or weaker depending on the timing correlation between their spikes. Furthermore, in recent years new in-vitro and in-vivo experiments have given evidences towards the existence of a Hebbian alike rule involving the growth of neurons’ spines based on their activity (Feldman and Brecht, 2005). This process is often also named structural plasticity.

The structural plasticity in the mammalian brain is dominated by different mechanisms that operate in different time scales (from minutes to weeks and months), and that span from single spines and axon boutons to entire dendritic arbors. This process is believed to be tightly related to the formation of memories and the process of learning (Holtmaat and Svoboda, 2009).

The long structural processes involving the creation and shaping of axons and complex dendritic arbors are related to the activity that takes place during neurogenesis (the process of creation of neurons in animals) and injury recovery. Despite their obvious relevance for understanding the operation of the brain, the rapid experience-dependent structural plasticity is more relevant for establishment of growth rules applicable to Spiking Neural Networks. These rapid dynamics typically correspond to dendritic spines growing towards axon boutons (Holtmaat and Svoboda, 2009).

2.1.2 Hebbian learning and evolution of the synapse strength

Fundamentals of Hebbian learning in spiking neurons:

There are two main different models for the neuron assessment of incoming pulses: rate-based or time-based. The first one establishes that the connection between two neurons will be stronger when the spiking rate of both is similar, whereas the second one calculates the average during time slots an evaluates each spike during that time slot.

The first model is easier to describe and implement in a computational system. However, it is considered to have some limitations that makes it insufficient for explaining learning in some situations e.g. the visual response time for many animals is less than 200ms, which makes it incompatible with rate-based input analysis, as neurons would not have enough amount of information for evaluating the input spike rate.

Literature supports the idea of biological neurons being more similar to time-based model, as they transfer information through spikes and the strength of the connections between depending on the time correlations between their spikes (Kempter et al., 1999). Eq. 2.1 is proposed for calculating the efficacy of a synapse of a neuron with another presynaptic neuron i,


where is the firing time of neuron i, is the firing time of the actual neuron, and W is the learning-window function. is a very small parameter that makes the learning evolve much slower than the actual network dynamics. The parameters and depend on .

Opposite to what is implemented in a typical ANN, this equation models neurons with binary outputs i.e. as long as the neurons are connected, the individual incoming pulses have always the same effect on the postsynaptic neuron.

Spike-timing-dependent synaptic modification induced by spike trains:

In (Kempter et al., 1999) it is assumed that the contribution of each pulse in a train of pulses in a neuron will be independent, and based on a general strength variation rule i.e. a formal representation of the Hebbian rule.

However, the authors in (Froemke and Dan, 2002) prove with in-vivo experiments that the influence of a pulse is strongly determined by the presence of previous pulses in the same neuron.

Hence, if two pulses happen within few difference in time, the effect of the second one gets considerably diminished, and the final change in the synapse strength will be mostly determined by the first one. This is formally presented by applying an exponential decrease to the effect of a pulse depending on the time difference with the preceding pulse. It is therefore introduced the neuron efficacy , calculated by using eq. 2.2,


where and are the timings of the and pulses of the neuron i, and is a time constant.

The efficacy is applied to the general equation eq. 2.3 for calculating the variation of the synapse strength between neurons i and j,


where is the function for calculating the spike-timing-dependent plasticity (STDP) between both neurons, and calculated by using eq. 2.1.2,

where A is a scaling factor, and a time constant.

2.1.3 An overview of the forces driving neurotaxis

In (Dickson, 2002) is offered a review of the known mechanisms of the axon guidance at that time (2000). On the one hand, it offers a review of the different guidance chemical cues present in the axon environment, stressing the fact that single cues may have more than one role in the axon guidance. There are four known families of guidance cues:

  • Netrins: They have the ability to both attract and repel the axon, depending on the receptors present on it. Their effect can span from the short range to millimetres.

  • Slits: They are large proteins that act as a repellent for certain receptors (Roundabout receptors, or Robo), hence they establish the borders of the axon growth. They also work as stimulants of the axon sensory axon branching and elongation

  • Sempahorins: These molecules are divided in 8 classes. They work as a short-range inhibitor for the growth, although it seems they may also work as attraction cues for some receptors. Their main role seems to be the avoidance of inappropriate cells contacts

  • Ephrins: They form molecular gradients that lead to the topographic order of an axon, though not its precise end. Thus, they signal the direction of the growth

The growth cone is formed by actin parallel oriented filaments (fillopodia), and an intervening networks of filaments, whereas the growth is directed by the extension and contraction of the microtubules. The turning in the growth cone can be explained by the signaling produced by different proteins contained in the cone structure. Moreover, it has been observed gradients of Calcium ions in the fillopodia that can lead to the turning of the cone.

One of the key features of the growth cone is the plasticity of its properties, which allows to react in a different way to the chemical cues, depending on the stage of the growth. This plasticity is obtained by at least three mechanisms:

  • Modulation by cyclic nucleotides: Inhibiting, lowering or rising the levels of certain proteins (cAMP or cGMP levels, or PKA or PKG) provokes the attraction or repulsion to certain cues.

  • Local translation in the growth cone: The translation of certain molecules through the axon provokes the synthesis of certain proteins. Blocking this translation, and hence inhibit the synthesis, inhibits as well the turning of the cone, although not the growth.

  • Switching responses at the midline: The axons may change their sensitivity to the cues after reaching and passing through them

“The ultimate challenge, after all, is to find out how a comparatively small number of guidance molecules generate such astonishingly complex patterns of neuronal wiring.” (Dickson, 2002)

2.1.4 Mechanics of the growth cone

The trip of the tip: Understanding the growth cone machinery

In (Lowery and Van Vactor, 2009) is offered a review of the chemical mechanisms and factors that influence the movement of the growth cone, which is the tip of the axon in a neuron. Using the metaphor of a road trip for explaining the whole system (See Fig. 2.1), the growth cone is considered a vehicle that has to drive through a roadway (adhesive substrate-bound cues), delimited by guard rails (repellent substrate-bound cues), and follows road signs for deciding the path (diffusible chemotropic cues). They also distinguish two main functions of the system i.e. the vehicle, which deals with the motion mechanisms of the growth cone for keep growing; and the navigator, which is responsible for deciding which path to take.

Figure 2.1: Representation of the main chemical cues interacting with the growth cone, and their function. Image obtained from (Lowery and Van Vactor, 2009).

Regarding the mechanisms behind the growth cone (the vehicle), there are three stages that are repeated continuously and that make the growth cone to progress on its trip. At first, during the protrusion stage, the filopodia and lamellipodia extend forward. Secondly, during the engorgement stage, the main body of the growth cone moves forward following the filopodia. At last, during the consolidation stage, the shape of the axon shaft is formed again. This all is achieved by a molecular clutch model, that allows the cytoskeleton to get anchored to the adhesive substrate (This is achieved thanks to the properties of a family of proteins called actin, which comprise the cytoskeleton of the cells). During this process, filopodia are considered to work both as exploration sensors and points of attachment. Fully understanding the clutch mechanism is paramount for being able to understand the overall logic governing the progression of the growth cone (Lowery and Van Vactor, 2009).

Regarding the navigator, there are two main components ruling the process of steering the growth cone and finding the adequate path. First of all, the aforementioned interaction of the actin structures with the chemical cues working as road signs. Depending on the chemical properties of the axon (for example, the presence of different kinds of neurotransmitters), the actin will be attracted or repelled by the different kinds of chemotropic cues (the road signs). Secondly, the microtubules present along the growth cone seem to have a very important role in the steering of the axon. Their polymer structure makes them show a dynamic instability. Thanks to this property, they act as sensors during the protrusion stage, interacting with actin cues and steering the growth cone towards the correct direction; but they act as well as inhibitors, granting stabilization against the guidance cues and acting as a scaffold for guidance of cue signaling.

The growth cone structure and mechanisms offer several alternatives to be considered when developing the rules for a model of the evolution of neuron topology. From the perspective of a multi-agent system, and given the importance of different elements in the growth and steering of the axon, some of these elements may be modeled as separate agents: A multi-agent model of the neuron may include soma agents (maybe another type for the dendrites), one or several axon segment agents, a growth cone agent, several filopodium agents, several microtubule agents, and several actin bundle agents.

Moreover, it is clear that the environment has a very important effect on the growth cone, and some of its properties could be included in the model of the MAS map. At least 4 different types of chemotropic cues are mentioned, which are crucial for determining the route the growth cone will take, due to their interaction with the growth cone elements.

However, several questions without an easy answer appear, mostly related to the how. How is determined the location and intensity of the chemical cues? How are the attraction and repulsion rules between the different chemical elements quantitatively decided? According to (Lowery and Van Vactor, 2009), the overall logic that governs this process is still emerging.

2.2 Spiking neural networks

Spiking Neural Networks (SNNs), which are often called 3rd generation neural networks, were developed under the idea of creating a system able to reproduce the behaviour of neural circuits. Although the ANNs were created as well as an attempt of replicating some of their behaviour, advances in neuroscience soon proved that the assumptions made for ANNs were far from the reality of biological neural circuits (Paugam-Moisy and Bohte, 2012).

The most innovative idea behind SNNs is taking into account the time evolution of input spikes for calculating the state of neurons. This means that SNNs are based on the evolution of the inputs in the time dimension, and they get excited when enough spikes have been recently received in their inputs. This idea also implies that the neurons in SNNs are event-driven computing units i.e. the computation is performed when certain events occur, opposed to time-driven processing, were the computation is performed at constant time intervals.

Despite there are several models for spiking neurons, most of them make use of the concept of the membrane potential. When spikes arrive to a neuron, itsmembrane potential increases and, depending on the model, the neuron will reach an excitation state after some time as a function of this potential. Finally, this excitation state entails the generation of a spike in the proper neuron, which is transmitted by its axon to other neurons connected to it. Moreover, this process involves time delays, and many models include stochastic processes as well. They are a way of representing a plethora of phenomena occurring in the neuron and its environment, such as the amount of and sensitivity to neurotransmitters, as well as the amount of them which can travel to the synapses; or the presence of chemical cues around the neuron, which can boost or inhibit the synapse. In any case, the modification of the membrane potential after an incoming spike is normally referred as postsynaptic potential, which can be excitatory (EPSP), or inhibitory (IPSP) (Paugam-Moisy and Bohte, 2012).

In any case, biological neurons get electrically charged after receiving spikes, and produce new spikes asynchronously based on the amount of received spikes. Neurons in an SNN work in an analogue way, where the timing between spikes is the most important way of transmitting information.

2.2.1 Firing models

The firing model is a crucial part of an SNN design, as it defines what is required for neurons to spike, and how they will behave after a spike occurs. The two main features that are evaluated in a firing model are their similarity with biological neurons and their computational simplicity. In the rest of this section the most popular firing models are introduced. In any case, many more have been designed, most of which derive from the following ones or combine some of their features.

Hodgkin-Huxley model:

The Hodgkin-Huxley model represents the spiking of the neurons by an electrical circuit consisting in a capacitor representing the capacitance of the neuron membrane, and three parallel conductances representing the different ion channels (potasium, calcium, etc.). Moreover, they are in series with batteries representing the equilibrium potentials of these channels. It was initially developed for modeling the behaviour of the squid nervous system, after a series of experiments during the 1940s. A generalized differential equation proved to be highly accurate for describing the properties of neurons’ action potential (Nelson, 2004).

Despite being very faithful to the biological neuron, this model is very complex and rather complicated for being implemented in an SNN with a big number of neurons. In Fig. 2.2 is depicted the electrical schematic of this model.

An important feature of this model is the existence of a temporary refractory time after a spike, where the occurrence of a second spike is very unlikely to happen. This time typically consists of few milliseconds.

Figure 2.2: Electrical schematic of the Hodgkin-Huxley model. Image obtained from (Nelson, 2004).
Leaky integrate and fire model:

The Leaky I-F model consists in the representation of the soma membrane as an R-C electrical circuit (See Fig. 2.3). It is derived from the Hodgkin-Huxley model, where the incoming electric pulses contribute to charge the capacitor in the model. Furthermore, it slowly discharges over time until reaching the rest voltage.

Figure 2.3: Electrical circuit of the Leaky I-F model and shape of the intensity and voltage plots when a pulse is received. Image obtained from (Gerstner et al., 2014).

The occurrence of a spike is determined by a threshold value. Whenever the membrane potential raises over that value, the neuron produces a spike. Immediately after that event happens, the voltage drops back to the rest voltage value. The existence of a refractory time is optional in this model.

Izhikevich model:

The Izhikevich model (Izhikevich, 2003) is one of the most suitable models for neuromorphic engineering, as it grants a remarkable trade-off between computational complexity and resemblance to actual neurons.

By using the differential equations shown in eqs. 2.2.1 and 2.2.1, this model is able to reproduce several different spiking patterns, which makes it possible to adjust to the different dynamic behaviours present in biological neurons. This is achieved by changing the values of the variables a, b, and c

This model adds complexity that is not needed in this project and some of its most relevant features would not be used for testing the neuron growth, such as the adaptation to several different types of neuron spike trains. Fig. 2.4 contains four different spiking patterns that can be obtained with eqs. 2.2.1 and 2.2.1.

Figure 2.4: Example of four different spiking patterns that can be obtained with the Izhikevich model, by varying the variables a, b, and c. Image obtained from (Izhikevich, 2003).

3.1 The neuron MAS model

A MAS model has been designed in order to build up a system that implements the behaviours and the designed set of rules that are explained in detail in the following sections. Therefore, the computation is distributed in several agents, and the final goal with this system is to distribute the intelligence between several computational units.

It has been developed a neuron model with dynamic dendrite spines that grow towards neighbouring neurons based on the correlation between their firing times. This model has been implemented and tested by using a multi-agent approach, where the entity representing a neuron is formed by 3 different types of agents.

Soma agent:

There is one soma agent per neuron, and it is the top agent in the neuron’s hierarchy tree i.e. it is the parent of the rest of the neuron’s agents. The application of the spiking rule is implemented in this agent i.e. the Leaky I-F model, as well as the stochastic process for deciding whether a spike occurs or not.

The spiking algorithm requires to gather the data of the incoming input pulses, which are received through the neuron dendrites. Both the intensity and the time at which these events happen are needed. They are fed to the Leaky I-F algorithm in order to obtain the soma membrane potential, which also takes into account the noise (See eq. 3.18

). This value is used for calculating the probability of a spike event happening, by following the next steps:

  1. Normalization of the membrane potential, being 0 the rest voltage , and 1 the threshold voltage (eq. 3.13).

  2. Calculation of the spiking probability, by feeding the normalized voltage to the Sigmoid function (eq.


  3. Decision of whether the neuron spikes or not, by using the previously calculated probability by the Bernouilli algorithm (eq. 3.4.2).

Moreover, this agent spawns a growth cone agent at the start of the simulation, and also after the current growth cone connects to another neuron. However, in some networks this functionality has been limited, for the sake of avoiding undesired complex structures i.e. some neurons can only spawn a limited amount of growth cones and thus dendrite spines.

This agent communicates certain information to other agents by emitting asynchronous events:

  • excited_neuron: The targets of this event are the growth cone agents of the same neuron. It is emitted when a spike happens in the soma, and it is required by the Hebb’s rule implemented in those agents.

  • electric_pulse: The targets of this event are the growth cone agents from different neurons in the neighbourhood. This event represents the propagation of an electric pulse through the environment when the current neuron gets excited.

  • assign_group: This event is used for informing the children agents of the soma the identity of their parent soma.

  • firing_time and stop_growth: These are auxiliary events used for recording the times at which the firings occur for later analysis, and to prevent the neuron’s cone to keep growing.

Growth cone agent:

This agent represents the tip of a neuron spine, and it contains the logic dealing with the Hebbian based spine growth. In order to apply Hebbian learning, the next steps are followed:

  1. Record the time at which the neuron’s soma is excited. Only the last excited_neuron is relevant for the algorithm.

  2. Detect electric_pulse events occurring in the environment, and record the time when they happen.

  3. When each of the previous events happen, calculate the time difference between both.

  4. Calculate the STDP (eq. 3.2) in the neuron i.e. the increase in the membrane potential.

  5. If the growth cone is not connected, calculate its acceleration and velocity by using the designed kinematics algorithm (see eqs. 3.3.1 and 3.6).

Therefore, it is sensitive to electric pulses traveling in the environment which, altogether with the firing state of the current neuron, are paramount for establishing the spine growth direction.

This agent communicates relevant information by emitting the following asynchronous events:

  • synapse: Its destination is the parent soma agent of the growth cone. It is used for transmitting received synaptic pulses to the soma for further processing.

  • cone_init: Event indicating that the growth cone has been correctly spawned and will start its normal operation.

  • cone_connected: Event for informing the parent soma agent that the cone has reached a destination. The soma will normally spawn a new growth cone after receiving this event.

  • cone_parent: and cone_kinematics: Auxiliary events used for recording the reached destination neuron and the historic velocity and acceleration values for further data analysis

Spine agent:

This agent represents a link in the spine of the neuron, and a new one is generated when the spine grows a certain length. Its purpose so far is merely graphical, as there are no functions associated with this agent.

However, further improvements of the model shall give intelligence to this agents, so they can react to electric pulses present in the environment, and eventually allow a dendrite to fork in more than one direction.

3.2 Spine growth rule

As explained in the introductory chapter of this document, the main goal of this project was to create intelligence by using a set of rules that make the structure of SNNs to evolve through time. Therefore, the first step for creating such model was to establish a rule for the growth of neuron spines.

In order to establish a rule for making neurons to grow towards each other, the Spike Time-Dependent Potential (STDP) process has been used, which is a reformulation of Hebbian learning,and it is described in (Froemke and Dan, 2002). According to them, is the function for calculating the STDP between 2 neurons, and it is calculated by using eq. 3.2,

where is the time difference between the spiking of both neurons involved, A is a scaling factor, and a time constant.

In Fig. 3.0(a) are depicted the results of an in-vivo experiment regarding the change in the excitatory postsynaptic potentials (EPSC) of a neuron in relation to the spiking time difference between the input pulse and the trigger of the neuron, and in Fig. 3.0(b) is depicted the plot of eq. 3.2. It can be observed that the shape of the implemented equation highly resembles the results obtained with real neurons.

Figure 3.1: a)Results of the in-vivo experiments obtained in (Bi and Poo, 1998), for the change in EPSC amplitude (in %), plotted against the time difference between spikes, and b) plot of the mathematical equation representing the STDP (eq. 3.2), for , and .

Furthermore, eq. 3.2 was modified so it takes into account the distance vector between neurons. This way they are more influenced by close neighbours than by distant ones. This modification has been inspired on the electric field equation, as electric forces are one of the main phenomena explaining the interaction between neurons. The resulting force vector is represented in eq. 3.2 for assessing the attraction force of one neuron towards a second one when both are triggering at similar times. Therefore, the attraction force of neuron j over neuron i is determined by


where is the EPSC amplitude calculated according to eq. 3.2, is the absolute distance between the growth cone of neuron i and the soma of neuron j, and is the unit vector joining both elements.

Implementing the negative side of the Hebbian rule

In the previous section the equation modeling Hebbian learning was introduced, which results in the variation of the synaptic strength between two neurons based on their firing times. This model is relatively simple to implement for modeling the positive growth of the neurons’ spines.

However, in order to have a full model of the Hebbian structure plasticity, the negative side of the Hebbian rule has to be implemented as well. Eq. 3.2 can be used for obtaining the negative STDP value corresponding to a presynaptic neuron firing after the postsynaptic one.

Obtaining the final repulsion force is not straightforward, and different approaches can be taken into consideration (see Fig. 3.2). The problem is, all of the shown approaches offer at least one major drawback:

  • The first approach would provoke the growth cone to travel in an unpredictable and unrealistic direction, moving away from both its soma and the second neuron

  • The second approach would create a very long spine doing a loop, although this could be fixed by, for example, making the spine to decay and ending up disappearing.

  • In both the second and third approaches the growth cone would not necessarily take a direction opposite to the second neuron. It could actually get closer to the neuron if it is located between the growth cone and the soma of the first one.

Figure 3.2: Basic approaches of the direction that the growth cone can take when it suffers a repulsion force from a second neuron. In a) it moves in the opposite direction to the second neuron, in b) it moves back towards the neuron’s soma, whereas in c) the cone goes backwards in the same direction it came from.

A fourth approach would be to implement a decay function which makes the attraction force between two neurons to get weaker. Therefore, for a system with two neurons, if there is a repulsion detected the force between them would start to decay until reaching zero. From that point, the spine would start to get weaker until it disappears.

Still, there would be gaps in this model, as the behaviour when the growth cone is influenced by more than one neuron would not be realistic and thus it needs to be further developed.

In any case, this was left as future work, and hence the negative side of Hebb’s rule has not been included in the design.

3.3 Spine growth dynamics

In order to be able to implement complex network topologies, the neuron model needed changes in some important aspects. Otherwise, it would not be possible to create the necessary connections for generating many networks.

First of all, the model so far describes neurons whose spines follow irregular trajectories, as they are considerably affected by random noise and stranger incoming pulses (experiments showing these effects are shown in chapter 4). Therefore, it is paramount to implement more stable dynamics for the spines in order to get robust connections between neurons. This problem is addressed in the subsection 3.3.1, where the concept of the spine acceleration is introduced for creating inertia in the agents.

Second of all, when the velocity is obtained from a constant acceleration, the obtained velocity has to be bounded, otherwise it would grow to the infinity, being a source of instability in the system. Hence, a drag force opposing the movement of the spine is implemented in the model. It is described in subsection 3.3.2.

3.3.1 Getting smoother trajectories

A drawback observed in the experiments was the high sensitivity of the spine growth to incoming pulses. Despite incoming white noise would not considerably affect the success of reaching the final destination, some disturbances provoked big spikes in the trajectory of the spine e.g. stranger points provoked sudden changes in the path of the spine.

In other words, it is desired for the spine to follow smooth trajectories, and to be robust against noise and undesired incoming pulses. The following are some alternatives that could achieve this:

  • Use a correction function that reduces the resulting attraction force if its direction deviates from the current neuron velocity. For example, using the cosine of half of the angle between the velocity and the force would reduce the effect of the force the more the force deviates from existing trajectory.

  • Using the attraction force for getting the second derivative of the position, instead of the first derivative. Said in other words, using the force for calculating the acceleration instead of the velocity.

  • Use the second order derivative equation typical of spring-mass systems for calculating the new position.

The first option presents the shortcoming that forces would tend to be ignored the more deviated they are from the current velocity.

A solution based on the second option has been included in the design, and is explained in the following subsection.

Movement based on the second derivative

By using this approach, the spine acceleration is obtained from the incoming electric force. This makes sense from the physics point of view, as the acceleration of an object is directly affected by the force that is applied to it according to Newton’s second law of motion (eq. 3.4).


Thus, last equation was used for calculating the spine acceleration from the electric force value obtained from eq. 3.2.

As it was explained in the previous section, the current velocity is used for calculating the new position of the agent by applying simple kinematics. The model has been modified by calculating first the acceleration with the previous equations and from there the instantaneous velocity of the agent (See eq. 3.3.1).

These equations were tested on a delay line (for details, see section 4.3). Fig. 3.3 shows the resulting layout. It can be seen that it is satisfactory, as the spines follow smooth trajectories. The experiment also shows that this model is robust against the presence of a stranger point emitting electric pulses with a given frequency.

Figure 3.3: RANA simulation of a delay line initialized with 8 neuron somas (red dots) using the neuron model defined by eqs. 3.3.1. The experiment was performed without environment noise, neither trigger Poisson noise, but a stranger pulse generator was included in the bottom right corner. The kinematics of the marked neuron are plotted in Fig. 3.4. The experiment is pseudo-deterministic, and thus the result is always the same.

In Fig. 3.4 are shown the plots of the acceleration and velocity of one of the neurons in the experiment. It can be observed that the pulse generator creates single pulse disturbances in the spine acceleration at 2 different points in the plot (Around iteration 800 and iteration 5800). However, the effect of these disturbances can not be observed in the plot of the velocity of the agent.

Furthermore, the result also shows that, in the sole presence of one attractive source, a spine will tend to constantly accelerate towards it, which will provoke the velocity to grow towards infinity. This is an undesired behaviour, which was dealt with by introducing a drag force. This is addressed in the subsection 3.3.2.

Figure 3.4: X and Y components of the velocity(red) and acceleration (blue) of the agent marked in Fig. 3.3. The plotted time spans from the simulation start until the dendrite reaches the destination soma.
Robustness against noise:

The new kinematic model was tested under different noise parameters in order to assess its robustness. Environment noise was added in the system, by following the same approach as in the former model i.e. The noise has an intensity that follows a Gaussian distribution with zero mean, and the direction of the noise follows a uniform distribution between 0 and 360 degrees.

The spines are still able to reach the desired target, even with a considerable amount of noise (see Fig. 3.5). The velocity and acceleration plots of one of the agents show that the noise is considerably big compared to the attraction of the other neurons, but it is still capable of reaching the target without being observed any undesired deviations from its path i.e. the spine has a certain inertia that prevents it from deviating from the current path.

Figure 3.5: Outcome of the bottom spine agent in the delay detector whown in Fig. 3.3, with an environment Gaussian noise of , and . Plot of a) the instantaneous acceleration of the agent in blue, and its average value in a window of 70 neighbouring samples in green; and b) the value of the instantaneous velocity in red, together with the filtered value of the acceleration. The acceleration was filtered only for visualization purposes.

3.3.2 Drag force

As mentioned in the previous section, the spine agents in the model show a constant and limitless increase in their velocities. This leads to cases where the spines travel too fast and can not manage to correct their trajectory and thus can not reach in some cases their desired goals. This can be observed in Fig. 3.6, where half of the spines grow passed their target somas.

Figure 3.6: RANA simulation of a delay line where four out of the eight neurons overgrow and do not reach their target goals due to the high inertia they have when they approach their destination.

In order to minimize this issue, the calculation of the resulting attraction force over a spine growth cone includes now a new component opposing the current movement i.e. a drag force, whose modulus is proportional to the current velocity and has opposite direction.

The drag force is a concept used in physics for describing the opposition of a fluid to the movement of an object through it. Its value is proportional to the object velocity, the fluid density, the cross section area, and the drag coefficient. For the sake of simplicity, the last 3 terms have been grouped under one single term, as all of them are considered constant in the current context. It will be named drag coefficient, , from now on. Thus, the drag force is obtained by using eq. 3.6,


Due to the lack of resemblance with typical fluid mechanics problems, the value of has been tuned empirically based on the desired kinematics of the spines, instead of using typical values used in fluid mechanics. In fact, this model is a simplification of the biological neuron, so the used used values are the ones that fit the best for achieving the desired networks.

Biological neurons do not freely navigate through a fluid inside the brain. Actually, their movement depends on the mechanic adhesion to the substrate, the forces created by the attraction and repulsion of surrounding chemical cues, and the mechanisms of the spines for growing. These are based on the contraction and extension of the fillopodia inside the neuron, as well as the alteration of the proteins distribution inside it (Lowery and Van Vactor, 2009).

Based on the previous experiments, it has been tried to stabilize the speed of the spine growth around when it is attracted by another neuron. Looking at the results plotted in Fig. 3.4, when the neuron has a velocity of the acceleration has an average value of . Thus, in order to stabilize the acceleration around this point, the drag force has to cancel out the attraction force. It is then obtained a drag coefficient of (Setting and in eq. 3.4, and solving eq. 3.6).

To sum up, the force acting on a spine growth is calculated using eq. 3.7, and entails the addition of 3 different components, which are the attraction forces from the rest of the neurons in the environment, the force provoked by the different sources of noise, and the drag force.


The previous equation was implemented for the same experiment (the delay) line. In Fig. 3.7 are depicted the plots of the obtained velocity and acceleration after implementing the drag force. It can be observed that the magnitude of the velocity gets limited, because the acceleration gets reduced when the velocity increases. The acceleration plot has been zoomed in compared to the previous experiment due to its smaller values. It can be clearly observed the stranger electric pulses in the acceleration, and how small their effect on the velocity is. In general, it can be observed that the introduced changes make the velocity more robust towards fast changes in the acceleration.

Figure 3.7: Result of the delay detector network with the second-derivative-based model, and the presence of a drag force, calculated with eq. 3.6 for a drag coefficient of . The experiment was performed with a stranger pulse generator, but without environment noise nor trigger Poisson noise. In b) is plotted the velocity and acceleration of the marked agent. The result of the simulation is deterministic.

The introduced changes have increased the robustness of the model against noise. In order to prove this, an experiment with the same delay line with an environment noise with variance

has been conducted. In Fig. 3.8 is depicted the result of una of the simulation runs, where the spine agents perform very irregular trajectories along the map.

Figure 3.8: Simulation of the delay line with a white Gaussian noise with a variance of . The drag force was calculated using eq. 3.6 for a drag coefficient of . The success rate after 1000 seconds of simulation is 0.57.

The success rate was 0.57 (i.e. 57% of the spines connected to the desired destination soma). It can also be observed that an unconnected spine was very close to reach its destination, which would have increased the success rate to 0.71. From the previous experiment, it is obtained that the average acceleration of a spine for reaching another soma oscillates around . Therefore, the noise variance is times bigger than the attraction force for reaching the desired soma. It is noticeable that the network obtained a high success rate given such big difference in the orders of magnitude of the noise and the attraction force.

3.3.3 Neuron agents with multiple spines

The design so far deals with neurons that consist in one soma agent from which only one spine agent can grow. This is a big limitation for the development of meaningful networks, as the implementation of logic functions and complex tasks involve as well a complex tree of connections between neurons.

The model described until this point has been applied for the delay line and the coincidence detector introduced in sections 4.3 and 4.4 in the next chapter (in Fig. 3.9 the result is depicted). It can be observed that the rightmost neuron - which corresponds to the output neuron - connects to one of the input neurons, whereas the spines of the 2 input neurons wander around, as they are basically ruled by the environment noise.

As the current model does not allow neurons to develop more than one spine, the network can not be completed.

Figure 3.9: Result of the coincidence detector network with the second-order derivative model, and the presence of a drag force, calculated with eq. 3.6 for a drag coefficient of . The experiment was performed with environment noise with , but without trigger Poisson noise. In b) is plotted the velocity and acceleration of the rightmost spine. The acceleration has been filtered by a mean filter with a kernel size of 2000, only for visualization purposes.

Hence, the model was modified in order to allow neurons to grow more than one spine, so they can have multiple input signals. The inception and growth of more than one spine has been modeled by the following rules:

  • A soma can generate a new spine if and only if all of the existing spines are already connected

  • Two spines of the same neuron can not grow towards the same destination.

  • A soma can only generate a finite number of spines. Although a function for limiting its number shall be included, this goal can be achieved as well naturally due to the properties of the environment.

These rules have been implemented in the neuron existing model, excepting the last rule. This is due to the fact that the spawning of new spines gets naturally limited by the amount of neighbouring neurons that are actually affecting the current neuron. This shall be revisited in future work, but it proved to be good enough for the current problem.

In Fig. 3.10 is depicted the result of this changes in the model for the coincidence detector model. It can be observed that the result is satisfactory, as the output neuron connects two spines to both input neurons.

Figure 3.10: Result of the coincidence detector network after using the rule for spawning more than one spine per neuron, with an environment Gaussian noise of . Thus, the rightmost neuron spawns a second spine after the first one gets connected. The input neurons’ spines wander around due to the environment noise

3.4 Neuron firing model

3.4.1 The Leaky integrate and fire spiking model

A fundamental topic regarding SNNs is determining which is the rule deciding when neurons spike. This has been a topic of interest for neuroscience since the beginning of the 20th century, and a plethora of models have been proposed so far (Gerstner et al., 2014). It is important to realize that despite biology prioritizes the resemblance with real neurons, computer science and neuromorphic engineering place much emphasis in the computational cost of implementing it.

Contrary to the Izhikevich model, the Leaky integrate and fire (I-F) model is a simple approach to the behaviour of the neuron membrane potential, where it is represented as an R-C electrical circuit (See Fig. 3.11). Therefore, the neuron membrane behaves as an electrical capacitance with a leak resistance towards its body, which has a characteristic rest voltage. External electric pulses are represented as intensity functions that charge the capacitor. Once the system voltage reaches a certain threshold, the neuron triggers a spike and the voltage is reset.

Figure 3.11: Leaky I-F model, based on an electric RC circuit modeling the neuron membrane. On the right is shown the shape of the membrane intensity and voltage. Image obtained from (Gerstner et al., 2014)

The instantaneous voltage potential of the neuron membrane is calculated by using eq. 3.8, which corresponds to the Leaky integrate and fire rule (Gerstner et al., 2014).


In the previous equations, is the instantaneous membrane potential, is its steady-state potential in the absence of external pulses, is the potential applied to the membrane due to an external pulse, the time since the pulse arrived, and is the RC circuit time constant, which is calculated with eq. 3.10.


The values of the model parameters have been initially chosen based on measurements done on real neurons. Namely, for cortical regular spiking pyramidal neuron, , , , and (Liu and Wang, 2001). From the values of and a value of is obtained by using eq. 3.10.

In Fig. 3.12 is shown the result of a test where the Leaky I-F model is applied. The central neuron stops getting externaly excited once it gets wired to the 2 input neurons. Moreover, the rightmost neuron is inactive until such event happens. Therefore, the central and right neurons get connected solely due to the spiking rule that is represented in eq. 3.8.

For an incoming electric pulse, the voltage increment in the membrane has been arbitrarily decided to be . Therefore, for reaching the threshold voltage, the neuron would have to receive at least 2 incoming pulses not very separate in time, allowing thus to replicate the behaviour of a coincidence detector.

Figure 3.12: Result of a coincidence detector network when using the Leaky I-F model for making central neuron to spike once it is connected to the input neurons. Contrary to the other connections, the connection between the central and output neuron is achieved by using the I-F model, as their excitation is not forced by external events.

3.4.2 Stochastic model for triggering the neuron

The Leaky integrate and fire model has been introduced in the previous section. It models the neuron as an RC circuit, where the instantaneous value of the membrane potential determines whether the neuron will trigger or not.

However, literature in neuroscience supports the idea that the triggering of the neurons is actually led by a stochastic process (Dayan and Abbott, 2001). According to this idea, the triggering times of the neurons are non-deterministic, and the membrane potential only contributes to increase the probability of spiking.

Therefore, a probabilistic function dependent of the neuron membrane potential has been introduced into the triggering model. As the outcome of this probabilistic function can only be 1 or 0, it has been used the Bernouilli distribution for modeling this behaviour (eq. 3.4.2),

where F is the existence of the firing event, and the expected value is calculated with the neuron’s membrane potential.

When the membrane potential has values closer to the rest potential , the probabilities of the neuron firing are low. On the other hand, when the membrane potential reaches values close to the threshold potential , the probability of the neuron firing is high, being the spike very likely to happen within few iterations. In order to represent this behaviour, it has been used the sigmoid function. It is a continuous function, which asymptotically grows towards 0 and 1 without the need of introducing artificial boundaries, and it covers the intermediate points by describing an S shape. It is represented by eq. 3.12, and its plot can be observed in Fig. 3.13,


where k is the growth rate of the Sigmoid function, and is its middle point. By modifying these 2 values, the function can be drifted towards one of the sides, and it can be made to grow faster towards the asymptotes.

Figure 3.13: Plot of the Sigmoid function S(x), compared to the membrane potential values. S(x) was plotted after normalizing between and , it was centered around the middle point of both values, and its growth rate was set at . The red vertical lines represent the voltage potentials and .

Moreover, the function was fed with normalized values of the membrane potential, by using the conversion formula shown in eq. 3.13. The values used for and are the membrane rest potential and the membrane threshold.


The result of the Sigmoid function (eq. 3.12) is the probability of the neuron to spike. Hence, it is the probability of obtaining a 1 in the Bernouilli process represented by eq. 3.4.2.

Tuning the stochastic model:

In order to implement the aforementioned stochastic process, a value for the probability

of the neuron to trigger at a discrete time step is required. For tuning this value, it has been used the better known value of the natural frequency of the neuron, altogether with the properties of the geometric distributions for setting


The geometric distribution can be used for representing the number of Bernouilli trials needed for obtaining one success for a constant probability. It is formally represented by eq. 3.14, whose result is the probability of obtaining one success after trials, given the success probability of one trial.


Given the natural frequency of the neuron is known, it can be obtained the typical number of iterations until a spike occurs by using eq. 3.15,


where is the time step of the simulation environment.

Assuming the membrane potential is constant and equal to , the Bernouilli probability

can be calculated by using the cumulative distribution function of the geometric distribution (eq.



Therefore, in order to centre the firing rate around the neuron’s natural rate when no external pulses are received (i.e. the membrane potential is always ), the previous equation is solved for , and , resulting .

Eq. 3.12 can be transformed into eq. 3.17 for getting the value of that satisfies the previous condition. This is, , , and . The obtained result is


The previous stochastic model for determining the firing time of the neuron proves to be realistic when the membrane potential is close to the rest potential, and also for values close to the threshold potential. However, intermediate values provoke the neuron to trigger after 1 to 3 iterations, which makes is a similar behaviour to the neuron at high voltage values.

This issue is addressed in the next section, and a modification for solving it is introduced.

3.4.3 Charge of a neuron by random noise

The introduced stochastic model consists in consecutive Bernouilli events for deciding whther a neuron spikes or not on each iteration. Its probability is determined by the membrane potential. Therefore, if a neuron is isolated from the surrounding ones, its membrane potential is constant and then the probability of spiking is constant too. Then, the process can be modeled by a geometric distribution, and hence the probability density function will be monotonic and decreasing. This behaviour can be observed in Fig.


Figure 3.14: Normalized histogram of the firing times of an isolated neuron, with a constant membrane potential . The mean of the triggering times is , and the median is 1002 [ms]. It was calculated over a set of 5000 samples. Note that only the range from 0 to 3000 [ms] is shown, although theoretically the values can span until infinity.

Although in these results the expected firing time of the neuron is close to its natural period, real neurons tend to show a behaviour resembling a Poisson distribution, where the firing is centered around the expected value, instead as being spread out like in the Fig.

3.14. As the process consists in a series of identical and independent Bernouilli trials, the distribution is drifted to the left side and hence trials will be more likely to occur close to the starting time.

In order to obtain a behaviour closer to the one of real neurons, white noise has been added to the membrane potential on each iteration, so the voltage tends to grow and the probability of spiking in a single trial increases with time. Therefore the distribution of the firing times will be squeezed to a more narrow area, following a bell-like shape.

Since the added noise follows a Gaussian distribution with mean , eq. 3.8 is transformed into eq. 3.18.


In Fig. 3.15 the evolution of the membrane potential is depicted after adding white Gaussian noise to the membrane potential on each iteration. In the second figure is depicted the histogram of the firing times. The neuron spikes at a fast pace, and this is due to the fact that the parameters of the Sigmoid function were kept unaltered from the previous model i.e. it was tuned in a way that the neuron would spike on average in the first 1000 ms under the assumption that it would stay with a membrane potential of , which is not the case anymore.

Figure 3.15: a) Evolution of the membrane potential during 2000 [ms]. The average firing time is . b) Normalized histogram of the firing times of an isolated neuron getting excited by Gaussian noise with . The Sigmoid function has and . It was calculated over a set of 10000 samples.

Since there is a voltage addition at each iteration due to the noise, the Leaky time difference parameter in the previous equation is the same at each iteration, and equal to the time step. In absence of other disturbances, the membrane potential grows monotonically towards an asymptote that is calculated as follows:


Therefore, with eq. 3.19 can be calculated the steady-state voltage i.e. the voltage at which the neuron will stabilize in absence of inputs.

In order to slow down the firing rate of the neuron, the Sigmoid parameters were modified, following a graphical rule of thumb: As a reference voltage, it is taken the voltage of the membrane after 75% of the natural period (62 mV). That voltage is used for solving the Sigmoid function (eq. 3.12), setting to the probability of the Bernouilli process for having 50% chances of triggering in the following 1000 trials. Moreover, was set to .

In Fig. 3.16 is depicted the ideal and simulated behaviour of a neuron with . Then, . From eq. 3.19, it is obtained that the steady-state membrane potential is . It is also depicted the histogram of the firing times of a neuron when Gaussian noise is added to it on each iteration. It can be observed that it resembles the shape of a Poisson distribution, and its sample mean is near the desired natural period.

Figure 3.16: a) Behaviour of the membrane potential for , when Gaussian noise with is added on each iteration. In blue, ideal evolution of the membrane potential. In green, behaviour of a neuron during a simulation. b) Normalized histogram of the firing times of an isolated neuron. The average firing time is . It was calculated over a set of 10000 samples.

4.1 Selection of a multi-agents system simulation tool: The RANA environment

In order to implement the designed system it was used a simulation tool for multi-agent systems, namely the RANA software framework (Jørgensen et al., 2015). This tool is a software project aimed at executing multi-agent simulations able to replicate behaviours in real-time.

In terms of software structure, RANA is divided in two main blocks:

  • The simulation core, written in C++, contains the code for running the graphical interface, iterating and rendering the simulation following real-time constraints, parsing the agents’ behaviours, and processing the communication between them.

  • The agent scope, where the behaviours of the different agents are defined. It is written in Lua, and requires a master agent as the entry point to the agent system description. In the code defining this master agents it can be specified the spawning of new agents with different behaviours.

The agents in the simulation environment can communicate with each other by using asynchronous events, which can be addressed to specific agents or broadcasted to all of them, and they can travel at a specific speed through the media. Moreover, they can include data that can be shared between the agents.

This project was preceded by a self-study project titled “Learning RANA”, covering the basic principles of the tool. That self-study also provided some improvements to the RANA project, mainly a new method for generating the movement of the agents based on their velocity. The findings and experiments done within that project are open and free to access at a GitHub repository222https://github.com/jlrandulfe/learning_RANA/wiki.

4.1.1 Implementation of the MAS design in RANA. The master agent

Beside the division of the tasks into the agents explained in section 3.1, the implementation in the RANA framework required the creation of a master agent. It is the entry point of the simulation, so RANA will spawn this agent and execute its routines as soon as it is created. This is due to the fact that RANA only allows to specify one type of agent for creating a simulation. Therefore, its main purpose is to create the agents structure and their distribution in a 2-D space, as well as feeding them with some crucial data. This main goal can be subdivided in the following tasks:

  • Spawning the neuron somas, according to a specified initial layout.

  • Emitting electric pulses to the different neurons, in order to force their spiking during the learning period. The time sequence of the firings follow a pattern based on the type of function that has to be implemented.

  • Commanding the neurons to not grow, for those whose growth is not desired.

  • Gathering data that will be used afterwards for analyzing the performance of the system.

4.1.2 Enhancements in the RANA framework

In order to implement the described MAS, two new functionalities were added to the RANA core. Namely, two functions were implemented in the core: One function for getting numbers following a Gaussian distribution, and one function for getting numbers following a Poisson distribution. The next steps were followed:

  1. Calculate values following Gaussian and Poisson distributions in two different class functions of the Phys class in the RANA core. They are calculated by using the corresponding functions of the C++ standard library 333https://en.cppreference.com/w/cpp/numeric/random/normal_distribution444https://en.cppreference.com/w/cpp/numeric/random/poisson_distribution.

  2. In the AgentLuaInterface class, two functions were added for interfacing each of the functions created in the previous step with the Lua environment.

  3. In the Lua modules statistic library, a call to each of the functions was added, so they are accessible from the Lua environment.

Auxiliary tools used for data analysis and experiment set-up

During this project, the next two RANA functionalities were used in order to ease the debugging and data analysis:

Log results into a .csv file:

In order to centralize the data collection, the master agent handles now certain incoming events in order to gather the desired data. Also, on the clean up function of the master agent all the collected data is written down to a local .csv file. This is performed by using the Lua IO library 555http://lua-users.org/wiki/IoLibraryTutorial.

Further work focused on the data analysis entailed developing scripts for parsing the obtained data and extracting meaningful results. The data analysis was performed in an independent environment using the Python programming language.

Set automatic experiments:

This second functionality creates automatic experiments by specifying in a Lua script the desired configuration of the experiment. So far, this configuration file commands to do experiments where the master agent is spawned and run for a specific amount of seconds. Moreover, it is also specified how many runs of this configuration will be executed.

In order to use this feature, the experiments are run from the command line, and the GUI is not executed. This is done with the following command:

    /<RANA-path>/Rana_qt --nogui -f /<project-path>/experiment.lua

This is a powerful tool that, altogether with the previous feature, allows to automatically execute and record the results of a given number of experiments. Furthermore, different parameters of the agents can be configured, so the influence of them in the result of the simulation can be easily quantitatively analyzed. The RANA wiki offers a good introduction and baseline of the topic 666https://github.com/sojoe02/RANA/wiki/Setup-of-mass-experiments

Network name Tested features Network layout
Neuron pair Neuron’s individual firing histogram. Hebbian learning based dendrite growth.
Delay line Hebbian learning based dendrite growth. Noise tolerance. Spine growth dynamics
Coincidence detector Dendrites growing towards multiple goals. Spine growth dynamics
Extended coincidence
2-D firing histogram. Timing behaviour. Leaky I-F model
Neuron reservoir Capacity to generalize the application. Ability to generate multiple functions.
Table 4.1: Summary of the experiments that have been performed for testing the designed model, as well as the different features that have been tested on each of the experiments.

4.2 Neuron pair

This is the simplest network that has been implemented. It consists in 2 neurons , which receive input pulses following different functions dependent on time i.e. .

First of all, this layout was very useful for testing Hebb’s rule i.e. the time difference between the spikes in both neurons was varied, so it could be easily analyzed how sensitive are the implemented equations to that parameter. Moreover, the effect of other parameters could be tested, such as the noise and the distance between neurons.

Secondly, this was the chosen layout for studying the behaviour of single neurons. For example, by blocking the spiking of one of the neurons, the probabilistic histogram of a neuron when it is isolated was tested.

4.3 Delay line

This layout was used for the initial stages of the design, aimed at testing the Hebbian rule as a simplified way of modeling the mechanism of the neuron growth and synapse. It was chosen because it is a simple network which could be used for implementing certain features in the first versions of the model.

The equation representing the Hebbian learning (eq. 2.1) was implemented in the neurons’ growth cones, and it defined the growth behaviour. Basically, once a neuron triggers, it looks for incoming pulses and gets the corresponding excitation level from each of them. Moreover, Hebb’s rule was applied into an euclidean space by getting a vector representing the growth cone velocity in a 2-D space.

Finally, a stochastic model was added to the environment and to the firing times of the neurons. The final results show the influence that the noise has in the final structure.

4.3.1 Implemented network

This experiment was intended to be a baseline for further research on the application of the Hebbian rule for making SNNs to dynamically grow. Therefore, the choice prioritized a network that is both well-known and simple. Hence, the delay lines used by the mammalian brain for locating the source of sounds were replicated. Initially proposed by L. Jefress (Jeffress, 1948), this theory is widely accepted, and in-vivo experiments have given empirical evidence of its existence (Konishi, 1993).

In Fig. 4.1 is depicted the schematic of the biological neural network that has been partially reproduced. When a sound reaches both ears of the animal, its information travels through 2 delay lines in opposite directions composed by several shared coincidence detectors, which are sub-networks of 3 neurons. When the sound is registered by both inputs of a coincidence detector, its output spikes and further layers of the brain can determine the spatial location of the perceived sound.

Figure 4.1: Schematic of the neural model for the location of sound sources in mammals and other animals. Obtained from (Konishi, 1993).

One important difference between the implemented network and the ones described in literature (Konishi, 1993) is that the line delay of the biological network is believed to be achieved by a unique wire with a specific time delay for the signals that it transmits. On the contrary, the designed network achieves the delay line by connecting neurons sequentially and making use of their trigger delay.

Despite not being a faithful reproduction of the biological model, it mocks the same behaviour. Moreover, the main purpose of the experiment was to test the Hebbian rule, where the evolution of the network is known and simple.

Environment set-up

In order to simplify the aforementioned network, the initial experiment was focused on reproducing the delay line connected to only one of the ears. This delay line is composed by several neurons that trigger with a certain time delay from each other. The network starts with the neurons unconnected, and if the implementation is successful, the neurons would end up connected to each other in a sequential shape, so when the first neuron registers an incoming pulse, the following ones would trigger after their corresponding delay. In Fig. 4.2 is depicted the distribution of the neurons in RANA for this initial layout, before and after the growth process takes place.

Figure 4.2: Initial layout of the first experiment, before and after the growth process. The red dots are the somas of the neurons, the blue dots are the axon/dendrites links, and the green dots represent their growth cones.

The shown layout has the main inconvenient that a neuron can only reach any of its two adjacent neurons, as it would have to jump over any of them in order to reach any other neuron. This would introduce a limitation to the neurons, as they could only grow towards their immediate neighbours.

Therefore, the layout was modified into a circular pattern in order to minimize this drawback. In Fig. 4.3 is depicted this layout both before and after the growth process takes place.

Figure 4.3: Initial layout of the first experiment after modifying the layout from a linear to a circular distribution of the agents, before and after the growth process. The red dots are the somas of the neurons, the blue dots are the axon/dendrites links, and the green dots represent their growth cones.

Once the geometrical layout of the network was set, the following assumptions and specifications were applied to the designed model:

  • the excitation of a neuron only needs an incoming pulse for happening. Therefore, the triggering is deterministic, and the stochastic model used for the spiking is not tested in this experiment.

  • The electric pulses of all neurons have the same intensity. Hence, the attraction force solely depends on the distance between agents, and on the spiking time difference.

  • The learning process is much slower than the neuron dynamics i.e. , where .

  • Once a neuron’s input is connected to a second neuron, the learning process for that neuron stops.

  • Neuron’s soma is the agent getting excited through the RANA events called synapse. Once excited, it propagates the excitation state to the rest of the neuron’s agents.

  • Neuron’s growth cone is the agent managing the incoming electric pulses, through the RANA events called electric_pulse.Hebb’s rule (eq. 3.2) is applied by merging this information with the timing of the soma’s excitation, .

  • Signal propagation speed is neglected.

  • Once a neuron is excited, it emits an electric pulse after a certain time delay.

Environment and process noise

In order to make the model more realistic, noise was added to the system. Thus, the robustness of the network could be tested. Moreover, learning is highly influenced by stochastic processes taking place in the neurons and neural circuits according to literature (Dayan and Abbott, 2001).

Namely, 2 different sources of noise have been added to the model:

  • Environment noise: Corresponding to disturbances caused by a plethora of unknown sources, or sources that can not be controlled or monitored, such as other neurons in the brain that are unconnected to the current system. By applying the central limit theorem, the intensity of this noise can be modeled as a normal distribution. Moreover, the direction of the noise follows a uniform distribution, from 0 to 360 degrees.

  • Process noise: Corresponding to the uncertainty of the triggering time of the neuron. Although it is likely that a neuron triggers right after being excited, the precise time of its triggering follows as well a stochastic model. It has been applied the Poisson distribution for modeling this noise.

Moreover, the experiment also tested the robustness of the system towards an intense source of noise at a specific location i.e. there may be a group of very active neurons close to the current neurons that may induce attraction to them. However, the current neurons should not grow towards the others in case that they belong to an independent process. This is assessed true if the triggering frequencies of both groups of neurons are different.

This second group of neurons has been modeled as a single electric pulse_generator agent, which sends an electric signal with a frequency uncoupled of the current group of neurons. Should this external group of neurons have a very intense field, the number of coincidences in time inside the Hebbian equation curve would be small enough to avoid the attraction between them.

4.3.2 Simulation results

Environment Gaussian noise:

In Fig. 4.4 are depicted the results of two simulations when the growth cones are affected by additive Gaussian noise, corresponding to the environment noise. In both cases the noise follows a normal distribution . Whereas the first image shows the pure effect of the noise, the second one shows the effect of the noise when the neurons are also following the Hebbian rule (eq. 3.2).

It is remarkable that even though the noise level is very high () compared to the value obtained by eq. 3.2 (), most of the growth cones are still able to make their way to the required destination somas. For the shown test, 5 out of 7 neurons managed to reach their destination somas. The results for lower levels of noise are better, and most of the times the neurons were able to connect successfully to the desired neighbours (For ).

The simulations showed a very low degree of repeatability in , in terms of the paths the dendrites would take. This was due to the high noise introduced to the growth. In the first case, the result always resembled the output in Fig. 4.3(a), with the dendrites showing a very chaotic trajectory around their somas. In the second case, most of the neurons managed to connect to the goal soma in all the experiments, although the degree of success varied for each repetition of the simulation.

Figure 4.4: Effect of additive white noise on the neuron growth, with , and . a) depicts the behaviour of the growth cone when it is only driven by the Gaussian pattern, and b) shows the addition of the noise over the Hebbian rule behaviour.
Process noise:

It corresponds to the delay in the triggering of the neurons after the input synapse gets excited. In order to model it, three different alternatives have been tried:

  • Add Poisson patterns with the same mean to all of the neurons. Even though the obtained noise values will be randomly determined for each neuron, and hence different, they will all have the same average values.

  • The mean of the Poisson pattern of each neuron is randomly decided at the beginning of the execution i.e. When a neuron is created, the Poisson distribution determining its process noise is characterized by , where is a uniform distribution between time values and .

  • The value of the noise level of each neuron is manually determined, based on its ID number. Thus, the effect of the process noise in a specific neuron can be more easily assessed.

It has been observed that if all of the neurons are disturbed with the same Poisson distribution, the effect will be almost neglected. This is due to the fact that for a given Poisson distribution , the relative delay in the triggering between neurons will be kept static, as all of them will tend to suffer the same amount of delay. Thus, the aforementioned second and third alternatives have been tried out on the neuron model.

The behaviour of the system when the 3 different approaches are applied is depicted in Fig. 4.5. It can be observed that the effect is almost none when every neuron is affected by the same distribution. When neurons triggering time follow Poisson distributions with randomly selected expected times, the results are unpredictable, as depending on the randomly obtained parameters, they will grow towards different directions. Finally, as it was expected, adding a trigger delay to only one neuron affects the direction of that neuron and its closest neighbours.

Figure 4.5: Effect of Poisson noise on the soma triggering time, when a) the triggering time of each soma is affected by a Poisson distribution with , b) the triggering time of each soma is affected by a Poisson distribution with , where is a uniform distribution between and , and c) when only the triggering time of one soma is affected by a Poisson distribution with .
Disturbances at specific locations:

As mentioned before, the neurons may be affected by two different noise sources from the environment. The first one, already introduced, relates to white noise due to all of the processes that are taking place in the brain. The second one, however, corresponds to intense neural activity taking place at a close location from the current neurons, so it has a specific direction and intensity. Moreover, it spikes at a characteristic frequency, which can not be the same as the current network - Otherwise they should end up being connected, according to Hebb’s rule.

Therefore, it was implemented a pulse_generator agent in the simulation, for replicating a neighbouring cluster of neurons uncoupled with the current set of neurons, but generating a very intense field.

The simulation is shown in Fig. 4.6. It can be observed that the growth cones tend to follow the desired path. However, they experience sudden changes corresponding to the spikes of the pulse_generator that happen when the neurons are excited.

Figure 4.6: Effect on the growth of an agent (near the bottom right corner) emitting pulses with a period of 1020 ms, whereas the rest of the neurons have a period of 1000 ms. The intensity of the emitted pulses is 10 (the intensity of the other neurons oscillate between 0.1 and 1). The environment noise was set to , and the triggering delay was not included. Neither the 2nd derivative model and the drag force were included in this simulation.

The effect of this noise is highly dependent on how the Hebbian rule is applied. Namely, it is very relevant how the influence of a neuron on a second one decays during time when they do not fire at similar times. The implemented model considers that the incoming intensity is reduced by a factor of when the incoming pulse is beyond the time limits of eq. 3.2. Changing the decay factor when a neuron triggers out of its range have a big impact on the reaction of a neuron towards incoming uncoupled pulses.

4.4 Coincidence detector

The delay detector described in the previous section served as an initial experiment for testing the ground knowledge of the Hebbian learning. In this section, a different neural network topology is introduced in order to delve deeper into the growth of the networks, and observe how they get affected by other factors.

Applying again the initial versions of the neuron model, a new topology has been implemented. The implemented function is a coincidence detector, which means that the output signal spikes when there is a spike detected in both of the input signals. In Fig. 4.7 is depicted the time evolution of an ideal network implementing a coincidence detector.

Furthermore, more features of the designed model were tested with this layout. As before, all the simulations and experiments have been performed using the RANA framework, where the neurons are represented by agents modeled in Lua.

4.4.1 Implemented network

The goal of the developed neuron model is to be able to reproduce a coincidence detector with two inputs. Therefore, a successful implementation of this network will generate spikes in the output signal whenever two spikes are detected at the same time in both inputs. In Fig. 4.7 is depicted the value of the output signal in relation with the input signals.

Figure 4.7: Schematic of the signals of an ideal coincidence detector. The two top graphs represent the pulses trains received in the inputs of the neurons, whereas the bottom one represents the pulses emitted by the neuron at its output. A neuron acting as a coincidence detector should only emit spikes when spikes arrive simultaneously to both of its inputs.
Assessment of the temporal values of input neuron spikes

In the model introduced in the previous section, a neuron gets excited every time an electric_pulse event reaches its soma. Furthermore, the excitation of the neuron is handled by the master agent, which sends electric pulses to the different neurons according to a given pattern.

Therefore, that experiment does not take into account the spiking of interconnected neurons. In fact, their outputs solely depend on the pattern of electric pulses sent by the master agent. Despite that experiment proved useful for testing a rule based on the Hebbian learning for making neurons’ spines grow towards other neurons, the experiment does not assess the firing rule of the neurons, which means that the effect of connected neurons is not analysed. The main purpose of this experiment was to implement and test the growth rule of the dendrite spines of the neurons.

Initial conditions

It is important to take into account that Hebbian learning and rapid structure growth is only one of the multiple mechanisms behind learning in neural circuits. Actually, the rapid growth of spines for connecting different neurons is a fine tuning of the network, and the main structure of the network is built up by longer and more complex processes taking place during the brain development (Feldman, 2009). Therefore, it is assumed that the initial setup of a network is formed by neurons getting excited by previously connected sources, which include many other neuronal layers in the brain, as well as previously established connections with external sensory signals. Hence, the modeled network receives inputs from the environment, which is assumed to be a black box, so the implemented model reacts to incoming signals which follow unknown rules.

This experiment has been set up with 3 neurons, where the 2 leftmost ones represent the input neurons, and the rightmost one represents the output neuron.

4.4.2 Simulation results

A network with the shape of a coincidence detector was obtained after simulating in RANA the described setup. The result for the basic experiment is shown in Fig. 4.8. The neurons start completely unwired, and after running the simulation for a certain amount of time, the two input neurons end up connected to the output neuron. Moreover, a stranger pulse was added in the bottom right corner of the map, which slightly leans the dendrites but is not able to alter the final result.

Figure 4.8: Evolution of a set of 3 agents that trigger following the behaviour of a coincidence detector i.e. the leftmost agents fire at the same time, whereas the rightmost agent fires 20ms later. The environment noise was set to , an stranger pulse generator with a rate of 1020ms and intensity of 10 was added at the bottom right corner, and the triggering delay was not included. The result is pseudo-deterministic (although the trajectories are subject to random deviations, the result is always the same).
Final consideration about the geometrical layout

If a set of three neurons must behave in certain way, their behaviour should be robust enough to keep working even if they are part of a larger network where other neurons are triggering with undetermined frequency and intensity.

Furthermore, the specific neurons that would need to be wired together in order to implement the desired function should not need to be initially determined. The way the network is connected is irrelevant, as long as it accomplishes its function. Therefore, for a given network there may be an undetermined number of layouts that can be considered successful.

4.5 Extended coincidence detector: Temporal pattern recognition

In the previous sections two different logic functions taking into account the time dimension have been introduced, namely a delay line and a coincidence detector. The former was used for testing the feasibility of using Hebb’s rule for generating the growth of neurons’ spines, whereas the latter explored the performance of the spine dynamics designed in this project and introduced in section 3.3. In order to test the designed rules on a more complicated scenario, a network that implements the features of the two aforementioned networks has been developed i.e. a coincidence detector where at least one of the inputs go through a delay line before reaching the output neuron.

The main purpose of using this layout was being able to test the firing rule for the neurons (let us recall that in the previous experiments the firing of the neurons was forced by the master agent). Moreover, this layout could be used for testing the whole integration of the designed rules, as it requires all of them to work in order to obtain the desired outcome. Therefore, it could be used as a golden standard for testing further modifications of the system and being able to have comparable results.

To sum up, with this layout it is possible to have a network that starts with an unconnected set of neurons that get wired by using Hebb’s rule and fire following the proposed modification of the Leaky I-F model (see section 3.4). Furthermore, the growth of the spines through the network environment is achieved by implementing the designed dynamics for the spines.

With all this in mind, an extended version of the coincidence detector introduced in section 4.4 has been designed.

4.5.1 Implemented network

In Fig. 4.9 is shown the schematic of the desired network after the learning process. The two input signals are fed to and respectively, and once they get excited they propagate electric pulses to the output neuron. It can be observed that there is an additional neuron between and , being its purpose the addition of a delay in the spiking of with regards to the spiking time of .

Figure 4.9: 2-neurons coincidence detector with one of them delayed by using an intermediate neuron. The tipycal delay of each neuron is 20ms, and the spiking period is 200ms
Learning process

When the execution of the systems starts, the neurons are not connected to each other. During a learning period of 500 s, the inputs of the neurons are forced by triggering 20 ms on average later than and , which fire simultaneously. Furthermore, the excitation of is forced 40 ms later than and 20 ms later than the other two neurons.

4.5.2 Simulation results

The described network has been implemented and simulated in RANA. The simulation result is depicted in Fig. 4.10. The initial layout consists in two input neurons in the left side, and output neuron in the right side, and an intermediate neuron for provoking a time delay in the bottom branch.

The neurons start unconnected, and during the simulation the first input neuron gets connected to the output neuron, whereas the second one gets connected to the intermediate neuron, which afterwards gets connected to the output neuron too.

Despite the randomness in the process, the default values for the set of parameters makes the network to always result in the same topology. The required learning window oscillates between 380 ms and 450 ms.

Figure 4.10: Simulation result in RANA of the 4-neurons extended coincidence detector. The output soma is the rightmost agent (red dot), and the input neurons are two leftmost somas. There is an extra neuron in the bottom branch for delaying the signal of one of the inputs. The delay in the firing of neurons was set to 20 ms.

4.6 Neuron reservoir

This network consists in a neuron reservoir with random initial connections, where the input signals are applied to some of the neurons, and another one is connected to the desired output signal during the learning period (see Fig. 4.11).

Figure 4.11: Sketch of a possible initial state of a neuron reservoir with the inputs and output of a coincidence detector.

As neurons are highly dominated by stochastic processes, the given reservoir may end up connected in a plethora of different ways. If the model is robust enough, the great majority of the times the final layout should reproduce the desired function in a satisfactory way. Moreover, the number of neurons in the reservoir may probably affect the success rate of the final layout.

If this layout were to be successful, further versions could deal with more than one function at the same time i.e. If the network implements 2 coincidence detectors at the same time, it would have their 4 inputs and 2 outputs in the same reservoir, or even sharing the inputs.

  • The output neuron is initially isolated i.e. it starts without any input connections.

  • The rest of the neurons are initially randomly connected (or disconnected) to each other.

  • System inputs are initially connected to random neurons (at least 1 neuron per input).

  • There is present in the system environment noise (Gaussian + specific disturbances), and process noise (Poisson based delay).

Model requirements:
  • The LTP is modeled by the positive side of the Hebbian rule (eq. 3.2)

  • The LTD is modeled by the negative side of the Hebbian rule (eq. 3.2)

  • There are needed 3 noise models for the environment white noise, the neuron trigger delay, and the presence of local disturbances.

  • The synapse between neurons decay over time

  • The output of the neurons is correlated with the connected inputs ( )

Controller parameters:

In order to do a quantitative assessment of the model, the influence of the following parameters in the final outcome was evaluated:

  • Number of neurons in the reservoir

  • Number of neurons connected to the inputs

  • Number of initial random connections

  • Noise intensity

4.6.1 Simulation results

The neuron reservoir network has been simulated in RANA, for a network with 21 neurons. In Fig. 4.12 is depicted the network layout at the beginning and end of the simulation. The network was fed with a function equivalent to a coincidence detector of three inputs i.e. there is an output neuron whose spiking is forced a given time after the simultaneous spikes of three input neurons are produced.

It is observed that the neurons in the system are able to grow their spines and connect to each other. What is more important, both the input and output neurons are connected to the main body of the network, and therefore the input spikes can be propagated to the rest of the network.

Doing a quantitative and meaningful analysis of this network is not trivial, as there are many signals and intermediate neurons involved. There is also a high degree of randomness and was difficult to predict the topology that was created.

Figure 4.12: Initial and final states of a neuron reservoir network simulated in RANA. The network contains 21 neurons, and the learning took place for 400 s. The three leftmost neurons are selected as input neurons. Their spiking are forced, and the growth of their dendrites is inhibited. The rightmost neuron, which is considered the output neuron of the system, and spikes with a certain delay relative to the input neurons.

4.7 Final results and discussion

In the previous sections of this chapter, a couple of layouts implementing the designed model have been introduced. Each of them has been tested with the RANA environment, and it has been proven in most cases that the result is satisfactory. However, when running the simulations the design parameters have been set to values that lead to the desired results, without doing a thorough assessment on their impact.

In this section are shown the results of some experiments that have been done in order to analyse the effect of some of the parameters of the designed model on its performance.

When not stated the contrary, the set of parameters used in the following experiments are fixed to the values shown in Table 4.2. Moreover, the initial neurons layout is crucial for determining the final outcome of the network.

Parameter name Value
-70 [mV]
-54 [mV]
Pulse amplitude 10 [mV]
Leaky 520 [ms]
Natural period 200 [ms]
Neuron delay 20 [ms]
0.02 [mV]
1 [N]
Table 4.2: Default set of parameters used in the RANA simulations.

The following results offer a quantitative measure of the performance of the model in 2 different aspects: The spiking probability of a neuron when it is isolated from other influences, and the correlation in the spiking of 3 neurons, when they fire at different time intervals. For the former, the neuron pair layout (see section 4.2) has been used, whereas the extended coincidence detector (see section 4.5) was the chosen layout for the latter.

4.7.1 Evaluation of the performance of a single neuron isolated from neighbours

In order to assess a specific set of parameters, it is relevant to take a look to their performance on a single neuron, when it does not have any sort of interaction with other neurons. As the introduced design consists in neurons following a stochastic model, they trigger even without the presence of external pulses.

If this first assessment results unfruitful, it will be unlikely that those parameters will be able to be used for producing useful networks. In order to perform the analysis, it was used the neuron pair network introduced in section 4.2, as it only contains two neurons and it is straightforward to remove the interaction between the two of them.

Characterization of the 1-D firing histogram of a neuron

Fig. 4.13 shows the histograms of a neuron which follows all the rules introduced in the chapter 3 for two different values for the characteristic time constant of the Leaky I-F model. Namely, this results were used in section 3.4 for validating the designed probabilistic model for making neurons to spike.

It is observed that the bigger the time constant, the more the histogram gets shifted to the right, keeping a shape similar to a Gaussian bell. This makes sense, as that parameter is the one defining how slow the charge of the membrane capacitance will go back to its original value ().

Figure 4.13: Normalized histogram of the firing times of an isolated neuron getting excited by white noise with when a) (the average firing time is ), and b) when The parameters of the Sigmoid function are and (the average firing time is ). The histograms were calculated over a set of 10000 samples.

4.7.2 Correlation between firing times in the extended coincidence detector

By using the extended coincidence detector introduced in section 4.5 it was evaluated how different parameters affect to the correlation between the time differences of the firing of the output neuron relative to the two input neurons.

Fig. 4.14 shows the correlation of the firing times for a successful set of parameters. It can be observed that most of the firings concentrate at , and .

Figure 4.14: Histogram of the firing times of the output neuron in the extended coincidence detector in relation with the input neurons. The electric pulses have an amplitude of 3 mV. Times in the X and Y axis represent the time difference between the triggering of the output neuron () and the input neurons ( and , respectively). The white area implies no triggers at all happening with those time difference values. The relative firing times of the input and output neurons are recorded for 10000 spikes of the output neuron.
Effect of the voltage intensity of the incoming pulses

In the next figure can be observed how the amplitude of the incoming pulses affect the timing performance of the network. In Fig. 4.15 is depicted the simulation when the electric pulses have an amplitude of 1 mV (whereas in the simulation shown in Fig. 4.14, they had an amplitude of 6 mV).

It is observed that the dispersion of the histogram is bigger, and it tends to get repeated at time intervals which correspond to the period of the neurons (200 ms).

Figure 4.15: Histogram of the firing times of the output neuron in the extended coincidence detector after the learning process, when the pulse intensity is 1 mV. The network learns during 500 [ms]. After that, the relative firing times of the input and output neurons are recorded for 10000 spikes of the output neuron.

It is also interesting to observe the results when one of the connections in the network was not successful i.e. when neurons and did not get connected during the simulation. This phenomena is depicted in Fig. 4.16, for an electric impulse amplitude of 1 mV. The dispersion is greater than in the same experiment when all the connections are correctly achieved.

Figure 4.16: Histogram of the firing times of the output neuron in the extended coincidence detector, when the pulse intensity is 1 mV, and the neurons and did not connect during the learning process. The relative firing times of the input and output neurons are recorded for 10000 spikes of the output neuron.

In general, it can be assumed that the correlation in the firing times tends to grow to more precise results when the impulse amplitude is bigger. For a specific circuit and function, this will result in more reliable results, giving always the correct output when the corresponding inputs spike.

However, this is actually an undesired behaviour, as it would imply an ad-hoc network, which could only work for a very specific function. In the extreme situation when the amplitude gets extremely high values, the neuron will tend do always spike within the same time difference, which would be a behaviour closer to that of an ANN.

4.8 Discussion

In the previous sections of this chapter the different simulations that have been run for testing the proposed model are summarized. Moreover, analytic results of the simulations show the non-deterministic nature of the model.

The results also show that the multi-agent approach and the implemented growing rules are able to evolve the topology of a small set of neurons and reproduce the function that was initially fed to the system.

The 2D histogram of the extended coincidence detector shows the relationship between the relative times of the input and output neurons firing. The network starts fully unconnected, and after a learning period, the connections between the neurons evolve until reaching certain topology. The graph shows a distribution in the firing times of the neurons centered around the time delays considered during the learning process.

In any case, the proposed simulations are simple, and the resulting networks require a deeper analysis in order to have a more solid model.

Usage of a 2D map

Due to the nature of the RANA simulation environment, the tests have been conducted for geometrical spaces with two dimensions. This is different to biological neural circuits, where the neurons are distributed in a 3D domain.

Despite the simulations were implemented in 2D and the obtained results are satisfactory, the flexibility of the model is limited due to this condition. Future networks will probably require complex structures which will not be possible to reach without a third dimension being implemented. Otherwise, it is likely that the dendrite spines of different neurons collide and block the growth of each other.

Crossing of dendrites in the geometric space

The larger a network is, the bigger the amount of dendrites, which occupy physical space (as well as the neuron somas). In this project, it has been omitted the collisions that could happen between dendrites, and they were allowed to ignore and cross other agents in their way. In the situation of a dendrite colliding with a soma that was not its goal, the implemented rules force the establishment of a connection between the dendrite and the soma. This is particularly undesired in the case of three aligned neurons, where it is necessary to get a connection between the two in the borders. Although the noise may let the connection be reached, there are high chances of it happening, no matter if it is a 2D or 3D space.

Some modifications of the model are proposed for shortcoming this issue:

  • Adding a repulsion force in the model that keeps the growth cone away from other somas when their firing times are not correlated i.e. when Hebb’s rule do not apply to them. Moreover, this would be a short range function, so it would only have a noticeable effect when the growth cone is close to the soma. The ideal behaviour would be the growth cone describing a circular trajectory around the undesired soma, until it is able to reach again its trajectory to the target one.

  • The studies covering the function of chemical cues in biological neural circuits suggest that some proteins present in the neural circuit accomplish this function. Therefore, a solution could be adding more agents to the MAS design which perform the function of such proteins. This would be a long line of work, and it is not clear that its results would present a high performance.

  • The combination of growth noise and connection pruning can probably make the network converge to the desired structure i.e. undesired connections would appear in many occasions, but the presence of noise would allow the dendrite to avoid the undesired soma and continue its path towards the destination in some of the tries. Moreover, adding a pruning process in the network would make the undesired connections to disappear.

Computational cost of the MAS model and RANA

An issue observed during the execution of the neuron reservoir experiment was that the simulations were very computationally expensive. Moreover, the speed of the simulation decreased exponentially, so an improvement in the processing hardware would still lead to very long simulation times after certain time of execution.

A thorough computing analysis of the designed MAS has not been done. However, it may be crucial for exploring the alternatives of this simulation framework, specially for implementing networks with more complexity.


  • (1)
  • Bi and Poo (1998) Bi, G.-q. and Poo, M.-m. (1998), ‘Synaptic modifications in cultured hippocampal neurons: dependence on spike timing, synaptic strength, and postsynaptic cell type’, Journal of neuroscience 18(24), 10464–10472.
  • Dayan and Abbott (2001) Dayan, P. and Abbott, L. F. (2001), ‘Theoretical neuroscience: computational and mathematical modeling of neural systems’.
  • Dickson (2002) Dickson, B. J. (2002), ‘Molecular mechanisms of axon guidance’, Science 298(5600), 1959–1964.
  • Feldman (2009) Feldman, D. E. (2009), ‘Synaptic mechanisms for plasticity in neocortex’, Annual review of neuroscience 32, 33–55.
  • Feldman and Brecht (2005) Feldman, D. E. and Brecht, M. (2005), ‘Map plasticity in somatosensory cortex’, Science 310(5749), 810–815.
  • Ferber and Weiss (1999) Ferber, J. and Weiss, G. (1999), Multi-agent systems: an introduction to distributed artificial intelligence, Vol. 1, Addison-Wesley Reading.
  • Froemke and Dan (2002) Froemke, R. C. and Dan, Y. (2002), ‘Spike-timing-dependent synaptic modification induced by natural spike trains’, Nature 416(6879), 433.
  • Furber et al. (2014) Furber, S. B., Galluppi, F., Temple, S. and Plana, L. A. (2014), ‘The spinnaker project’, Proceedings of the IEEE 102(5), 652–665.
  • Gerstner et al. (2014) Gerstner, W., Kistler, W. M., Naud, R. and Paninski, L. (2014), Neuronal dynamics: From single neurons to networks and models of cognition, Cambridge University Press.
  • Ghosh-Dastidar and Adeli (2009) Ghosh-Dastidar, S. and Adeli, H. (2009), ‘A new supervised learning algorithm for multiple spiking neural networks with application in epilepsy and seizure detection’, Neural networks 22(10), 1419–1431.
  • Graves et al. (2008) Graves, A., Liwicki, M., Fernández, S., Bertolami, R., Bunke, H. and Schmidhuber, J. (2008), ‘A novel connectionist system for unconstrained handwriting recognition’, IEEE transactions on pattern analysis and machine intelligence 31(5), 855–868.
  • Graves et al. (2013) Graves, A., Mohamed, A.-r. and Hinton, G. (2013), Speech recognition with deep recurrent neural networks, in ‘2013 IEEE international conference on acoustics, speech and signal processing’, IEEE, pp. 6645–6649.
  • Holtmaat and Svoboda (2009) Holtmaat, A. and Svoboda, K. (2009), ‘Experience-dependent structural synaptic plasticity in the mammalian brain’, Nature Reviews Neuroscience 10(9), 647.
  • Indiveri and Horiuchi (2011) Indiveri, G. and Horiuchi, T. K. (2011), ‘Frontiers in neuromorphic engineering’, Frontiers in neuroscience 5, 118.
  • Izhikevich (2003) Izhikevich, E. M. (2003), ‘Simple model of spiking neurons’, IEEE Transactions on neural networks 14(6), 1569–1572.
  • Jeffress (1948) Jeffress, L. A. (1948), ‘A place theory of sound localization.’, Journal of comparative and physiological psychology 41(1), 35.
  • Jørgensen et al. (2015) Jørgensen, S. V., Demazeau, Y. and Hallam, J. (2015), Rana, a real-time multi-agent system simulator, in ‘2015 IEEE/WIC/ACM International Conference on Web Intelligence and Intelligent Agent Technology (WI-IAT)’, Vol. 2, IEEE, pp. 92–95.
  • Kempter et al. (1999) Kempter, R., Gerstner, W. and Van Hemmen, J. L. (1999), ‘Hebbian learning and spiking neurons’, Physical Review E 59(4), 4498.
  • Kohonen (1990) Kohonen, T. (1990), ‘The self-organizing map’, Proceedings of the IEEE 78(9), 1464–1480.
  • Konishi (1993) Konishi, M. (1993), ‘Listening with two ears’, Scientific American 268(4), 66–73.
  • Krizhevsky et al. (2012)

    Krizhevsky, A., Sutskever, I. and Hinton, G. E. (2012), Imagenet classification with deep convolutional neural networks,

    in ‘Advances in neural information processing systems’, pp. 1097–1105.
  • Liu and Wang (2001) Liu, Y.-H. and Wang, X.-J. (2001), ‘Spike-frequency adaptation of a generalized leaky integrate-and-fire model neuron’, Journal of computational neuroscience 10(1), 25–45.
  • Lowery and Van Vactor (2009) Lowery, L. A. and Van Vactor, D. (2009), ‘The trip of the tip: understanding the growth cone machinery’, Nature reviews Molecular cell biology 10(5), 332.
  • Merolla et al. (2014) Merolla, P. A., Arthur, J. V., Alvarez-Icaza, R., Cassidy, A. S., Sawada, J., Akopyan, F., Jackson, B. L., Imam, N., Guo, C., Nakamura, Y. et al. (2014), ‘A million spiking-neuron integrated circuit with a scalable communication network and interface’, Science 345(6197), 668–673.
  • Nelson (2004) Nelson, M. E. (2004), ‘Electrophysiological models’, Databasing the brain: from data to knowledge pp. 285–301.
  • Patel and Poo (1982) Patel, N. and Poo, M.-M. (1982), ‘Orientation of neurite growth by extracellular electric fields’, Journal of Neuroscience 2(4), 483–496.
  • Paugam-Moisy and Bohte (2012) Paugam-Moisy, H. and Bohte, S. (2012), ‘Computing with spiking neuron networks’, Handbook of natural computing pp. 335–376.
  • Soman et al. (2016) Soman, S., Suri, M. et al. (2016), ‘Recent trends in neuromorphic engineering’, Big Data Analytics 1(1), 15.
  • Thurau et al. (2003)

    Thurau, C., Bauckhage, C. and Sagerer, G. (2003), Combining self organizing maps and multilayer perceptrons to learn bot-behaviour for a commercial game.,


    ‘GAME-ON’, Citeseer, p. 119.


  • (1)
  • Bi and Poo (1998) Bi, G.-q. and Poo, M.-m. (1998), ‘Synaptic modifications in cultured hippocampal neurons: dependence on spike timing, synaptic strength, and postsynaptic cell type’, Journal of neuroscience 18(24), 10464–10472.
  • Dayan and Abbott (2001) Dayan, P. and Abbott, L. F. (2001), ‘Theoretical neuroscience: computational and mathematical modeling of neural systems’.
  • Dickson (2002) Dickson, B. J. (2002), ‘Molecular mechanisms of axon guidance’, Science 298(5600), 1959–1964.
  • Feldman (2009) Feldman, D. E. (2009), ‘Synaptic mechanisms for plasticity in neocortex’, Annual review of neuroscience 32, 33–55.
  • Feldman and Brecht (2005) Feldman, D. E. and Brecht, M. (2005), ‘Map plasticity in somatosensory cortex’, Science 310(5749), 810–815.
  • Ferber and Weiss (1999) Ferber, J. and Weiss, G. (1999), Multi-agent systems: an introduction to distributed artificial intelligence, Vol. 1, Addison-Wesley Reading.
  • Froemke and Dan (2002) Froemke, R. C. and Dan, Y. (2002), ‘Spike-timing-dependent synaptic modification induced by natural spike trains’, Nature 416(6879), 433.
  • Furber et al. (2014) Furber, S. B., Galluppi, F., Temple, S. and Plana, L. A. (2014), ‘The spinnaker project’, Proceedings of the IEEE 102(5), 652–665.
  • Gerstner et al. (2014) Gerstner, W., Kistler, W. M., Naud, R. and Paninski, L. (2014), Neuronal dynamics: From single neurons to networks and models of cognition, Cambridge University Press.
  • Ghosh-Dastidar and Adeli (2009) Ghosh-Dastidar, S. and Adeli, H. (2009), ‘A new supervised learning algorithm for multiple spiking neural networks with application in epilepsy and seizure detection’, Neural networks 22(10), 1419–1431.
  • Graves et al. (2008) Graves, A., Liwicki, M., Fernández, S., Bertolami, R., Bunke, H. and Schmidhuber, J. (2008), ‘A novel connectionist system for unconstrained handwriting recognition’, IEEE transactions on pattern analysis and machine intelligence 31(5), 855–868.
  • Graves et al. (2013) Graves, A., Mohamed, A.-r. and Hinton, G. (2013), Speech recognition with deep recurrent neural networks, in ‘2013 IEEE international conference on acoustics, speech and signal processing’, IEEE, pp. 6645–6649.
  • Holtmaat and Svoboda (2009) Holtmaat, A. and Svoboda, K. (2009), ‘Experience-dependent structural synaptic plasticity in the mammalian brain’, Nature Reviews Neuroscience 10(9), 647.
  • Indiveri and Horiuchi (2011) Indiveri, G. and Horiuchi, T. K. (2011), ‘Frontiers in neuromorphic engineering’, Frontiers in neuroscience 5, 118.
  • Izhikevich (2003) Izhikevich, E. M. (2003), ‘Simple model of spiking neurons’, IEEE Transactions on neural networks 14(6), 1569–1572.
  • Jeffress (1948) Jeffress, L. A. (1948), ‘A place theory of sound localization.’, Journal of comparative and physiological psychology 41(1), 35.
  • Jørgensen et al. (2015) Jørgensen, S. V., Demazeau, Y. and Hallam, J. (2015), Rana, a real-time multi-agent system simulator, in ‘2015 IEEE/WIC/ACM International Conference on Web Intelligence and Intelligent Agent Technology (WI-IAT)’, Vol. 2, IEEE, pp. 92–95.
  • Kempter et al. (1999) Kempter, R., Gerstner, W. and Van Hemmen, J. L. (1999), ‘Hebbian learning and spiking neurons’, Physical Review E 59(4), 4498.
  • Kohonen (1990) Kohonen, T. (1990), ‘The self-organizing map’, Proceedings of the IEEE 78(9), 1464–1480.
  • Konishi (1993) Konishi, M. (1993), ‘Listening with two ears’, Scientific American 268(4), 66–73.
  • Krizhevsky et al. (2012)

    Krizhevsky, A., Sutskever, I. and Hinton, G. E. (2012), Imagenet classification with deep convolutional neural networks,

    in ‘Advances in neural information processing systems’, pp. 1097–1105.
  • Liu and Wang (2001) Liu, Y.-H. and Wang, X.-J. (2001), ‘Spike-frequency adaptation of a generalized leaky integrate-and-fire model neuron’, Journal of computational neuroscience 10(1), 25–45.
  • Lowery and Van Vactor (2009) Lowery, L. A. and Van Vactor, D. (2009), ‘The trip of the tip: understanding the growth cone machinery’, Nature reviews Molecular cell biology 10(5), 332.
  • Merolla et al. (2014) Merolla, P. A., Arthur, J. V., Alvarez-Icaza, R., Cassidy, A. S., Sawada, J., Akopyan, F., Jackson, B. L., Imam, N., Guo, C., Nakamura, Y. et al. (2014), ‘A million spiking-neuron integrated circuit with a scalable communication network and interface’, Science 345(6197), 668–673.
  • Nelson (2004) Nelson, M. E. (2004), ‘Electrophysiological models’, Databasing the brain: from data to knowledge pp. 285–301.
  • Patel and Poo (1982) Patel, N. and Poo, M.-M. (1982), ‘Orientation of neurite growth by extracellular electric fields’, Journal of Neuroscience 2(4), 483–496.
  • Paugam-Moisy and Bohte (2012) Paugam-Moisy, H. and Bohte, S. (2012), ‘Computing with spiking neuron networks’, Handbook of natural computing pp. 335–376.
  • Soman et al. (2016) Soman, S., Suri, M. et al. (2016), ‘Recent trends in neuromorphic engineering’, Big Data Analytics 1(1), 15.
  • Thurau et al. (2003)

    Thurau, C., Bauckhage, C. and Sagerer, G. (2003), Combining self organizing maps and multilayer perceptrons to learn bot-behaviour for a commercial game.,

    in ‘GAME-ON’, Citeseer, p. 119.