format: text

text:
# Dynamical mean-field theory of noisy spiking neuron ensembles: Application to the Hodgkin-Huxley model

## I INTRODUCTION

It is well known that a small cluster of cortex may contain thousands of similar neurons. Each neuron which receives spikes from hundreds of other neurons, generates spikes propagating along the axon towards synapses exciting neurons in the next stage. Dynamics of an individual neuron with voltage-dependent ionic channels can be described by Hodgkin-Huxley-type (HH) model , or by reduced, simplified neuron models such as integrate-and-fire (IF), FitzHugh-Nagumo (FN)  and Hindmarsh-Rose (HR) models . Although the response of a single neuron in vitro is rather accurate, that in vivo is not reliable . This is due to noisy environment in living brains, where various kinds of noises are reported to be ubiquitous (for a review see ). In recent years, the population of neuron ensembles has been recognized to play important roles in the information transmission (pooling effect) -. Then it is necessary for us to theoretically investigate high-dimensional, stochastic differential equations (DEs) describing the large-scale noisy neuron ensemble. In order to make our discussion concrete, let us consider ensembles consisting of $`N`$-unit neurons, each of which is described by $`K`$-dimensional coupled DEs: for example, $`K`$ =1, 2, 3 and 4 for IF, FN, HR and HH neuron models, respectively. Dynamics of such neuron ensembles, expressed by $`KN`$-dimensional stochastic DEs, has been so far investigated with the use of the two approaches: (i) direct simulations and (ii) analytical methods such as Fokker-Planck equation (FPE) and the moment method. Simulations have been made for large-scale networks mostly consisting of IF neurons. Since the CPU time to simulate networks by conventional methods is proportional to $`N^2`$, it is rather difficult to simulate realistic neuron clusters in spite of recent computer development. In FPE dynamics of neuron ensembles is described by the population activity. Although FPE is the powerful method formally applicable to the case of arbitrary $`K`$ and $`N`$ , actual calculations have been made mostly for $`N=\mathrm{}`$ ensembles of a $`K=1`$ model with the use of the mean-field and/or diffusion approximations . Similar population density approaches have been recently developed for a large-scale neuronal clusters . The moment method initiated by Rodriguez and Tuckwell (RT) has been applied to single FN  and HH neurons . When the moment method is applied to a single neuron model with $`K`$ variables, $`K`$-dimensional stochastic DEs are replaced by $`(1/2)K(K+3)`$-dimensional deterministic DEs. When the moment method is applied to $`N`$-unit neuron ensembles under consideration, $`KN`$-dimensional stochastic DEs are replaced by $`N_{eq}`$-dimensional deterministic DEs where $`N_{eq}=(1/2)KN(KN+3)`$ . For example, in the case of $`K=2`$ (FN model), the number of equations is $`N_{eq}=`$ 230, 20 300 and 2 003 000 for $`N=`$ 10, 100 and 1000, respectively. In the case of $`K=4`$ (HH model), we get $`N_{eq}=`$ 860, 80 600 and 8 006 000 for $`N=`$ 10, 100 and 1000, respectively. These figures are too large for us to make simulations for realistic neuron clusters. In their subsequent paper , RT transplanted the result of the moment method for HH neuron ensembles to FPE-type equation which has not been solved yet.

In a previous study (Ref. is hereafter referred to as I), the present author proposed a semi-analytical dynamical mean-field approximation (DMA), in which equations of motions for means, variances and covariances of local and global variables were derived for $`N`$-unit FN neuron ensemble. The original $`2N`$-dimensional stochastic DEs are replaced by eight-dimensional deterministic DEs: $`N_{eq}=8`$ is much smaller than corresponding figures in the moment method mentioned above. DMA calculations in I on the spiking-time precision and the synchronization in FN neuron ensembles are in good agreement with direct simulations. The feasibility of DMA has been demonstrated in I.

The purpose of the present paper is two folds. The first purpose is to extend DMA of I to general neuron ensembles subject to white noises described by $`KN`$-dimensional stochastic DEs, which will be replaced by $`K(K+2)`$-dimensional deterministic DEs. The second purpose of the present paper is to apply the generalized DMA to an ensemble of HH neurons, which is more realistic than FN neuron model previously studied in I. Since Hodgkin and Huxley proposed the HH model in 1952 , much studies have been intensively made on properties of the HH model. Responses of a single, pairs and ensembles HH neurons mostly to direct and sinusoidal currents have been investigated. In recent years, responses of HH neurons to spike-train inputs have been studied -. The stochastic resonance (SR) of HH neurons for sinusoidal and spike inputs with various kinds of added noises has been investigated -. These studies have shown that noise can play a constructive role in signal transmission against our conventional wisdom. In most studies on SR, however, noises added to ensemble neurons are considered to be independent of each other. Quite recently effects of spatially correlated noises on SR have been investigated , which shows that although common noises work to enhance the synchronization in neuron ensembles, they are not effective for SR, in contrast to independent noises. We will adopt in this study, spatially correlated white noises in order to clarify respective effects of common and independent noises on the response of ensemble neurons.

The paper is organized as follows: In Sec. II, we extend a DMA theory to general neuron ensembles described by $`KN`$ stochastic DEs. Our DMA theory is applied to HH neuron ensembles in Sec. III. Some numerical results on HH neuron ensembles are presented in Sec. IV. Conclusions and discussions are given in Sec. V.

## II DMA for A General Neuron Ensemble

### A Equation of motions

We assume an ensemble of $`N`$-unit neurons ($`N2`$), each of which is described by $`K`$-dimensional non-linear differential equations (DEs). Dynamics of a given neuron ensemble is expressed by

$`{\displaystyle \frac{dv_i}{dt}}`$ $`=`$ $`F^{(1)}(\{u_{qi}\})+\left({\displaystyle \frac{w}{N1}}\right){\displaystyle \underset{j(i)}{}}G(v_j(t))+K^{(e)}(t)+\xi _i(t),`$ (1)
$`{\displaystyle \frac{du_{pi}}{dt}}`$ $`=`$ $`F^{(p)}(\{u_{qi}\}),\text{(}p=2\text{ to }K\text{)}`$ (2)
where $`v_i=u_{pi}`$ with $`p=1`$ denotes the membrane potential of a neuron $`i`$ ($`=1`$ to $`N`$), $`u_{pi}`$ with $`p=2`$ to $`K`$ stands for auxiliary variables and $`F^{(p)}`$ is functions of ($`\{u_{qi}\}`$). The synaptic-coupling strength $`w`$ is assumed to be constant, $`G(v)=1/[1+\mathrm{exp}[(v\theta )/ϵ]]`$ is the sigmoid function with the threshold $`\theta `$ and the width $`ϵ`$ , and $`K^{(e)}`$ stands for an applied external input whose explicit form will be given later \[Eq. (56)\]. The last term of Eq. (1) expresses the spatially correlated white noises, $`\xi _i(t)`$, given by

$`<\xi _i(t)>`$ $`=`$ $`0,`$ (3)
$`<\xi _i(t)\xi _j(t^{})>`$ $`=`$ $`[\beta _0^2\delta _{ij}+\beta _1^2(1\delta _{ij})]\delta (tt^{})`$ (4)
$`=`$ $`(\beta _C^2+\beta _I^2\delta _{ij})\delta (tt^{}),`$ (5)
where $`\beta _C=\beta _1`$ and $`\beta _I=\sqrt{\beta _0^2\beta _1^2}`$ denote the magnitudes of common and independent noises, respectively, and the bracket $`<>`$ expresses the stochastic average ; the case of $`\beta _1=0`$ ($`\beta _1=\beta _0`$) stands for independent (common) noises only.

In order to derive DEs in DMA theory, we first define the global variables for the ensemble by

$`U_p(t)`$ $`=`$ $`(1/N){\displaystyle \underset{i}{}}u_{pi}(t),`$ (6)
and their averages by

$`\mu _p(t)`$ $`=`$ $`\mu _{u_p}(t)=<U_p(t)>.`$ (7)
Deviations from these averages of local variables are given by

$`\delta u_{pi}(t)`$ $`=`$ $`u_{pi}(t)\mu _{u_p}(t),`$ (8)
and those of global variables given by

$`\delta U_p(t)`$ $`=`$ $`U_p(t)\mu _{u_p}(t).`$ (9)
Next we define the variances and co-variances between local variables given by (argument $`t`$ is neglected hereafter)

$`\gamma _{p,q}`$ $`=`$ $`\gamma _{u_p,u_q}={\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}<\delta u_{pi}\delta u_{qi}>,`$ (10)
and those between global variables given by

$`\rho _{p,q}`$ $`=`$ $`\rho _{u_p,u_q}=<\delta U_p\delta U_q>,`$ (11)
It is noted that $`\gamma _{u_p,u_q}`$ expresses fluctuations in local variables while $`\rho _{u_p,u_q}`$ those in global variables.

We assume that the noise intensity is weak and that the distribution function $`p(𝐳)`$ for $`KN`$-dimensional random variables of $`𝐳`$ = ($`\{u_{pi}\}`$) is given by the Gaussian distribution concentrated near the mean point of $`𝝁`$=$`(\{\mu _{u_p}\})`$ . Numerical simulations have shown that for weak noises, the distribution of $`v(t)`$ of the membrane potential of a single HH neuron nearly obeys the Gaussian distribution, although for strong noises, the distribution of $`v(t)`$ deviates from the Gaussian, taking a bimodal form . Similar behavior of the membrane-potential distribution has been reported also in a FN neuron model . By using Eq. (7), we express Eqs. (1) and (2) in a Taylor expansion of $`\delta u_{pi}`$ up to the fourth-order terms. The average yields DEs for the means of $`d\mu _{u_p}/dt`$ \[Eq. (16)\]. DEs of variances and covariances may be obtained by using the equations of motions of $`\delta u_{pi}`$. For example, DE for $`d\gamma _{u_p,u_q}/dt`$ is given by

$$\frac{d\gamma _{u_p,u_q}}{dt}=\frac{1}{N}\underset{i}{}<\left(\frac{\delta u_{pi}}{t}\right)\delta u_{qi}+\delta u_{pi}\left(\frac{\delta u_{qi}}{t}\right)>,$$
(12)
with

$`{\displaystyle \frac{\delta u_{pi}}{t}}`$ $`=`$ $`{\displaystyle \underset{q}{}}F_{u_q}^{(u_p)}\delta u_{qi}+{\displaystyle \frac{1}{2}}{\displaystyle \underset{q}{}}{\displaystyle \underset{r}{}}F_{u_qu_r}^{(u_p)}(\delta u_{qi}\delta u_{ri}\gamma _{u_q,u_r})+{\displaystyle \frac{1}{6}}{\displaystyle \underset{q}{}}{\displaystyle \underset{r}{}}{\displaystyle \underset{s}{}}F_{u_qu_ru_s}^{(u_p)}\delta u_{qi}\delta u_{ri}\delta u_{ti}`$ (13)
$`+`$ $`\delta _{p1}\{\xi _i+{\displaystyle \frac{1}{N1}}{\displaystyle \underset{k(i)}{}}[G_{u_1}\delta u_{1k}+{\displaystyle \frac{1}{2}}G_{u_1u_1}(\delta u_{1k}^2\gamma _{1,1})+{\displaystyle \frac{1}{6}}G_{u_1u_1u_1}\delta u_{1k}^3]+K^{(e)}\},`$ (14)
where $`q,r`$ and $`s`$ run from 1 to $`K`$, $`F^{(u_p)}=F^{(p)}`$, $`F_{u_q}^{(u_p)}=F^{(p)}/u_q`$, $`F_{u_qu_r}^{(u_p)}=^2F^{(p)}/u_qu_r`$ and $`F_{u_qu_ru_s}^{(u_p)}=^3F^{(p)}/u_qu_ru_s`$ evaluated at the means of ({$`\mu _{u_p}`$}), and similar derivatives for $`G`$. In the process of calculations of means, variances and covariances, we have taken into account the fourth-order moment contributions with the use of the Gaussian approximation, as given by

$`<\delta u_{pi}\delta u_{qi}\delta u_{ri}\delta u_{si}>`$ $``$ $`<\delta u_{pi}\delta u_{qi}><\delta u_{ri}\delta u_{si}>+<\delta u_{pi}\delta u_{ri}><\delta u_{qi}\delta u_{si}>`$ (15)
$`+`$ $`<\delta u_{pi}\delta u_{si}><\delta u_{qi}\delta u_{ri}>,`$ (16)
$`{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}<\delta u_{pi}\delta u_{qi}\delta u_{ri}\delta u_{si}>`$ $``$ $`\gamma _{u_p,u_q}\gamma _{u_r,u_s}+\gamma _{u_p,u_r}\gamma _{u_q,u_s}+\gamma _{u_p,u_s}\gamma _{u_q,u_r},`$ (17)
$`{\displaystyle \frac{1}{N^2}}{\displaystyle \underset{i}{}}{\displaystyle \underset{j}{}}<\delta u_{pi}\delta u_{qj}\delta u_{rj}\delta u_{sj}>`$ $``$ $`\rho _{u_p,u_q}\gamma _{u_r,u_s}+\rho _{u_p,u_r}\gamma _{u_q,u_s}+\rho _{u_p,u_s}\gamma _{u_q,u_r},`$ (18)
The importance of including the fourth-order term has been pointed out by Tanabe and Pakdaman  in the improved moment method for a noisy FN neuron.

After some manipulations, we get DEs for means, variances and covariances given by (details being given in Appendix A of I):

$`{\displaystyle \frac{d\mu _{u_p}}{dt}}`$ $`=`$ $`F^{(u_p)}+{\displaystyle \frac{1}{2}}{\displaystyle \underset{q}{}}{\displaystyle \underset{r}{}}F_{u_q,u_r}^{(u_p)}\gamma _{u_q,u_r}+\delta _{p1}[wU_0+K^{(e)}],`$ (19)
$`{\displaystyle \frac{d\gamma _{u_p,u_q}}{dt}}`$ $`=`$ $`{\displaystyle \underset{r}{}}[F_{u_r}^{(u_p)}\gamma _{u_q,u_r}+F_r^{(u_q)}\gamma _{u_p,u_r}]+\beta _0^2\delta _{p1}\delta _{q1}+wU_1[\delta _{p1}\zeta _{u_q,u_1}+\delta _{q1}\zeta _{u_p,u_1}]`$ (22)
$`+{\displaystyle \frac{1}{6}}{\displaystyle \underset{r}{}}{\displaystyle \underset{s}{}}{\displaystyle \underset{w}{}}[F_{u_ru_su_w}^{(u_p)}(\gamma _{u_q,u_r}\gamma _{u_s,u_w}+\gamma _{u_q,u_s}\gamma _{u_r,u_w}+\gamma _{u_q,u_w}\gamma _{u_r,u_s})`$
$`+F_{u_ru_su_w}^{(u_q)}(\gamma _{u_p,u_r}\gamma _{u_s,u_w}+\gamma _{u_p,u_s}\gamma _{u_r,u_w}+\gamma _{u_p,u_w}\gamma _{u_r,u_s})],`$
$`{\displaystyle \frac{d\rho _{u_p,u_q}}{dt}}`$ $`=`$ $`{\displaystyle \underset{r}{}}[F_{u_r}^{(u_p)}\rho _{u_q,u_r}+F_{u_r}^{(u_q)}\rho _{u_p,u_r}]+[{\displaystyle \frac{1}{N}}\beta _0^2+\left(1{\displaystyle \frac{1}{N}}\right)\beta _1^2]\delta _{p1}\delta _{q1}`$ (26)
$`+wU_1[\delta _{p1}\rho _{u_q,u_1}+\delta _{q1}\rho _{u_p,u_1}]`$
$`+{\displaystyle \frac{1}{6}}{\displaystyle \underset{r}{}}{\displaystyle \underset{s}{}}{\displaystyle \underset{t}{}}[F_{u_ru_su_t}^{(u_p)}(\rho _{u_q,u_r}\gamma _{u_s,u_t}+\rho _{u_q,u_s}\gamma _{u_r,u_t}+\rho _{u_q,u_t}\gamma _{u_r,u_s})`$
$`+F_{u_ru_su_t}^{(u_q)}(\rho _{u_p,u_r}\gamma _{u_s,u_t}+\rho _{u_p,u_s}\gamma _{u_r,u_t}+\rho _{u_p,u_t}\gamma _{u_r,u_s})],`$
with

$`\zeta _{u_p,u_q}`$ $`=`$ $`\left({\displaystyle \frac{1}{N1}}\right)(N\rho _{u_p,u_q}\gamma _{u_p,u_q}),`$ (27)
$`U_0`$ $`=`$ $`{\displaystyle \frac{1}{N}}{\displaystyle \underset{j}{}}<G(v_j)>=G+{\displaystyle \frac{1}{2}}G_{vv}\gamma _{v,v},`$ (28)
$`U_1`$ $`=`$ $`G_v+{\displaystyle \frac{1}{2}}G_{vvv}\gamma _{v,v},`$ (29)
where $`U_0`$ expresses output spikes of the ensemble, $`v_j=u_{1j}`$, and arguments of $`r`$, $`s`$ and $`w`$ in the sums run from 1 to $`K`$. The original $`KN`$-dimensional stochastic DEs are transformed to $`N_{eq}`$-dimensional deterministic DEs where $`N_{eq}=K+K(K+1)=K(K+2)`$.

### B Property of DMA

In previous Sec. IIA, DMA has been derived with the use of equations of motions for moments. It is, however, possible to alternatively derive DMA from the conventional moment method with a reduction in numbers of variables, as was shown in I for FN neuron ensembles. In Appendix A, we present a derivation of DMA from the moment method for a general neuron ensemble under consideration.

We should note that the noise contribution is $`\beta _0^2`$ in Eq. (17) while that is $`[(1/N)\beta _0^2+(11/N)\beta _1^2]`$ in Eq.(18). When model parameters of $`\beta _0`$, $`\beta _1`$, $`w`$ and $`N`$ are varied, the ratio of $`\rho _{v,v}/\gamma _{v,v}`$ changes. In particular, in the case of $`w=0`$, we get

$`{\displaystyle \frac{\rho _{v,v}}{\gamma _{v,v}}}`$ $`=`$ $`{\displaystyle \frac{1}{N}}+\left(1{\displaystyle \frac{1}{N}}\right)\left({\displaystyle \frac{\beta _1}{\beta _0}}\right)^2,`$ (30)
$`=`$ $`{\displaystyle \frac{1}{N}},\text{for }\beta _1=0`$ (31)
$`=`$ $`1.\text{for }\beta _1=\beta _0`$ (32)
Equation (23) agrees with the central-limit theorem for independent noises while Eq. (24) expresses the result for common noises. On the other hand, in the opposite limit of $`w\mathrm{}`$, we get $`\rho _{v,v}/\gamma _{v,v}1`$. The change in the ratio of $`\rho _{v,v}/\gamma _{v,v}`$ reflects on the the firing time distributions and the degree of synchronization in neuron ensembles, as will be discussed in the followings.

Firing Time Distributions

The $`n`$th firing time of a given neuron $`i`$ in the ensemble is defined as the time when the membrane potential $`v_i(t)`$ crosses the threshold $`\theta `$ from below:

$$t_{oin}=\{tv_i(t)=\theta ;\dot{v_i}>0\}.$$
(33)
The distribution of firing times of $`t_{oin}`$ of a given neuron $`i`$ is given by

$$Z_{\mathrm{}}(t)\varphi \left(\frac{tt_o^{}}{\delta t_o\mathrm{}}\right)\frac{d}{dt}\left(\frac{\mu _v}{\sigma _{\mathrm{}}}\right)\mathrm{\Theta }(\dot{\mu _v}),$$
(34)
with the normal distribution function given by

$$\varphi (x)=\frac{1}{\sqrt{2\pi }}\mathrm{exp}\left(\frac{x^2}{2}\right),$$
(35)
and

$$\delta t_o\mathrm{}=\frac{\sigma _{\mathrm{}}}{\dot{\mu _v}},$$
(36)
where $`\sigma _{\mathrm{}}=\sqrt{\gamma _{v,v}}`$ and $`\dot{\mu _v}=d\mu _v/dt`$ evaluated at $`t=t_o^{}`$ where $`\mu _v(t_o^{})=\theta `$. In the limit of vanishing $`\beta `$, Eq. (26) reduces to

$$Z_{\mathrm{}}(t)=\delta (tt_o^{}).$$
(37)
Similarly we may define the $`m`$th firing time relevant to the global variable $`V(t)=(1/N)_iv_i(t)`$ as

$$t_{gm}=\{tV(t)=\theta ;\dot{V}(t)>0\}.$$
(38)
The distribution of firing times of $`t_{gm}`$ is given by

$$Z_g(t)=\varphi \left(\frac{tt_o^{}}{\delta t_{og}}\right)\frac{d}{dt}\left(\frac{\mu _1}{\sigma _g}\right)\mathrm{\Theta }(\dot{\mu _v}),$$
(39)
with

$$\delta t_{og}=\frac{\sigma _g}{\dot{\mu _v}}$$
(40)
where $`\sigma _g=\sqrt{\rho _{v,v}}`$. In particular, in the case of no couplings, we get

$`{\displaystyle \frac{\delta t_{og}}{\delta t_o\mathrm{}}}`$ $`=`$ $`\sqrt{{\displaystyle \frac{1}{N}}+(1{\displaystyle \frac{1}{N}})\left({\displaystyle \frac{\beta _1}{\beta _0}}\right)^2}.\text{(for }w=0\text{)}`$ (41)
Synchronous Response

The synchronization ratio is defined by

$$S(t)=\frac{(\rho _{v,v}/\gamma _{v,v}1/N)}{(11/N)}=\frac{\zeta _{v,v}}{\gamma _{v,v}},$$
(42)
with

$$\zeta _{v,v}=\left(\frac{1}{N1}\right)(N\rho _{v,v}\gamma _{v,v})=\frac{1}{N(N1)}\underset{i}{}\underset{j(i)}{}<\delta v_i\delta v_j>,$$
(43)
expressing the averaged covariance for the variable of ($`\{\delta v_i\}`$). $`S(t)`$ changes as the model parameters of $`\beta _0`$, $`\beta _1`$, $`w`$ and $`N`$ are varied. It is easy to see from Eqs. (23) and (24) that $`S=0`$ (the asynchronous state) for $`w=0`$ and $`\beta _1\beta _0`$, while $`S=1`$ (the completely synchronous state) for $`w\beta _0^2`$ or $`\beta _1=\beta _0`$. In particular, for $`w=0`$, we get

$$S(t)=\left(\frac{\beta _1}{\beta _0}\right)^2,\text{for }w=0$$
(44)
which implies that the synchronization is induced by common noises.

## III DMA for HH Neuron Ensembles

### A Equation of Motions

For the HH neuron model ($`K=4`$), $`F^{(p)}`$ in Eq. (1) is given by

$`F^{(1)}=F^{(v)}(v_i,m_i,h_i,n_i)`$ $`=`$ $`{\displaystyle \frac{1}{C}}[g_{\mathrm{Na}}m_i^3h_i(v_iv_{\mathrm{Na}})+g_\mathrm{K}n_i^4(v_iv_\mathrm{K})+g_\mathrm{L}(v_iv_\mathrm{L})],`$ (45)
$`F^{(p)}=F^{(u_p)}(v_i,u_{pi})`$ $`=`$ $`[a_{u_p}(v_i)+b_{u_p}(v_i)]u_{pi}+a_{u_p}(v_i),\text{(}p=2\text{ to 4)}`$ (46)
In Eqs. (37) and (38) $`u_{1i}=v_i`$ expresses the membrane potential of a neuron $`i`$, and $`u_{2i}=m_i`$, $`u_{3i}=h_i`$ and $`u_{4i}=n_i`$ denote gate variables of Na and K channels for which $`a_{u_p}(v)`$ and $`b_{u_p}(v)`$ ($`p=2`$ to 4) are given by

$$a_m(v)=\frac{0.1(v+40)}{[1e^{(v+40)/10}]},$$
(47)
$$b_m(v)=4e^{(v+65)/18},$$
(48)
$$a_h(v)=0.07e^{(v+65)/20},$$
(49)
$$b_h(v)=\frac{1}{[1+e^{(v+35)/10}]}.$$
(50)
$$a_n(v)=\frac{0.01(v+55)}{[1e^{(v+55)/10}]},$$
(51)
$$b_n(v)=0.125e^{(v+65)/80}.$$
(52)
In Eq. (37), the reversal potentials of Na, K channels and leakage are $`v_{\mathrm{Na}}=50`$ mV, $`v_\mathrm{K}=77`$ mV and $`v_\mathrm{L}=54.5`$ mV: the maximum values of corresponding conductances are $`g_{\mathrm{Na}}=120`$ $`\mathrm{mS}/cm^2`$, $`g_\mathrm{K}=36`$ $`\mathrm{mS}/cm^2`$ and $`g_\mathrm{L}=0.3`$ $`\mathrm{mS}/cm^2`$: the capacitance of the membrane is $`C=1`$ $`\mu \mathrm{F}/cm^2`$. From functional forms for $`F^{(v)}`$ and $`F^{(u_p)}`$ given by Eqs. (37)-(44), we get $`F_{v,v}=0,F_{u_q}^{(u_p)}=F_{u_p}^{(u_p)}\delta _{pq},F_{v,u_q}^{(u_p)}=F_{v,u_p}^{(u_p)}\delta _{pq}`$ and $`F_{u_p,u_q}^{(u_p)}=0`$. Numbers of non-vanishing third-order derivatives are six for $`F^{(v)}`$ \[$`F_{vmm}^{(v)}`$, $`F_{vmh}^{(v)}`$, $`F_{vnn}^{(v)}`$, $`F_{mmm}^{(v)}`$, $`F_{nnn}^{(v)}`$ and $`F_{mmh}^{(v)}`$\] and two for each $`F^{(u_p)}`$ $`(p=2`$ to 4) \[$`F_{vvv}^{(u_p)}`$ and $`F_{vvu_p}^{(u_p)}`$\].

After some manipulations with Eqs. (16)-(18), we get DEs for means, variances and covariances given by ($`p,q=2`$ to 4)

$`{\displaystyle \frac{d\mu _v}{dt}}`$ $`=`$ $`F^{(v)}+{\displaystyle \frac{1}{2}}{\displaystyle \underset{p=2}{\overset{4}{}}}{\displaystyle \underset{q=2}{\overset{4}{}}}F_{u_pu_q}^{(v)}\gamma _{u_p,u_q}+{\displaystyle \underset{p=2}{\overset{4}{}}}F_{vu_p}^{(v)}\gamma _{v,u_p}+wU_0+K^{(e)},`$ (53)
$`{\displaystyle \frac{d\mu _{u_p}}{dt}}`$ $`=`$ $`F^{(u_p)}+{\displaystyle \frac{1}{2}}F_{v,v}^{(u_p)}\gamma _{v,v}+F_{v,u_p}^{(u_p)}\gamma _{v,u_p},`$ (54)
$`{\displaystyle \frac{d\gamma _{v,v}}{dt}}`$ $`=`$ $`2[F_v^{(v)}\gamma _{v,v}+{\displaystyle \underset{p=2}{\overset{4}{}}}F_{u_p}^{(v)}\gamma _{v,u_p}]+\beta _0^2+2wU_1\zeta _{v,v}+X_{v,v},`$ (55)
$`{\displaystyle \frac{d\gamma _{v,u_p}}{dt}}`$ $`=`$ $`(F_v^{(v)}+F_{u_p}^{(u_p)})\gamma _{v,u_p}+{\displaystyle \underset{q=2}{\overset{4}{}}}F_{u_q}^{(v)}\gamma _{u_q,u_p}+F_v^{(u_p)}\gamma _{v,v}+w\zeta _{v,u_p}+X_{v,u_p},`$ (56)
$`{\displaystyle \frac{d\gamma _{u_p,u_q}}{dt}}`$ $`=`$ $`(F_{u_p}^{(u_p)}+F_{u_q}^{(u_q)})\gamma _{u_p,u_q}+F_v^{(u_p)}\gamma _{v,u_q}+F_v^{(u_q)}\gamma _{v,u_p}+X_{u_p,u_q},`$ (57)
$`{\displaystyle \frac{d\rho _{v,v}}{dt}}`$ $`=`$ $`2[F_v^{(v)}\rho _{v,v}+{\displaystyle \underset{p=2}{\overset{4}{}}}F_{u_p}^{(v)}\rho _{v,u_p}]+[{\displaystyle \frac{1}{N}}\beta _0^2+(1{\displaystyle \frac{1}{N}})\beta _1^2]+2wU_1\rho _{v,v}+Y_{v,v},`$ (58)
$`{\displaystyle \frac{d\rho _{v,u_p}}{dt}}`$ $`=`$ $`(F_v^{(v)}+F_{u_p}^{(u_p)})\rho _{v,u_p}+{\displaystyle \underset{q=2}{\overset{4}{}}}F_{u_q}^{(v)}\rho _{u_q,u_p}+F_v^{(u_p)}\rho _{v,v}+wU_1\rho _{v,u_p}+Y_{v,u_p},`$ (59)
$`{\displaystyle \frac{d\rho _{u_p,u_q}}{dt}}`$ $`=`$ $`(F_{u_p}^{(u_p)}+F_{u_q}^{(u_q)})\rho _{u_p,u_q}+F_v^{(u_p)}\rho _{v,u_q}+F_v^{(u_q)}\rho _{v,u_p}+Y_{u_p,u_q},`$ (60)
with

$`\zeta _{u_p,u_q}`$ $`=`$ $`\left({\displaystyle \frac{1}{N1}}\right)(N\rho _{u_p,u_q}\gamma _{u_p,u_q}),`$ (61)
$`U_0`$ $`=`$ $`{\displaystyle \frac{1}{N}}{\displaystyle \underset{j}{}}<G(v_j)>=G+{\displaystyle \frac{1}{2}}G_{vv}\gamma _{v,v},`$ (62)
$`U_1`$ $`=`$ $`G_v+{\displaystyle \frac{1}{2}}G_{vvv}\gamma _{v,v},`$ (63)
where $`F^{(v)}`$, $`F_v^{(v)}=F^{(v)}/v`$ et al. evaluated at means of $`(\mu _v,\mu _m,\mu _h,\mu _n)`$. In Eqs. (45)-(52), $`X_{v,v}`$ and $`Y_{v,v}`$ et al. denote the contributions from the fourth-order terms, whose explicit expressions are given by Eqs. (B1)-(B6) in Appendix B becuase they are rather lengthy. Although calculations of the fourth-order terms are rather tedious, they play important roles in stabilizing DEs. This is numerically demonstrated in Appendix B for the case of $`N=1`$.

The original $`4N`$-dimensional stochastic DEs given by Eqs. (37) and (38) are transformed to 24-dimensional deterministic DEs given by Eqs. (45)-(52) with Eqs. (B1)-(B6): four means, ($`\mu _v`$, $`\mu _m`$, $`\mu _h`$, $`\mu _n`$), ten moments for local variables ($`\gamma _{v,v}`$, $`\gamma _{m,m}`$, $`\gamma _{h,h}`$, $`\gamma _{n,n}`$, $`\gamma _{v,m}`$, $`\gamma _{v,h}`$, $`\gamma _{v,n}`$, $`\gamma _{m,h}`$, $`\gamma _{h,n}`$, $`\gamma _{m,n}`$), and ten moments for global variables ($`\rho _{v,v}`$, $`\rho _{m,m}`$, $`\rho _{h,h}`$, $`\rho _{n,n}`$, $`\rho _{v,m}`$, $`\rho _{v,h}`$, $`\rho _{v,n}`$, $`\rho _{m,h}`$, $`\rho _{h,n}`$, $`\rho _{m,n}`$).

In this subsection, DMA for HH model has been obtained by the method of equations of motions of means, variances and covariances of local and global variables. We may, however, derive it from the moment method, as mentioned before. In Appendix C, DEs in the moment method are presented for HH model.

We expect that our DMA equations given by Eqs.(45)-(52) and (B1)-(B6) may show much variety depending on model parameters such as the strength of white noise ($`\beta _0`$, $`\beta _1`$), couplings ($`w`$) and the ensemble size ($`N`$). In the next Sec. IV, we will present some numerical DMA calculations, which are compared with simulation results. DMA equations have been solved by the fourth-order Runge-Kutta method with a time step of 0.01 ms for the initial conditions of $`\mu _v=65.0`$, $`\mu _m=0.0528`$, $`\mu _h=0.597`$, $`\mu _n=0.317`$, and $`\gamma _{u_p,u_q}=\rho _{u_p,u_q}=0`$ ($`u_p`$, $`u_q`$=$`v`$, $`m`$, $`h`$ and $`n`$). Direct simulations have been performed by solving $`4N`$-dimensional DEs given by Eqs.(37) and (38) by using also the fourth-order Runge-Kutta method with a time step of 0.01. Simulation results are the average of 100 trials otherwise noticed.

## IV Calculated Results of HH neuron ensembles

### A Firing time distribution

In the present study, we pay our attention to the response of the HH neuron ensembles to a single spike input applied to all neurons in the ensemble, given by

$$K^{(e)}(t)=\left(\frac{I_i}{C}\right)\alpha (tt_i),$$
(64)
with the alpha function:

$$\alpha (t)=\left(\frac{t}{\tau _s}\right)e^{(1t/\tau _s)}\mathrm{\Theta }(t),$$
(65)
where $`\mathrm{\Theta }(x)=1`$ for $`x0`$ and 0 otherwise, $`I_i`$ stands for the magnitude of an input spike, $`C`$ the membrane capacitance \[Eq. (37)\], $`t_i`$ the input time of a spike, and $`\tau _s`$ (=1 ms) the time constant of synapses. We get the critical magnitude of $`I_{ic}`$ = 3.62 $`\mu `$A/$`\mathrm{cm}^2`$, below which firings of neuron cannot take place without noises ($`\beta _0=\beta _1=0`$). We have adopted the value of $`I_i`$ =5 $`\mu `$A/$`\mathrm{cm}^2`$ for a study of the response to a supra-threshold input. We express the coupling constant $`w`$ by $`w=J/C`$ with $`J`$ in units of $`\mu `$A/$`\mathrm{cm}^2`$. The time, voltage, current and noise intensity are hereafter expressed in units of ms, mV, $`\mu `$A/$`\mathrm{cm}^2`$ and V/s, respectively, though they are sometimes omitted for a simplicity of our explanation. We have adopted parameters of $`\theta =0`$ mV and $`ϵ=10`$ mV in the sigmoid function $`G(v)`$ such that output $`U_0`$ is similar to the result given by the alpha function \[see Fig. 1(a)\]. Adopted parameter values of $`\beta _0`$, $`\beta _1`$, $`J`$ and $`N`$ will be explained shortly.

Figures 1(a), 1(b) and 1(c) show the time courses of $`\mu _v`$, $`\sigma _{\mathrm{}}(=\sqrt{\gamma _{v,v}})`$ and $`\sigma _g(=\sqrt{\rho _{v,v}})`$, respectively, when a single spike is applied at $`t=100`$ ms. Solid and dashed curves express the results of DMA and direct simulations, respectively, which are calculated with parameters of $`\beta _0=0.1`$, $`\beta _1=0`$, $`J=0`$ and $`N=100`$. States of neurons in an ensemble when an input spike is injected at $`t=100`$ ms, are randomized because noises have been already added since $`t=0`$. We note that $`\mu _v`$ obtained by DMA is in very good agreement with that obtained by simulations as shown in Fig. 1(a), where an external input of $`K^{(e)}(t)`$ and an output of $`U_0(t)`$ are also plotted. Figures 1(b) and 1(c) show that $`\sigma _{\mathrm{}}`$ and $`\sigma _g`$ calculated by DMA are again in good agreement with those of simulations. We note that the relation given by Eq.(22): $`\sigma _g=\sigma _{\mathrm{}}/N`$ valid for $`w=J/C=0`$, is supported by our calculations.

Figure 2(a) shows $`Z_{\mathrm{}}`$, the firing probability of local variable, which is calculated for $`\beta _0=0.1`$, $`\beta _1=0`$, $`J=0`$ and $`N=100`$. Firings occur at $`t103.6`$ ms with a delay of about 3.6 ms. Fluctuations of firing times of local variable, $`\delta t_o\mathrm{}`$, are 0.066 ms in DMA while it is 0.069 ms in simulations which is the root-mean-square (RMS) value of firing times defined by Eq. (25). In contrast, Fig. 2(b) shows $`Z_g`$, the firing probability of global variables. Fluctuations of firing times of global variable $`\delta t_{og}`$ are 0.0066 ms in DMA and it is 0.0083 ms in simulations, respectively. We note that $`\delta t_{og}`$ is much smaller than $`\delta t_o\mathrm{}`$ \[Eq. (33)\].

Noise-strength dependence

When the noise strength is increased, the distribution of membrane potentials is widen and fluctuations of firing times are increased, as was discussed in Sec. IIB. Filled squares in Fig. 3(a) show the $`\beta _0`$ dependence of $`\delta t_o\mathrm{}`$ obtained by DMA theory with $`\beta _1=0`$, $`J=0`$ and $`N=100`$, while open squares express the RMS value of firing times obtained by simulations. The agreement between the two methods is in fairly good for $`\beta _0<0.1`$ but becomes worse for $`\beta >0.1`$. In contrast, filled circles in Fig. 3(a) show the $`\beta _0`$ dependence of $`\delta t_{og}`$ relevant to the global variable obtained by DMA theory and open circles stand for RMS values of firing times in simulations. We note that $`\delta t_{og}`$ is much smaller than $`\delta t_o\mathrm{}`$ because $`\delta t_{og}=\delta t_o\mathrm{}/\sqrt{N}`$ \[Eq. (33)\].

As $`\beta _1`$ is increased for a fixed $`\beta _0`$, the contribution from common noises increases while that from independent noises decreases ($`\beta _C=\beta _1`$, $`\beta _I=\sqrt{\beta _0^2\beta _1^2}`$). The $`\beta _1`$ dependence of firing-time fluctuations is shown in Fig. 3(b). Filled squares and circles denote the results of $`t_o\mathrm{}`$ and $`t_{og}`$, respectively, obtained by DMA, and open squares and circles those by simulations. Figure 3(b) shows that $`\delta t_{og}`$ is almost linearly increased as $`\beta _1`$ is increased, while $`\delta t_o\mathrm{}`$ remains constant. In the limit of $`\beta _1=\beta _0=0.1`$, for which only common noises are applied ($`\beta _C=0.1`$ and $`\beta _I=0`$), we get $`\delta t_{og}=\delta t_o\mathrm{}`$, which shows that common noises do not work to reduce global fluctuations.

Ensemble-size dependence

Filled squares in Fig. 4(a) show the $`N`$ dependence of $`\delta t_o\mathrm{}`$ relevant to local fluctuations for $`\beta _0=0.1`$, $`\beta _1=0`$ and $`J=0`$, obtained by DMA theory, while open squares express that obtained by simulations. We note that $`\delta t_o\mathrm{}`$ is independent of $`N`$ because of no couplings ($`J=0`$). In contrast, $`\delta t_{og}`$ relevant to global fluctuations inversely decreases when the size $`N`$ is increased, as shown by filled and open circles which are obtained by DMA theory and simulations, respectively. The relation: $`\delta t_{og}(1/\sqrt{N})`$, holds as given by Eq. (33) for $`\beta _1=0`$. Figure 4(b) shows a similar plot for a finite value of $`\beta _1=0.05`$ with $`\beta _0=0.1`$, $`J=0`$ and $`N=100`$. In the limit of $`N\mathrm{}`$, the ratio of $`\delta t_{og}/\delta t_o\mathrm{}`$ approaches a finite value of $`\beta _1/\beta _0=0.5`$ \[Eq. (33)\].

Coupling-strength dependence

So far we have neglected the coupling of $`J`$, which is now introduced. Filled squares in Fig. 5(a) show the $`J`$ dependence of $`\delta t_o\mathrm{}`$ calculated by DMA theory for $`\beta _0=0.1`$, $`\beta _1=0`$ and $`N=100`$, while open squares that obtained by simulations. Filled and open circles express $`\delta t_{og}`$ in the DMA theory and simulations, respectively. We note that $`\delta t_o\mathrm{}`$ is much reduced as $`J`$ is increased although there is little change in $`\delta t_{og}`$. Figure 5(b) shows a similar plot of the $`J`$ dependence of firing time accuracy for finite $`\beta _1=0.05`$ with $`\beta _0=0.1`$ and $`N=100`$. Again a reduction in $`\delta t_o\mathrm{}`$ as increasing $`J`$ is more significant than that of $`\delta _{og}`$.

### B Synchronization ratio

One of important effects of the couplings is to yield synchronous firings in ensemble neurons. Figures 6(a) and 6(b) show the time course of the synchronization ratio $`S(t)`$ for $`J=100`$ and 200 $`\mu `$A/$`\mathrm{cm}^2`$, respectively, with $`\beta _0=0.1`$, $`\beta _1=0`$ and $`N=100`$: solid and dashed curves denote the results of DMA and simulations, respectively. Fairly large fluctuations in simulation results are due to a lack of trial number of one hundred, which is a limit of our computer facility. A comparison between Figs. 6(a) and 6(b) shows that $`S(t)`$ is increased as $`J`$ is increased: the maximum value of $`S(t)`$ in Fig. 6(b) is $`S_{max}=0.019`$ which is larger than $`S_{max}=0.007`$ in Fig. 6(b). Figure 6(c) shows the time course of $`S(t)`$ for a finite $`\beta _1=0.05`$ with $`\beta _0=0.1`$, $`J=100`$ and $`N=100`$. A significant increase in $`S`$ is realized at $`100\stackrel{<}{}t\stackrel{<}{}120`$ ms which is induced by an applied spike \[note the difference in vertical scales of Figs. 6(a), 6(b) and 6(c)\]. We note a fairly large value of $`S=0.25`$ even without an applied input spike at $`t\stackrel{<}{}100`$ or $`t\stackrel{>}{}120`$. This expresses the synchronization among the membrane potentials of ensemble neurons induced by added noises although they do not induce firings. In order to distinguish the synchronization with firings from that without firings, we define the firing-induced synchronization ratio, $`S^{}(t)`$, given by

$$S^{}(t)=S(t)S_b,$$
(66)
where $`S_b=(\beta _1/\beta _0)^2`$ denotes the background synchronization induced by noises only \[Eq. (36)\]. We get $`S_{max}=0.369`$, $`S_{max}^{}=0.119`$ and $`S_b=0.25`$ in Fig. 6(c). From a comparison of Fig. 6(c) with Fig. 6(a), we note that $`S^{}(t)`$ is also much increased by common noises.

An increase in $`S(t)`$ by an increase of $`J`$ is clearly shown in Fig. 7(a), where the maximum of $`S(t)`$ ($`S_{max}`$) is plotted as a function of $`J`$. A disagreement between results of DMA and simulations for $`J<50`$ is due to fluctuations in simulations because of insufficient trial number as mentioned above. The dependence of $`S_{max}`$ on the size $`N`$ is shown in Fig. 7(b) where $`\beta _0=0.1`$, $`\beta _1=0`$ and $`J=100`$. $`S_{max}`$ is decreased as increasing $`N`$. Figure 7(c) expresses the $`\beta _0`$ dependence of $`S_{max}`$ for $`\beta _1=0.05`$, $`J=100`$ and $`N=100`$. At $`\beta _0=\beta _1=0.05`$, we get $`S_{max}=1`$, which is decreased as increasing $`\beta _0`$. Filled squares in Fig. 7( c) denote $`S_{max}^{}`$, which shows the maximum around $`\beta _00.08`$. In contrast, Fig. 7(d) show the $`\beta _1`$ dependence of $`S_{max}`$ for $`\beta _0=0.1`$, $`J=100`$ and $`N=100`$. $`S_{max}`$ is increased as increasing $`\beta _1`$, and approaches unity as $`\beta _1\beta _0`$ (= 0.1). We note that $`S_{max}^{}`$ has the maximum at $`\beta _10.07`$.

## V Conclusion and Discussion

In previous Sec. IV, we have reported DMA calculations for a single spike input to HH neuron ensembles. DMA calculations and simulations have shown that

(a) $`\delta t_o\mathrm{}`$ increases as increasing $`\beta _0`$, or decreasing $`J`$, independently of $`\beta _1`$ and $`N`$,

(b) $`\delta t_{og}`$ increases as increasing $`\beta _0`$ or $`\beta _1`$, or decreasing $`N`$, independently $`J`$, and

(c) $`S_{max}`$ increases as increasing $`\beta _1`$ or $`J`$, or decreasing $`\beta _0`$ or $`N`$.

In order to understand these behaviors, we have tried to obtain phenomenological, analytical expressions for $`\delta t_o\mathrm{}`$, $`\delta t_{og}`$ and $`S_{max}`$ as functions of $`\beta _0`$, $`\beta _1`$, $`J`$ and $`N`$. For small $`J`$, we express $`\gamma _{v,v}`$ and $`\rho _{v,v}`$ in power series of $`J`$ at $`t=t_o^{}`$ where neurons fire, given by (see Appendix E of I, )

$`\gamma _{v,v}`$ $``$ $`\beta _0^2[1(a_1J+a_2J^2+)],`$ (67)
$`\rho _{v,v}`$ $``$ $`\beta _0^2\left[{\displaystyle \frac{1}{N}}+\left(1{\displaystyle \frac{1}{N}}\right)\left({\displaystyle \frac{\beta _1}{\beta _0}}\right)^2\right].`$ (68)
It is noted that in the limit of $`J=0`$ ($`w=0`$), Eqs. (59) and (60) reduce to Eq. (22). Substituting Eqs. (59) and (60) to Eqs. (28) and (32), we get

$`\delta t_o\mathrm{}`$ $``$ $`\beta _0[1{\displaystyle \frac{1}{2}}(a_1J+a_2J^2+)],`$ (69)
$`\delta t_{og}`$ $``$ $`\beta _0\left[{\displaystyle \frac{1}{N}}+\left(1{\displaystyle \frac{1}{N}}\right)\left({\displaystyle \frac{\beta _1}{\beta _0}}\right)^2\right]^{1/2}.`$ (70)
Equations (61) and (62) may explain the behavior of $`\delta t_o\mathrm{}`$ and $`\delta t_{og}`$ in items (a) and (b) mentioned above.

Next we will obtain the analytical expression for $`S_{max}`$. For small $`J`$, we get $`\gamma _{v,v}`$ and $`\rho _{v,v}`$ in power series of $`J`$ at $`t=t_o^{(m)}`$ where $`S(t)`$ takes the maximum value, given by (see Appendix E of I, )

$`\gamma _{v,v}`$ $``$ $`\beta _0^2[1(b_1J+b_2J^2+)],`$ (71)
$`\rho _{v,v}`$ $``$ $`\beta _0^2\left[{\displaystyle \frac{1}{N}}+\left(1{\displaystyle \frac{1}{N}}\right)\left({\displaystyle \frac{\beta _1}{\beta _0}}\right)^2\right].`$ (72)
Substituting Eqs. (63) and (64) to Eq. (34), we get

$`S_{\mathrm{max}}`$ $`=`$ $`\left({\displaystyle \frac{\beta _1}{\beta _0}}\right)^2+[{\displaystyle \frac{1}{N}}+(1{\displaystyle \frac{1}{N}})\left({\displaystyle \frac{\beta _1}{\beta _0}}\right)^2][b_1J+(b_2+b_1^2)J^2+].`$ (73)
Equation (65) is consistent with the item (c): $`S_{max}`$ increases as increasing $`\beta _1`$ and $`J`$, or decreasing $`\beta _0`$ and $`N`$. In the case of $`\beta _1=0`$, Eq. (65) shows that $`S_{max}`$ is independent of $`\beta _0`$, which is supported in DMA calculation and simulations (not shown). Expressions given by Eqs. (61), (62) and (65) are useful in a phenomenological sense. In principle, expressions as given by Eqs. (59) and (60) may be derived from DMA equations given by Eqs. (45)-(52) although we have not unfortunately succeeded in getting them because of their complexity.

Numerical calculations in Sec. IV have been made for the response to a single spike input. DMA is, however, applicable to arbitrary inputs. This will be demonstrated by adding spike trains to HH neuron ensembles, given by

$$K^{(e)}(t)=(\frac{I_i}{C})\underset{n}{}\alpha (tt_{in}),$$
(74)
where $`t_{in}`$ expresses the $`n`$th input time. Figures 8(a) and 8(b) show the time courses of $`\mu _v`$ and $`\sigma _{\mathrm{}}`$ (=$`\sqrt{\gamma _{v,v}}`$), respectively, for Poisson spike trains with the average interspike interval (ISI) of 25 ms; solid and dashed curves express results DMA and simulations, respectively, for $`\beta _0=0.1`$, $`\beta _1=0`$, $`J=0`$ and $`N=100`$. The time course of $`\mu _v`$ of DMA is in good agreement with that of simulations. A comparison between the input $`K^{(e)}`$ and output $`U_0`$ shows that when ISI of input is shorter than about 10 ms, HH neurons cannot respond because of the refractory period. Figures 8(b) shows that $`\sigma _{\mathrm{}}`$ of DMA is also in good agreement with that of simulations.

To summarize, the DMA theory previously proposed for FN neuron ensemble in I, has been generalized to an ensemble described by $`KN`$-dimensional stochastic DEs, which has been be replaced by $`K(K+2)`$-dimensional deterministic DEs expressed by means and second-order moments: contributions from the fourth-order moments are taken in account by the Gaussian decoupling approximation. DMA has been applied to HH neuron ensembles, for which we get 24-dimensional deterministic DEs. We have studied effects of noise, the coupling strength and the ensemble size on the firing time precision and the firing synchronization for single-spike inputs, obtaining the following results:

(i) the firing-time accuracy of the order of one-tenth ms is possible in a large-scale HH neuron ensemble, even without couplings,

(ii) the spike transmission is improved with the synchronous response by increasing the coupling strength, and

(iii) the synchronization is increased by common noises but decreased by independent noises.

The item (i) and (ii) are consistent with the SR results in HH neuron ensembles -. Although they are quite similar to the case of FN model discussed in I, their quantitative discussions are possible with the use of the realistic HH model. The item (iii) agrees with the result of Ref.  for SR in HH neuron ensembles subject of common and independent noises.

Our calculations have demonstrated the feasibility of DMA, whose advantages may be summarized as follows:

(1) because of the semi-analytical nature of DMA, some results may be derived without numerical calculations,

(2) DMA is free from the weak-coupling constraint although it assumes weak noises,

(3) a tractable small number of DEs makes calculations feasible for large-scale neuron ensembles with a fairly short computational time,

(4) DMA may be applicable to ensembles with fluctuations not only due to noises but also due to some inhomogeneities in model parameters, and

(5) DMA can be applied to more general stochastic systems besides neuron models.

As for the item (3), we may point out that for example, the CPU time of DMA calculations for a 200 ms time course of a $`N=100`$ HH neuron ensemble with the use of 1.8 GHz PC is 2 s, which is about 2500 times faster than the CPU time of 85 min ($``$ 5000 s) for direct simulations with 100 trials. It is necessary to stress the importance of the fourth-order contributions in stabilizing solutions of DMA, which is numerically demonstrated in Appendix B. Although expressions for the fourth-order contributions are lengthy, we are much benefited from them once they are derived and planted into computer programs . This paper is the first report on our DMA calculations of HH neuron ensembles. We are now under consideration to incorporate the time delay in the coupling terms in Eq. (1), with which the HH neuron ensemble may show the intrigue behavior like chaos. Such calculations will be reported in a separate paper.

## Acknowledgments

The author would like to express his sincere thanks to Professor Hideo Nitta for critical reading of the manuscript. This work is partly supported by a Grant-in-Aid for Scientific Research from the Japanese Ministry of Education, Culture, Sports, Science and Technology.

## A Derivation of DMA from the moment method for general neuron ensembles

In the moment method, we define the means, variances and covariances given by

$`m_{u_p}^i`$ $`=`$ $`<u_{pi}>,`$ (A1)
$`C_{u_p,u_q}^{i,j}`$ $`=`$ $`<\mathrm{\Delta }u_{pi}\mathrm{\Delta }u_{qj}>,`$ (A2)
where $`\mathrm{\Delta }u_{pi}=u_{pi}m_{u_p}^i`$. Assuming the weak couplings and adopting the Gaussian decoupling approximations for the fourth-order moments, we get DEs for general neuron ensembles described by Eqs. (1) and (2):

$`{\displaystyle \frac{dm_{u_p}^i}{dt}}`$ $`=`$ $`F^{(u_p)}+{\displaystyle \frac{1}{2}}{\displaystyle \underset{q}{}}{\displaystyle \underset{r}{}}F_{u_{qi},u_{ri}}^{(u_p)}C_{u_q,u_r}^{(i,i)}`$ (A4)
$`+\delta _{p1}[\left({\displaystyle \frac{w}{N1}}\right){\displaystyle \underset{k(i)}{}}(G+{\displaystyle \frac{1}{2}}G_{u_{1k}u_{1k}}C_{u_1,u_1}^{k,k})+K^{(e)}],`$
$`{\displaystyle \frac{dC_{u_p,u_q}^{i,j}}{dt}}`$ $`=`$ $`{\displaystyle \underset{r}{}}[F_{u_{ri}}^{(u_p)}C_{u_q,u_r}^{i,j}+F_{u_{rj}}^{(u_q)}C_{u_p,u_r}^{i,j}]+[\beta _0^2\delta _{ij}+\beta _1^2(1\delta _{ij})]\delta _{p1}\delta _{q1}`$ (A8)
$`+\delta _{p1}\left({\displaystyle \frac{w}{N1}}\right){\displaystyle \underset{k(i)}{}}G_{u_{1k}}C_{u_q,u_1}^{j,k}+\delta _{q1}\left({\displaystyle \frac{w}{N1}}\right){\displaystyle \underset{k(j)}{}}G_{u_{1k}}C_{u_p,u_1}^{i,k}`$
$`+{\displaystyle \frac{1}{6}}{\displaystyle \underset{r}{}}{\displaystyle \underset{s}{}}{\displaystyle \underset{w}{}}[F_{u_{ri}u_{si}u_{wi}}^{(u_p)}(C_{u_r,u_s}^{i,i}C_{u_w,u_q}^{i,j}+C_{u_r,u_w}^{i,i}C_{u_s,u_q}^{i,j}+C_{u_s,u_w}^{i,i}C_{u_r,u_q}^{i,j})`$
$`+F_{u_{rj}u_{sj}u_{wj}}^{(u_q)}(C_{u_r,u_s}^{j,j}C_{u_w,u_p}^{j,i}+C_{u_r,u_w}^{j,j}C_{u_s,u_p}^{j,i}+C_{u_s,u_w}^{j,j}C_{u_r,u_p}^{j,i})],`$
where $`F^{(u_p)}=F^{(p)}`$, $`F_{u_{ri}}^{(u_p)}=F^{(p)}/u_{ri}`$ and $`F_{u_{ri}u_{si}u_{wi}}^{(p)}=^{(3)}F^{(p)}/u_{ri}u_{si}u_{wi}`$ evaluted for the means of $`(\{m_{u_p}^i\})`$, and the last term in Eq. (A4) denotes the fourth-order contribution. The number of DEs is $`N_{eq}=KN+(1/2)KN(KN+1)=(1/2)KN(KN+3)`$.

In order to derive DMA from the moment method, we define the quantities given by

$`\overline{\mu }_\kappa `$ $`=`$ $`{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}m_\kappa ^i,`$ (A9)
$`\overline{\gamma }_{\kappa ,\lambda }`$ $`=`$ $`{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}C_{\kappa ,\lambda }^{i,i}+d_{\kappa ,\lambda },`$ (A10)
$`\overline{\rho }_{\kappa ,\lambda }`$ $`=`$ $`{\displaystyle \frac{1}{N^2}}{\displaystyle \underset{i}{}}{\displaystyle \underset{j}{}}C_{\kappa ,\lambda }^{i,j},\text{(}\kappa ,\lambda =u_p,u_q\text{)}`$ (A11)
where

$`d_{\kappa ,\lambda }`$ $`=`$ $`{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}\delta m_\kappa ^i\delta m_\lambda ^i,`$ (A12)
$`\delta m_\kappa ^i`$ $`=`$ $`m_\kappa ^i\mu _\kappa .`$ (A13)
We may show that Eqs. (A3) and (A4) with Eqs. (A3)-(A7) yield Eqs. (16)-(18): $`\overline{\mu }_\kappa =\mu _\kappa `$, $`\overline{\gamma }_{\kappa ,\lambda }=\gamma _{\kappa ,\lambda }`$ and $`\overline{\rho }_{\kappa ,\lambda }=\rho _{\kappa ,\lambda }`$. Then the moment method yields the same results as DMA as far as the averaged quantities are concerned (see also Appendix B of I).

## B The fourth-order contributions in DMA for HH neuron ensembles

The fourth-order contributions given by $`X_{\kappa ,\lambda }`$ and $`Y_{\kappa ,\lambda }`$ ($`\kappa ,\lambda =v,u_p`$) in Eqs. (45)-(52) are expressed by

$`X_{v,v}`$ $`=`$ $`F_{vmm}^{(v)}(\gamma _{v,v}\gamma _{m,m}+2\gamma _{v,m}\gamma _{v,m})+F_{vmh}^{(v)}(\gamma _{v,v}\gamma _{m,h}+2\gamma _{v,m}\gamma _{v,h})`$ (B3)
$`+F_{vnn}^{(v)}(\gamma _{v,v}\gamma _{n,n}+2\gamma _{v,n}\gamma _{v,n})+F_{mmm}^{(v)}\gamma _{v,m}\gamma _{m,m}+F_{nnn}^{(v)}\gamma _{v,n}\gamma _{n,n}`$
$`+F_{mmh}^{(v)}(\gamma _{v,h}\gamma _{m,m}+2\gamma _{v,m}\gamma _{m,h}),`$
$`X_{v,u_p}`$ $`=`$ $`{\displaystyle \frac{1}{2}}[F_{vmm}^{(v)}(\gamma _{v,u_p}\gamma _{m,m}+2\gamma _{u_p,m}\gamma _{v,m})+F_{vmh}^{(v)}(\gamma _{v,u_p}\gamma _{m,h}+\gamma _{u_p,m}\gamma _{v,h}+\gamma _{u_p,h}\gamma _{v,m})`$ (B7)
$`+F_{vnn}^{(v)}(\gamma _{v,u_p}\gamma _{n,n}+2\gamma _{u_p,n}\gamma _{v,n})+F_{mmm}^{(v)}\gamma _{u_p,m}\gamma _{m,m}+F_{nnn}^{(v)}\gamma _{u_p,n}\gamma _{n,n}`$
$`+F_{mmh}^{(v)}(\gamma _{u_p,h}\gamma _{m,m}+2\gamma _{u_p,m}\gamma _{m,h})+F_{vvv}^{(u_p)}\gamma _{v,v}\gamma _{v,v}`$
$`+F_{vvu_p}^{(u_p)}(\gamma _{v,u_p}\gamma _{v,v}+2\rho _{v,v}\gamma _{v,u_p})],`$
$`X_{u_p,u_q}`$ $`=`$ $`{\displaystyle \frac{1}{2}}[F_{vvv}^{(u_p)}\gamma _{v,u_q}\gamma _{v,v}+F_{vvv}^{(u_q)}\gamma _{v,u_p}\gamma _{v,v}+(F_{vvu_p}^{(u_p)}+F_{vvu_q}^{(u_q)})(\gamma _{u_p,u_q}\gamma _{v,v}+2\gamma _{v,u_q}\gamma _{v,u_p})],`$ (B8)
$`Y_{v,v}`$ $`=`$ $`F_{vmm}^{(v)}(\rho _{v,v}\gamma _{m,m}+2\rho _{v,m}\gamma _{v,m})+F_{vmh}^{(v)}(\rho _{v,v}\gamma _{m,h}+\rho _{v,m}\gamma _{v,h}+\rho _{v,h}\gamma _{v,m})`$ (B11)
$`+F_{vnn}^{(v)}(\rho _{v,v}\gamma _{n,n}+2\rho _{v,n}\gamma _{v,n})+F_{mmm}^{(v)}\rho _{v,m}\gamma _{m,m}+F_{nnn}^{(v)}\rho _{v,n}\gamma _{n,n}`$
$`+F_{mmh}^{(v)}(\rho _{v,h}\gamma _{m,m}+2\rho _{v,m}\gamma _{m,h}),`$
$`Y_{v,u_p}`$ $`=`$ $`{\displaystyle \frac{1}{2}}[F_{vmm}^{(v)}(\rho _{v,u_p}\gamma _{m,m}+2\rho _{u_p,m}\gamma _{v,m})+F_{vmh}^{(v)}(\rho _{v,u_p}\gamma _{m,h}+\rho _{u_p,m}\gamma _{v,h}+\rho _{u_p,h}\gamma _{v,m})`$ (B15)
$`+F_{vnn}^{(v)}(\rho _{v,u_p}\gamma _{n,n}+2\rho _{u_p,n}\gamma _{v,n})+F_{mmm}^{(v)}\rho _{u_p,m}\gamma _{m,m}+F_{nnn}^{(v)}\rho _{u_p,n}\gamma _{n,n}`$
$`+F_{mmh}^{(v)}(\rho _{u_p,h}\gamma _{m,m}+2\rho _{u_p,m}\gamma _{m,h})+F_{vvv}^{(u_p)}\rho _{v,v}\gamma _{v,v}`$
$`+F_{vvu_p}^{(u_p)}(\gamma _{v,u_p}\gamma _{v,v}+2\rho _{v,v}\gamma _{v,u_p})],`$
$`Y_{u_p,u_q}`$ $`=`$ $`{\displaystyle \frac{1}{2}}[F_{vvv}^{(u_p)}\rho _{v,u_q}\gamma _{v,v}+F_{vvu_p}^{(u_p)}(\rho _{u_p,u_q}\gamma _{v,v}+2\rho _{v,u_q}\gamma _{v,u_p})+F_{vvv}^{(u_q)}\rho _{v,u_p}\gamma _{v,v}`$ (B17)
$`+F_{vvu_q}^{(u_q)}(\rho _{u_p,u_q}\gamma _{v,v}+2\rho _{v,u_p}\gamma _{v,u_q})]`$
where $`F_{vmh}^{(v)}=^3F^{(v)}/vmh`$ et al. Although calculations and computer programming of fourth-order contributions given by Eqs. (B1)-(B6) are rather tedious, they play important roles in stabilizing the solution of DEs .

Here we demonstrate the importance of the fourth-order contributions in the case of a single HH neuron ($`N=1`$) for which $`w=0`$, $`\gamma _{\kappa ,\lambda }=\rho _{\kappa ,\lambda }`$ and $`X_{\kappa ,\lambda }=Y_{\kappa ,\lambda }`$ in Eqs. (45)-(52) and (B1)-(B6). Figure 9(a) shows the time course of $`\mu _v`$ for $`\beta _0=0.1`$ and $`\beta _1=0`$ when the constant input of $`I_i=10`$ $`\mu `$A/$`\mathrm{cm}^2`$ is applied at $`t=0`$ ms. The solid and dashed curves express the results of DMA and the simulation (100 trials), respectively. The dotted curve denote the result of DMA2 (the second-order DMA) in which the fourth-order contributions are neglected ($`X_{\kappa ,\lambda }=Y_{\kappa ,\lambda }=0`$). For $`t<60`$ ms, all results seem to be in good agreement. At $`t\stackrel{>}{}60`$, however, the solution of DMA2 becomes unstable and significantly deviates from those of DMA and the simulation. From the time course of $`\sigma _{\mathrm{}}=\sqrt{\gamma _{v,v}}`$ shown in Fig. 9(b), we note that such deviation of DMA2 already starts from $`t30`$ ms. The solution of DMA2 is stable at $`\beta 0.037`$ for the constant current of $`I_i=10`$ $`\mu `$A/$`\mathrm{cm}^2`$.

Figures 10(a) and 10(b) show the time courses of $`\mu _v`$ and $`\sigma _{\mathrm{}}`$ for $`\beta _0=0.2`$ and $`\beta _1=0`$ when we apply the periodic spike train input given by Eq. (66) with $`I_i=5`$ $`\mu `$A/$`\mathrm{cm}^2`$ and a constant ISI of 25 ms. Figure 10(b) clearly shows that the result of DMA2 deviates from those of DMA and simulations from the first spike input and that the result of DMA2 diverges at the second spike input. The solution of DMA2 is stable only at $`\beta 0.178`$ for this periodic spike.

## C The moment method for HH neuron ensembles

We will derive DEs in the moment method for HH neuron ensembles, defining the means, variances and covariances given by

$`m_v^i`$ $`=`$ $`<v_i>,`$ (C1)
$`m_{u_p}^i`$ $`=`$ $`<u_{pi}>,`$ (C2)
$`C_{v,v}^{i,j}`$ $`=`$ $`<\mathrm{\Delta }v_i\mathrm{\Delta }v_j>,`$ (C3)
$`C_{v,u_p}^{i,j}`$ $`=`$ $`<\mathrm{\Delta }v_i\mathrm{\Delta }u_{pi}>,`$ (C4)
$`C_{u_p,u_q}^{i,j}`$ $`=`$ $`<\mathrm{\Delta }u_{pi}\mathrm{\Delta }u_{qj}>,`$ (C5)
where $`\mathrm{\Delta }v_i=v_im_v^i`$ and $`\mathrm{\Delta }u_{pi}=u_{pi}m_{u_p}^i`$. Adopting the weak-coupling approximation and the Gaussian approximation for the fourth-order terms, we get DEs given by

$`{\displaystyle \frac{dm_v^i}{dt}}`$ $`=`$ $`F^{(v_i)}+{\displaystyle \frac{1}{2}}{\displaystyle \underset{p=2}{\overset{4}{}}}{\displaystyle \underset{q=2}{\overset{4}{}}}F_{u_p,u_q}^{(v_i)}C_{u_p,u_q}^{i,i}+{\displaystyle \underset{p=2}{\overset{4}{}}}F_{v,u}^{(v_i)}C_{v,u_p}^{i,i}`$ (C7)
$`+\left({\displaystyle \frac{w}{N1}}\right){\displaystyle \underset{k(i)}{}}(G^{(v_k)}+{\displaystyle \frac{1}{2}}G_{vv}^{(v_k)}C_{v,v}^{k,k})+K^{(e)}(t),`$
$`{\displaystyle \frac{dm_{u_p}^i}{dt}}`$ $`=`$ $`F^{(u_{pi})}+{\displaystyle \frac{1}{2}}F_{v,v}^{(u_{pi})}C_{v,v}^{i,i}+F_{v,u_p}^{(u_{pi})}C_{v,u_p}^{i,i},`$ (C8)
$`{\displaystyle \frac{dC_{v,v}^{i,j}}{dt}}`$ $`=`$ $`2F_v^{(v_i)}C_{v,v}^{i,j}+{\displaystyle \underset{p=2}{\overset{4}{}}}F_{u_p}^{(v_i)}(C_{v,u_p}^{i,j}+C_{v,u_p}^{j,i})+[\beta _0^2\delta _{ij}+\beta _1^2(1\delta _{ij})]`$ (C10)
$`+\left({\displaystyle \frac{w}{N1}}\right)[{\displaystyle \underset{k(i)}{}}G_v^{(v_k)}C_{v,v}^{j,k}+{\displaystyle \underset{k(j)}{}}G_v^{(v_k)}C_{v,v}^{i,k}]+Z_{v,v}^{i,j},`$
$`{\displaystyle \frac{dC_{v,u_p}^{i,j}}{dt}}`$ $`=`$ $`(F_v^{(v_i)}+F_{u_p}^{(u_{pi})})C_{v,u_p}^{i,j}+{\displaystyle \underset{q=2}{\overset{4}{}}}F_{u_q}^{(v_i)}C_{u_q,u_p}^{i,j}+F_v^{(u_{pi})}C_{v,v}^{i,j}`$ (C12)
$`+\left({\displaystyle \frac{w}{N1}}\right){\displaystyle \underset{k(i)}{}}G_v^{(v_k)}C_{v,u_p}^{k,j}+Z_{v,u_p}^{i,j},`$
$`{\displaystyle \frac{dC_{u_p,u_q}^{i,j}}{dt}}`$ $`=`$ $`(F_{u_p}^{(u_{pi})}C_{u_q,u_p}^{i,j}+F_{u_q}^{(u_{qi})}C_{u_p,u_q}^{i,j})+F_v^{(u_{pi})}C_{v,u_q}^{i,j}+F_v^{(u_{qi})}C_{v,u_p}^{j,i}+Z_{u_p,u_q}^{i,j},`$ (C13)
with

$`Z_{v,v}^{i,j}`$ $`=`$ $`{\displaystyle \frac{1}{2}}[F_{vmm}^{(v)}(C_{v,v}^{i,j}C_{m,m}^{j,j}+2C_{v,m}^{i,j}C_{v,m}^{j,j}+C_{v,v}^{i,j}C_{m,m}^{i,i}+2C_{v,m}^{j,i}C_{v,m}^{i,j})`$ (C18)
$`+F_{vmh}^{(v)}(C_{v,v}^{i,j}C_{m,h}^{j,j}+C_{v,m}^{i,j}C_{v,h}^{j,j}+C_{v,h}^{i,j}C_{v,m}^{j,j}+C_{v,v}^{i,j}C_{m,h}^{i,i}+C_{v,m}^{j,i}C_{v,h}^{i,i}+C_{v,h}^{j,i}C_{v,m}^{i,i})`$
$`+F_{vnn}^{(v)}(C_{v,v}^{i,j}C_{n,n}^{j,j}+2C_{v,n}^{i,j}C_{v,n}^{j,j}+C_{v,v}^{i,j}C_{n,n}^{i,i}+2C_{v,n}^{j,i}C_{v,n}^{i,i})`$
$`+F_{mmm}^{(v)}(C_{v,m}^{i,j}C_{m,m}^{j,j}+C_{v,m}^{j,i}C_{m,m}^{i,i})+F_{nnn}^{(v)}(C_{v,n}^{i,j}C_{n,n}^{j,j}+C_{v,n}^{j,i}C_{n,n}^{i,i})`$
$`+F_{mmh}^{(v)}(C_{v,h}^{i,j}C_{m,m}^{j,j}+2C_{v,m}^{i,j}C_{m,h}^{j,j}+C_{v,h}^{j,i}C_{m,m}^{i,i}+2C_{v,m}^{j,i}C_{m,h}^{i,i})],`$
$`Z_{v,u_p}^{i,j}`$ $`=`$ $`{\displaystyle \frac{1}{2}}[F_{vmm}^{(v)}(C_{v,u_p}^{i,j}C_{m,m}^{i,i}+2C_{m,u_p}^{i,j}C_{v,m}^{i,i})+F_{vmh}^{(v)}(C_{v,u_p}^{i,j}C_{m,h}^{i,i}+C_{m,u_p}^{i,j}C_{v,h}^{i,i}+C_{h,u_p}^{i,j}C_{v,m}^{i,i})`$ (C22)
$`+F_{vnn}^{(v)}(C_{v,u_p}^{i,j}C_{n,n}^{i,i}+2C_{n,u_p}^{i,j}C_{v,n}^{i,i})+F_{mmm}^{(v)}C_{m,u_p}^{i,j}C_{m,m}^{i,i}+F_{nnn}^{(v)}C_{n,u_p}^{i,j}C_{n,n}^{i,i}`$
$`+F_{mmh}^{(v)}(C_{h,u_p}^{i,j}C_{m,m}^{i,i}+2C_{m,u_p}^{i,j}C_{m,h}^{i,i})+F_{vvv}^{(u_p)}C_{v,v}^{i,j}C_{v,v}^{j,j}`$
$`+F_{vvu_p}^{(u_p)}(C_{v,u_p}^{i,j}C_{v,v}^{j,j}+2C_{v,v}^{i,j}C_{v,u_p}^{j,j})],`$
$`Z_{u_p,u_q}^{i,j}`$ $`=`$ $`{\displaystyle \frac{1}{2}}[F_{vvv}^{(u_p)}C_{v,u_q}^{i,j}V_{v,v}^{i,i}+F_{vvu_p}^{(u_p)}(C_{u_p,u_q}^{i,j}C_{v,v}^{i,i}+2C_{v,u_q}^{i,j}C_{v,u_p}^{i,i})+F_{vvv}^{(u_q)}C_{u_p,v}^{i,j}C_{v,v}^{j,j}`$ (C24)
$`+F_{vvu_q}^{(u_q)}(C_{u_p,u_q}^{i,j}C_{v,v}^{j,j}+2C_{u_p,v}^{i,j}C_{v,u_q}^{j,j})],`$
where $`F^{(v_i)}`$, $`F_v^{(v_i)}=F/v_i`$ and $`F_{v,u_p}^{(v_i)}=^2F/v_iu_{pi}`$ evaluated for the averages of $`(m_v^i,\{m_{u_p}^i\})`$, and similar derivatives for $`F^{(u_{pi})}`$ and $`G^{(v_k)}`$. Equations given by Eqs. (C6)-(C13) denote the result of the fourth-order moment method. The second-order moment method was applied to a single HH neuron by RT , whose result is given by Eqs. (C6)-(C10) when we set $`i=j=1`$, $`w=0`$, $`\beta _1=0`$ and $`Z_{v,v}=Z_{v,u_p}=Z_{u_p,u_q}=0`$. Equations (C6)-(C13) lead to DMA equation given by Eqs. (45)-(52) and (B1)-(B6) if we adopt the relations as given by Eqs. (A5)-(A7). In particular, in the case of $`N=1`$, Eqs. (C6)-(C13) are identical with Eqs. (45)-(52) and (B1)-(B6) of DMA if we read $`m_\kappa ^1=\mu _\kappa `$, $`C_{\kappa ,\lambda }^{1,1}=\gamma _{\kappa ,\lambda }=\rho _{\kappa ,\lambda }`$, and $`Z_{\kappa ,\lambda }^{1,1}=X_{\kappa ,\lambda }=Y_{\kappa ,\lambda }`$.