In their seminal work on systemic risk in financial markets, Eisenberg and Noe [Larry Eisenberg and Thomas Noe, 2001] proposed and studied a model with n firms embedded into a network of debt relations. We analyze this model from a game-theoretic point of view. Every firm is a rational agent in a directed graph that has an incentive to allocate payments in order to clear as much of its debt as possible. Each edge is weighted and describes a liability between the firms. We consider several variants of the game that differ in the permissible payment strategies. We study the existence and computational complexity of pure Nash and strong equilibria, and we provide bounds on the (strong) prices of anarchy and stability for a natural notion of social welfare. Our results highlight the power of financial regulation - if payments of insolvent firms can be centrally assigned, a socially optimal strong equilibrium can be found in polynomial time. In contrast, worst-case strong equilibria can be a factor of Ω(n) away from optimal, and, in general, computing a best response is an NP-hard problem. For less permissible sets of strategies, we show that pure equilibria might not exist, and deciding their existence as well as computing them if they exist constitute NP-hard problems.

We consider probabilistic PCA and related factor models from a Bayesian perspective. These models are in general not identifiable as the likelihood has a rotational symmetry. This gives rise to complicated posterior distributions with continuous subspaces of equal density and thus hinders efficiency of inference as well as interpretation of obtained parameters. In particular, posterior averages over factor loadings become meaningless and only model predictions are unambiguous. Here, we propose a parameterization based on Householder transformations, which remove the rotational symmetry of the posterior. Furthermore, by relying on results from random matrix theory, we establish the parameter distribution which leaves the model unchanged compared to the original rotationally symmetric formulation. In particular, we avoid the need to compute the Jacobian determinant of the parameter transformation. This allows us to efficiently implement probabilistic PCA in a rotation invariant fashion in any state of the art toolbox. Here, we implemented our model in the probabilistic programming language Stan and illustrate it on several examples.

Blackwell's theorem shows the equivalence of two preorders on the set of information channels. Here, we restate, and slightly generalize, his result in terms of random variables. Furthermore, we prove that the corresponding partial order is not a lattice; that is, least upper bounds and greatest lower bounds do not exist.

The information that two random variables Y, Z contain about a third random variable X can have aspects of shared information (contained in both Y and Z), of complementary information (only available from (Y, Z) together) and of unique information (contained exclusively in either Y or Z). Here, we study measures SĨ of shared, UĨ unique and CĨ complementary information introduced by Bertschinger et al. [1] which are motivated from a decision theoretic perspective. We find that in most cases the intuitive rule that more variables contain more information applies, with the exception that SĨ and CĨ information are not monotone in the target variable X. Additionally, we show that it is not possible to extend the bivariate information decomposition into SĨ, UĨ and CĨ to a non-negative decomposition on the partial information lattice of Williams and Beer [2]. Nevertheless, the quantities UĨ, SĨ and CĨ have a well-defined interpretation, even in the multivariate setting.

In this paper we analyze the relationship between the computational ca- pabilities of randomly connected networks of threshold gates in the time- series domain and their dynamical properties. In particular we propose a complexity measure which we find to assume its highest values near the edge of chaos, i.e. the transition from ordered to chaotic dynamics. Furthermore we show that the proposed complexity measure predicts the computational capabilities very well: only near the edge of chaos are such networks able to perform complex computations on time series. Ad- ditionally a simple synaptic scaling rule for self-organized criticality is presented and analyzed. 1 Introduction It has been proposed that extensive computational capabilities are achieved by systems whose dynamics is neither chaotic nor ordered but somewhere in between order and chaos. This has led to the idea of "computation at the edge of chaos". Early evidence for this hypothesis has been reported e.g. in [1]. The results of numerous computer simulations carried out in these studies suggested that there is a sharp transition between ordered and chaotic dynamics. Later on this was confirmed by Derrida and others [2]. They used ideas from statistical physics to develop an accurate mean-field theory which allowed to determine the critical parameters analytically. Because of the physical background, this theory focused on the autonomous dynamics of the system, i.e. its relaxation from an initial state (the input) to some terminal state (the output) without any external influences. In contrast to such "off-line" computations, we will focus in this article on time-series computations, i.e. mappings, also called filters, from a time-varying input signal to a time- varying output signal. Such "online" or real-time computations describe more adequately the input to output relation of systems like animals or autonomous robots which must react in real-time to a continuously changing stream of sensory input. The purpose of this paper is to analyze how the computational capabilities of randomly connected recurrent neural networks in the domain of real-time processing and the type of dynamics induced by the underlying distribution of synaptic weights are related to each other. In particular, we will show that for the types of neural networks considered in this pa- per (defined in Sec. 2) there also exists a transition from ordered to chaotic dynamics. This phase transition is determined using an extension of the mean-field approach described in [3] and [4] (Sec. 3). As the next step we propose a novel complexity measure (Sec. 4) which input network activity 10 20 neuron # 30 0 20 40 0 20 40 0 20 40 timesteps timesteps timesteps 0.4 0.2 0.8 0 0.6 -0.2 critical -0.4 ordered chaotic 0.4 mean activity -0.6 0.1 1 10 2 Figure 1: Networks of randomly connected threshold gates can exhibit ordered, critical and chaotic dynamics. In the upper row examples of the temporal evolution of the network state xt are shown (black: xi,t = 1, white: xi,t = 0, input as indicated above) for three different networks with parameters taken from the ordered, critical and chaotic regime, respectively. Parameters: K = 5, N = 500, u = -0.5, r = 0.3 and and 2 as indicated in the phase plot below. The background of the phase plot shows the mean activity a (see Sec. 3) of the networks depending on the parameters and 2. can be calculated using the mean-field theory developed in Sec. 3 and serves as a predic- tor for the computational capability of a network in the time-series domain. Employing a recently developed framework for analyzing real-time computations [5, 6] we investigate in Sec. 5 the relationship between network dynamics and the computational capabilities in the time-series domain. In Sec. 6 of this paper we propose and analyze a synaptic scaling rule for self-organized criticality (SOC) for the types of networks considered here. In con- trast to previous work [7], we do not only check that the proposed rule shows adaptation towards critical dynamics, but also show that the computational capabilities of the network are actually increased if the rule is applied. Relation to previous work: In [5], the so-called liquid state machine (LSM) approach was proposed and used do analyze the computational capabilities in the time-series domain of randomly connected networks of biologically inspired network models (composed of leaky integrate-and-fire neurons). We will use that approach to demonstrate that only near the edge of chaos, complex computations can be performed (see Sec. 5). A similar analysis for a restricted case (zero mean of synaptic weights) of the network model considered in this paper can be found in [4]. 2 The Network Model and its Dynamics We consider input driven recurrent networks consisting of N threshold gates with states xi {0, 1}. Each node i receives nonzero incoming weights wij from exactly K randomly chosen nodes j. Each nonzero connection weight wij is randomly drawn from a Gaussian distribution with mean and variance 2. Furthermore, the network is driven by an exter- nal input signal u() which is injected into each node. Hence, in summary, the update of the network state xt = (x1,t, . . . , xN,t) is given by xi,t = ( N w j=1 ij xj,t-1 + ut-1) which is applied to all neurons in parallel and where (h) = 1 if h 0 and (h) = 0 otherwise. In the following we consider a randomly drawn binary input signal u(): at each time step ut assumes the value u + 1 with probability r and the value u with probability 1 - r. This network model is similar to the one we have considered in [4]. However it differs in two important aspects: a) By using states xi {0, 1} we emphasis the asymmet- ric information encoding by spikes prevalent in biological neural systems and b) it is more general in the sense that the Gaussian distribution from which the non-zero weights are drawn is allowed to have an arbitrary mean R. This implies that the network activity a N t = 1 x N i=1 i,t can vary considerably for different parameters (compare Fig. 1) and enters all the calculations discussed in the rest of the paper. The top row of Fig. 1 shows typical examples of ordered, critical and chaotic dynamics (see the next section for a definition of order and chaos). The system parameters corresponding to each type of dynamics are indicated in the lower panel (phase plot). We refer to the (phase) transition from the ordered to the chaotic regime as the critical line (shown as the solid line in the phase plot). Note that increasing the variance 2 of the weights consistently leads to chaotic behavior. 3 The Critical Line: Order and Fading Memory versus Chaos To define the chaotic and ordered phase of an input driven network we use an approach which is similar to that proposed by Derrida and Pomeau [2] for autonomous systems: consider two (initial) network states with a certain (normalized) Hamming distance. These states are mapped to their corresponding successor states (using the same weight matrix) with the same input in each case and the change in the Hamming distance is observed. If small distances tend to grow this is a sign of chaos whereas if the distance tends to decrease this is a signature of order. Following closely the arguments in [4, 3] we developed a mean-field theory (see [8] for all details) which allows to calculate the update dt+1 = f (dt, at, ut) of the normalized Hamming distance dt = |{i : xi,t = ~ xi,t}|/N between two states xt and ~ xt as well as the update at+1 = A(at, ut) of the network activity in one time step. Note that dt+1 depends on the input ut (in contrast to [3]) and also on the activity at (in contrast to [4]). Hence the two-dimensional map Fu(dt, at) := (dt+1, at+1) = (f (dt, at, ut), A(at, ut)) describes the time evolution of dt and at given the input times series u(). Let us consider the steady state of the averaged Hamming distance f as well as the steady state of the averaged network activity a, i.e. (f , a) = limt F tu .1 If f = 0 we know that any state differences will eventually die out and the network is in the ordered phase. If on the other hand a small difference is amplified and never dies out we have f = 0 and the network is in the chaotic phase. Whether f = 0 or f = 0 can be decided by looking at the slope of the function f (, , ) at its fixed point f = 0. Since at does not depend on dt we calculate the averaged steady state activity a and determine the slope of the map rf (d, a, u + 1) + (1 - r)f (d, a, u) at the point (d, a) = (0, a). Accordingly we say that the network is in the ordered, critical or chaotic regime if < 1, = 1 or > 1 respectively. In [8] it is shown that the so called critical line = 1 where the phase transition from ordered to chaotic behavior occurs is given by K-1 K - 1 1 Pbf = an(1 - a)K-1-n(rQ(1, n, u + 1) + (1 - r)Q(1, n, u)) = (1) n K n=0 where Pbf denotes the probability (averaged over the inputs and the network activity) that a node will change its output if a single out of its K input bits is flipped.2 Examples of 1F tu denotes t-fold composition of the map Fu(, ) where in the k-th iteration the input uk is applied and denotes the average over all possible initial conditions and all input signals with a given statistics determined by u and r. 2The actual single bit-flip probability Q depends on the number n of inputs which are 1 and the K = 5 K = 10 0.4 0.4 0.1 0.3 0.3 0.07 0.2 0.2 0.08 0.06 0.1 0.1 0.05 0 0 0.06 -0.1 -0.1 0.04 -0.2 0.04 -0.2 0.03 NM-Separation NM-Separation -0.3 -0.3 0.02 -0.4 0.02 -0.4 0.01 -0.5 -0.5 -0.6 0 -0.6 0 0.1 1 10 0.1 1 10 2 2 Figure 2: N M -separation assumes high values on the critical line. The gray coded image shows the N M -separation in dependence on and 2 for K denoted in the panels, r = 0.3, u = -0.5 and b = 0.1. The solid line marks the critical values for and 2. critical lines that were calculated from this formula (marked by the solid lines) can be seen in Fig. 2 for K = 5 and K = 10.3 The ordered phase can also be described by using the notion of fading memory (see [5] and the references therein). Intuitively speaking in a network with fading memory a state xt is fully determined by a finite history ut-T , ut-T +1, . . . , ut-1, ut of the input u(). A slight reformulation of this property (see [6] and the references therein) shows that it is equivalent to the requirement that all state differences vanish, i.e. being in the ordered phase. Fading memory plays an important role in the "liquid state machine" framework [5] since together with the separation property (see below) it would in principle allow an appropriate readout function to deduce the recent input, or any function of it, from the network state. If on the other hand the network does not have fading memory (i.e. is in the chaotic regime) a given network state xt also contains "spurious" information about the initial conditions and hence it is hard or even impossible to deduce any features of the recent input. 4 NM-Separation as a Predictor for Computational Power The already mentioned separation property [5] is especially important if a network is to be useful for computations on input time-series: only if different input signals separate the network state, i.e. different inputs result in different states, it is possible for a readout function to respond differently. Hence it is necessary that any two different input time series for which the readout function should produce different outputs drive the recurrent network into two sufficiently different states. The mean field theory we have developed (see [8]) can be extended to describe the update dt+1 = s(dt, ...) of the Hamming distance that result from applying differ- ent inputs u() and ~ u() with a mean distance of b := Pr {ut = ~ ut}, i.e. the separa- tion. In summary the three-dimensional map Su,~u(dt, at, ~ at) := (dt+1, at+1, ~at+1) = (s(dt, at, ~at, ut, ~ ut), A(at, ut), A(~at, ~ ut)) fully describes the time evolution of the Ham- ming distance and the network activities. Again we consider the steady state of the averaged Hamming distance s and the network activities a, ~ a, i.e. (s, a, ~ a) = limt St . u,~ u The overall separation for a given input statistics (determined by u, r, and b) is then given by s. However, this overall separation measure can not be directly related to the computa- external input u and is given by Q(1, n, u) = -u (, n, n2) 1 - (-u - , , 2) d + - (, n, n2)(-u - , , 2)d where , denote the Gaussian density and cumulative den- -u sity respectively (see [8] for a detailed explanation). 3For each value of = -0.6 + k 0.01, k = 0 . . . 100 a search was conducted to find the value for 2 such that = 1. Numerical iterations of the function A were used to determine a. A B C 3bit parity (K = 5) 3bit parity (K = 10) 5bit random boolean functions 0.4 5 0.8 0.2 4 4 0.6 0 3 3 -0.2 MC (MI) 2 MC (MI) 0.4 mean MI 2 -0.4 0.2 1 1 -0.6 0 0 0 0.01 0.1 1 10 100 0.01 0.1 1 10 100 0.01 0.1 1 10 100 2 2 2 Figure 3: Real-time computation at the edge of chaos. A The gray coded image (an in- terpolation between the data points marked with open diamonds) shows the performance of trained networks in dependence of the parameters and 2 for the delayed 3-bit par- ity task. Performance is measured as the memory capacity M C = I(v, y()) where I(v, y()) is the mutual information between the classifier output v() and the target func- tion y() = PARITY(u t t- , ut- -1, ut- -2) measured on a test set. B Same as panel A but for K = 10. C Same as panel A but for an average over 50 randomly drawn Boolean functions f of the last 5 time steps, i.e. yt = f (ut, ut-1, ..., ut-4). tional power since chaotic networks separate even minor differences in the input to a very high degree. The part of this separation that is caused by the input distance b and not by the distance of some initial state is therefore given by s - f because f measures the state distance that is caused by differences in the initial states and remains even after long runs with the same inputs (see Sec. 3). Note that f is always zero in the ordered phase and non-zero in the chaotic phase. Since we want the complexity measure, which we will call N M -separation, to be a predictor for computational power we correct s - f by a term which accounts for the separation due to an all-dominant input drive. A suitable measure for this "imme- diate separation" i is the average increase in the Hamming distance if the system is run for a long time (t ) with equal inputs u() = ~ u() and then a single step with an input pair (v, ~ v) with an average difference of b = Pr {v, = ~ v} is applied: i = lim 1 t rv(1-r)1-vb|v-~v|(1-b)1-|v-~v| s(, , , v, ~ v) St -f . Hence v,~ v=0 u,u a measure of the network mediated separation N Msep due to input differences is given by N Msep = s - f - i (2) In Fig. 2 the N M -separation resulting from an input difference of b = 0.1 is shown in dependence of the network parameters and 2.4 Note that the N M -separation peaks very close to the critical line. Because of the computational importance of the separation property this also suggests that the computational capabilities of the networks will peak at the onset of chaos, which is confirmed in the next section. 5 Real-Time Computations at the Edge of Chaos To access the computational power of a network we make use of the so called "liquid state machine" framework which was proposed by Maass et.al. [5] and independently by Jaeger [6]. They put forward the idea that any complex time-series computation can be imple- mented by composing a system which consists of two conceptually different parts: a) a 4For each value of = -0.6 + k 0.05, k = 0 . . . 20, 10 values for 2 where chosen near the critical line and 10 other values where equally spaced (on a logarithmic scale) over the interval [0.02,50]. For each such pair (, 2) extensive numerical iterations of the map S where performed to obtain accurate estimates of s, f and i. Hopefully these numerical estimates can be replaced by analytic results in the future. properly chosen general-purpose recurrent network with "rich" dynamics and b) a read- out function that is trained to map the network state to the desired outputs (see [5, 6, 4] for more details). This approach is potentially successful if the general-purpose network encodes the relevant features of the input signal in the network state in such a way that the readout function can easily extract it. We will show that near the critical line the net- works considered in this paper encode the input in such a way that a simple linear classifier C(xt) = (w xt + w0) suffices to implement a broad range of complex nonlinear fil- ters. Note that in order to train the network for a given task only the parameters w RN , w0 R of the linear classifier are adjusted such that the actual network output vt = C(xt) is as close as possible to the target values yt. To access the computational power in a principled way networks with different parameters were tested on a delayed 3-bit parity task for increasing delays and on randomly drawn Boolean functions of the last 5 input bits. Note that these tasks are quite complex for the networks considered here since most of them are not linear separable (i.e. the parity function) and require memory. Hence to achieve good performance it is necessary that a state xt contains information about several input bits ut , t < t in a nonlinear transformed form such that a linear classifier C is sufficient to perform the nonlinear computations. The results are summarized in Fig. 3 where the performance (measured in terms of mutual information) on a test set between the network output and the target signal is shown for various parameter settings (for details see [4]). The highest performance is clearly achieved for parameter values close to the critical line where the phase transition occurs. This has been noted before [1]. In contrast to these previous results the networks used here are not optimized for any specific task but their computational capabilities are assessed by evaluating them for many different tasks. Therefore a network that is specifically designed for a single task will not show a good performance in this setup. These considerations suggest the following hypotheses regarding the computational function of generic recurrent neural circuits: to serve as a general-purpose temporal integrator, and simultaneously as a kernel (i.e., nonlinear projection into a higher dimensional space) to facilitate subsequent (linear) readout of information whenever it is needed. 6 Self-Organized Criticality via Synaptic Scaling Since the computational capabilities of a network depend crucially on having almost critical dynamics an adaptive system should be able to adjust its dynamics accordingly. Equ. (1) states that critical dynamics are achieved if the probability Pbf that a single bit- flip in the input shows up in the output should on average (over the external and internal input statistics given by u, r and a respectively) be equal to 1 . To allow for a rule that K can adjust the weights of each node a local estimate of Pbf must be available. This can be accomplished by estimating Pbf from the margin of each node, i.e. the distance of the internal activation from the firing threshold. Intuitively a node with an activation that is much higher or lower than its firing threshold is rather unlikely to change its output if a single bit in its input is flipped. Formally P i of node i is given by the average (over the bf internal and external input statistics) of the following quantity: 1 N (w K ij (1 - 2xj,t-1)(1 - 2xi,t) - mi,t) (3) j=1,wij =0 where mi,t = N w j=1 ij xj,t-1 + ut-1 denotes the margin of node i (see [8] for details). Each node now applies synaptic scaling to adjust itself towards the critical line. Accord- ingly we arrive at the following SOC-rule: 1 w (t) > 1 w 1+ ij if P esti bf K ij (t + 1) = (4) (1 + ) wij(t) if P esti(t) < 1 bf K A 50 100 neuron # 150 200 100 200 300 400 500 600 700 timesteps B C 1.5 1.5 1 1 bf bf K*P K*P 0.5 0.5 0 0 0 100 200 300 400 500 600 700 0 100 200 300 400 500 600 700 timesteps timesteps Figure 4: Self-organized criticality. A Time evolution of the network state xt starting in a chaotic regime while the SOC-rule (4) is active (black: xi,t = 1, white: xi,t = 0). Parameters: N = 500, K = 5, u = -0.5, r = 0.3, = 0 and initial 2 = 100. B Estimated Pbf . The dotted line shows how the node averaged estimate of Pbf evolves over time for the network shown in A. The running average of this estimate (thick black line) as used by the SOC-rule clearly shows that Pbf approaches its critical value (dashed line). C Same as B but for K = 10 and initial 2 = 0.01 in the ordered regime. where 0 < 1 is the learning rate and P esti(t) is a running average of the formula bf in Equ. (3) to estimate P i . Applying this rule in parallel to all nodes of the network is bf then able to adjust the network dynamics towards criticality as shown in Fig. 45. The upper row shows the time evolution of the network states xt while the SOC-rule (4) is running. It is clearly visible how the network dynamics changes from chaotic (the initial network had the parameters K = 5, = 0 and 2 = 100) to critical dynamics that respect the input signal. The lower row of Fig. 4 shows how the averaged estimated bit-flip probability 1 N P esti(t) approaches its critical value for the case of the above network and one N i=1 bf that started in the ordered regime (K = 10, = 0, 2 = 0.01). Since critical dynamics are better suited for information processing (see Fig. 3) it is ex- pected that the performance on the 3-bit parity task improves due to SOC. This is con- firmed in Fig. 5 which shows how the memory capacity M C (defined in Fig. 3) grows for networks that were initialized in the chaotic and ordered regime respectively. Note that the performance reached by these networks using the SOC-rule (4) is as high as for networks where the critical value for 2 is chosen apriori and stays at this level. This shows that rule (4) is stable in the sense that it keeps the dynamics critical and does not destroy the computational capabilities. 7 Discussion We developed a mean-field theory for input-driven networks which allows to determine the position of the transition line between ordered and chaotic dynamics with respect to the 5Here a learning rate of = 0.01 and an exponentially weighted running average with a time constant of 15 time steps were used. A B K = 5, start 2 = 100 (chaotic) K = 10, start 2 = 0.01 (ordered) 5 4 4 3 3 MC [bits] 2 2 MC [bits] 1 1 0 0 0 500 1000 1500 2000 0 500 1000 1500 2000 SOC steps SOC steps Figure 5: Time evolution of the performance with activated SOC-rule. A The plot shows the memory capacity M C (see Fig. 3) on the 3-bit parity task averaged over 25 networks ( standard deviation as error-bars) evaluated at the indicated time steps. At each evaluation time step the network weights were fixed and the M C was measured as in Fig. 3 by training the corresponding readouts from scratch. The networks were initialized in the chaotic regime. B Same as in A but for K = 10 and networks initialized in the ordered regime. parameters controlling the network connectivity and input statistics. Based on this theory we proposed a complexity measure (called N M -separation) which assumes its highest values at the critical line and shows a clear correlation with the computational power for real-time time-series processing. These results provide further evidence for the idea of "computation at the edge of chaos" [1] and support the hypothesis that dynamics near the critical line are expected to be a general property of input driven dynamical systems which support complex real-time computations. Therefore our analysis and the proposed complexity measure provide a new approach towards discovering dynamical principles that enable biological systems to do sophisticated information processing. Furthermore we have shown that a local rule for synaptic scaling is able to adjust the weights of a network towards critical dynamics. Additionally networks adjusting them- selves by this rule have been found to exhibit enhanced computational capabilities. Thereby systems can combine task-specific optimization provided by (supervised) learning rules with self-organization of its dynamics towards criticality. This provides an explanation how specific information can be processed while still being able to react to incoming sig- nals in a flexible way. Acknowledgement This work was supported in part by the PASCAL project #IST-2002-506778 of the European Community.

What makes a neural microcircuit computationally powerful? Or more precisely, which measurable quantities could explain why one microcir- cuit C is better suited for a particular family of computational tasks than another microcircuit C ? We propose in this article quantitative measures for evaluating the computational power and generalization capability of a neural microcircuit, and apply them to generic neural microcircuit mod- els drawn from different distributions. We validate the proposed mea- sures by comparing their prediction with direct evaluations of the com- putational performance of these microcircuit models. This procedure is applied first to microcircuit models that differ with regard to the spatial range of synaptic connections and with regard to the scale of synaptic efficacies in the circuit, and then to microcircuit models that differ with regard to the level of background input currents and the level of noise on the membrane potential of neurons. In this case the proposed method allows us to quantify differences in the computational power and gen- eralization capability of circuits in different dynamic regimes (UP- and DOWN-states) that have been demonstrated through intracellular record- ings in vivo. 1 Introduction Rather than constructing particular microcircuit models that carry out particular computa- tions, we pursue in this article a different strategy, which is based on the assumption that the computational function of cortical microcircuits is not fully genetically encoded, but rather emerges through various forms of plasticity ("learning") in response to the actual distribution of signals that the neural microcircuit receives from its environment. From this perspective the question about the computational function of cortical microcircuits C turns into the questions: a) What functions (i.e. maps from circuit inputs to circuit outputs) can the circuit C learn to compute. b) How well can the circuit C generalize a specific learned computational function to new inputs? We propose in this article a conceptual framework and quantitative measures for the in- vestigation of these two questions. In order to make this approach feasible, in spite of numerous unknowns regarding synaptic plasticity and the distribution of electrical and bio- chemical signals impinging on a cortical microcircuit, we make in the present first step of this approach the following simplifying assumptions: Particular neurons ("readout neurons") learn via synaptic plasticity to extract specific information encoded in the spiking activity of neurons in the circuit. We assume that the cortical microcircuit itself is highly recurrent, but that the impact of feedback that a readout neuron might send back into this circuit can be neglected.1 We assume that synaptic plasticity of readout neurons enables them to learn arbitrary linear transformations. More precisely, we assume that the input to such readout neuron can be approximated by a term n-1 w i=1 ixi(t), where n - 1 is the number of presynaptic neurons, xi(t) results from the output spike train of the ith presynaptic neuron by filtering it according to the low-pass filtering property of the membrane of the readout neuron,2 and wi is the efficacy of the synaptic connection. Thus wixi(t) models the time course of the contribution of previous spikes from the ith presynaptic neuron to the membrane potential at the soma of this readout neuron. We will refer to the vector x(t) as the circuit state at time t. Under these unpleasant but apparently unavoidable simplifying assumptions we propose new quantitative criteria based on rigorous mathematical principles for evaluating a neural microcircuit C with regard to questions a) and b). We will compare in sections 4 and 5 the predictions of these quantitative measures with the actual computational performance achieved by 132 different types of neural microcircuit models, for a fairly large number of different computational tasks. All microcircuit models that we consider are based on bio- logical data for generic cortical microcircuits (as described in section 3), but have different settings of their parameters. 2 Measures for the kernel-quality and generalization capability of neural microcircuits One interesting measure for probing the computational power of a neural circuit is the pair- wise separation property considered in [Maass et al., 2002]. This measure tells us to what extent the current circuit state x(t) reflects details of the input stream that occurred some time back in the past (see Fig. 1). Both circuit 2 and circuit 3 could be described as being chaotic since state differences resulting from earlier input differences persist. The "edge-of- chaos" [Langton, 1990] lies somewhere between points 1 and 2 according to Fig. 1c). But the best computational performance occurs between points 2 and 3 (see Fig. 2b)). Hence the "edge-of-chaos" is not a reliable predictor of computational power for circuits of spik- ing neurons. In addition, most real-world computational tasks require that the circuit gives a desired output not just for 2, but for a fairly large number m of significantly different inputs. One could of course test whether a circuit C can separate each of the m pairs of 2 1This assumption is best justified if such readout neuron is located for example in another brain area that receives massive input from many neurons in this microcircuit and only has diffuse back- wards projection. But it is certainly problematic and should be addressed in future elaborations of the present approach. 2One can be even more realistic and filter it also by a model for the short term dynamics of the synapse into the readout neuron, but this turns out to make no difference for the analysis proposed in this article. 8 a 4 b state separation 0.25 c 4 7 2 circuit 3 0.2 0 2 6 3 0 1 2 3 5 0.2 1 0.15 0.7 scale 2 4 0.1 0.5 W circuit 2 0.1 0.3 3 1 00 1 2 3 state separation 2 0.1 0.05 0.1 1 0.05 0.05 circuit 1 0 0.5 1 1.4 2 3 4 6 8 0 1.4 1.6 1.8 2 2.2 0 1 2 3 t [s] Figure 1: Pointwise separation property for different types of neural microcircuit models as specified in section 3. Each circuit C was tested for two arrays u and v of 4 input spike trains at 20 Hz over 3 s that differed only during the first second. a) Euclidean differences between resulting circuit states xu(t) and xv(t) for t = 3 s, averaged over 20 circuits C and 20 pairs u, v for each indicated value of and Wscale (see section 3). b) Temporal evolution of xu(t) - xv(t) for 3 different circuits with values of , Wscale according to the 3 points marked in panel a) ( = 1.4, 2, 3 and Wscale = 0.3, 0.7, 2 for circuit 1, 2, and 3 respectively). c) Pointwise separation along a straight line between point 1 and point 2 of panel a). such inputs. But even if the circuit can do this, we do not know whether a neural readout from such circuit would be able to produce given target outputs for these m inputs. Therefore we propose here the linear separation property as a more suitable quantitative measure for evaluating the computational power of a neural microcircuit (or more precisely: the kernel-quality of a circuit; see below). To evaluate the linear separation property of a circuit C for m different inputs u1, . . . , um (which are in this article always functions of time, i.e. input streams such as for example multiple spike trains) we compute the rank of the n m matrix M whose columns are the circuit states xu (t i 0 ) resulting at some fixed time t0 for the preceding input stream ui. If this matrix has rank m, then it is guaranteed that any given assignment of target outputs yi R at time t0 for the inputs ui can be implemented by this circuit C (in combination with a linear readout). In particular, each of the 2m possible binary classifications of these m inputs can then be carried out by a linear readout from this fixed circuit C. Obviously such insight is much more informative than a demonstration that some particular classification task can be carried out by such circuit C. If the rank of this matrix M has a value r < m, then this value r can still be viewed as a measure for the computational power of this circuit C, since r is the number of "degrees of freedom" that a linear readout has in assigning target outputs yi to these inputs ui (in a way which can be made mathematically precise with concepts of linear algebra). Note that this rank-measure for the linear separation property of a circuit C may be viewed as an empirical measure for its kernel-quality, i.e. for the complexity and diversity of nonlinear operations carried out by C on its input stream in order to boost the classification power of a subsequent linear decision-hyperplane (see [Vapnik, 1998]). Obviously the preceding measure addresses only one component of the computational per- formance of a neural circuit C. Another component is its capability to generalize a learnt computational function to new inputs. Mathematical criteria for generalization capability are derived in [Vapnik, 1998] (see ch. 4 of [Cherkassky and Mulier, 1998] for a compact ac- count of results relevant for our arguments). According to this mathematical theory one can quantify the generalization capability of any learning device in terms of the VC-dimension of the class H of hypotheses that are potentially used by that learning device.3 More pre- 3The VC-dimension (of a class H of maps H from some universe Suniv of inputs into {0, 1}) is defined as the size of the largest subset S Suniv which can be shattered by H. One says that S Suniv is shattered by H if for every map f : S {0, 1} there exists a map H in H such that H(u) = f (u) for all u S (this means that every possible binary classification of the inputs u S cisely: if VC-dimension (H) is substantially smaller than the size of the training set Strain, one can prove that this learning device generalizes well, in the sense that the hypothesis (or input-output map) produced by this learning device is likely to have for new examples an error rate which is not much higher than its error rate on Strain, provided that the new examples are drawn from the same distribution as the training examples (see equ. 4.22 in [Cherkassky and Mulier, 1998]). We apply this mathematical framework to the class HC of all maps from a set Suniv of inputs u into {0, 1} which can be implemented by a circuit C. More precisely: HC consists of all maps from Suniv into {0, 1} that a linear readout from circuit C with fixed internal parameters (weights etc.) but arbitrary weights w Rn of the readout (that classifies the circuit input u as belonging to class 1 if w xu(t0) 0, and to class 0 if w xu(t0) < 0) could possibly implement. Whereas it is very difficult to achieve tight theoretical bounds for the VC-dimension of even much simpler neural circuits, see [Bartlett and Maass, 2003], one can efficiently estimate the VC-dimension of the class HC that arises in our context for some finite ensemble Suniv of inputs (that contains all examples used for training or testing) by using the following mathematical result (which can be proved with the help of Radon's Theorem): Theorem 2.1 Let r be the rank of the n s matrix consisting of the s vectors xu(t0) for all inputs u in Suniv (we assume that Suniv is finite and contains s inputs). Then r VC-dimension(HC) r + 1. We propose to use the rank r defined in Theorem 2.1 as an estimate of VC-dimension(HC ), and hence as a measure that informs us about the generalization capability of a neural microcircuit C. It is assumed here that the set Suniv contains many noisy variations of the same input signal, since otherwise learning with a randomly drawn training set Strain Suniv has no chance to generalize to new noisy variations. Note that each family of computational tasks induces a particular notion of what aspects of the input are viewed as noise, and what input features are viewed as signals that carry information which is rel- evant for the target output for at least one of these computational tasks. For example for computations on spike patterns some small jitter in the spike timing is viewed as noise. For computations on firing rates even the sequence of interspike intervals and temporal rela- tions between spikes that arrive from different input sources are viewed as noise, as long as these input spike trains represent the same firing rates. Examples for both families of computational tasks will be discussed in this article. 3 Models for generic cortical microcircuits We test the validity of the proposed measures by comparing their predictions with direct evaluations of the computational performance for a large variety of models for generic cor- tical microcircuits consisting of 540 neurons. We used leaky-integrate-and-fire neurons4 and biologically quite realistic models for dynamic synapses.5 Neurons (20 % of which were randomly chosen to be inhibitory) were located on the grid points of a 3D grid of dimensions 6 6 15 with edges of unit length. The probability of a synaptic connection can be carried out by some hypothesis H in H). 4Membrane voltage V dVm m modeled by m = -(V dt m -Vresting )+Rm (Isyn(t)+Ibackground + Inoise), where m = 30 ms is the membrane time constant, Isyn models synaptic inputs from other neurons in the circuits, Ibackground models a constant unspecific background input and Inoise models noise in the input. 5Short term synaptic dynamics was modeled according to [Markram et al., 1998], with distribu- tions of synaptic parameters U (initial release probability), D (time constant for depression), F (time constant for facilitation) chosen to reflect empirical data (see [Maass et al., 2002] for details). from neuron a to neuron b was proportional to exp(-D2(a, b)/2), where D(a, b) is the Euclidean distance between a and b, and regulates the spatial scaling of synaptic connec- tivity. Synaptic efficacies w were chosen randomly from distributions that reflect biological data (as in [Maass et al., 2002]), with a common scaling factor Wscale. 8 0.7 b 4 a 2 3 0.65 1 0.7 scale 2 0.5 W 0.3 1 0.6 0 50 100 150 200 0 50 100 150 200 0.1 t [ms] t [ms] 0.05 0.5 1 1.4 2 3 4 6 8 Figure 2: Performance of different types of neural microcircuit models for classification of spike patterns. a) In the top row are two examples of the 80 spike patterns that were used (each consisting of 4 Poisson spike trains at 20 Hz over 200 ms), and in the bottom row are examples of noisy variations (Gaussian jitter with SD 10 ms) of these spike patterns which were used as circuit inputs. b) Fraction of examples (for 200 test examples) that were correctly classified by a linear readout (trained by linear regression with 500 training examples). Results are shown for 90 different types of neural microcircuits C with varying on the x-axis and Wscale on the y-axis (20 randomly drawn circuits and 20 target classification functions randomly drawn from the set of 280 possible classification functions were tested for each of the 90 different circuit types, and resulting correctness-rates were averaged. The mean SD of the results is 0.028.). Points 1, 2, 3 defined as in Fig. 1. Linear readouts from circuits with n - 1 neurons were assumed to compute a weighted sum n-1 w i=1 ixi(t) + w0 (see section 1). In order to simplify notation we assume that the vector x(t) contains an additional constant component x0(t) = 1, so that one can write w x(t) instead of n-1 w i=1 ixi(t) + w0. In the case of classification tasks we assume that the readout outputs 1 if w x(t) 0, and 0 otherwise. 4 Evaluating the influence of synaptic connectivity on computational performance Neural microcircuits were drawn from the distribution described in section 3 for 10 differ- ent values of (which scales the number and average distance of synaptically connected neurons) and 9 different values of Wscale (which scales the efficacy of all synaptic connec- tions). 20 microcircuit models C were drawn for each of these 90 different assignments of values to and Wscale. For each circuit a linear readout was trained to perform one (randomly chosen) out of 280 possible classification tasks on noisy variations u of 80 fixed spike patterns as circuit inputs u. The target performance of any such circuit input was to output at time t = 100 ms the class (0 or 1) of the spike pattern from which the preceding circuit input had been generated (for some arbitrary partition of the 80 fixed spike patterns into two classes. Each spike pattern u consisted of 4 Poisson spike trains over 200 ms. Per- formance results are shown in Fig. 2b for 90 different types of neural microcircuit models. We now test the predictive quality of the two proposed measures for the computational power of a microcircuit on spike patterns. One should keep in mind that the proposed measures do not attempt to test the computational capability of a circuit for one particu- lar computational task, but for any distribution on Suniv and for a very large (in general infinitely large) family of computational tasks that only have in common a particular bias regarding which aspects of the incoming spike trains may carry information that is relevant for the target output of computations, and which aspects should be viewed as noise. Fig. 3a explains why the lower left part of the parameter map in Fig. 2b is less suitable for any 8 8 8 a b c 20 4 450 4 450 4 2 400 2 400 2 3 15 1 350 1 350 1 0.7 0.7 0.7 scale 2 0.5 W 0.5 0.5 10 300 300 0.3 0.3 0.3 1 250 250 5 0.1 0.1 200 0.1 200 0.05 0.05 0.05 0 0.5 1 1.4 2 3 4 6 8 0.5 1 1.4 2 3 4 6 8 0.5 1 1.4 2 3 4 6 8 Figure 3: Values of the proposed measures for computations on spike patterns. a) Kernel-quality for spike patterns of 90 different circuit types (average over 20 circuits, mean SD = 13; For each circuit, the average over 5 different sets of spike patterns was used).6 b) Generalization capability for spike patterns: estimated VC-dimension of HC (for a set Suniv of inputs u consisting of 500 jittered versions of 4 spike patterns), for 90 different circuit types (average over 20 circuits, mean SD = 14; For each circuit, the average over 5 different sets of spike patterns was used). c) Difference of both measures (mean SD = 5.3). This should be compared with actual computational performance plotted in Fig. 2b. Points 1, 2, 3 defined as in Fig. 1. such computation, since there the kernel-quality of the circuits is too low. Fig. 3b explains why the upper right part of the parameter map in Fig. 2b is less suitable, since a higher VC-dimension (for a training set of fixed size) entails poorer generalization capability. We are not aware of a theoretically founded way of combining both measures into a single value that predicts overall computational performance. But if one just takes the difference of both measures then the resulting number (see Fig. 3c) predicts quite well which types of neural microcircuit models perform well for the particular computational tasks considered in Fig. 2b. 5 Evaluating the computational power of neural microcircuit models in UP- and DOWN-states Data from numerous intracellular recordings suggest that neural circuits in vivo switch be- tween two different dynamic regimes that are commonly referred to as UP- and DOWN states. UP-states are characterized by a bombardment with synaptic inputs from recurrent activity in the circuit, resulting in a membrane potential whose average value is signifi- cantly closer to the firing threshold, but also has larger variance. We have simulated these different dynamic regimes by varying the background current Ibackground and the noise current Inoise. Fig. 4a shows that one can simulate in this way different dynamic regimes of the same circuit where the time course of the membrane potential qualitatively matches data from intracellular recordings in UP- and DOWN-states (see e.g. [Shu et al., 2003]). We have tested the computational performance of circuits in 42 different dynamic regimes (for 7 values of Ibackground and 6 values of Inoise) with 3 complex nonlinear computations on firing rates of circuit inputs.7 Inputs u consisted of 4 Poisson spike trains with time- varying rates (drawn independently every 30 ms from the interval of 0 to 80 Hz for the first two and the second two of 4 input spike trains, see middle row of Fig. 4a for a sample). Let f1(t) (f2(t)) be the actual sum of rates normalized to the interval [0, 1] for the first 6The rank of the matrix consisting of 500 circuit states xu(t) for t = 200 ms was computed for 500 spike patterns over 200 ms as described in section 2, see Fig. 2a. 7Computations on firing rates were chosen as benchmark tasks both because UP states were con- jectured to enhance the performance for such tasks, and because we want to show that the proposed measures are applicable to other types of computational tasks than those considered in section 4. 16 a 100 UP-state [mV] 14 m 50 V 12 0 16 100 DOWN-state [mV] 14 m 50 V 12 0 300 350 400 450 500 350 400 450 500 t [ms] t [ms] b c d 10 10 10 120 6 70 6 6 0.2 4.5 UP 4.5 100 4.5 60 0.15 3.2 3.2 3.2 80 I noise 50 0.1 1.9 DOWN 1.9 60 1.9 40 0.05 40 30 0 0.6 0.6 20 0.6 11.5 12 12.5 13.5 14.3 11.5 12 12.5 13.5 14.3 11.5 12 12.5 13.5 14.3 e f g 10 10 10 0.25 6 6 6 0.3 4.5 0.7 4.5 4.5 0.2 3.2 3.2 3.2 I noise 0.6 1.9 1.9 0.15 1.9 0.25 0.5 0.1 0.2 0.6 0.6 0.6 11.5 12 12.5 13.5 14.3 11.5 12 12.5 13.5 14.3 11.5 12 12.5 13.5 14.3 I I I background background background Figure 4: Analysis of the computational power of simulated neural microcircuits in different dy- namic regimes. a) Membrane potential (for a firing threshold of 15 mV) of two randomly selected neurons from circuits in the two parameter regimes marked in panel b), as well as spike rasters for the same two parameter regimes (with the actual circuit inputs shown between the two rows). b) Estimates of the kernel-quality for input streams u with 34 different combinations of firing rates from 0, 20, 40 Hz in the 4 input spike trains (mean SD = 12). c) Estimate of the VC-dimension for a set Suniv of inputs consisting of 200 different spike trains u that represent 2 different combinations of firing rates (mean SD = 4.6). d) Difference of measures from panels b and c (after scaling each lin- early into a common range [0,1]). e), f), g): Evaluation of the computational performance (correlation coefficient; all for test data; mean SD is 0.06, 0.04, and 0.03 for panels e), f), and g) respectively.) of the same circuits in different dynamic regimes for computations involving multiplication and abso- lute value of differences of firing rates (see text). The theoretically predicted parameter regime with good computational performance for any computations on firing rates (see panel d) agrees quite well with the intersection of areas with good computational performance in panels e, f, g. two (second two) input spike trains computed from the time interval [t - 30ms, t]. The computational tasks considered in Fig. 4 were to compute online (and in real-time) every 30 ms the functions f1(t) f2(t) (see panel e), to decide whether the value of the product f1(t) f2(t) lies in the interval [0.1, 0.3] or lies outside of this interval (see panel f), and to decide whether the absolute value of the difference f1(t) - f2(t) is greater than 0.25 (see panel g). We wanted to test whether the proposed measures for computational power and general- ization capability were able to make reasonable predictions for this completely different parameter map, and for computations on firing rates instead of spike patterns. It turns out that also in this case the kernel-quality (Fig. 4b) explains why circuits in the dynamic regime corresponding to the left-hand side of the parameter map have inferior computa- tional power for all three computations on firing rates (see Fig. 4 e,f,g). The VC-dimension (Fig. 4c) explains the decline of computational performance in the right part of the pa- rameter map. The difference of both measures (Fig. 4d) predicts quite well the dynamic regime where high performance is achieved for all three computational tasks considered in Fig. 4 e,f,g. Note that Fig. 4e has high performance in the upper right corner, in spite of a very high VC-dimension. This could be explained by the inherent bias of linear readouts to compute smooth functions on firing rates, which fits particularly well to this particular target output. If one estimates kernel-quality and VC-dimension for the same circuits, but for computa- tions on sparse spike patterns (for an input ensemble Suniv similarly as in section 4), one finds that circuits at the lower left corner of this parameter map (corresponding to DOWN- states) are predicted to have better computational performance for these computations on sparse input. This agrees quite well with direct evaluations of computational performance (not shown). Hence the proposed quantitative measures may provide a theoretical founda- tion for understanding the computational function of different states of neural activity.

Depending on the connectivity, recurrent networks of simple computational units can show very different types of dynamics, ranging from totally ordered to chaotic. We analyze how the type of dynamics (ordered or chaotic) exhibited by randomly connected networks of threshold gates driven by a time-varying input signal depends on the parameters describing the distribution of the connectivity matrix. In particular, we calculate the critical boundary in parameter space where the transition from ordered to chaotic dynamics takes place. Employing a recently developed framework for analyzing real-time computations, we show that only near the critical boundary can such networks perform complex computations on time series. Hence, this result strongly supports conjectures that dynamical systems that are capable of doing complex computational tasks should operate near the edge of chaos, that is, the transition from ordered to chaotic dynamics.