Dynamic analyses of information encoding in neural ensembles

Neural spike train decoding algorithms and techniques to compute Shannon mutual information are important methods for analyzing how neural systems represent biological signals. Decoding algorithms are also one of several strategies being used to design controls for brain-machine interfaces. Developing optimal strategies to design decoding algorithms and compute mutual information are therefore important problems in computational neuroscience. We present a general recursive filter decoding algorithm based on a point process model of individual neuron spiking activity and a linear stochastic state-space model of the biological signal. We derive from the algorithm new instantaneous estimates of the entropy, entropy rate, and the mutual information between the signal and the ensemble spiking activity. We assess the accuracy of the algorithm by computing, along with the decoding error, the true coverage probability of the approximate 0.95 confidence regions for the individual signal estimates. We illustrate the new algorithm by reanalyzing the position and ensemble neural spiking activity of CA1 hippocampal neurons from two rats foraging in an open circular environment. We compare the performance of this algorithm with a linear filter constructed by the widely used reverse correlation method. The median decoding error for Animal 1 (2) during 10 minutes of open foraging was 5.9 (5.5) cm, the median entropy was 6.9 (7.0) bits, the median information was 9.4 (9.4) bits, and the true coverage probability for 0.95 confidence regions was 0.67 (0.75) using 34 (32) neurons. These findings improve significantly on our previous results and suggest an integrated approach to dynamically reading neural codes, measuring their properties, and quantifying the accuracy with which encoded information is extracted.     

Nonparametric modeling of neural point processes via stochastic gradient boosting regression

Statistical nonparametric modeling tools that enable the discovery and approximation of functional forms (e.g., tuning functions) relating neural spiking activity to relevant covariates are desirable tools in neuroscience. In this article, we show how stochastic gradient boosting regression can be successfully extended to the modeling of spiking activity data while preserving their point process nature, thus providing a robust nonparametric modeling tool. We formulate stochastic gradient boosting in terms of approximating the conditional intensity function of a point process in discrete time and use the standard likelihood of the process to derive the loss function for the approximation problem. To illustrate the approach, we apply the algorithm to the modeling of primary motor and parietal spiking activity as a function of spiking history and kinematics during a two-dimensional reaching task. Model selection, goodness of fit via the time rescaling theorem, model interpretation via partial dependence plots, ranking of covariates according to their relative importance, and prediction of peri-event time histograms are illustrated and discussed. Additionally, we use the tenfold cross-validated log likelihood of the modeled neural processes (67 cells) to compare the performance of gradient boosting regression to two alternative approaches: standard generalized linear models (GLMs) and Bayesian P-splines with Markov chain Monte Carlo (MCMC) sampling. In our data set, gradient boosting outperformed both Bayesian P-splines (in approximately 90% of the cells) and GLMs (100%). Because of its good performance and computational efficiency, we propose stochastic gradient boosting regression as an off-the-shelf nonparametric tool for initial analyses of large neural data sets (e.g., more than 50 cells; more than 10(5) samples per cell) with corresponding multidimensional covariate spaces (e.g., more than four covariates). In the cases where a functional form might be amenable to a more compact representation, gradient boosting might also lead to the discovery of simpler, parametric models.     

The time-rescaling theorem and its application to neural spike train data analysis.

Measuring agreement between a statistical model and a spike train data series, that is, evaluating goodness of fit, is crucial for establishing the model's validity prior to using it to make inferences about a particular neural system. Assessing goodness-of-fit is a challenging problem for point process neural spike train models, especially for histogram-based models such as perstimulus time histograms (PSTH) and rate functions estimated by spike train smoothing. The time-rescaling theorem is a well-known result in probability theory, which states that any point process with an integrable conditional intensity function may be transformed into a Poisson process with unit rate. We describe how the theorem may be used to develop goodness-of-fit tests for both parametric and histogram-based point process models of neural spike trains. We apply these tests in two examples: a comparison of PSTH, inhomogeneous Poisson, and inhomogeneous Markov interval models of neural spike trains from the supplementary eye field of a macque monkey and a comparison of temporal and spatial smoothers, inhomogeneous Poisson, inhomogeneous gamma, and inhomogeneous inverse gaussian models of rat hippocampal place cell spiking activity. To help make the logic behind the time-rescaling theorem more accessible to researchers in neuroscience, we present a proof using only elementary probability theory arguments. We also show how the theorem may be used to simulate a general point process model of a spike train. Our paradigm makes it possible to compare parametric and histogram-based neural spike train models directly. These results suggest that the time-rescaling theorem can be a valuable tool for neural spike train data analysis.     

Performance analysis of two EM-based measurement bias estimation processes for tracking systems

In target tracking, the measurements collected by sensors can be biased in some real scenarios, e.g., due to systematic error. To accurately estimate the target trajectory, it is essential that the measurement bias be identified in the first place. We investigate the iterative bias estimation process based on the expectation-maximization (EM) algorithm, for cases where sufficiently large numbers of measurements are at hand. With the assistance of extended Kalman filtering and smoothing, we derive two EM estimation processes to estimate the measurement bias which is formulated as a random variable in one state-space model and a constant value in another. More importantly, we theoretically derive the global convergence result of the EM-based measurement bias estimation and reveal the link between the two proposed EM estimation processes in the respective state-space models. It is found that the bias estimate in the second state-space model is more accurate and of less complexity. Furthermore, the EM-based iterative estimation converges faster in the second state-space model than in the first one. As a byproduct, the target trajectory can be simultaneously estimated with the measurement bias, after processing a batch of measurements. These results are confirmed by our simulations.     

Approximate conditional least squares estimation of a nonlinear state-space model via an unscented Kalman filter

The problem of estimating a nonlinear state-space model whose state process is driven by an ordinary differential equation (ODE) or a stochastic differential equation (SDE), with discrete-time data is studied. A new estimation method is proposed based on minimizing the conditional least squares (CLS) with the conditional mean function computed approximately via the unscented Kalman filter (UKF). Conditions are derived for the UKF–CLS estimator to preserve the limiting properties of the exact CLS estimator, namely, consistency and asymptotic normality, under the framework of infill asymptotics, i.e. sampling is increasingly dense over a fixed domain. The efficacy of the proposed method is demonstrated by simulation and a real application.     

A latent variable regression model for capture-recapture data

Capture-recapture methods are used to estimate the prevalence of diseases in the field of epidemiology. The information used for estimation purposes are available from multiple lists, whereby giving rise to the problems of list dependence and heterogeneity. In this paper, modelling is focused on the heterogeneity part. We present a new binomial latent class model which takes into account both the observed and unobserved heterogeneity within capture-recapture data. We adopt the conditional likelihood approach and perform estimation via the EM algorithm. We also derive the mathematical expressions for the computation of the standard error of the unknown population size. An application to data on diabetes patients in a town in northern Italy is discussed.     

Local influence analysis of nonlinear structural equation models with nonignorable missing outcomes from reproductive dispersion models

Nonlinear structural equation models with nonignorable missing outcomes from reproductive dispersion models are proposed to identify the relationship between manifest variables and latent variables in modern educational, medical, social and psychological studies. The nonignorable missing mechanism is specified by a logistic regression model. An EM algorithm is developed to obtain the maximum likelihood estimates of the structural parameters and parameters in the logistic regression model. Assessment of local influence is investigated in nonlinear structural equation models with nonignorable missing outcomes from reproductive dispersion models on the basis of the conditional expectation of the complete-data log-likelihood function. Some local influence diagnostics are obtained via observations of missing data and latent variables that are generated by the Gibbs sampler and Metropolis-Hastings algorithm on the basis of the conformal normal curvature. A simulation study and a real example are used to illustrate the application of the proposed methodologies.     

Maximum likelihood estimation of an extended latent Markov model for clustered binary panel data

Computational aspects concerning a model for clustered binary panel data are analyzed. The model is based on the representation of the behavior of a subject (individual panel member) in a given cluster by means of a latent process. This latent process is decomposed into a cluster-specific component and an individual-specific component. The first component follows a first-order Markov chain, whereas the second is time-invariant and is represented by a discrete random variable. An algorithm for computing the joint distribution of the response variables is introduced. The algorithm may be used even in the presence of a large number of subjects in the same cluster. An Expectation-Maximization (EM) scheme for the maximum likelihood estimation of the model is also described together with the estimation of the Fisher information matrix on the basis of the numerical derivative of the score vector. The estimate of this matrix is used to obtain standard errors for the parameter estimates and to check the identifiability of the model and the convergence of the EM algorithm. The approach is illustrated by means of an application to a data set concerning Italian employees’ illness benefits.     

Dynamic analysis of neural encoding by point process adaptive filtering

Neural receptive fields are dynamic in that with experience, neurons change their spiking responses to relevant stimuli. To understand how neural systems adapt their representations of biological information, analyses of receptive field plasticity from experimental measurements are crucial. Adaptive signal processing, the well-established engineering discipline for characterizing the temporal evolution of system parameters, suggests a framework for studying the plasticity of receptive fields. We use the Bayes' rule Chapman-Kolmogorov paradigm with a linear state equation and point process observation models to derive adaptive filters appropriate for estimation from neural spike trains. We derive point process filter analogues of the Kalman filter, recursive least squares, and steepest-descent algorithms and describe the properties of these new filters. We illustrate our algorithms in two simulated data examples. The first is a study of slow and rapid evolution of spatial receptive fields in hippocampal neurons. The second is an adaptive decoding study in which a signal is decoded from ensemble neural spiking activity as the receptive fields of the neurons in the ensemble evolve. Our results provide a paradigm for adaptive estimation for point process observations and suggest a practical approach for constructing filtering algorithms to track neural receptive field dynamics on a millisecond timescale.     

Analyzing functional connectivity using a network likelihood model of ensemble neural spiking activity

Analyzing the dependencies between spike trains is an important step in understanding how neurons work in concert to represent biological signals. Usually this is done for pairs of neurons at a time using correlation-based techniques. Chornoboy, Schramm, and Karr (1988) proposed maximum likelihood methods for the simultaneous analysis of multiple pair-wise interactions among an ensemble of neurons. One of these methods is an iterative, continuous-time estimation algorithm for a network likelihood model formulated in terms of multiplicative conditional intensity functions. We devised a discrete-time version of this algorithm that includes a new, efficient computational strategy, a principled method to compute starting values, and a principled stopping criterion. In an analysis of simulated neural spike trains from ensembles of interacting neurons, the algorithm recovered the correct connectivity matrices and interaction parameters. In the analysis of spike trains from an ensemble of rat hippocampal place cells, the algorithm identified a connectivity matrix and interaction parameters consistent with the pattern of conjoined firing predicted by the overlap of the neurons' spatial receptive fields. These results suggest that the network likelihood model can be an efficient tool for the analysis of ensemble spiking activity.     

Detection of hidden structures in nonstationary spike trains

We propose an algorithm for simultaneously estimating state transitions among neural states and nonstationary firing rates using a switching state-space model (SSSM). This algorithm enables us to detect state transitions on the basis of not only discontinuous changes in mean firing rates but also discontinuous changes in the temporal profiles of firing rates (e.g., temporal correlation). We construct estimation and learning algorithms for a nongaussian SSSM, whose nongaussian property is caused by binary spike events. Local variational methods can transform the binary observation process into a quadratic form. The transformed observation process enables us to construct a variational Bayes algorithm that can determine the number of neural states based on automatic relevance determination. Additionally, our algorithm can estimate model parameters from single-trial data using a priori knowledge about state transitions and firing rates. Synthetic data analysis reveals that our algorithm has higher performance for estimating nonstationary firing rates than previous methods. The analysis also confirms that our algorithm can detect state transitions on the basis of discontinuous changes in temporal correlation, which are transitions that previous hidden Markov models could not detect. We also analyze neural data recorded from the medial temporal area. The statistically detected neural states probably coincide with transient and sustained states that have been detected heuristically. Estimated parameters suggest that our algorithm detects the state transitions on the basis of discontinuous changes in the temporal correlation of firing rates. These results suggest that our algorithm is advantageous in real-data analysis.     

基于混合滤波最大期望算法的高速列车建模

针对高速列车非线性单质点模型的特殊结构及含有隐含变量问题, 提出一种基于混合滤波的最大期望辨识方法. 借助递阶辨识理论, 将高铁列车状态空间模型分解为线性子系统模型和非线性子系统模型. 进而, 分别利用卡尔曼滤波和粒子滤波对速度和位移状态进行联合估计. 最后, 使用最大期望方法辨识高铁列车子系统模型参数, 解决了隐含变量辨识问题. 和传统方法相比, 本文所提出方法计算量小, 且具有较高的辨识精度. 仿真对比实验结果验证了该方法的有效性.     

Computing confidence intervals for point process models

Characterizing neural spiking activity as a function of intrinsic and extrinsic factors is important in neuroscience. Point process models are valuable for capturing such information; however, the process of fully applying these models is not always obvious. A complete model application has four broad steps: specification of the model, estimation of model parameters given observed data, verification of the model using goodness of fit, and characterization of the model using confidence bounds. Of these steps, only the first three have been applied widely in the literature, suggesting the need to dedicate a discussion to how the time-rescaling theorem, in combination with parametric bootstrap sampling, can be generally used to compute confidence bounds of point process models. In our first example, we use a generalized linear model of spiking propensity to demonstrate that confidence bounds derived from bootstrap simulations are consistent with those computed from closed-form analytic solutions. In our second example, we consider an adaptive point process model of hippocampal place field plasticity for which no analytical confidence bounds can be derived. We demonstrate how to simulate bootstrap samples from adaptive point process models, how to use these samples to generate confidence bounds, and how to statistically test the hypothesis that neural representations at two time points are significantly different. These examples have been designed as useful guides for performing scientific inference based on point process models.     

Pattern filtering for detection of neural activity, with examples from HVc activity during sleep in zebra finches

The detection of patterned spiking activity is important in the study of neural coding. A pattern filtering approach is developed for pattern detection under the framework of point processes, which offers flexibility in combining temporal details and firing rates. The detection combines multiple steps of filtering in a coarse-to-fine manner. Under some conditional Poisson assumptions on the spiking activity, each filtering step is equivalent to classifying by likelihood ratios all the data segments as targets or as background sequences. Unlike previous studies, where global surrogate data were used to evaluate the statistical significance of the detected patterns, a localized p-test procedure is developed, which better accounts for firing modulation and nonstationarity in spiking activity. Common temporal structures of patterned activity are learned using an entropy-based alignment procedure, without relying on metrics or pair-wise alignment. Applications of pattern filtering to single, presumptive interneurons recorded in the nucleus HVc of zebra finch are illustrated. These demonstrate a match between the auditory-evoked response to playback of the individual bird's own song and spontaneous activity during sleep. Small temporal compression or expansion, or both, is required for optimal matching of spontaneous patterns to stimulus-evoked activity.     

Fast ML Estimation for the Mixture of Factor Analyzers via an ECM Algorithm

In this brief, we propose a fast expectation conditional maximization (ECM) algorithm for maximum-likelihood (ML) estimation of mixtures of factor analyzers (MFA). Unlike the existing expectation–maximization (EM) algorithms such as the EM in Ghahramani and Hinton, 1996, and the alternating ECM (AECM) in McLachlan and Peel, 2003, where the missing data contains component–indicator vectors as well as latent factors, the missing data in our ECM consists of component–indicator vectors only. The novelty of our algorithm is that closed-form expressions in all conditional maximization (CM) steps are obtained explicitly, instead of resorting to numerical optimization methods. As revealed by experiments, the convergence of our ECM is substantially faster than EM and AECM regardless of whether assessed by central processing unit (CPU) time or number of iterations.      

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.     

Dual information representation with stable firing rates and chaotic spatiotemporal spike patterns in a neural network model.

Although various means of information representation in the cortex have been considered, the fundamental mechanism for such representation is not well understood. The relation between neural network dynamics and properties of information representation needs to be examined. We examined spatial pattern properties of mean firing rates and spatiotemporal spikes in an interconnected spiking neural network model. We found that whereas the spatiotemporal spike patterns are chaotic and unstable, the spatial patterns of mean firing rates (SPMFR) are steady and stable, reflecting the internal structure of synaptic weights. Interestingly, the chaotic instability contributes to fast stabilization of the SPMFR. Findings suggest that there are two types of network dynamics behind neuronal spiking: internally-driven dynamics and externally driven dynamics. When the internally driven dynamics dominate, spikes are relatively more chaotic and independent of external inputs; the SPMFR are steady and stable. When the externally driven dynamics dominate, the spiking patterns are relatively more dependent on the spatiotemporal structure of external inputs. These emergent properties of information representation imply that the brain may adopt a dual coding system. Recent experimental data suggest that internally driven and externally driven dynamics coexist and work together in the cortex.     

Efficient Markov chain Monte Carlo methods for decoding neural spike trains

Stimulus reconstruction or decoding methods provide an important tool for understanding how sensory and motor information is represented in neural activity. We discuss Bayesian decoding methods based on an encoding generalized linear model (GLM) that accurately describes how stimuli are transformed into the spike trains of a group of neurons. The form of the GLM likelihood ensures that the posterior distribution over the stimuli that caused an observed set of spike trains is log concave so long as the prior is. This allows the maximum a posteriori (MAP) stimulus estimate to be obtained using efficient optimization algorithms. Unfortunately, the MAP estimate can have a relatively large average error when the posterior is highly nongaussian. Here we compare several Markov chain Monte Carlo (MCMC) algorithms that allow for the calculation of general Bayesian estimators involving posterior expectations (conditional on model parameters). An efficient version of the hybrid Monte Carlo (HMC) algorithm was significantly superior to other MCMC methods for gaussian priors. When the prior distribution has sharp edges and corners, on the other hand, the "hit-and-run" algorithm performed better than other MCMC methods. Using these algorithms, we show that for this latter class of priors, the posterior mean estimate can have a considerably lower average error than MAP, whereas for gaussian priors, the two estimators have roughly equal efficiency. We also address the application of MCMC methods for extracting nonmarginal properties of the posterior distribution. For example, by using MCMC to calculate the mutual information between the stimulus and response, we verify the validity of a computationally efficient Laplace approximation to this quantity for gaussian priors in a wide range of model parameters; this makes direct model-based computation of the mutual information tractable even in the case of large observed neural populations, where methods based on binning the spike train fail. Finally, we consider the effect of uncertainty in the GLM parameters on the posterior estimators.     

Nonlinear system identification under missing observations: The case of unknown model structure

This article presents an algorithm for identification of nonlinear state-space models when the “true” model structure of a process is unknown. In order to estimate the parameters in a state-space model, one needs to know the model structure and have an estimate of states. An approximation of the model structure is obtained using radial basis functions centered around a maximum a posteriori estimate of the state trajectory. A particle filter approximation of smoothed states is then used in conjunction with expectation maximization algorithm for estimating the parameters. The proposed approach is extended to handle missing observations and illustrated through a real application.     

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.