A sequential Monte Carlo approach to estimate biophysical neural models from spikes

A biophysical and statistical modeling paradigm for connecting neural physiology and function

To understand single neuron computation, it is necessary to know how specific physiological parameters affect neural spiking patterns that emerge in response to specific stimuli. Here we present a computational pipeline combining biophysical and statistical models that provides a link between variation in functional ion channel expression and changes in single neuron stimulus encoding. More specifically, we create a mapping from biophysical model parameters to stimulus encoding statistical model parameters. Biophysical models provide mechanistic insight, whereas statistical models can identify associations between spiking patterns and the stimuli they encode. We used public biophysical models of two morphologically and functionally distinct projection neuron cell types: mitral cells (MCs) of the main olfactory bulb, and layer V cortical pyramidal cells (PCs). We first simulated sequences of action potentials according to certain stimuli while scaling individual ion channel conductances. We then fitted point process generalized linear models (PP-GLMs), and we constructed a mapping between the parameters in the two types of models. This framework lets us detect effects on stimulus encoding of changing an ion channel conductance. The computational pipeline combines models across scales and can be applied as a screen of channels, in any cell type of interest, to identify ways that channel properties influence single neuron computation.

Estimating Distributions of Parameters in Nonlinear State Space Models with Replica Exchange Particle Marginal Metropolis–Hastings Method

Extracting latent nonlinear dynamics from observed time-series data is important for understanding a dynamic system against the background of the observed data. A state space model is a probabilistic graphical model for time-series data, which describes the probabilistic dependence between latent variables at subsequent times and between latent variables and observations. Since, in many situations, the values of the parameters in the state space model are unknown, estimating the parameters from observations is an important task. The particle marginal Metropolis–Hastings (PMMH) method is a method for estimating the marginal posterior distribution of parameters obtained by marginalization over the distribution of latent variables in the state space model. Although, in principle, we can estimate the marginal posterior distribution of parameters by iterating this method infinitely, the estimated result depends on the initial values for a finite number of times in practice. In this paper, we propose a replica exchange particle marginal Metropolis–Hastings (REPMMH) method as a method to improve this problem by combining the PMMH method with the replica exchange method. By using the proposed method, we simultaneously realize a global search at a high temperature and a local fine search at a low temperature. We evaluate the proposed method using simulated data obtained from the Izhikevich neuron model and Lévy-driven stochastic volatility model, and we show that the proposed REPMMH method improves the problem of the initial value dependence in the PMMH method, and realizes efficient sampling of parameters in the state space models compared with existing methods.

Parameter Estimation in Multiple Dynamic Synaptic Coupling Model Using Bayesian Point Process State-Space Modeling Framework

It is of great interest to characterize the spiking activity of individual neurons in a cell ensemble. Many different mechanisms, such as synaptic coupling and the spiking activity of itself and its neighbors, drive a cell's firing properties. Though this is a widely studied modeling problem, there is still room to develop modeling solutions by simplifications embedded in previous models. The first shortcut is that synaptic coupling mechanisms in previous models do not replicate the complex dynamics of the synaptic response. The second is that the number of synaptic connections in these models is an order of magnitude smaller than in an actual neuron. In this research, we push this barrier by incorporating a more accurate model of the synapse and propose a system identification solution that can scale to a network incorporating hundreds of synaptic connections. Although a neuron has hundreds of synaptic connections, only a subset of these connections significantly contributes to its spiking activity. As a result, we assume the synaptic connections are sparse, and to characterize these dynamics, we propose a Bayesian point-process state-space model that lets us incorporate the sparsity of synaptic connections within the regularization technique into our framework. We develop an extended expectation-maximization. algorithm to estimate the free parameters of the proposed model and demonstrate the application of this methodology to the problem of estimating the parameters of many dynamic synaptic connections. We then go through a simulation example consisting of the dynamic synapses across a range of parameter values and show that the model parameters can be estimated using our method. We also show the application of the proposed algorithm in the intracellular data that contains 96 presynaptic connections and assess the estimation accuracy of our method using a combination of goodness-of-fit measures.

Capturing Multiple Timescales of Adaptation to Second-Order Statistics With Generalized Linear Models: Gain Scaling and Fractional Differentiation

Single neurons can dynamically change the gain of their spiking responses to take into account shifts in stimulus variance. Moreover, gain adaptation can occur across multiple timescales. Here, we examine the ability of a simple statistical model of spike trains, the generalized linear model (GLM), to account for these adaptive effects. The GLM describes spiking as a Poisson process whose rate depends on a linear combination of the stimulus and recent spike history. The GLM successfully replicates gain scaling observed in Hodgkin-Huxley simulations of cortical neurons that occurs when the ratio of spike-generating potassium and sodium conductances approaches one. Gain scaling in the GLM depends on the length and shape of the spike history filter. Additionally, the GLM captures adaptation that occurs over multiple timescales as a fractional derivative of the stimulus envelope, which has been observed in neurons that include long timescale afterhyperpolarization conductances. Fractional differentiation in GLMs requires long spike history that span several seconds. Together, these results demonstrate that the GLM provides a tractable statistical approach for examining single-neuron adaptive computations in response to changes in stimulus variance.

Continuous Prediction of Cognitive State Using A Marked-Point Process Modeling Framework

Behavioral outcomes in many cognitive tasks are often recorded in a trial structure at discrete times. To adapt to this structure, neural encoder and decoder models have been built to take into account the trial organization to characterize the connection between brain dynamics and behavior, e.g. through latent dynamical models. The challenge of these models is that they are limited to discrete trial times while neural data is continuous. Here, we propose a marked-point process framework to characterize multivariate behavioral outcomes recorded during a trial-structured cognitive task, to build an estimation of cognitive state at a fine time resolution. We propose a state-space marked-point process modeling framework to characterize the relationship between observed behavior and underlying dynamical cognitive processes. We define the framework for a class of behavioral readouts by a response time and a discrete mark signifying an observed binary decision, and develop the state estimation and system identification steps. We define the filter and smoother for the marked-point process observation and develop an EM algorithm to estimate the model's free parameters. We demonstrate this modeling approach in a behavioral readout captured while participants perform an emotional conflict resolution task (ECR) and show that we can estimate underlying cognitive processes at a fine temporal resolution beyond the trial by trial approach.

Inferring synaptic inputs from spikes with a conductance-based neural encoding model

Descriptive statistical models of neural responses generally aim to characterize the mapping from stimuli to spike responses while ignoring biophysical details of the encoding process. Here, we introduce an alternative approach, the conductance-based encoding model (CBEM), which describes a mapping from stimuli to excitatory and inhibitory synaptic conductances governing the dynamics of sub-threshold membrane potential. Remarkably, we show that the CBEM can be fit to extracellular spike train data and then used to predict excitatory and inhibitory synaptic currents. We validate these predictions with intracellular recordings from macaque retinal ganglion cells. Moreover, we offer a novel quasi-biophysical interpretation of the Poisson generalized linear model (GLM) as a special case of the CBEM in which excitation and inhibition are perfectly balanced. This work forges a new link between statistical and biophysical models of neural encoding and sheds new light on the biophysical variables that underlie spiking in the early visual pathway.

Parameter Estimation in Synaptic Coupling Model Using a Point Process Modeling Framework

Biophysical models are widely used to characterize temporal dynamics of the brain networks on different topological and spatial scales. In parallel, the state-space modeling framework with point process observations has been successfully applied in characterizing spiking activity of neuronal ensembles in response to different dynamical covariates. Parameter estimation in biophysical models is generally done heuristically, which hampers their applicability and interpretability. Heuristic parameter estimation becomes an intractable problem when the number of model parameters grows. Here, we propose an algorithm for estimating biophysical model parameters using point-process models and a state-space framework. The framework provides methods for parameter estimation as well as model validation. We demonstrate the application of this methodology to the problem of estimating the parameters of a dynamic synapse model. We generate simulation data for the dynamic synapse across a range of parameters values and assess the estimation accuracy of our method using a combination of goodness-of-fit measures. The proposed methodology can be applied broadly to parameter estimation problems across a broad range of biophysical models, including Hodgkin-Huxley models and network models.

Real-time particle filtering and smoothing algorithms for detecting abrupt changes in neural ensemble spike activity

Sequential change-point detection from time series data is a common problem in many neuroscience applications, such as seizure detection, anomaly detection, and pain detection. In our previous work (Chen Z, Zhang Q, Tong AP, Manders TR, Wang J. J Neural Eng 14: 036023, 2017), we developed a latent state-space model, known as the Poisson linear dynamical system, for detecting abrupt changes in neuronal ensemble spike activity. In online brain-machine interface (BMI) applications, a recursive filtering algorithm is used to track the changes in the latent variable. However, previous methods have been restricted to Gaussian dynamical noise and have used Gaussian approximation for the Poisson likelihood. To improve the detection speed, we introduce non-Gaussian dynamical noise for modeling a stochastic jump process in the latent state space. To efficiently estimate the state posterior that accommodates non-Gaussian noise and non-Gaussian likelihood, we propose particle filtering and smoothing algorithms for the change-point detection problem. To speed up the computation, we implement the proposed particle filtering algorithms using advanced graphics processing unit computing technology. We validate our algorithms, using both computer simulations and experimental data for acute pain detection. Finally, we discuss several important practical issues in the context of real-time closed-loop BMI applications. NEW & NOTEWORTHY Sequential change-point detection is an important problem in closed-loop neuroscience experiments. This study proposes novel sequential Monte Carlo methods to quickly detect the onset and offset of a stochastic jump process that drives the population spike activity. This new approach is robust with respect to spike sorting noise and varying levels of signal-to-noise ratio. The GPU implementation of the computational algorithm allows for parallel processing in real time.

New insights into olivo-cerebellar circuits for learning from a small training sample

Artificial intelligence such as deep neural networks exhibited remarkable performance in simulated video games and 'Go'. In contrast, most humanoid robots in the DARPA Robotics Challenge fell down to ground. The dramatic contrast in performance is mainly due to differences in the amount of training data, which is huge and small, respectively. Animals are not allowed with millions of the failed trials, which lead to injury and death. Humans fall only several thousand times before they balance and walk. We hypothesize that a unique closed-loop neural circuit formed by the Purkinje cells, the cerebellar deep nucleus and the inferior olive in and around the cerebellum and the highest density of gap junctions, which regulate synchronous activities of the inferior olive nucleus, are computational machinery for learning from a small sample. We discuss recent experimental and computational advances associated with this hypothesis.

Segmental Bayesian estimation of gap-junctional and inhibitory conductance of inferior olive neurons from spike trains with complicated dynamics

The inverse problem for estimating model parameters from brain spike data is an ill-posed problem because of a huge mismatch in the system complexity between the model and the brain as well as its non-stationary dynamics, and needs a stochastic approach that finds the most likely solution among many possible solutions. In the present study, we developed a segmental Bayesian method to estimate the two parameters of interest, the gap-junctional (gc ) and inhibitory conductance (gi ) from inferior olive spike data. Feature vectors were estimated for the spike data in a segment-wise fashion to compensate for the non-stationary firing dynamics. Hierarchical Bayesian estimation was conducted to estimate the gc and gi for every spike segment using a forward model constructed in the principal component analysis (PCA) space of the feature vectors, and to merge the segmental estimates into single estimates for every neuron. The segmental Bayesian estimation gave smaller fitting errors than the conventional Bayesian inference, which finds the estimates once across the entire spike data, or the minimum error method, which directly finds the closest match in the PCA space. The segmental Bayesian inference has the potential to overcome the problem of non-stationary dynamics and resolve the ill-posedness of the inverse problem because of the mismatch between the model and the brain under the constraints based, and it is a useful tool to evaluate parameters of interest for neuroscience from experimental spike train data.

Clusterless Decoding of Position from Multiunit Activity Using a Marked Point Process Filter

Point process filters have been applied successfully to decode neural signals and track neural dynamics. Traditionally these methods assume that multiunit spiking activity has already been correctly spike-sorted. As a result, these methods are not appropriate for situations where sorting cannot be performed with high precision, such as real-time decoding for brain-computer interfaces. Because the unsupervised spike-sorting problem remains unsolved, we took an alternative approach that takes advantage of recent insights into clusterless decoding. Here we present a new point process decoding algorithm that does not require multiunit signals to be sorted into individual units. We use the theory of marked point processes to construct a function that characterizes the relationship between a covariate of interest (in this case, the location of a rat on a track) and features of the spike waveforms. In our example, we use tetrode recordings, and the marks represent a four-dimensional vector of the maximum amplitudes of the spike waveform on each of the four electrodes. In general, the marks may represent any features of the spike waveform. We then use Bayes's rule to estimate spatial location from hippocampal neural activity. We validate our approach with a simulation study and experimental data recorded in the hippocampus of a rat moving through a linear environment. Our decoding algorithm accurately reconstructs the rat's position from unsorted multiunit spiking activity. We then compare the quality of our decoding algorithm to that of a traditional spike-sorting and decoding algorithm. Our analyses show that the proposed decoding algorithm performs equivalent to or better than algorithms based on sorted single-unit activity. These results provide a path toward accurate real-time decoding of spiking patterns that could be used to carry out content-specific manipulations of population activity in hippocampus or elsewhere in the brain.

Identifying and tracking simulated synaptic inputs from neuronal firing: insights from in vitro experiments

Accurately describing synaptic interactions between neurons and how interactions change over time are key challenges for systems neuroscience. Although intracellular electrophysiology is a powerful tool for studying synaptic integration and plasticity, it is limited by the small number of neurons that can be recorded simultaneously in vitro and by the technical difficulty of intracellular recording in vivo. One way around these difficulties may be to use large-scale extracellular recording of spike trains and apply statistical methods to model and infer functional connections between neurons. These techniques have the potential to reveal large-scale connectivity structure based on the spike timing alone. However, the interpretation of functional connectivity is often approximate, since only a small fraction of presynaptic inputs are typically observed. Here we use in vitro current injection in layer 2/3 pyramidal neurons to validate methods for inferring functional connectivity in a setting where input to the neuron is controlled. In experiments with partially-defined input, we inject a single simulated input with known amplitude on a background of fluctuating noise. In a fully-defined input paradigm, we then control the synaptic weights and timing of many simulated presynaptic neurons. By analyzing the firing of neurons in response to these artificial inputs, we ask 1) How does functional connectivity inferred from spikes relate to simulated synaptic input? and 2) What are the limitations of connectivity inference? We find that individual current-based synaptic inputs are detectable over a broad range of amplitudes and conditions. Detectability depends on input amplitude and output firing rate, and excitatory inputs are detected more readily than inhibitory. Moreover, as we model increasing numbers of presynaptic inputs, we are able to estimate connection strengths more accurately and detect the presence of connections more quickly. These results illustrate the possibilities and outline the limits of inferring synaptic input from spikes.

Likelihood Methods for Point Processes with Refractoriness

Likelihood-based encoding models founded on point processes have received significant attention in the literature because of their ability to reveal the information encoded by spiking neural populations. We propose an approximation to the likelihood of a point-process model of neurons that holds under assumptions about the continuous time process that are physiologically reasonable for neural spike trains: the presence of a refractory period, the predictability of the conditional intensity function, and its integrability. These are properties that apply to a large class of point processes arising in applications other than neuroscience. The proposed approach has several advantages over conventional ones. In particular, one can use standard fitting procedures for generalized linear models based on iteratively reweighted least squares while improving the accuracy of the approximation to the likelihood and reducing bias in the estimation of the parameters of the underlying continuous-time model. As a result, the proposed approach can use a larger bin size to achieve the same accuracy as conventional approaches would with a smaller bin size. This is particularly important when analyzing neural data with high mean and instantaneous firing rates. We demonstrate these claims on simulated and real neural spiking activity. By allowing a substantive increase in the required bin size, our algorithm has the potential to lower the barrier to the use of point-process methods in an increasing number of applications.

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.

Control strategies for underactuated neural ensembles driven by optogenetic stimulation

Motivated by experiments employing optogenetic stimulation of cortical regions, we consider spike control strategies for ensembles of uncoupled integrate and fire neurons with a common conductance input. We construct strategies for control of spike patterns, that is, multineuron trains of action potentials, up to some maximal spike rate determined by the neural biophysics. We emphasize a constructive role for parameter heterogeneity, and find a simple rule for controllability in pairs of neurons. In particular, we determine parameters for which common drive is not limited to inducing synchronous spiking. For large ensembles, we determine how the number of controllable neurons varies with the number of observed (recorded) neurons, and what collateral spiking occurs in the full ensemble during control of the subensemble. While complete control of spiking in every neuron is not possible with a single input, we find that a degree of subensemble control is made possible by exploiting dynamical heterogeneity. As most available technologies for neural stimulation are underactuated, in the sense that the number of target neurons far exceeds the number of independent channels of stimulation, these results suggest partial control strategies that may be important in the development of sensory neuroprosthetics and other neurocontrol applications.