Exact subthreshold integration with continuous spike times in discrete-time neural network simulations

Very large networks of spiking neurons can be simulated efficiently in parallel under the constraint that spike times are bound to an equidistant time grid. Within this scheme, the subthreshold dynamics of a wide class of integrate-and-fire-type neuron models can be integrated exactly from one grid point to the next. However, the loss in accuracy caused by restricting spike times to the grid can have undesirable consequences, which has led to interest in interpolating spike times between the grid points to retrieve an adequate representation of network dynamics. We demonstrate that the exact integration scheme can be combined naturally with off-grid spike events found by interpolation. We show that by exploiting the existence of a minimal synaptic propagation delay, the need for a central event queue is removed, so that the precision of event-driven simulation on the level of single neurons is combined with the efficiency of time-driven global scheduling. Further, for neuron models with linear subthreshold dynamics, even local event queuing can be avoided, resulting in much greater efficiency on the single-neuron level. These ideas are exemplified by two implementations of a widely used neuron model. We present a measure for the efficiency of network simulations in terms of their integration error and show that for a wide range of input spike rates, the novel techniques we present are both more accurate and faster than standard techniques.     

Exact simulation of integrate-and-fire models with exponential currents

Neural networks can be simulated exactly using event-driven strategies, in which the algorithm advances directly from one spike to the next spike. It applies to neuron models for which we have (1) an explicit expression for the evolution of the state variables between spikes and (2) an explicit test on the state variables that predicts whether and when a spike will be emitted. In a previous work, we proposed a method that allows exact simulation of an integrate-and-fire model with exponential conductances, with the constraint of a single synaptic time constant. In this note, we propose a method, based on polynomial root finding, that applies to integrate-and-fire models with exponential currents, with possibly many different synaptic time constants. Models can include biexponential synaptic currents and spike-triggered adaptation currents.     

On the simulation of nonlinear bidimensional spiking neuron models

Bidimensional spiking models are garnering a lot of attention for their simplicity and their ability to reproduce various spiking patterns of cortical neurons and are used particularly for large network simulations. These models describe the dynamics of the membrane potential by a nonlinear differential equation that blows up in finite time, coupled to a second equation for adaptation. Spikes are emitted when the membrane potential blows up or reaches a cutoff θ. The precise simulation of the spike times and of the adaptation variable is critical, for it governs the spike pattern produced and is hard to compute accurately because of the exploding nature of the system at the spike times. We thoroughly study the precision of fixed time-step integration schemes for this type of model and demonstrate that these methods produce systematic errors that are unbounded, as the cutoff value is increased, in the evaluation of the two crucial quantities: the spike time and the value of the adaptation variable at this time. Precise evaluation of these quantities therefore involves very small time steps and long simulation times. In order to achieve a fixed absolute precision in a reasonable computational time, we propose here a new algorithm to simulate these systems based on a variable integration step method that either integrates the original ordinary differential equation or the equation of the orbits in the phase plane, and compare this algorithm with fixed time-step Euler scheme and other more accurate simulation algorithms.     

Model-based decoding, information estimation, and change-point detection techniques for multineuron spike trains

One of the central problems in systems neuroscience is to understand how neural spike trains convey sensory information. Decoding methods, which provide an explicit means for reading out the information contained in neural spike responses, offer a powerful set of tools for studying the neural coding problem. Here we develop several decoding methods based on point-process neural encoding models, or forward models that predict spike responses to stimuli. These models have concave log-likelihood functions, which allow efficient maximum-likelihood model fitting and stimulus decoding. We present several applications of the encoding model framework to the problem of decoding stimulus information from population spike responses: (1) a tractable algorithm for computing the maximum a posteriori (MAP) estimate of the stimulus, the most probable stimulus to have generated an observed single- or multiple-neuron spike train response, given some prior distribution over the stimulus; (2) a gaussian approximation to the posterior stimulus distribution that can be used to quantify the fidelity with which various stimulus features are encoded; (3) an efficient method for estimating the mutual information between the stimulus and the spike trains emitted by a neural population; and (4) a framework for the detection of change-point times (the time at which the stimulus undergoes a change in mean or variance) by marginalizing over the posterior stimulus distribution. We provide several examples illustrating the performance of these estimators with simulated and real neural data.     

Determination of firing times for the stochastic Fitzhugh-Nagumo neuronal model

We present for the first time an analytical approach for determining the time of firing of multicomponent nonlinear stochastic neuronal models. We apply the theory of first exit times for Markov processes to the Fitzhugh-Nagumo system with a constant mean gaussian white noise input, representing stochastic excitation and inhibition. Partial differential equations are obtained for the moments of the time to first spike. The observation that the recovery variable barely changes in the prespike trajectory leads to an accurate one-dimensional approximation. For the moments of the time to reach threshold, this leads to ordinary differential equations that may be easily solved. Several analytical approaches are explored that involve perturbation expansions for large and small values of the noise parameter. For ranges of the parameters appropriate for these asymptotic methods, the perturbation solutions are used to establish the validity of the one-dimensional approximation for both small and large values of the noise parameter. Additional verification is obtained with the excellent agreement between the mean and variance of the firing time found by numerical solution of the differential equations for the one-dimensional approximation and those obtained by simulation of the solutions of the model stochastic differential equations. Such agreement extends to intermediate values of the noise parameter. For the mean time to threshold, we find maxima at small noise values that constitute a form of stochastic resonance. We also investigate the dependence of the mean firing time on the initial values of the voltage and recovery variables when the input current has zero mean.     

Analytical integrate-and-fire neuron models with conductance-based dynamics and realistic postsynaptic potential time course for event-driven simulation strategies

In a previous paper (Rudolph & Destexhe, 2006), we proposed various models, the gIF neuron models, of analytical integrate-and-fire (IF) neurons with conductance-based (COBA) dynamics for use in event-driven simulations. These models are based on an analytical approximation of the differential equation describing the IF neuron with exponential synaptic conductances and were successfully tested with respect to their response to random and oscillating inputs. Because they are analytical and mathematically simple, the gIF models are best suited for fast event-driven simulation strategies. However, the drawback of such models is they rely on a nonrealistic postsynaptic potential (PSP) time course, consisting of a discontinuous jump followed by a decay governed by the membrane time constant. Here, we address this limitation by conceiving an analytical approximation of the COBA IF neuron model with the full PSP time course. The subthreshold and suprathreshold response of this gIF4 model reproduces remarkably well the postsynaptic responses of the numerically solved passive membrane equation subject to conductance noise, while gaining at least two orders of magnitude in computational performance. Although the analytical structure of the gIF4 model is more complex than that of its predecessors due to the necessity of calculating future spike times, a simple and fast algorithmic implementation for use in large-scale neural network simulations is proposed.     

Analysis of a ring-based loss network

Loss networks are a class of resource allocation models which have proved useful in the study and design of communication networks. A major aim when analysing these models is to obtain measures of performance such as congestion probabilities and throughput. Even for moderately sized networks it is often not possible to calculate the quantities of interest exactly in reasonable time. Accordingly there is motivation to investigate methods of accurate approximation. This paper is concerned with estimating the performance measures of a ring-shaped network with branches. Two methods are presented: an asymptotic analysis made possible by the network's cyclic structure, and a technique utilising the network's symmetry which could be applied more generally.     

Higher-order distributional properties in closed queueing networks

In many real-life computer and networking applications, the distributions of service times, or times between arrivals of requests, or both, can deviate significantly from the memoryless negative exponential distribution that underpins the product-form solution for queueing networks. Frequently, the coefficient of variation of the distributions encountered is well in excess of one, which would be its value for the exponential. For closed queueing networks with non-exponential servers there is no known general exact solution, and most, if not all, approximation methods attempt to account for the general service time distributions through their first two moments.We consider two simple closed queueing networks which we solve exactly using semi-numerical methods. These networks depart from the structure leading to a product-form solution only to the extent that the service time at a single node is non-exponential. We show that not only the coefficients of variation but also higher-order distributional properties can have an important effect on such customary steady-state performance measures as the mean number of customers at a resource or the resource utilization level in a closed network.Additionally, we examine the state that a request finds upon its arrival at a server, which is directly tied to the resulting quality of service. Although the well-known Arrival Theorem holds exactly only for product-form networks of queues, some approximation methods assume that it can be applied to a reasonable degree also in other closed queueing networks. We investigate the validity of this assumption in the two closed queueing models considered. Our results show that, even in the case when there is a single non-exponential server in the network, the state found upon arrival may be highly sensitive to higher-order properties of the service time distribution, beyond its mean and coefficient of variation.This dependence of mean numbers of customers at a server on higher-order distributional properties is in stark contrast with the situation in the familiar open M/G/1 queue. Thus, our results put into question virtually all traditional approximate solutions, which concentrate on the first two moments of service time distributions.     

Capacity Gain of Mixed Multicast/Unicast Transport Schemes in a TV Distribution Network

This paper presents three approaches to estimate the required resources in an infrastructure where digital TV channels can be delivered in unicast or multicast (broadcast) mode. Such situations arise for example in Cable TV, IPTV distribution networks or in (future) hybrid mobile TV networks. The three approaches presented are an exact calculation, a Gaussian approximation and a simulation tool. We investigate two scenarios that allow saving bandwidth resources. In a static scenario, the most popular channels are multicast and the less popular channels rely on unicast. In a dynamic scenario, the list of multicast channels is dynamic and governed by the users' behavior. We prove that the dynamic scenario always outperforms the static scenario. We demonstrate the robustness, versatility and the limits of our three approaches. The exact calculation application is limited because it is computationally expensive for cases with large numbers of users and channels, while the Gaussian approximation is good exactly for such systems. The simulation tool takes long to yield results for small blocking probabilities. We explore the capacity gain regions under varying model parameters. Finally, we illustrate our methods by discussing some realistic network scenarios using channel popularities based on measurement data as much as possible.     

Advances in Design and Application of Spiking Neural Networks

This paper presents new findings in the design and application of biologically plausible neural networks based on spiking neuron models, which represent a more plausible model of real biological neurons where time is considered as an important feature for information encoding and processing in the brain. The design approach consists of an evolutionary strategy based supervised training algorithm, newly developed by the authors, and the use of different biologically plausible neuronal models. A dynamic synapse (DS) based neuron model, a biologically more detailed model, and the spike response model (SRM) are investigated in order to demonstrate the efficacy of the proposed approach and to further our understanding of the computing capabilities of the nervous system. Unlike the conventional synapse, represented as a static entity with a fixed weight, employed in conventional and SRM-based neural networks, a DS is weightless and its strength changes upon the arrival of incoming input spikes. Therefore its efficacy depends on the temporal structure of the impinging spike trains. In the proposed approach, the training of the network free parameters is achieved using an evolutionary strategy where, instead of binary encoding, real values are used to encode the static and DS parameters which underlie the learning process. The results show that spiking neural networks based on both types of synapse are capable of learning non-linearly separable data by means of spatio-temporal encoding. Furthermore, a comparison of the obtained performance with classical neural networks (multi-layer perceptrons) is presented.     

Standard clock simulation and ordinal optimization applied to admission control in integrated communication networks

In this paper we apply the ideas of ordinal optimization and the technique of Standard Clock (SC) simulation to the voice-call admission-control problem in integrated voice/data multihop radio networks. This is an important problem in networking that is not amenable to exact analysis by means of the usual network modeling techniques. We first describe the use of the SC approach on sequential machines, and quantify the speedup in simulation time that is achieved by its use in a number of queueing examples. We then develop an efficient simulation model for wireless integrated networks based on the use of the SC approach, which permits the parallel simulation of a large number of admission-control policies, thereby reducing computation time significantly. This model is an extension of the basic SC approach in that it incorporates fixed-length data packets, whereas SC simulation is normally limited to systems with exponentially distributed interevent times. Using this model, we demonstrate the effectiveness of ordinal-optimization techniques, which provide a remarkably good ranking of admission-control policies after relatively short simulation runs, thereby facilitating the rapid determination of good policies. Moreover, we demonstrate that the use of crude, inaccurate analytical and simulation models can provide highly accurate policy rankings that can be used in conjunction with ordinal-optimization methods, provided that they incorporate the key aspects of system operation.     

Analytical integrate-and-fire neuron models with conductance-based dynamics for event-driven simulation strategies

Event-driven simulation strategies were proposed recently to simulate integrate-and-fire (IF) type neuronal models. These strategies can lead to computationally efficient algorithms for simulating large-scale networks of neurons; most important, such approaches are more precise than traditional clock-driven numerical integration approaches because the timing of spikes is treated exactly. The drawback of such event-driven methods is that in order to be efficient, the membrane equations must be solvable analytically, or at least provide simple analytic approximations for the state variables describing the system. This requirement prevents, in general, the use of conductance-based synaptic interactions within the framework of event-driven simulations and, thus, the investigation of network paradigms where synaptic conductances are important. We propose here a number of extensions of the classical leaky IF neuron model involving approximations of the membrane equation with conductance-based synaptic current, which lead to simple analytic expressions for the membrane state, and therefore can be used in the event-driven framework. These conductance-based IF (gIF) models are compared to commonly used models, such as the leaky IF model or biophysical models in which conductances are explicitly integrated. All models are compared with respect to various spiking response properties in the presence of synaptic activity, such as the spontaneous discharge statistics, the temporal precision in resolving synaptic inputs, and gain modulation under in vivo-like synaptic bombardment. Being based on the passive membrane equation with fixed-threshold spike generation, the proposed gIF models are situated in between leaky IF and biophysical models but are much closer to the latter with respect to their dynamic behavior and response characteristics, while still being nearly as computationally efficient as simple IF neuron models. gIF models should therefore provide a useful tool for efficient and precise simulation of large-scale neuronal networks with realistic, conductance-based synaptic interactions.     

Analytical modeling of multistage, multipath networks

Because of their ability to tolerate faults, multipath, multistage networks provide useful interconnection schemes for large-scale parallel computers. However, the analytical models that have been used to analyze the performance of Banyan networks cannot be used to evaluate the performance of multipath networks. We present here what we believe to be the first analytical model that allows calculation of the bandwidth of the general class of unbuffered, packet-switched, multipath, multistage networks. The equations yielded by the model can be solved either exactly or by Monte Carlo approximation. The model agrees well with the results of a more complex simulation and provides a first step towards solution of the open problem of modeling of buffered, packet-switched, multipath, multistage networks     

On Taylor series expansions for waiting times in tandem queues: an algorithm for calculating the coefficients and an investigation of the approximation error

Recently, a Taylor series expansion was developed for expected stationary waiting times in open (max,+)-linear stochastic systems with Poisson input process; these systems cover various instances of queueing networks.As an application, we present an algorithm for calculating the coefficients for infinite capacity tandem queueing networks with discrete service-time distributions. The algorithm works quite efficiently if the random vector of the service times of all servers is concentrated at a small number of atoms. We investigate the relative error of the Taylor approximation by simulation; in many cases, it follows very well a simple expression which holds exactly for independent, exponentially distributed servers.     

Spike train probability models for stimulus-driven leaky integrate-and-fire neurons

Mathematical models of neurons are widely used to improve understanding of neuronal spiking behavior. These models can produce artificial spike trains that resemble actual spike train data in important ways, but they are not very easy to apply to the analysis of spike train data. Instead, statistical methods based on point process models of spike trains provide a wide range of data-analytical techniques. Two simplified point process models have been introduced in the literature: the time-rescaled renewal process (TRRP) and the multiplicative inhomogeneous Markov interval (m-IMI) model. In this letter we investigate the extent to which the TRRP and m-IMI models are able to fit spike trains produced by stimulus-driven leaky integrate-and-fire (LIF) neurons. With a constant stimulus, the LIF spike train is a renewal process, and the m-IMI and TRRP models will describe accurately the LIF spike train variability. With a time-varying stimulus, the probability of spiking under all three of these models depends on both the experimental clock time relative to the stimulus and the time since the previous spike, but it does so differently for the LIF, m-IMI, and TRRP models. We assessed the distance between the LIF model and each of the two empirical models in the presence of a time-varying stimulus. We found that while lack of fit of a Poisson model to LIF spike train data can be evident even in small samples, the m-IMI and TRRP models tend to fit well, and much larger samples are required before there is statistical evidence of lack of fit of the m-IMI or TRRP models. We also found that when the mean of the stimulus varies across time, the m-IMI model provides a better fit to the LIF data than the TRRP, and when the variance of the stimulus varies across time, the TRRP provides the better fit.     

Approximate transient analysis of large stochastic models with WinPEPSY-QNS

In this paper we present a new algorithm for the approximate transient analysis of large stochastic models. The new algorithm is based on the self-correcting analysis principle for continuous-time Markov chains (CTMC). The approach uses different time dependent aggregations of the CTMC of a stochastic model. With the method of uniformization the transient state probabilities of each aggregated CTMC for a time step are calculated. The derived probabilities are used for the construction of stronger aggregations, which are applied for the correction of the transition rates of the previous aggregations. This is done step by step, until the final time is reached. High aggregations of the original continuous-time Markov chain lead to a time and space efficient computational effort. Therefore the approximate transient analysis method based on the self-correcting aggregation can be used for models with large state spaces. For queuing networks with phase-type distributions of the service times this newly developed algorithm is implemented in WinPEPSY-QNS, a tool for performance evaluation and prediction of stochastic models based on queuing networks. It consists of a graphical editor for the construction of queuing networks and an easy-to-use evaluation component, which offers suitable analysis methods. The newly implemented algorithm is used for the analysis of several examples, and the results are compared to the results of simulation runs where exact values cannot be achieved.     

Response time analysis of parallel computer and storage systems

Fork-join structures have gained increased importance in recent years as a means of modeling parallelism in computer and storage systems. The basic fork-join model is one in which a job arriving at a parallel system splits into K independent tasks that are assigned to K unique, homogeneous servers. In the paper, a simple response time approximation is derived for parallel systems with exponential service time distributions. The approximation holds for networks modeling several devices, both parallel and nonparallel. (In the case of closed networks containing a stand-alone parallel system, a mean response time bound is derived.) In addition, the response time approximation is extended to cover the more realistic case wherein a job splits into an arbitrary number of tasks upon arrival at a parallel system. Simulation results for closed networks with stand-alone parallel subsystems and exponential service time distributions indicate that the response time approximation is, on average, within 3 percent of the seeded response times. Similarly, simulation results with nonexponential distributions also indicate that the response time approximation is close to the seeded values. Potential applications of our results include the modeling of data placement in disk arrays and the execution of parallel programs in multiprocessor and distributed systems     

Learning the structure of dynamic Bayesian networks from time series and steady state measurements

Dynamic Bayesian networks (DBN) are a class of graphical models that has become a standard tool for modeling various stochastic time-varying phenomena. In many applications, the primary goal is to infer the network structure from measurement data. Several efficient learning methods have been introduced for the inference of DBNs from time series measurements. Sometimes, however, it is either impossible or impractical to collect time series data, in which case, a common practice is to model the non-time series observations using static Bayesian networks (BN). Such an approach is obviously sub-optimal if the goal is to gain insight into the underlying dynamical model. Here, we introduce Bayesian methods for the inference of DBNs from steady state measurements. We also consider learning the structure of DBNs from a combination of time series and steady state measurements. We introduce two different methods: one that is based on an approximation and another one that provides exact computation. Simulation results demonstrate that dynamic network structures can be learned to an extent from steady state measurements alone and that inference from a combination of steady state and time series data has the potential to improve learning performance relative to the inference from time series data alone.