Open access peer-reviewed chapter

Simulation of Neural Behavior

By Tatsuo Kitajima, Zonggang Feng and Azran Azhim

Submitted: October 23rd 2015Reviewed: April 29th 2016Published: August 24th 2016

DOI: 10.5772/64028

Downloaded: 760

Abstract

The brain is an organ that takes the central role in advanced information processing. There exist great many neurons in our brain, which build complicated neural networks. All information processing in the brain is accomplished by neural activity in the form of neural oscillations. In order to understand the mechanisms of information processing, it is necessary to clarify functions of neurons and neural networks. Although the current progress of experiment technology is remarkable, only experiments by themselves cannot uncover the behavior of only a single neuron. Computational neuroscience is a research field, which fills up the deficiency in experiments. By modeling the essential features of a neuron or a neural network, we can analyze their fundamental properties by computer simulation. In this chapter, one aspect of computational neuroscience is described. At the first, the cell membrane and a neuron can be modeled by using an RC circuit. Next, the Hodgkin-Huxley model is introduced, which has the function of generation of action potentials. Furthermore, many neurons show the subthreshold resonance phenomena, and the cell membrane is necessary to be modeled by an RLC circuit. Finally, some simulation results are shown, and properties of such neuronal behaviors are discussed.

Keywords

  • cell membrane
  • action potential
  • neural oscillation
  • subthreshold resonance phenomenon
  • RLC circuit

1. Introduction

Our brain is an extraordinary microsome and has been completely shrouded in mystery. However, its mystery has been just a little bit by bit solved owing to recent advances in experimental technologies and tremendous development of computers. Many people can simply say “brain,” but it is a general term for a collection of six main regions, that is, cerebrum, diencephalon, midbrain, cerebellum, pons, and medulla oblongata. The brain is an organ that takes the central role in advanced information processing, such as visual, auditory, speech or language faculties, motion control, recognition, emotion, and so on. According to advances of experiment and computer technology, the research of brain science or neuroscience has been made not only in the fields of medicine, biology, biochemistry, pharmacology, and psychology but also in the field of engineering.

The present-day computers have outstanding processing capacity. For example, they can find the data that satisfy some requirements among huge quantities of data (database) or can calculate over five trillion figures of pi(π). Therefore, many people are inclined to think that our brain will be able to be replaced by computer in near future. Surely, computers excel at processing of digitized data and processing by following a standard algorithm. However, it can hardly execute processing, such as recognition of ambiguity figures (such as illusionism) or inference based on imperfect information, which our brain can instantaneously carry out. Reason for this comes from differences in ways of information processing of the computer and our brain. The current computers, called von Neumann computer, are grounded in sequential processing by using central processing units (CPUs) and memory storages, while on the other hand, our brain bases on parallel and distributed processing through neural networks whose components are neurons.

There exist tens of billions of neurons in our brain, which build neural networks in complicated arrangement. All information processing in the brain is accomplished by neural activity in the form of neural oscillations that cause cortical oscillations (delta, theta, alpha, beta, or gamma oscillation). In order to clarify the mechanisms of advanced information processing in the brain, such as learning and memory, it is necessary to understand functions and features of neurons and neural networks. Although the current progress in experiment technology and measuring system is remarkable, only experiments by themselves cannot uncover the behavior of only a single neuron, because even a single neuron has complex biophysical characteristics and never stops growth. Computational neuroscience is a research field which fills up such a deficiency in experiments. By modeling the essential features of a neuron or a neural network at multiple spatial-temporal scales, we can capture and analyze the fundamental properties of a neuron or a neural network by computer simulation. Moreover, we can even offer some suggestions to experimental study by taking into account the probable results obtained from the simulation.

2. Neuron model with a low-pass filter property

2.1. Electrical circuit model of the cell membrane

Neurons play a key role in almost all brain functions. Fundamental function of neurons is to generate action potentials when they received sufficient stimuli from the environment. Once action potentials are generated, they are transmitted to other neurons so as to communicate information from one neuron to another. There exist many types of neurons in the brain, such as pyramidal neurons in the hippocampus and neocortex (Figure 1(a)), motor neurons in motor cortex (Figure 1(b)), or Purkinje cells in the cerebellum (Figure 1(c)) [1]. Although their shapes are different, they have basically the same structure. As shown in Figure 2(a), a neuron is composed of three parts, that is, the soma (cell body) where action potentials are generated, the dendrite that receives inputs from other neurons, and the axon along which action potentials are transmitted to axon terminals. One thing especially worth mentioning, the dendrite of a neuron has hundreds to thousands of spines, on which axon terminals of other neurons connect. This junction is called a synapse, through which information are transmitted from one neuron to another (Figure 2(b)). Actually, there exist two kinds of synapses, one of which is an electric synapse and the other is a chemical synapse [1]. The former is a junction where neurons are directly contacted each other and information are electrically conducted from one neuron to another. This junction is also called a gap junction. On the other hand, the latter one is a junction with a cleft, called a synaptic cleft, into which neurochemical transmitters are released from the axon terminal and they bind to receptors on the spine head. Electrical synapses are found at the sites that require the fastest possible response, such as nociceptive reflex, whereas chemical synapses are found in almost all neurons of the brain. Figure 2(b) shows an example of a chemical synapse. Both synapses have a very important role in signal processing between neurons.

Figure 1.

Various types of neurons. (a) Pyramidal neuron (cortex), (b) motor neuron (spinal cord),and (c) Purkinje cell (cerebellum).

Figure 2.

(a) Schematic neuron (Structure of neuron). (b) Synaptic connection at the synapse.

Surfaces of a neuron are covered with the cell membrane, which separates the interior of cell from the exterior environment. The cell membrane is composed of protein, lipid, and carbohydrate [1]. As shown in Figure 3(a), it is composed of two layers of phospholipid molecules, each of which has a hydrophilic head (circle) and hydrophobic tail (two waved lined), and both hydrophobic tails face each other inside the cell membrane. This structure is called lipid bilayer. Furthermore, the concentration of the extracellular ions, such as Na+, Cl-, or Ca2+, is higher than the intracellular one. Contrarily, the concentration of intracellular ion, such as K+, is higher than the extracellular one. In addition, many types of ion channels (protein) are penetrating the cell membrane. Those ion channels are normally closed. If neurochemical transmitters released from the presynaptic axon terminal bind to receptors of the corresponding ion channels on the spine head, those ion channels are activated and open. Subsequently, specific ion flow occurs according to their ionic gradients. At the resting state, those channels are closed and no ionic flows occur except for small leakage.

Figure 3.

Cell membrane. (a) Cross section of a cell membrane lipid bilayer and (b) equivalent RC circuit model of cell membrane.

Based on the above properties, the cell membrane has the following electrical properties:

  1. It is lipid bilayer, that is, it is composed of two parallel plates. Thus, the cell membrane has characteristics similar to “capacitance,” C.

  2. The difference between intracellular and extracellular ion concentrations corresponds to “power source,” Ei (i = Na+, K+, Cl or Ca2+).

  3. The ionic flowability of opening ion channels is thought of as “resistance” Ri or “conductance” 1/Ri.

  4. The corresponding flows of Na+, K+, Cl, or Ca2+ through ion channels are “current,” INa, IK, ICl, or ICa.

With these points in mind, the cell membrane can be modeled by an equivalent RC circuit, which is shown in Figure 3(b). Once an RC circuit is obtained, we can obtain its dynamics by using Ohm’s law, Kirchhoff’s law, or other knowledge of electrical circuit theory. From Figure 3(b), the following equation is obtained:

CdVdt=1RL(VEL)1Ri(VEi)+Iinp,E1

where V is the membrane potential of the cell membrane, C is the membrane capacitance, RL is the leakage resistance, EL is the reversal potential, Ri is the flowability of ion i(i = Na+, K+, Cl, or Ca2+), Ei is the corresponding ionic equilibrium potential, and Iinp is the specific input current given to the cell membrane. As a neuron is covered with the cell membrane, a synapse or a soma can be also expressed by using an RC circuit. Accordingly, we can study the synaptic properties or neuronal characteristics by using computer simulations.

2.2. Generation of action potentials (Hodgkin-Huxley model)

In this section, we give one model that can generate an action potential, which is the basic function of a neuron. When a dendritic spine receives stimuli from an axon terminal of other neuron, the membrane potential of a spine head changes depending on that stimulus. Those potential changes are transmitted to the soma (strictly speaking, the axon hillock in the neighborhood of the soma) through dendrites and integrated there. If the accumulated potential of the soma exceeds the threshold, an action potential is generated. Generated action potentials are transmitted to axon terminals along the axon. Based on this knowledge, McCulloch and Pitts expressed a neuron as a product-sum threshold element in 1943 [2]. Their model is a formal neuron model, called McCulloch-Pitts model, and is shown in Figure 4.

Figure 4.

Formal neuron model (The McCulloch-Pitts model).

The McCulloch-Pitts model is expressed as follows:

u(t)=i=1Nwi(t)xi(t),E2
y(t)=φ(u(t))={1:u(t)θ0:u(t)<θ,E3

where xi(t) is an input from ith neuron, wi(t) is a weight from a neuron i, u(t) is a state (potential) of a neuron, y(t) is its output, and θ is a threshold. In this model, if a state u(t) exceeds a threshold θ, output 1 is send to other neurons. Notice that, however, McCulloch-Pitts model does not consider a refractory period, during which neurons cannot or find it hard to generate the next action potential.

As the McCulloch-Pitts model was a very easy model for engineers to understand the mechanism of generation of action potentials, many engineers have applied this model to study basic neuronal behaviors. The most prominent example is the application to the perceptron, which was known as one of the powerful tools for some kinds of pattern recognition problems. Although there exist many variations of the McCulloch-Pitts model, one of them uses a sigmoid function instead of a step function expressed by Eq. (3). This kind of model is applied to the back propagation algorithm and recently the deep learning method, because a sigmoid function is a differentiable function. However, the practical neurons are not so simple as the McCulloch-Pitts model and the back propagation algorithm. Therefore, more profound considerations were necessary to describe complicated neuronal behaviors.

In 1952, Hodgkin and Huxley developed one mathematical model that explains the generation of an action potential (impulse or spike) based on physiological experiments for a squid giant axon [3]. As described in the previous section, the extracellular concentration of Na+ is higher than the intracellular one, and the intracellular concentration K+ is higher than the extracellular one, and the cell membrane has both Na+ permeable channel (Na channel) and K+ permeable channel (K channel). Hodgkin and Huxley found that both Na and K channels are voltage-dependently activated, that is, the activation and inactivation of these channels are affected by the membrane potential of the cell membrane. They also elucidated that action potentials are generated by increased or decreased activation and inactivation of Na channel and increased or decreased activation of K channel. Based on the results of physiological experiments for a squid giant axon, they showed that an action potential is generated whenever the cell membrane is depolarized over the threshold. They proposed a schematic electrical circuit model that can explain the mechanism for generation of action potentials, called the Hodgkin-Huxley model (HH model).

Figure 5.

The Hodgkin-Huxley model. (a) Conductance-based electrical circuit of the Hodgkin-Huxley model and (b) simulation result.

Figure 5(a) shows the HH model and its dynamics is given as follows:

Iinp=CdVdt+IL+IN+IK,E4

where V is the membrane potential of the cell, C is the capacitance, IL is the leak current, INa and IK are currents through Na channel and K channel, respectively, and Iinp is the input current. They proposed the empirical formulae, which appropriately indicate activation and inactivation properties of Na channel and activation properties of K channel, that is, the change of ionic permeability of Na+ and K+. Currents IL, INa, and IK are given as follows [3]:

IL=g¯L(VEL),E5
IN=g¯Nam(V,t)3h(V,t)(VEN),E6
IK=g¯Kn(V,t)4(VEK),E7

where g¯Lis the leakage conductance, g¯Naand g¯Kare the amplitude of Na channel conductance and K channel conductance, respectively, EL is the resting potential, and ENa and EK are equilibrium potentials of Na channel and K channel. m(V, t) and h(V, t) are activation and inactivation variables of Na channel, and n(V,t) is an activation variable of K channel. They gave the following empirical formula:

dx(V,t)dt=αx(V)(1x(V,t))βx(V)x(V,t),(x=m,h,n)E8
αm(V)=0.1(25V)exp(25V10)1,βm(V)=4exp(V18),E9
αh(V)=0.07exp(V20),βh(V)=1exp(30V10)+1,E10
αn(V)=0.01(10V)exp(10V10)1,βm(V)=0.125exp(V80).E11

Figure 5(b) shows one example of computer simulation results for the HH model. When a continuous DC input is given to the HH model, action potentials can be generated at certain interval, that is, with a refractory period. Regardless of the strength of inputs, action potentials have the same shape and size. All differential equations were solved by the fourth-order Rung-Kata method by using C++. Parameters used here were as C = 1μF/cm2, g¯L=0.3mS/cm2,g¯Na=0.3mS/cm2,g¯K=0.3mS/cm2, ENa= 115 mV, EK = −12 mV, Iinp = 10 mA/cm2, and the resting potential = 0 mV.

3. Neuron model with a band-pass filter property

3.1. Subthreshold resonance phenomenon

As described in Section 2.2, neurons can generate action potentials depending on the strength of DC input stimuli. As shown in Figure 6(a), for small DC inputs (dark blue, light blue, dashed red lines), action potentials are not generated by reason that membrane potentials do not exceed the threshold. However, if a larger input (red line) is given to a neuron, the membrane potential can exceed the threshold and as a result, an action potential is generated. On the contrary, when AC inputs are given to a neuron, outputs of a neuron are unlike the cases of DC inputs, apart from whether the membrane potential exceeds the threshold or not. We consider three AC inputs (blue, red, and green in Figure 6(b)), whose amplitudes are equal but their frequencies are different (f1<f2<f3). By using an AC input with frequency f1(blue), the membrane potential (blue line) is assumed to be obtained under the threshold level, that is, in a subthreshold level. If the input frequency increases from f1 to f2 (red), the membrane potential (red line) is still in a subthreshold level, but its amplitude becomes larger than that of frequency f1 (blue line). However, if the input frequency further increases to f3 (green), the amplitude of the membrane potential (green line) reduces and becomes smaller than that of frequency f2 (red line). Instead of AC inputs with a single frequency, let an AC input whose frequency increases with time be given to this neuron. This kind of AC is called a chirp current. Then, its membrane potential has the shape with an expanded center section as shown in Figure 6(c), that is, the membrane potential takes the maximum at a specific frequency, however, remains at a subthreshold level. As its FFT shows, this neuron has a band-pass property, that is, frequency selectivity. These kinds of oscillatory phenomena in a subthreshold level are called the subthreshold resonance phenomena.

Figure 6.

Subthreshold resonance phenomena. (a) DC inputs with different amplitudes, (b) AC inputs with different frequencies, and (c) a chirp current input whose frequency increases as time increases.

Subthreshold resonance oscillations have been found in many excitatory and/or inhibitory neurons in the whole brain. Mauro et al. first reported a subthreshold resonance oscillation in squid giant axon [4]. Koch discussed these resonance oscillations in relation to the cable theory [5]. Since then, these resonance phenomena have been observed in many neurons in various regions of the brain, such as trigeminal root ganglion [6], inferior olive [7], and thalamus [8, 9]. These subthreshold resonance phenomena have been also reported in cortical neurons [1014], and in the 2000s, also in hippocampal neurons in CA1 [15, 16]. Although it is suspected that frequency selectivity of neurons should play an important role in behavioral or perceptual functions in animals, their practical roles have still been unclear. Recently, Narayanan and Johnston [17] reported that subthreshold resonance oscillations in hippocampal CA1 neurons are closely related to the long-term synaptic plasticity, which is currently considered as one of possible foundations of learning and memory [18, 19]. So, it is very interesting and attractive to study those resonance oscillatory features, in order to clarify the mechanisms of higher information processing functions in the brain, such as learning, short-term memory, or working memory.

As already described, the cell membrane is usually modeled by an RC circuit. However, if a chirp current is given to an RC circuit, a membrane potential shows only a property of low-pass filter shown in Figure 7(a). On electrical circuit theory, resonance circuit must contain inductive elements, that is, inductance L. Indeed, if a chirp current is given to an RLC circuit, a membrane potential shows a band-pass property shown in Figure 7(b). Having many neurons, the subthreshold resonance phenomena indicate that those neurons must have some kind of inductive factor. So exactly, what is a distinguishing major role of such inductive characteristics in the cell membrane? By advances of experimental technique, it has been reported that many kinds of voltage-dependent ion channels have an important role in subthreshold phenomena. Such ion channels involved in the subthreshold resonance phenomena are different from neuron to neuron that belongs to brain regions. Among them, slow non-inactivating K+ channel (Krs channels) [10], hyperpolarization-activated cationic channel (h channel) [12], and persistent Na+ channels (NaP channel) [13] are well known. In addition to these channels, voltage-dependent Ca2+ channels in neurons and/or dendritic spines [8] are also concerned in the subthreshold resonance oscillation.

Figure 7.

Calculated membrane potential for a chirp current and its magnitude of FFT. (a) RC circuit and (b) RLC circuit.

3.2. Inductive property of voltage-dependent ion channels

A hyperpolarization-activated cation channel (h channel) and a persistent sodium channel (NaP channel) are known to mediate the subthreshold resonance oscillation observed in entorhinal cortical neurons [12, 14]. In this section, we show how such voltage-dependent ion channels have inductive properties. We consider a compartment neuron model with h channel and NaP channel as shown in Figure 8. Its dynamics are expressed by the following conductance-based equations [13]:

CdVdt=ILIhINaPIinp,E12

Figure 8.

A compartmental neuron model with h channel and NaP channel.

where V is the membrane potential, IL is the leak current, Ih and INaP are the currents through h channel and NaP channel, respectively, and Iinp is the input current. The leak current IL is given by

IL=g¯L(VEL).E13

where g¯Lis the leak conductance and EL is the resting potential. Currents Ih and INaP are given as follows:

Ih=gh(V)(VEh)=g¯h{0.65mhf(V)+0.35mhs(V)}(VEh).E14
INaP=gNaP(V)(VEN)=g¯NaPmNaP(V)(VEN).E15

where gh(V) and gNaP(V) are, respectively, the h channel conductance and the NaP channel conductance, g¯hand g¯NaPare, respectively, the maximum amplitude of gh(V) and gNaP(V), and Eh and EN are the equilibrium potentials for K+ through h channel and Na+ for NaP channel, respectively. mhf(V) and mhs(V) are, respectively, the fast activation and slow activation variables of h channel, and mNaP(V) is an activation variable of NaP channel. They satisfy the following equations:

dmxdt=1τx(mxmx),(x=hf,hs,NaP)E16
mhf=11+exp(V+79.29.78),τhf=1+0.51exp(V1.710)+exp(V+34052),E17
mhs=11+exp(V+71.37.9),τhs=1+5.6exp(V1.714)+exp(V+26043),E18
mNaP=11+exp(V+386.5),τNaP=0.15ms.E19

3.2.1. Equivalent admittance (impedance) of h channel

Let V* be the equilibrium potential, Ih*be the h-current at V*. From Eq. (14), Ih*satisfies the following relation:

Ih*=g¯h{0.65mhf(V*)+0.35mhs(V*)}(V*Eh).E20

When the membrane potential V(t) changes from V* to V* + δV(t), where δV(t) is a small variation of the membrane potential from V*, the current Ih(t) also changes from Ih*to Ih*+δIh*(t), where δIh(t) is a small variation of h current caused by δV(t). Ih*+δIh(t)satisfies the following relation:

Ih*+δIh(t)=g¯h{0.65mhf(V*+δV(t))+0.35mhf(V*+δV(t))}(V*+δV(t)Eh).E21

Let mhf(V* + δV) approximate by mhf (V*) + δmhf and mhs(V* + δV) by mhs(V*) + δmhs for a small variation δV. Then, Eq. (21) can be expressed by the following equation:

Ih*+δIh(t)=g¯h{0.65mhf(V*)+0.35mhs(V*)}(V*Eh)+g¯h{0.65δmhf+0.35δmhs}(V*Eh)+g¯h{0.65mhf(V*)+0.35mhs(V*)}δV(t)+g¯h{0.65δmhf+0.35δmhs}δV(t).E22

By subtracting Eq. (20) from Eq. (22) and dropping the higher-order variation terms than the second, which appeared on the right-hand side of Eq. (22), the following equation is obtained:

δIhg¯h{0.65δmhf+0.35δmhs}(V*EK)+g¯h{0.65mhf*+0.35mhs*}δV.E23

As an activation variable mhf(V) satisfies Eq. (16), mhf(V*) and mhf(V* + δV) must satisfy the following equations:

dmhf(V*)dt=1τhf(V*){mhf(V*)mhf},E24
dmhf(V*+δV)dt=1τhf(V*+δV){mhf(V*+δV)mhf(V*+δV)},E25

where left terms of Eqs. (24) and (25), dmhf(V*)/dtand dmhf(V*+δV)/dt, represent the quantity dmhf/dtevaluated at V* and V* + δV, respectively. By approximating also τhf(V*+δV)by τhf(V*)+δτhf, where δτhfis a small variation caused by δV, Eq. (25) is written as follows:

d[mhf(V*)+δmhf]dt1τhf(V*)+δτhf{mhf(V*)+δmhf[mhf(V*)+δmhf]}1τhf(V*)[1δτhfτhf(V*)+...]{mhf(V*)+δmhf(mhf(V*)+δmhf)}.E26

By dropping the higher-order variation terms than the second and subtracting Eq. (24) from Eq. (26), the following equation is obtained:

mhfdt1τhf(V*){δmhfδmhf}.E27

Furthermore, as a small variation δmhfmay be approximately expressed by [dmhf(V*)/dV]δV, Eq. (27) becomes as follows:

mhfdt1τhf(V*){dmhf(V*)dVδVδmhf}.E28

Using the differential operator p instead of time derivative (d/dt), Eq. (28) can be written as follows:

(p+1τhf(V*))δmhf=1τhf(V*)(dmhf(V*)dV)δV.E29

From this relation, a small variation δmhf is definitely expressed by δV as follows:

δmhf=1τhf(V*)dmhf(V*)dVp+1τhf(V*)δV.E30

Notice that dmhf(V*)/dV in the numerator of the right-hand side can be directly calculated from Eq. (17), that is, it is given by

dmhf(V*)dV=19.78exp(V*+79.29.78){1+exp(V*+79.29.78)}2.E31

Exactly in the same way, a small variation δmhs is expressed by δV as follows:

δmhs=1τhs(V*)dmhs(V*)dVp+1τhs(V*)δV.E32

By substituting Eqs. (30) and (32) into Eq. (23), δIh is finally expressed by δV(t) as follows:

δIhδVg¯h{0.65mhf*+0.35mhs*}+g¯h{0.651τhf(V*)dmhf(V*)dVp+1τhf(V*)+0.351τhs(V*)dmhs(V*)dVp+1τhs(V*)}(V*Eh).E33

As δIh(t)/δV(t)represents admittance, Eq. (33) shows an equivalent admittance of h channel for a small variation δV and its admittance can be expressed by parallel coupling circuits of one conductance and two admittances. That is, the first term of Eq. (33) represents a conductance of h channel, which is expressed by the inverse of a pure resistance Rh; the second term is an admittance of a fast activation variable Yhf, which can be expressed by the inverse of series coupling of an inductance Lhf and a resistance Rhf, that is, Yhf = 1/(Rhf + pLhf); and the third term is an admittance element of a slow activation variable Yhs, which is also expressed by the inverse of series coupling of an inductance Lhs and a resistance Rhs, that is, Yhs = 1/(Rhs+ pLhs). Figure 9 shows an equivalent RLC circuit of h channel, where Rh, Rhf, Rhf, Rhs, and Lhs are given as follows:

Rh=1g¯h{0.65mhf*+0.35mhs*},E34
Rhf=10.65g¯hdmhf(V*)dV(V*Eh),Lhf=τhf*Rhf,E35
Rhs=10.35g¯hdmhs(V*)dV(V*Eh),Lhs=τhs*Rhs.E36

Figure 9.

An equivalent RLC circuit for h channel.

3.2.2. Equivalent admittance (impedance) of NaP channel

As in the case of with h channel, let V* be the equilibrium potential and INaP*be the NaP current at V*. From Eq. (15), INaP*satisfies the following relation:

INaP*=g¯NaPmNaP(V*)(V*EN).E37

When the membrane potential V(t) changes from V* to V* + δV(t), the current INaP (t) also changes from INaP*to INaP*+δINaP(t), where δINaP (t) is a small variation of NaP current caused by δV(t). Then, INaP*+δINaP(t)satisfies the following relation:

INaP*+δINaP(t)=g¯NaPmNaP(V*+δV(t))(V*+δV(t)EN).E38

Let mNaP(V* + δV) approximate by mNaP(V*) + δmNaP for a small variation δV. Then, Eq. (38) can be expressed by the following equation:

INaP*+δINaP(t)=g¯NaPmNaP(V*)(V*EN)+g¯NaPδmNaP(V*EN)+g¯NaPmNaP(V*)δV(t)+g¯NaPδmNaPδV(t).E39

By subtracting Eq. (37) from Eq. (39) and dropping the higher-order variation terms than the second, which appeared on the right-hand side of Eq. (39), the following equation is obtained:

δINaPg¯NaPδmNaP(V*EN)+g¯NaPmNaP*δV.E40

By following the same procedure from Eq. (24) to Eq. (26) except that τNaP is constant (0.15 ms), a small variation δmNaP can be expressed as follows:

mNaPdt10.15{δmNaPδmNaP}.E41

By approximating a small variation δmNaP by [dmNaP (V*)/dV] ⋅ δV, Eq. (41) becomes

mNaPdt10.15{dmNaP(V*)dVδVδmNaP}.E42

By using the differential operator p, a small variation δmNaP can be expressed by δV as follows:

δmNaP=10.15dmNaP(V*)dVp+10.15δV.E43

By substituting Eq. (43) into Eq. (40), δINaP is finally expressed by δV(t) as follows:

δINaPδVg¯NaPmNaP(V*)+g¯NaP10.15dmNaP(V*)dVp+10.15(V*EN).E44

Eq. (44) shows an equivalent admittance of NaP channel for a small variation δV, and it can be expressed by parallel coupling circuits of one conductance and one admittance. That is, the first term of Eq. (44) represents a conductance of NaP channel, which is expressed by the inverse of a pure resistance RNaP, and the second term is an admittance of an activation variable YNaP, which is expressed by the inverse of series coupling of an inductance LP and a resistance RP, that is, YNaP = 1/(RP+ pLP). Figure 10 shows an equivalent RLC circuit of NaP channel, where RNaP, RP, and LP are given as follows:

RNaP=1g¯NaPmNaP(V*),E45
RP=110.15g¯NaPdmNaP(V*)dV(V*EN),LP=0.15RP,E46

Figure 10.

An equivalent RLC circuit for NaP channel.

Figure 11.

An equivalent RLC circuit for a compartment model with h channel and NaP channel.

3.2.3. Equivalent RLC circuit of a neuron with h channel and NaP channel

By combining the results of Sections 3.2.1 and 3.2.2, an equivalent RLC circuit for a neuron model with h channel and NaP channel is obtained. Figure 11 shows its equivalent RLC circuit. In this section, we show some simulation results for a compartment neuron model (Figure 8) and its equivalent RLC circuit(Figure 11). A chirp current given to both a compartment neuron model and its equivalent RLC circuit is described as follows:

Iinp=Ainpsin(ω(t)t)+id,ω(t)=2πftT,E47

where the angular frequency ω(t) increases from 0 to 2πf over the period [0, T]. id is a DC bias current, which is set to zero in this subsection. The following parameter values were used in simulations; C = 1.5 μF/cm2, g¯L=0.15mS/cm2, EL = −65 mV, g¯NaP=0.5mS/cm2, g¯h=1.5mS/cm2, EN = 55 mV, EK = −90 mV, Eh = −20mV.

Figure 12.

Membrane potential and the magnitude and the phase of its FFT. Simulation result for (a) a compartment model (Figure 8) and (b) its equivalent RLC circuit (Figure 11).

Figure 12(a) shows one simulation result for a compartment model with h channel and NaP channel. The membrane potential V and the amplitude and the phase of its FFT are shown. As the magnitude of FFT shows, this neuron model has a band-pass property. Figure 12(b) shows the simulation result for its equivalent RLC circuit. Comparing Figure 12(a) and (b), the membrane potentials and their FFTs are exactly similar. This fact indicates that the derived RLC circuit represents almost the same properties of a compartment neuron model within a small variation of the membrane potential. In other words, voltage-dependent h channel and NaP channel may surely have inductive properties and contribute to the subthreshold resonance phenomena.

3.3. Properties of voltage responses to oscillatory current inputs

A neuron can generate action potentials whenever its membrane potential exceeds the threshold except during a refractory period. If the membrane potential of a neuron stays in a subthreshold level, a neuron cannot generate action potentials. However, as described in the previous section, many neurons in the brain have the subthreshold resonance properties. This fact indicates that a neuron may be able to generate an action potential when AC inputs whose frequencies are close to the resonance frequency of a neuron are given, because the resulting membrane potential for that input has the potential to exceed the threshold by the effects of the subthreshold resonance property. Actually, this fact has been observed in neurons of the brain. For example, Hutcheon et al. studied subthreshold voltage responses to AC inputs in neurons from the sensorimotor cortex of rats [11]. Figure 13 shows one of their results. Cases 1–3 show effects of the input amplitude: (Case 1) if a chirp current with a small amplitude is given to a neuron, the membrane potential does not exceed the threshold and no action potentials are generated; (Case 2) if its input amplitude increases a little bit, the membrane potential becomes larger but it still stays in a subthreshold level and no action potentials are generated; and (Case 3) if its amplitude increases more, the membrane potential exceeds for AC inputs with frequencies close to a resonance frequency of a neuron. As a result, action potentials are generated around that frequency. Cases 4 and 5 show effects of a DC bias input: (Case 4) if a DC-bias in the input current increases a little bit from 0μA/cm2, the membrane potential cannot exceed the threshold and no action potentials are generated and (Case 5) if the value of a DC bias is raised more, the membrane potential can exceed the threshold and action potentials are generated.

Figure 13.

Experimental results: Subthreshold voltage responses for a chirp current input with different amplitudes and DC bias currents observed in sensory cortex of rats [11].

We consider here a compartment neuron model with h channel and NaP channel described in the previous section, into which the HH model is incorporated in order to generate action potentials. Figure 14 shows its integrated neuron model, and the following equation is obtained:

CdVdt=ILIhINaPINIK+Iinp.E48

Figure 14.

A compartment model with h channel and NaP channel into which the Hodgkin-Huxley model (the HH-model) is incorporated.

Eq. (48) is an extension form of Eq. (12), into which two currents IN and IK of the HH model are added. Detail dynamics of IL, Ih, and INaP are given by Eqs. (13)–(15) in Section 3.2. Dynamics of INa and IK are also given by Eqs. (6) and (7) in Section 2.2. In addition to parameter values shown in the previous section, the following values were used in simulations: g¯Na=52mS/cm2,g¯K=52mS/cm2, EN = 55 mV, and EK = −90 mV.

Figure 15(a) shows the membrane potential for a chirp current input with Ainp = 2.5 μA/cm2 and id = 0 μA/cm2. For this input, the membrane potential cannot exceed the threshold, and no action potentials are generated. However, if the amplitude of AC input (Ainp) increases, the situation changes. Figure 15(b) shows the membrane potentials for a chirp current input whose amplitude increases to 3.3 μA/cm2. In this case, the membrane potential exceeds the threshold for the input frequency close to the resonance frequency of this neuron model. As shown in Figure 12, it is 13 Hz in this neuron model. On the other hand, Figure 15(c) shows the membrane potential for a chirp input with Ainp = 2.5 μA/cm2 and increased id = 1 μA/cm2. Almost the same response as Figure 12(b) is obtained. However, resulting membrane potentials in Figures 13 and 15 are not completely identical, because neurons used in experiments by Hutcheon et al. and a neuron model used in this simulations are different. However, simulation results follow a similar tendency as experimental results by Hutcheon et al.[11].

Figure 15.

Simulated membrane potential for a chirp current input with different amplitude and DC bias current. (a) Response for Ainp = 2.5 μA/cm2 and id = 0 μA/cm2. (b) Response for Ainp = 3.3 μA/cm2 and id = 0 μA/cm2. (c) Response for Ainp = 2.5 μA/cm2 and id = 1 μA/cm2.

By computer simulations, effects of the amplitude of AC input, its frequency, and a DC-bias current on the membrane potential were studied [20]. Some results are shown here. Figure 16 shows the membrane potentials for the AC inputs with a single frequency, setting Ainp = 3.3 μA/cm2. Figure 16(a) shows the results for 2 Hz input frequency. In this case, no action potentials are generated, because the membrane potential cannot exceed the threshold at all. Figure 16(b) shows that action potentials are generated, because the input frequency is 13 Hz, which is close to the resonance frequency of this neuron model. Thus, the maximum amplitude of membrane potentials exceeds the threshold by the effect of the subthreshold resonance property. In the case of Figure 16(c), for 30 Hz input frequency, no action potentials are generated, because the membrane potential is less than the threshold. On the other hand, if a DC bias current input increases to 1μA/cm2, as Figure 16(d) shows, more action potentials are generated for 13 Hz input frequency input than the case of no DC bias input. This indicates the effect of a DC bias current on generation of action potentials. Comparing Figure 16(b) and (d), it is clarified that the more spikes are observed than the case of no DC bias current input. For the case of id = 1 μA/cm2,a neuron model can also generate an action potential for 30Hz input frequency, as shown in Figure 16(e). As no action potential is generated for the case of id = 0 μA/cm2, this is also caused by the effect of a DC bias current input.

Figure 16.

Simulated membrane potential for an AC input with a single frequency (Ainp is fixed to 2.5 μA/cm2). (a) Response for id = 0 μA/cm2 and f = 2 Hz. (b) Response for id = 0 μA/cm2 and f = 13 Hz. (c) Response for id = 0 μA/cm2 and f = 30 Hz. (d) Response for id = 1 μA/cm2 and f = 13 Hz. (e) Response for id = 1 μA/cm2 and f = 30 Hz.

4. Conclusion

Until now, it has been mainly proposed that the cell membrane and neurons are modeled by an RC circuit. However, from the fact that many neurons in various regions of the brain have band pass properties, a neuron should be modeled by an RLC circuit. In this chapter, an equivalent RLC circuit was developed for a neuron model with h and NaP channels, and it was clarified that the subthreshold resonance property of this neuron model comes from inductive properties of h and NaP channels, especially, h channel. Furthermore, by incorporating the Hodgkin-Huxley dynamics into this neuron model, we showed the relation between the subthreshold resonance oscillation and the generation of action potentials. By computer simulations, it was shown that the amplitude of an AC input and a DC bias current input strongly play a role in the generation of action potentials, coupled with AC input frequencies.

It is presumed that the subthreshold resonance phenomena may relate closely to various practical neuron activities and behaviors in our brain, such as the sensitivity to external noises in sensory system. So, it is very important and interesting to study firing patterns of a neuron or properties of burst oscillations by using computer simulations, in order to clarify the mechanisms of higher information processing in the brain.

Needless to say, mutual collaboration of experimental research and modeling research is of large significance, in order to create more sophisticated models based on the most recent findings from experiments, and in order to develop new experimental methods to verify the facts suggested from simulation results.

Acknowledgments

We thank Mr. Babak V. Ghaffari for preparing figures.

© 2016 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution 3.0 License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

How to cite and reference

Link to this chapter Copy to clipboard

Cite this chapter Copy to clipboard

Tatsuo Kitajima, Zonggang Feng and Azran Azhim (August 24th 2016). Simulation of Neural Behavior, Numerical Simulation - From Brain Imaging to Turbulent Flows, Ricardo Lopez-Ruiz, IntechOpen, DOI: 10.5772/64028. Available from:

chapter statistics

760total chapter downloads

More statistics for editors and authors

Login to your personal dashboard for more detailed statistics on your publications.

Access personal reporting

Related Content

This Book

Next chapter

Numerical Simulations of Dynamics Behaviour of the Action Potential of the Human Heart's Conduction System

By Beata Jackowska-Zduniak

Related Book

First chapter

Mechanical Models of Microtubules

By Slobodan Zdravković

We are IntechOpen, the world's leading publisher of Open Access books. Built by scientists, for scientists. Our readership spans scientists, professors, researchers, librarians, and students, as well as business professionals. We share our knowledge and peer-reveiwed research papers with libraries, scientific and engineering societies, and also work with corporate R&D departments and government entities.

More About Us