1.2. Mathematical modeling of stochastic neural activity
1.2.1. Definition of some important measures
The output of a neuron, its spike train, is completely determined by the spike times.
Mathematically, it is convenient to model a spike train as a superposition of Dirac delta functions,
x(t) =
∑
i
δ(t−ti), (1.1)
where{ti}are the spike times. Motivated by the observed variability of neural spiking, it is useful to consider a spike train a stochastic process (a point process). A particular set of spike times is then a realization of that process (out of an ensemble).
Theinstantaneous firing ratecan be written as the first moment of the spike train,
r(t) =hx(t)i, (1.2)
whereh·idenotes ensemble averaging.
1.2. Mathematical modeling of stochastic neural activity
One of the simplest models of neural spiking is theinhomogeneous Poisson process. This process is completely defined by prescribing a time-dependent firing rate r(t) and de-manding that the probability that a spike occurs at a given time is independent of whether or not spikes occur at other times. A way to approximate such a process is to discretize the time axis into bins of length ∆t, drawing a uniformly distributed random number Ri ∼ U(0, 1)for each bin and registering a spike wheneverRi <∆t·r(i·∆t). For∆t →0, this yields an inhomogeneous Poisson Process. A homogeneous Poisson process is ob-tained when the firing rate is constant in time,r(t)≡r0.
A stochastic process is calledstationaryif its moments do not depend on absolute time (but only time differences). Clearly, an inhomogeneous Poisson process that encodes a time-varying signal is not stationary. We can, however, consider the signal a stochas-tic process as well: Usually, we do not know which parstochas-ticular realization of a stimulus a neuron is bound to encounter; we can, at best, make statements about the statistics the stimulus will obey. If the signal process is stationary, then, by considering not only the spike-train ensemble (conditioned on a particular stimulus) but also the ensemble of stimuli, we are again dealing with a stationary process.
A stationary spike train with spike times{ti}can also be thought of as a sequence of interspike intervals(ISIs),
Ti =ti−ti−1, (1.3)
that are drawn from a certain distribution ρ(T). ISIs in a sequence will in general be correlated; if, by contrast, they are independent, then the spike train is arenewal process.
The mean ISI is related to thestationary firing rate r0,
hTi= 1
r0. (1.4)
A useful measure to quantify the irregularity of a spike train is thecoefficient of variation (CV)
CV =
ph∆T2i
hTi , (1.5)
with∆T= T− hTi.
Turning back to the example of a (homogeneous) Poisson process, it is straightforward to show that the ISI distribution is exponential,
ρ(T) =r0e−r0T, (1.6)
and that theCV =1. The Poisson process thus often serves as a reference when assessing the regularity of spiking – spike trains with aCV < 1 (CV > 1 ) are more (less) regular than a Poisson process.
Theauto-correlation functionof a stationary spike train is defined as
Kxx(τ):=hx(t)x(t+τ)i − hx(t)i2. (1.7) This can expressed via the joint probability densityP(t1,t2), whereP(t1,t2)dt2gives the
Chapter 1. Introduction
probability to observe a spike in the time interval(t1,t1+dt)and a spike in(t2,t2+dt) (the marginal probability density P(t1)and the conditional densityP(t1|t2)are defined analogously). One has
where m(τ) is the spike-triggered rate, i.e. the probability that after a reference spike at time t, there is a different spike at time t+τ (this does not need to be the first spike after the reference spike). For a (homogeneous) Poisson process, m(τ) = r0, and thus Kxx(τ) =r0δ(τ).
One can also study second-order statistics in the Fourier domain. For a stationary process, thepower spectrumis the Fourier transform of the correlation function,
Sxx(f) =
Z∞
−∞
dτe2πi fτKxx(τ). (1.9)
This relation is the so-called Wiener-Khinchin theorem (Gardiner, 1985; Risken, 1989). A different definition of the power spectrum, which is more easily calculated in simulations, is
δ(f − f0)Sxx(f) =xe(f)ex∗(f0), (1.10) where the tilde denotes the Fourier transform (defined as in eq. (1.9)), and the asterisk denotes complex conjugation. In simulations we have to use finite time windows; there, we use the Fourier transform
xeT(f) =
Z T
0 dt e2πi f tx(t) (1.11)
and approximate the power spectrum as
Sxx(f) = 1 T
xeT(f)xe∗T(f0). (1.12) For a renewal spike train, the power spectrum is related to the Fourier transform of the
1.2. Mathematical modeling of stochastic neural activity
ISI density (Stratonovich, 1967),
Sxx(f) =r01− |ρe(f)|2
|1−ρe(f)|2. (1.13)
From eq. (1.13), it can be shown that the limit of vanishing frequency of a renewal spike train is given by
limf→0Sxx(f) =r0C2V. (1.14) For spike-train auto-correlation functions which contain noδ-peak except the one atτ= 02, the high-frequency limit of the power spectrum is
flim→∞Sxx(f) =r0. (1.15)
Turning again to the example of a (homogeneous) Poisson spike train, the power spec-trum is
Sxx(f)≡r0 (1.16)
(the process has equal power at all frequencies, it iswhite), which is consistent with the two limits.
The definitions eq. (1.7), eq. (1.9), eq. (1.10), and eq. (1.12) generalize to the case where one considerstwoprocessess(t)andx(t). One has thecross-correlation,
Ksx(τ) =hs(t)x(t+τ)i − hs(t)i hx(t)i, (1.17) and the correspondingcross spectrum,
δ(f− f0)Ssx(f) =es(f)xe∗(f0). (1.18) When the influence of a signals(t)on the spike trainx(t)is weak, one may use linear response theory to calculate the cross spectrum: Conditioned on the signal, the first mo-ment ofx(t)(the firing rate) is time-dependent. Linear response theory assumes that the effect of the signal on this rate is captured by convolution with a kernelK(τ),
r(t) =hx(t)i ≈r0+
Z∞
−∞
dτK(τ)s(t−τ), (1.19)
whereK(τ)is causal, i.e. K(τ < 0) = 0. In Fourier space, the convolution turns into a multiplication, so that, for f 6=0,
er(f) =hxe(f)i=χ(f)es(f) (1.20) where χ(f)is the susceptibility. From eq. (1.18) and eq. (1.20), one finds that in linear
2In Chapter 3, we will encounter a spike-train correlation function that falls outside of this class.
Chapter 1. Introduction
response
Ssx(f)≈χ∗(f)Sss(f), (1.21) whereSss(f)is the signal power spectrum.
A measure that will be widely used in the present thesis is thecoherence function,
Csx(f) = |Ssx(f)|2
Sxx(f)Sss(f). (1.22) The coherence is a measure between zero and one; it can be thought of as a correlation coefficient in frequency space. In the neural context, it quantifies how well certain fre-quency components of the signal can be linearly reconstructed from the spike train. We will come back to this measure when introducing information-theoretic concepts.