Inference of intrinsic spiking irregularity based on the Kullback–Leibler information

Abstract

We have recently established an empirical Bayes method that extracts both the intrinsic irregularity and the time-dependent rate from a spike sequence [Koyama, S., Shinomoto, S., 2005. Empirical Bayes interpretations of random point events. J. Phys. A: Math. Gen. 38, L531–L537]. In the present paper, we examine an alternative method based on the more fundamental principle of minimizing the Kullback–Leibler information from the original distribution of spike sequences to a model distribution. Not only the empirical Bayes method but also the Kullback–Leibler information method exhibits a switch of the most plausible interpretation of the spikes between (I) being derived irregularly from a nearly constant rate, and (II) being derived rather regularly from a significantly fluctuating rate.The model distributions selected by both methods are similar for the same spike sequences derived from a given rate-fluctuating gamma process.

Introduction

Characterizing the neuronal firing is the first step to understand the function of the brain Dayan and Abbott, 2001, Rieke et al., 1997, Tuckwell, 1988. It has been revealed in recent studies that cortical neurons exhibit characteristic firing patterns specific to individual neurons.With a new coefficient termed the local variation of inter-spike intervals $L_{\text{V}}$, neurons in the medial motor cortical areas were found to be classified into two groups according to their specific firing patterns, “quasi-regular” and “likely random” (Shinomoto et al., 2003). Neurons in area TE were also classified according to specific firing patterns, one similar to the “likely random” type in the medial motor cortical areas, and the other exhibiting a “clumpy-bursty” firing pattern unique to TE (Shinomoto et al., 2005). Those neurons exhibiting distinctly different firing patterns were found to be distributed mainly in layers V–VI and layers II–III, respectively. This result suggests a possibility of identifying layer location of the neurons by their firing characteristics.

In this way, cortical neurons were successfully classified by means of the local variation, $L_{\text{V}}$, but not by means of the coefficient of variation, $C_{\text{V}}$. $L_{\text{V}}$ reflects only the local stepwise variability of inter-spike intervals, depends less on the firing rate, and represents the firing characteristics intrinsic and specific to individual neurons, while $C_{\text{V}}$ reflects the global variability of an entire spike sequence and is sensitive to rate fluctuations. It would be due to the fluctuating firing rate of the in vivo spike sequences that the classification by means of $C_{\text{V}}$ was unsatisfactory.

We recently proposed an empirical Bayes method that enables us to grasp not only such intrinsic irregularity, but also the fluctuating firing rate (Koyama and Shinomoto, 2005). In the first phase of the Bayes method, the fluctuating firing rate $\lambda (t)$ is estimated from a sequence of spikes $\{t_i\}=\{t_1\text{,}\dots \text{,}{t}_n\}$. For this purpose, we introduce a rate-fluctuating point process as a model of the spike sequence, $p_{\kappa }(\{t_i\}|\{\lambda (t)\})$, which includes a hyperparameter, $\kappa$, that specifies the intrinsic irregularity of the spike train. In order to avoid the possible over-fitting that takes place in estimating a continuous rate $\lambda (t)$ from a finite number of spikes, we introduce a prior distribution of $\lambda (t)$ as $p_{\beta }(\{\lambda (t)\})$, in which the hyperparameter $\beta$ indicates the tendency to flatness of the time-dependent rate. The posterior distribution of the rate $\lambda (t)$ is obtained using the Bayes theorem.

In the second phase, a set of the hyperparameters $\{\kappa \text{,}\beta \}$ is determined so that the “marginal likelihood function”, or the “evidence”:

$$p_{\kappa \text{,}\beta }(\{t_i\})=\int p_{\kappa }(\{t_i\}|\{\lambda (t)\})p_{\beta }(\{\lambda (t)\})\mathrm{d}\{\lambda (t)\}$$

is maximized for the given data $\{t_i\}$ according to the principle of the empirical Bayes method Carlin and Louis, 2000, MacKay, 1992.

We have demonstrated the effectiveness of this procedure and the suitability of the hyper parameters $\{{\hat{\kappa }}_{\text{EB}}\text{,}{\hat{\beta }}_{\text{EB}}\}$ selected accordingly, in the light of the simulation study (Koyama and Shinomoto, 2005). However, the empirical Bayes inferences could be evaluated as “overconfident”, since the same set of data is used twice; first in the prior to select $\beta$ and then in the posterior to select $\lambda (t)$.

In the present paper, we examine an alternative method based on the more fundamental principle of minimizing the Kullback–Leibler (KL) information from the true distribution of spike sequences to a model distribution. The KL information is a coefficient measuring a non-negative asymmetrical “distance” from one probability distribution to another, and the model distribution with a lower value of the KL information approximates the original probability distribution better (Cover and Thomas, 1991).

For a set of spike trains derived from a time dependent gamma process, the KL information from a true distribution to a model distribution may have multiple minima in the parameter space of the model distribution. This fact implies that there are multiple potentially reasonable parameter sets for the model distribution. We call the alternative choices for a set of parameters the “interpretations”. By increasing the fluctuation of the rate of occurrence, the plausible interpretation switches from (I) irregularly derived from a constant rate, to (II) regularly derived from a fluctuating rate. Second, we compare the distribution functions selected by the KL information method with the one selected by the empirical Bayes method (Koyama and Shinomoto, 2005). We show that the results obtained by both methods are similar not only qualitatively but also quantitatively.

This paper is organized as follows. In Section 2, we propose the method of inference based on the KL information. In Section 3, we apply the KL information method to spike sequences generated by numerical simulation to search for the optimal interpretations of given spike sequences. The results are then compared with those obtained by the empirical Bayes method. Section 4 is devoted to discussions.