Functional differentiations in evolutionary reservoir computing networks

We extended RCs Dominey, Arbib, and Joseph (1995); Maass, Natschläger, and Markram (2002); Jaeger and Haas (2004), adopting the optimization principle via evolutionary dynamics to the internal random recurrent networks, called ERCs. We employed a continuous-valued, discrete-timestep model of ERC. The network consists of $N$ neurons, whose dynamics is described by

where $x_{i}(t)$ is a state of the $i$-th neuron at time $t$, $w_{ij}$ and $w^{(\text{in})}_{k}$ are weights of recurrent connections and weights of connections to an input signal $I_{k}(t)$, respectively, and $w_{i0}$ is a constant term. In addition, $\alpha_{i}$ represents a decay constant, and $\xi_{i}(t)$ is a noise term. The initial state of $x_{i}$ was randomly taken from a uniform distribution over $[-0.5,0.5]$. The network consisting of these neurons was constructed with sparse connections, such that only 10% of synaptic connections on the network were nonzero.

Because of effective computations on differentiations, we used the following additional structure of internal networks. The internal networks were divided into two groups of neurons that construct an input and an output layer within internal networks, unlike a standard RC that is constructed by homogeneous random networks. Neurons in the input layer receive input signals but have no direct connections to output units, called readout units. In contrast, neurons in the output layer do not receive input signals directly but have direct connections to the readout units (Fig. 1).

The model was trained to perform the discrimination task of two different kinds of patterns: simultaneous discrimination of spatial and temporal patterns. Correspondingly, the model possesses two distinct output vectors, spatial $y^{\text{(sp)}}$ and temporal $y^{\text{(temp)}}$ vectors. Each output vector was determined by a linearly weighted sum of the internal states, as shown in the following equation.

where (*) indicates spatial or temporal patterns, $i=1,...,L$ are indices of readout units and $w^{\text{(*)}}_{ij}$ are weights for readout. In the present simulation, output weights are determined, using ridge regression. For input patterns denoted by an asterisk, the following loss function was minimized:

where the first term indicates the loss stemming from the differences between output and target patterns $p_{i}^{(*)}(t)$ during a training period of length $T_{\text{tr}}$ and the second term indicates the loss from the weight $L^{2}$ regularization.

For the presented tasks, reservoir computing networks with random internal networks, which belong to conventional RCs, did not show good performance, even after the output weights were learned using the ridge regression. This implies that the structure of internal networks must be changed for such tasks. This is a reason why we introduced a genetic algorithm for changing the internal networks and extended conventional RCs to evolutionary ones. For the first time, we prepared an initial ensemble of randomly connected networks, and second, we performed the following procedures (1) to (3) for each network.

1. (1)

   Collection of data for learning. Input a pattern sequence consisting of spatial and temporal patterns, and continue to update internal states according to Eq. (1), and discard the first transient data.

2. (2)

   Learning period. Perform supervised learning of the readout weights with ridge regression, using the internal states and teacher signals collected in (1).

3. (3)

   Test period. Update the states of both internal and readout units, and then calculate the root mean square of errors between teacher signals and outputs for both spatial and temporal patterns. Then, evaluate Eq.(4).

For all networks, the following procedure (4) was performed.

1. (4)

   Evolutionary period. Among all networks, select networks with smaller errors for which the next generation of the networks were produced using mutation and crossover.

In the present simulation, the targets for mutation and crossover are the weights of recurrent connections and decay constants, $\alpha_{i}$. Mutation of the weights of recurrent connections was performed, randomly reconnecting $4$% of the total connections and adding Gaussian noise, $N(0,\sigma_{w}^{2})$ in a randomly selected 40% of the weights of all internal connections. Decay constants $\alpha_{i}$ were also changed by randomly adding Gaussian noise, $N(0,\sigma_{\alpha}^{2})$, where the upper and lower cutoffs of $\alpha_{i}$ were commonly given as [$\alpha_{\text{low}},\alpha_{\text{high}}]$. As for crossover, we randomly selected two survived networks and constructed one new network, where a half connection of a new network stems from a randomly chosen half connection of one survived network and the other half connection from the other survived network. For each generation, we repeated (1) to (4), and the selected networks show high adaptability to the given task. As a constraint for evolutionary dynamics, the following loss function was adopted, which indicates the summation of mean squared errors of the activity of all readout units, based on the teacher signals for both spatial and temporal patterns. Here, an average was taken over both of all readout units and a test period of length $T_{\text{te}}$.

To clarify an indication of the functional differentiation of neurons, we adopted mutual information between neural activity in the internal network and teacher signals. Mutual information was calculated for spatial and temporal patterns, separately according to the following steps; Discretize neuronal states $x_{i}(t)$ by equally dividing the interval between its minimum and maximum values into $n_{d}=8$ states, which are represented by $X_{i}(t)\in\{1,...,n_{d}\}$. Then, calculate the joint-probability distributions $p(X_{i},l)=p(X_{i}(t),l(t-4))$ and $p(X_{i},m)=p(X_{i}(t),m(t-4))$, whose distributions provide the basis of the calculation of mutual information for spatial and temporal patterns defined in the following manner.

Computation results of evolutionary dynamics showed high accuracy for both spatial and temporal patterns, which were around 90$\%$ in the case that the number of neuron $N$ was $64$ (Fig. 2(b)). The process of evolution is shown in Fig. 2(a). Although the performance of the early networks was low, it rapidly improved during the first hundred generations, thereby converging to a certain equilibrium state.

In this experiment, the objective function depends on the discriminations independently performed for spatial and temporal patterns but independent of the performance of the combination task. Nevertheless, the network performance can be estimated through learning only the readout units for each generation, because a common feature of RCs is that multiple-readout units can learn in parallel. In Fig. 2(b), the performance of the network in the case that the combination task was learned by ridge regression for each network, and each generation was also indicated by the purple dotted line and dots. Because the network evolution was not optimized for the combination task, the network performance had an accuracy of around 70%.

We show the output time series for spatial and temporal patterns at initial and evolved networks in Fig. 3. The evolved network successfully outputs correct signals for the classification of spatial and temporal patterns, except for the first few steps after exchanging the inputs. Evolutions of spectral radii of the recurrent weight matrix and the distribution of decay constant $\alpha_{i}$ are discussed in Appendix B.

To investigate functional differentiation of each neuron through evolution, we examined mutual information between the state of neuron $x_{i}(t)$ and spatial patterns $l(t)$, and between $x_{i}(t)$ and temporal patterns $m(t)$. If two quantities are negatively correlated over neuron ensembles, this implies the existence of many differentiated neurons that primarily represent only one type of information. The correlation coefficient between the two mutual information values in the input layer was about $-0.1$, indicating little change with evolution, and variance between trials was large (Fig. 4(a)). In contrast, correlation in the output layer significantly decreased in early evolution stages, when classification performance was improving, finally reaching a trial average of about $-0.4$ (Fig. 4(a)). Focusing on the output layer that showed a negative correlation, we made histograms to approximate the joint probability distribution of the two mutual information values $I^{\text{(sp)}}$ and $I^{\text{(temp)}}$ in the output layer neurons (Fig. 4(b)). Compared with the initial networks, the distribution extends to the upper left and lower right after evolution. In particular, neurons with higher $I^{\text{(temp)}}$ increased in the evolved networks. This result implies the appearance of the functional differentiation of the neural specificity of sensory information. To test robustness of the emergence of functional differentiation, we conducted two additional experiments: one using spatiotemporal patterns in which every component has a different oscillatory phase and another using a different set of frequencies. Appendix B presents the results.

We also investigated how the network structure changed from a random network after evolutionary optimization was established. We found that the feedforward connections from the input layer to the output layer in the internal network were strengthened much more than the other connections (Fig. 5). Furthermore, the feedback connections from the output to the input layer were rather weakened, and were also weakened for internal feedback connections in both the input and output layers. We estimated the effective connectivity of the network type for the feedforward or feedback type of networks by such connection strength. The ratio of effective connectivity of the feedback to the feedforward networks varied depending on the system size $N$: about 2:3 at $N=32$, and 3:4– 6:7 at $N=64$. This tendency of the reduction of feedback connections relative to the feedforward connections does not contradict the structure of cortical local connections that was recently found by Seeman et al.Seeman et al. (2018), where the ratio of the number of connections was reported as 1:10 for mouse primary visual cortex and 1:5– 1:7 for human frontal and temporal cortices. Furthermore, the evolved structure shown in the present ERC is similar to the evolutionary change of the network structure of hippocampus. In the biological evolution from reptiles to mammals, the hippocampus evolved from rather randomly uniform networks consisting of a small-cell and a large-cell layer to differentiated feedforward networks including feedback connections Treves (2004). It is also known in mammals that some amount of feedback connections exist in CA3, whereas only a few feedback connections exist in CA1. Thus, the ratio of feedback connections to feedforward connections in the mammalian hippocampus looks consistent with our findings in ERCs. However, it is questionable if the constraint adopted here was used in biological evolution. In biological evolution, sexual and social constraints as well as natural selections would have been more significant than the information constraint adopted here.

The structural changes described may induce functional changes. It is known that uniformly random networks in reptiles give rise to a simple memory, which can be realized by a single association. In contrast, RNNs with inhibitory neurons in mammals can produce an episodic memory, which can be realized by a successive association of memories (see, for example, Tsuda, Koerner, and Shimizu (1987); Tsukada, Yamaguti, and Tsuda (2013)). The hippocampal CA3 plays a role in the formation of such sequences of memories, while the hippocampal CA1 receives such sequences, and encodes the sequential information in a form of fractal geometry, called Cantor coding (see, for example, Tsuda (2001); Fukushima et al. (2007); Tsuda and Kuroda (2001); Kuroda et al. (2009); Yamaguti et al. (2011)). Considering these findings, the evolution of the network structure to the feedforward network including feedback connections was a decisive milestone in the evolution of memory. Therefore, the present ERC possesses sufficient structures of neural networks to yield high performance for not only the discrimination of patterns but also memory formation.

For comparison with the presented task, called here the separation task, we investigated the development of networks in the case of the combination task, where the same network size and parameters were used. The combination task identifies specific combinations of the spatial and temporal patterns. Learning succeeded with a high accuracy (Fig. 6, 7), but a similar differentiation was not observed in the internal network structure, that is, $I^{(\text{sp})}$ and $I^{(\text{temp})}$ did not show negative correlations(Fig. 8(a)). The joint probability distribution of the evolved networks (Fig. 8(b)) shows a concentration along a diagonal line, which was not seen in the separation task. This implies that a different type of optimization occurred. Specifically, because many neurons possessing both spatial and temporal information were observed in the output layer, such neurons possibly do not separate sensory information, but rather unify sensations, thereby processing a combination of sensory inputs. The emergence of such neurons may be related to the concept of “mixed selectivity,” which has recently been intensively studied both experimentally and theoretically in the neocortex Rigotti et al. (2013).

We investigated two kinds of conditions on which differentiation depends: the size of the network and the strength of the weights. As for the network size, the computation results showed that a decreasing network size corresponds to an increasing number of neurons that specifically respond to either input signal of spatial or temporal patterns (Fig. 9). Classification accuracy increased with network size (Fig. 9(a)). However, the degree of functional differentiation in the output layer, as measured by the strength of negative correlation of mutual information, was higher at $N=32$ and $N=64$ than at $N=128$ (Fig. 9(b)). There was no negative correlation in any case in the combination task. As for network structure, we found that ratios of feedforward connections increased and those of feedback connections decreased as the network size decreased (Fig. 9(c)).

To investigate the effect of penalty for weight size, the following penalty term, $\lambda\sum_{ij}{w_{ij}^{2}/N^{2}}$, on the connections of the internal network was added to the loss function and evolved the network. Here, we show the computation results in the case of $N=128$. Because the network with $N=64$ sufficiently differentiates even without such an additional penalty term, the effect of the additional term was very weak. In contrast, in $N=128$, because the degree of differentiations was low in the case without the penalty term (Fig. 9(b), we investigated how differentiations changed by such an additional term (Fig. 10). As the strength of the penalty $\lambda$ increased, accuracy fell at $\lambda=50$ and rose again at $\lambda=100$ (Fig. 10(a)). A detailed investigation of the reasons for this decrease and increase is beyond the scope of this paper. Negative correlations between two mutual information values became stronger as $\lambda$ increased, implying an increase in differentiated neurons (Fig. 10(b)). The feedforward-to-feedback ratio did not strongly depend on the penalty strength (Fig. 10(c)). It looks the penalty strength did not affect the differentiation on feedforward/feedback structure of the network. In summary, there was an evident range of penalty strengths that enhances differentiation without adversely affecting classification accuracy.