Spike Train Data Research Articles

A semiparametric Bayesian model for detecting synchrony among multiple neurons.

We propose a scalable semiparametric Bayesian model to capture dependencies among multiple neurons by detecting their cofiring (possibly with some lag time) patterns over time. After discretizing time so there is at most one spike at each interval, the resulting sequence of 1s (spike) and 0s (silence) for each neuron is modeled using the logistic function of a continuous latent variable with a gaussian process prior. For multiple neurons, the corresponding marginal distributions are coupled to their joint probability distribution using a parametric copula model. The advantages of our approach are as follows. The nonparametric component (i.e., the gaussian process model) provides a flexible framework for modeling the underlying firing rates, and the parametric component (i.e., the copula model) allows us to make inferences regarding both contemporaneous and lagged relationships among neurons. Using the copula model, we construct multivariate probabilistic models by separating the modeling of univariate marginal distributions from the modeling of a dependence structure among variables. Our method is easy to implement using a computationally efficient sampling algorithm that can be easily extended to high-dimensional problems. Using simulated data, we show that our approach could correctly capture temporal dependencies in firing rates and identify synchronous neurons. We also apply our model to spike train data obtained from prefrontal cortical areas.

Discovering Functional Neuronal Connectivity from Serial Patterns in Spike Train Data

Repeating patterns of precisely timed activity across a group of neurons (called frequent episodes) are indicative of networks in the underlying neural tissue. This letter develops statistical methods to determine functional connectivity among neurons based on nonoverlapping occurrences of episodes. We study the distribution of episode counts and develop a two-phase strategy for identifying functional connections. For the first phase, we develop statistical procedures that are used to screen all two-node episodes and identify possible functional connections (edges). For the second phase, we develop additional statistical procedures to prune the two-node episodes and remove false edges that can be attributed to chains or fan-out structures. The restriction to nonoverlapping occurrences makes the counting of all two-node episodes in phase 1 computationally efficient. The second (pruning) phase is critical since phase 1 can yield a large number of false connections. The scalability of the two-phase approach is examined through simulation. The method is then used to reconstruct the graph structure of observed neuronal networks, first from simulated data and then from recordings of cultured cortical neurons.

Analyzing large-scale spiking neural data with HRLAnalysisâ„¢

The additional capabilities provided by high-performance neural simulation environments and modern computing hardware has allowed for the modeling of increasingly larger spiking neural networks. This is important for exploring more anatomically detailed networks but the corresponding accumulation in data can make analyzing the results of these simulations difficult. This is further compounded by the fact that many existing analysis packages were not developed with large spiking data sets in mind. Presented here is a software suite developed to not only process the increased amount of spike-train data in a reasonable amount of time, but also provide a user friendly Python interface. We describe the design considerations, implementation and features of the HRLAnalysis™ suite. In addition, performance benchmarks demonstrating the speedup of this design compared to a published Python implementation are also presented. The result is a high-performance analysis toolkit that is not only usable and readily extensible, but also straightforward to interface with existing Python modules.

A Unified Approach to Linking Experimental, Statistical and Computational Analysis of Spike Train Data

A fundamental issue in neuroscience is how to identify the multiple biophysical mechanisms through which neurons generate observed patterns of spiking activity. In previous work, we proposed a method for linking observed patterns of spiking activity to specific biophysical mechanisms based on a state space modeling framework and a sequential Monte Carlo, or particle filter, estimation algorithm. We have shown, in simulation, that this approach is able to identify a space of simple biophysical models that were consistent with observed spiking data (and included the model that generated the data), but have yet to demonstrate the application of the method to identify realistic currents from real spike train data. Here, we apply the particle filter to spiking data recorded from rat layer V cortical neurons, and correctly identify the dynamics of an slow, intrinsic current. The underlying intrinsic current is successfully identified in four distinct neurons, even though the cells exhibit two distinct classes of spiking activity: regular spiking and bursting. This approach – linking statistical, computational, and experimental neuroscience – provides an effective technique to constrain detailed biophysical models to specific mechanisms consistent with observed spike train data.

Introduction to neural spike train data for phase-amplitude analysis

Statistical analysis of spike trains is one of the central problems in neural coding, and can be pursued in several ways. One option is model-based, i.e. assume a parametric or semi-parametric model, such as the Poisson model, for spike train data and use it in decoding spike trains. The other option is metric-based, i.e. choose a metric for comparing the numbers and the placements of spikes in different trains, and does not need a model. A prominent idea in the latter approach is to derive metrics that are based on measurements of time-warpings of spike trains needed in the alignments of corresponding spikes. We propose the use of ideas developed in functional data analysis, namely the definition and separation of phase-amplitude components, as a novel tool for analyzing spike trains and decoding underlying neural signals. For concreteness, we introduce a real spike train dataset taken from experimental recordings of the primary motor cortex of a monkey while performing certain arm movements. To facilitate functional data analysis, one needs to smooth the observed discrete spike trains with Gaussian kernels.

Analysis of spike train data: Discussion of results

In this paper we summarize and comment on different techniques presented for analyzing spike train data.

Analysis of spike train data: An application of $k$-mean alignment

We analyze the spike train data by means of the $k$-mean alignment algorithm in a double perspective: data as non periodic and data as periodic. In the first analysis, we show that alignment is not needed to identify paths. Indeed, without allowing for warping, we detect four clusters strongly associated to the four possible paths. In the second analysis, by exploiting the circular nature of data and allowing for shifts, we detect two clusters distinguishing between spike trains presenting higher or lower neuronal activity during the bottom-left/bottom-right movement respectively. In this latter case, the alignment procedure is able to match the four movements across paths.

Analysis of spike train data: Alignment and comparisons using the extended Fisher-Rao metric

We present a metric-based framework for analyzing statistical variability of the neural spike train data that was introduced in an earlier paper on this section [14]. Treating the smoothed spike trains as functional data, we apply the extended Fisher-Rao Riemannian metric, first introduced in Srivastava et al. [9], to perform: (1) pairwise alignment of spike functions, (2) averaging of multiple functions, and (3) alignment of spike functions to the mean. The last item results in separation phase and amplitude components from the functional data. Further, we utilize proper metrics on these components for classification of activities represented by spike trains. This approach is based on the square-root slope function (SRSF) representation of functions that transforms the Fisher-Rao metric into the standard $\mathbb{L}^{2}$ metric and, thus, simplifies computations. We compare our registration results with some current methods and demonstrate an application of our approach in neural decoding to infer motor behaviors.

Analysis of spike train data: A multivariate mixed effects model for phase and amplitude

This discussion of the spike train data applies the approach of Hadjipantelis et al. developed for Linguistic data to jointly model both the phase and amplitude functions via a mixed effects model. The approach shows that care needs to be taken when assessing amplitude and phase functions, particularly if separate analysis is undertaken, as there can be high levels of correlation between the warping and amplitude component functions.

Analysis of spike train data: Classification and Bayesian alignment

We analyze a data set of spike trains obtained under four different experimental conditions. We model the data curves via mixtures of normal densities. The peak locations in the fitted curves are modeled via a non-homogeneous Poisson process and classification of the spike trains into groups may be done based on the estimated spacings between peaks. We employ a Bayesian, MCMC-based registration method to align the fitted curves and summarize the data using meaningful functional statistics and posterior intervals.

Analysis of spike train data: Comparison between the real and the simulated data

This paper compares the real experimental data with the data simulated by Wu and Srivastava, 2011. It turns out that the real data exhibit a different composition of the phase and the amplitude variation from the simulated data, where these two types of variation are separated using the Fisher Rao curve registration. As a result, for the real data the original functions are a better choice of data objects for path discrimination, while for the simulated data the domain warping functions are better.

Estimating synaptic connections from multiple spike trains based on a coupled escape rate model

An approach to study the mechanism of information processing in a neural circuit is to identify synaptic connectivity between neurons. As recent advances in experimental techniques facilitate simultaneous recording of a large population of neurons, estimation of synaptic connectivity from multiple neural spike train data has become a major goal in computational neuroscience. Although it is important to develop an accurate method, it is difficult to validate a method on the basis of experimental data because of the lack of knowledge regarding synaptic connectivity. In this contribution, we have developed a method to estimate synaptic connection on the basis of a coupled escape rate model (CERM), which is an extension of the escape rate model [1]. The strength of the synaptic connection was determined by maximizing the likelihood function. We applied this method as well as model-free methods (transfer entropy and cross-correlation [2]) to synthetic spike train data. The synthetic data were generated by a cortical network model, which consisted of thousands of Hodgkin-Huxley type neurons [3]. These methods are applied to the spike data generated by the cortical network model with different topologies of synaptic connectivity (regular, small-world, and random). It should be noted that the true connectivity is known and the estimation performance can be evaluated by the receiver operating characteristic (ROC) analysis. Our results demonstrated that all methods gave better results for nonrandom (regular and small-world) networks than for random networks. Furthermore, CERM is considered to be the optimal method because it does not involve binning and requires smallest observation length. Therefore, despite high computational costs, CERM is a most suitable method to estimate synaptic connections from multiple spike train data.

SPIKY: a graphical user interface for monitoring spike train synchrony

With the growing availability of multi-unit recordings there is increasing demand for methods which provide the possibility to study similarity patterns of activity across many neurons. Accordingly, a wide variety of approaches to quantify the similarity (or dissimilarity) between two or more spike trains has been suggested. Recently, the ISI- and the SPIKE-distance [1,2] have been proposed as parameter-free and time-scale independent measures of spike train synchrony. The key property of both measures is that they are time-resolved since they rely on instantaneous estimates of spike train dissimilarity. This makes it possible to track changes in instantaneous clustering, i.e., time-localized patterns of (dis)similarity among multiple spike trains. The SPIKE-distance also comes in a causal variant [2] which is defined such that the instantaneous values of dissimilarity are defined from past information only so that time-resolved spike train synchrony can be estimated in real-time. For both the regular and the real-time SPIKE-distance, there are several levels of information reduction [3]. The starting point is the most detailed representation in which one instantaneous value is obtained for each pair of spike trains. This results in a matrix of size 'number of sampled time instants' × 'squared number of spike trains' (i.e. #(tn)N2). By selecting a pair of spike trains one obtains a bivariate dissimilarity profile whereas the selection of a time instant yields an instantaneous matrix of pairwise spike train dissimilarities which can be used to divide the spike trains into instantaneous clusters, i.e., groups of spike trains with low intra-group and high inter-group dissimilarity. Another way to reduce the information is averaging. The spatial average over spike train pairs yields a dissimilarity profile for the selected (sub)population, whereas temporal averaging leads to a bivariate distance matrix for the selected interval or the selected trigger points. Finally, application of the remaining average results in one distance value which describes the overall level of synchrony for a group of spike trains over a given time interval. The Matlab source codes for calculating and visualizing both the ISI- and the SPIKE-distance have been made publicly available and have already been widely used in various contexts. However, the use of these codes is not very intuitive and their application requires some basic knowledge of Matlab. Thus it became desirable to provide a more user-friendly and accessible interface. Here we address this need and present the graphical user interface SPIKY [4,5]. This interactive program facilitates the application of the ISI- and the SPIKE-distance to both simulated and real data. SPIKY includes a spike train generator for testing purposes, as well as masks for selecting the analysis window and the neuronal subpopulation of interest. Once given a set of spike train data, it calculates the desired measure and allows visualization of all the different representations mentioned above (such as measure profiles and pairwise dissimilarity matrices). It even includes the possibility to generate movies which are very useful in order to track the varying patterns of (dis)similarity. Finally, we also have increased the high computation speed even further by transferring the most time-consuming parts of the original Matlab code to Matlab executables (MEX) with the new subroutines written in C.

Identification of sparse neural functional connectivity using penalized likelihood estimation and basis functions

One key problem in computational neuroscience and neural engineering is the identification and modeling of functional connectivity in the brain using spike train data. To reduce model complexity, alleviate overfitting, and thus facilitate model interpretation, sparse representation and estimation of functional connectivity is needed. Sparsities include global sparsity, which captures the sparse connectivities between neurons, and local sparsity, which reflects the active temporal ranges of the input-output dynamical interactions. In this paper, we formulate a generalized functional additive model (GFAM) and develop the associated penalized likelihood estimation methods for such a modeling problem. A GFAM consists of a set of basis functions convolving the input signals, and a link function generating the firing probability of the output neuron from the summation of the convolutions weighted by the sought model coefficients. Model sparsities are achieved by using various penalized likelihood estimations and basis functions. Specifically, we introduce two variations of the GFAM using a global basis (e.g., Laguerre basis) and group LASSO estimation, and a local basis (e.g., B-spline basis) and group bridge estimation, respectively. We further develop an optimization method based on quadratic approximation of the likelihood function for the estimation of these models. Simulation and experimental results show that both group-LASSO-Laguerre and group-bridge-B-spline can capture faithfully the global sparsities, while the latter can replicate accurately and simultaneously both global and local sparsities. The sparse models outperform the full models estimated with the standard maximum likelihood method in out-of-sample predictions.

Robust analysis of semiparametric renewal process models

A rate model is proposed for a modulated renewal process comprising a single long sequence, where the covariate process may not capture the dependencies in the sequence as in standard intensity models. We consider partial likelihood-based inferences under a semiparametric multiplicative rate model, which has been widely studied in the context of independent and identical data. Under an intensity model, gap times in a single long sequence may be used naively in the partial likelihood with variance estimation utilizing the observed information matrix. Under a rate model, the gap times cannot be treated as independent and studying the partial likelihood is much more challenging. We employ a mixing condition in the application of limit theory for stationary sequences to obtain consistency and asymptotic normality. The estimator's variance is quite complicated owing to the unknown gap times dependence structure. We adapt block bootstrapping and cluster variance estimators to the partial likelihood. Simulation studies and an analysis of a semiparametric extension of a popular model for neural spike train data demonstrate the practical utility of the rate approach in comparison with the intensity approach.

Method for stationarity-segmentation of spike train data with application to the Pearson cross-correlation

Correlations among neurons are supposed to play an important role in computation and information coding in the nervous system. Empirically, functional interactions between neurons are most commonly assessed by cross-correlation functions. Recent studies have suggested that pairwise correlations may indeed be sufficient to capture most of the information present in neural interactions. Many applications of correlation functions, however, implicitly tend to assume that the underlying processes are stationary. This assumption will usually fail for real neurons recorded in vivo since their activity during behavioral tasks is heavily influenced by stimulus-, movement-, or cognition-related processes as well as by more general processes like slow oscillations or changes in state of alertness. To address the problem of nonstationarity, we introduce a method for assessing stationarity empirically and then "slicing" spike trains into stationary segments according to the statistical definition of weak-sense stationarity. We examine pairwise Pearson cross-correlations (PCCs) under both stationary and nonstationary conditions and identify another source of covariance that can be differentiated from the covariance of the spike times and emerges as a consequence of residual nonstationarities after the slicing process: the covariance of the firing rates defined on each segment. Based on this, a correction of the PCC is introduced that accounts for the effect of segmentation. We probe these methods both on simulated data sets and on in vivo recordings from the prefrontal cortex of behaving rats. Rather than for removing nonstationarities, the present method may also be used for detecting significant events in spike trains.

Persistent activity and interactions in the brainstem respiratory network during hyperventilatory apnea

Brainstem respiratory neuron interactions during hyperventilatory apnea (HA) and the dissolution and reemergence of the motor pattern for breathing remain poorly understood. We used multi‐array recording, cross‐correlation, and gravity clustering methods to evaluate spike train data from the ventral respiratory column (VRC), pontine respiratory group, and medullary raphé nuclei of decerebrate vagotomized neuromuscularly‐blocked cats. Phrenic activity, blood pressure, end‐tidal CO2, and arterial blood gas measures of O2, CO2 and pH were monitored. During HA: (i) tonic activities persisted in subsets of VRC inspiratory (I) and expiratory (E) neurons; (ii) short‐time scale correlations indicated interactions among neurons that expressed E and non‐respiratory modulated discharge patterns during control conditions, including neurons with functional inhibitory actions on I neurons; (iii) selective stimulation of carotid chemoreceptors by close injection of CO2‐saturated saline transiently evoked rhythmic respiratory network activity and suppressed tonic E neuron activity; (iv) respiratory and slower (0.04 – 0.1 Hz) network oscillations were disrupted. We conclude that VRC tonic E neurons contribute to I drive inhibition during HA (cf. Sears et al. 1982) and document the loss of multiple rhythms during HA. Support: NIH grant NS19814

The potential of corticomuscular and intermuscular coherence for research on human motor control

The potential of corticomuscular and intermuscular coherence for research on human motor control

Estimating summary statistics in the spike-train space

Estimating sample averages and sample variability is important in analyzing neural spike trains data in computational neuroscience. Current approaches have focused on advancing the use of parametric or semiparametric probability models of the underlying stochastic process, where the probabilistic distribution is characterized at each time point with basic statistics such as mean and variance. To directly capture and analyze the average and variability in the observation space of the spike trains, we focus on a data-driven approach where statistics are defined and computed in a function space in which the spike trains are viewed as individual points. Based on the definition of a "Euclidean" metric, a recent paper introduced the notion of the mean of a set of spike trains and developed an efficient algorithm to compute it under some restrictive conditions. Here we extend this study by: (1) developing a novel algorithm for mean computation that is quite general, and (2) introducing a notion of covariance of a set of spike trains. Specifically, we estimate the covariance matrix using the geometry of the warping functions that map the mean spike train to each of the spike trains in the dataset. Results from simulations as well as a neural recording in primate motor cortex indicate that the proposed mean and covariance successfully capture the observed variability in spike trains. In addition, a "Gaussian-type" probability model (defined using the estimated mean and covariance) reasonably characterizes the distribution of the spike trains and achieves a desirable performance in the classification of the spike trains.

Detection of bursts and pauses in spike trains

Midbrain dopaminergic neurons in vivo exhibit a wide range of firing patterns. They normally fire constantly at a low rate, and speed up, firing a phasic burst when reward exceeds prediction, or pause when an expected reward does not occur. Therefore, the detection of bursts and pauses from spike train data is a critical problem when studying the role of phasic dopamine (DA) in reward related learning, and other DA dependent behaviors. However, few statistical methods have been developed that can identify bursts and pauses simultaneously. We propose a new statistical method, the Robust Gaussian Surprise (RGS) method, which performs an exhaustive search of bursts and pauses in spike trains simultaneously. We found that the RGS method is adaptable to various patterns of spike trains recorded in vivo, and is not influenced by baseline firing rate, making it applicable to all in vivo spike trains where baseline firing rates vary over time. We compare the performance of the RGS method to other methods of detecting bursts, such as the Poisson Surprise (PS), Rank Surprise (RS), and Template methods. Analysis of data using the RGS method reveals potential mechanisms underlying how bursts and pauses are controlled in DA neurons.