Modeling Excitable Cells with Memristors

: This paper presents an in-depth analysis of an excitable membrane of a biological system by proposing a novel approach that the cells of the excitable membrane can be modeled as the networks of memristors. We provide compelling evidence from the Chay neuron model that the state-independent mixed ion channel is a nonlinear resistor, while the state-dependent voltage-sensitive potassium ion channel and calcium-sensitive potassium ion channel function as generic memristors from the perspective of electrical circuit theory. The mechanisms that give rise to periodic oscillation, aperiodic (chaotic) oscillation, spikes, and bursting in an excitable cell are also analyzed via a small-signal model, a pole-zero diagram of admittance functions, local activity, the edge of chaos, and the Hopf bifurcation theorem. It is also proved that the zeros of the admittance functions are equivalent to the eigen values of the Jacobian matrix, and the presence of the positive real parts of the eigen values between the two bifurcation points lead to the generation of complicated electrical signals in an excitable membrane. The innovative concepts outlined in this paper pave the way for a deeper understanding of the dynamic behavior of excitable cells, offering potent tools for simulating and exploring the fundamental characteristics of biological neurons.


Introduction
The electrical activities of neurons are characterized by a diverse array of dynamic phenomena, such as action potential, oscillation, spike, chaos and bursting.Understanding these qualitative features are essential for unraveling the principles underlying neuronal excitability.The popular Hodgkin-Huxley (HH) model developed in 1952 consists of membrane voltage, potassium conductance, sodium conductance and leakage conductance, provide a framework for understanding the propagation of action potential based on the squid giant axon experiments [1].Recent analysis revealed that the potassium ionchannel and sodium ion-channel in the HH model, initially interpreted as time-varying conductances are in fact generic memristors from the perspective of electrical circuit theory [2][3][4][5].The HH model has spurred significant interests to design electrical circuit models and observe the experimental results in the wide varieties of complex system of the membrane potential, nervous system, barnacle giant muscle fiber, Purkinje fibers, solitary hair cells, auditory periphery, mini review of neuromorphic architectures and implementation, organic synapses for neuromorphic electronics, and photochromic compounds [6][7][8][9][10][11][12][13][14][15][16].Similarly, extensive researches have been conducted to observe the varieties of oscillations in pancreatic β-cells inspired by the HH model [17][18][19][20][21][22][23][24][25][26].The mathematical model of an excitable membrane in pancreatic β-cells consist of voltage-sensitive channels that allow the Na + and Ca 2+ to enter the cell and, voltage-sensitive K + channels and voltageinsensitive K+ channel which allow to leave K + ion and activate intracellular calcium ion respectively [27][28][29].Therefore, the outward current carried by K + ions pass through the voltage and calcium-sensitive channels, and inward current carried by Na + and Ca 2+ ions pass through the voltage-sensitive Na + and Ca 2+ channels.However, the above models consist of complicated nonlinear differential equations associated with membrane voltage.Later a modified model was presented by Chay [30], assuming the β-cells of the voltage-sensitive Na + conductance is almost inactive, and the inward current is almost exclusively carried by Ca 2+ ions through the voltage-sensitive Ca 2+ channel.Therefore, the assumption of a mixed effective conductance was formulated without affecting the results by expressing the total inward current in terms of a single mixed conductance gI, and reversal potential EI of the two functionally independent Na + and Ca 2+ channels.The model consists of three nonlinear differential equations in contrast to the other complicated models of an excitable membrane of pancreatic β-cell.Our studies in this paper typically focus on the simplified Chay neuron model of an excitable cell [30].
The scientific novelty of this study is to model the excitable cells using memristive theory.By characterizing the state-independent voltage-sensitive mixed ion channel gI as a nonlinear resistor, and the state-dependent voltage-sensitive potassium ion-channel gK,V and calcium-sensitive potassium ion-channel gK,Ca as time-invariant memristors in the Chay neuron model, this research introduces a novel approach to study the behavior of ion channels in excitable cells.This unique modeling framework extends memristive theory to the realm of neuroscience, opening up new avenues for investigating the complex dynamics of excitable cells and their role in neural information processing.Moreover, the study employs comprehensive analytical tools such as small signal equivalent circuit model, pole-zero diagrams, the local activity principle, edge of chaos theory, and Hopf bifurcation theorem with the goal of gaining deeper insights in to the dynamic behavior the excitable cells.By integrating these analytical tools, the study provides a comprehensive perspective on the dynamic behavior of excitable cells in the framework of memristive theory, potentially uncovering new insights and relationships that were previously unexplored.The contributions of the study have the potential to advance our understanding of excitable cell dynamics and their implications for neural function.
The paper is organized as follows.We introduce the Chay neuron model and its comparison analyses with HH, FitzHugh-Nagumo and Morris-Lecar(ML) models in section 2. Section 3 describes the pinched hysteresis fingerprints of ion-channel memristors.Section 4 presents Direct Current (DC) analysis of Memristive Chay neuron model.Section 5 provides the small-signal analysis.Section 6 explores the application of the local activity principle, edge of chaos theorem, and Hopf bifurcations in memristive Chay neuron.Finally, section 7 concludes the paper.

Chay Neuron Model of an Excitable Cell
Excitable cells are specialized cells in the body and neurons that are capable of generating electrical signals in response to certain stimuli.These cells are crucial for the functioning of various physiological processes, including nerve signaling, muscle contraction, and hormone release.Excitability in these cells is primarily due to the presence of specialized proteins called ion-channels in their cell membranes.These ion channels control the movement of ions such as sodium (Na + ), potassium (K + ), calcium (Ca 2+ ), and chloride (Cl -) across the membrane, leading to changes in the cell's membrane potential and the generation of electrical signals.The study of excitable cells encompasses a wide array of topics, and our primary aim is to present a unified model for both neuronal and secretory excitable membranes based on the Chay neuron model.The Chay neuron model, which focuses on a simplified representation of neuronal and secretory excitable membranes, aims to provide a unified framework for understanding the complex electrical activity observed in excitable cells.This model typically involves just three ordinary differential equations(ODEs) to capture the essential features of an excitable cell membrane.The model consists of (a) mixed ion-channel gI (b) the state-dependent voltagesensitive potassium ion-channel gK,V (c) calcium-sensitive potassium ion-channel gK,Ca and (d) leakage channel are described by the following differential equation: ,Ca L 1 V 50 20 0.1 20 Fig. 1(a) shows the typical circuit of Chay model with external current stimulus, denoted as I 1 .It consists of membrane potential V of capacitance Cm, potentials EI , EK and EL for mixed Na + -Ca 2+ ions, K + and leakage ions respectively.The conductances gI, gK,V, gK,Ca, and gL, represent the voltage-sensitive mixed ion-channel, voltage-sensitive potassium ion-channel, calcium-sensitive potassium ion-channel and leakage channel respectively.In the 1 Electrical model is not given in the Chay paper [30].We have designed a typical electrical circuit model following the differential equation of the membrane potential.The symbolic representation of the conductances and potentials are assumed in different notations compared to the original representation.Fig. 1(a) shows an electrical circuit model following the conventional assumption of HH model.upcoming session, we will provide a rigorous proof that the state independent mixed ionchannel functions as a nonlinear resistor.However, the commonly held belief regarding state-dependent ion-channels exhibiting time-varying conductances is found to be conceptually incorrect from the perspective of electrical circuit theory.Contrary to this conventional assumption, these ion-channels do not adhere to time-varying conductance principles.Instead, they align more accurately with the characteristics of time-invariant generic memristors from a circuit theoretic standpoint.A rigorous proof will be demonstrated in the subsequent section.The parameters for this model are summarized in Table 1 2 and list of abbreviations of the model parameters are illustrated in Appendix.
The comparison analyses of the HH model [1], FitzHugh-Nagumo model [31], ML model [8], and the Chay model [30] are sumarized in Table 2 along with their respective strengths and limitations.It is notable that each model possesses distinct advantages and drawbacks making them suitable for different research contexts and questions.The choice of model depends on the level of detail required, computational resources available, and specific phenomena under investigation.This study predominantly centers on the Chay neuron model of excitable cells. .As, we are assuming the value of membrane capacitance (Cm)=1mF/cm 2 , we use the unit of all the conductances of the ion-channels g=mS/cm 2 throughout this study, which is also the equivalent unit g* of original Chay model.Due to the periodic and dynamic nature of the conductance g * (g), it can also be considered as the "conductance periodic factor".
(b) Chay [30] We are proposing a framework that the cells of excitable membranes can be modeled as the networks of memristors.
Novel model of excitable cells to capture multiple neuronal states, such as action potentials, periodic oscillations, aperiodic oscillations, spikes and bursting patterns.

Limited validation in
experimental contexts and lack details for some applications.

Pinched Hysteresis Fingerprints of the Ion-Channel Memristor
A generic memristor driven by a current source or voltage source is a two-terminal electrical circuit element whose instantaneous current or voltage obeys a state-dependent Ohm's law.A generic memristor driven by a current source can be expressed as follows in terms of state n x : ( ) ( , ,..., ; ) where R(x) is the memristance of the memristor and depends on "n" (n≥1 ) states variables,.
Similarly, a voltage-controlled memristor is defined in terms of the memductance G(x) and the state variables 12 , ..., n x x x , as follows: ( ) 12 , ,..., Eqs. ( 2) and (3) play significant importances to distinguish the memristive and nonmemristive system [32][33].They provide evidence that the state independent voltagesensitive mixed ion-channel functions as a nonlinear resistor and, state dependent voltagesensitive potassium ion-channel and calcium-sensitive potassium ion-channel behaves as timeinvariant generic memristors.

Voltage-sensitive mixed ion-channel nonlinear resistor
The time varying voltage sensitive mixed ion-channel with input voltage vI and current iI in the second element (from left) in Fig. 1(a) is given by, ( , ) and the conductance of the voltage sensitive mixed ion channel is given by 3 ( , ) where  ∞  ℎ ∞ are computed using (1e) and (   is performed for the three different frequencies namely, f=100 Hz, 200 Hz, and 1 KHz, respectively.Fig. 2 shows the corresponding output nonlinear waveform on currente iI vs.
voltage vI plane for these frequencies, confirming that the mixed ion-channel exhibits the properties of a nonlinear resistor only.

Voltage-sensitive potassium ion-channel memristor
Let us define the voltage across the voltage-sensitive potassium ion-channel shown in third (from left) element in Fig. 1(a) is vK,V and current is iK,V , then and current entering to the channel is , , , () where the memductance is given by 4 , , () and the state equation describing the channel in terms of n can be simplified from 1(b) as, ( ) Note that (5b)-(5d) are identical to the voltage controlled generic memristor defined in (3a)-(3b) with first order differential equation.Hence, the time-varying conductance shown in Fig. 1(a) of voltage-sensitive potassium ion-channel is replaced with voltagesensitive potassium ion-channel memristor as shown in the third element (from left) in Fig.

1(b).
We observed the memristive fingerprint of the voltage-sensitive potassium ion-channel memristor by applying sinusoidal bipolar signal under different frequencies.This property asserts that beyond some frequency f*, the pinched hysteresis loops characterized by a memristor shrinks to a single-valued function through the origin as frequency f > f* tends to infinity.To verify this property, a sinusoidal voltage source vK,V(t)=100sin(2πft)mV is applied to the voltage-sensitive potassium ion-channel with frequencies f=100 KHz, 500 KHz, and 4 MHz respectively.As shown in Fig. 3, the zero crossing pinched hysteresis loops shrink as the frequencies increase and tend to a straight line at 4 MHz which confirms that the voltage-sensitive potassium ion-channel is a generic memristor.All of these pinched hysteresis loops exhibit the fingerprints of a memristor [33].

Calcium-sensitive potassium ion-channel memristor
Let us consider the input voltage of the calcium-sensitive potassium ion-channel, the fourth element (from left) in Fig. 1(a) is vK,Ca 4 and current is iK,Ca then the current entering to the channel is given by ,Ca ,Ca ,Ca (Ca) where and the memductance of the calcium-sensitive potassium channel is given by ,Ca ,Ca (Ca) 1 (6c) 4 Since the same potential EK is shared by the voltage-sensitive potassium ion-channel memristor and calcium-sensitive potassium ion-channel memristor, the voltage assumed V-EK = vK,V in (5a) and V-EK = vK,Ca in (6b) are identical.The voltages vK,V and vK,Ca are assumed to distinguish the input voltage applied to voltage-sensitive potassium ion-channel memristor and calcium-sensitive potassium ion-channel memristor, respectively.The state equation in terms of calcium concentration from 6(b) and (1c) is given by, ( ) ( )

DC analysis of Memristive Chay Model of an Excitable Cell
The primary objectives to analyze the DC behavior of the memristive Chay model is to identify its equilibrium points of the nonlinear equations.(mV) insights can be gained into the behavior of the excitable cell under different conditions, such as varying input stimuli or parameter values, and can be expressed as a function of current I as:

Small-Signal Circuit Model
The small-signal equivalent circuit is the linearized method to predict the response of electronic circuits when a small input signal is applied to an equilibrium point Q.The objective of this section is to analyze the small-signal response of voltage-sensitive mixed ion-channel nonlinear resistor, voltage-sensitive potassium ion-channel memristor and calcium-sensitive potassium ion-channel memristor using Taylor series expansion and Laplace transformation.

Small-signal circuit model of the mixed ion-channel nonlinear resistor
The small signal equivalent circuit of the mixed ion-channel nonlinear resistor at an equilibrium point QI 5 on the DC VI-II curve is derived as follows ( ) ( ) Applying Taylor series expansion to the voltage-sensitive mixed ion-channel nonlinear resistor defined in (8a)-(8b) at the DC operating point QI, we get 00 ( ) Where h.o.t denotes higher order terms and coefficeints can be computed as, () ( ) Linearize (8c) by neglecting the h.o.t.then,

( )
Taking the Laplace transform of (8f), we obtain 5 The equilibrium point QI at vi=VI is obtained by solving 4(b).
12 ( ) ( ) ( ) The admittance YI(s; QI) of the small-signal equivalent circuit of the voltage sensitive mixed ion-channel nonlinear resistor at the DC operating point QI is given by, where 1/ ( )    respectively.The explicit formulas for computing coefficient a12(QI) are given in Table 3 for readers' convenience.

Small-signal circuit model of the voltage-sensitive potassium ion-channel memristor
The small-signal circuit model of the voltage sensitive potassium ion-channel memristor at an equilibrium point QK,V 6 on the DC VK,V-IK,V curve is derived by defining 6 The equilibrium point QK,V at vK,V = VK,V is obtained from (5d) by solving f(n;VK,V) = 0 for n = nK,V.

Small-signal circuit model of the calcium-sensitive potassium ion-channel memristor
The small-signal circuit model of the calcium-sensitive potassium-channel memristor at an equilibrium point QK,Ca 7 in the DC VK,Ca-IK,Ca curve is derived by defining from (6a) in a Taylor series about the equilibrium point (Ca(QK.Ca), VCa(QK,Ca)), we obtain 7 The equilibrium point QK,Ca at vK,Ca = VK,Ca is obtained from (6d) by solving f(Ca;VK,Ca) = 0 for Ca = CaK,Ca.The explicit formula for Ca(VK,Ca) is given in Table 5.

Small-signal circuit model of the memristive Chay model
Let The corresponding range of local activity, edge of chaos 1 and edge of chaos 2 at equilibrium voltage Vm(Q) (resp.I) are also given in Fig. 13 for readers' convenience.We will cover the details of these regimes in the section on locally activity and edge of chaos.
The circuit element R1,I is obtained by calculating the small signal model of the voltagesensitive mixed ion-channel nonlinear resistor from Table 3 at equilibrium voltage , 50 20 0.1 20 1 0.07 1 0.07 1 ) ( ) ) Similarly, LK,V, R1K,V , and R2K,V are calculated from the small-signal equivalent circuit of the voltage sensitive potassium ion-channel memristor from Table 4 and LK,Ca, R1K,Ca , and R2K,Ca are calculated from the small signal equivalent circuit of the calcium-sensitive potassium ion channel memristor from Table 5 at equilibrium voltage Vm(Q) respectively.Note that VK,V+EK and VK,Ca+EK must be replaced by Vm(Q) in Table 4 and Table 5 by the small signal model of the voltage-sensitive potassium ion-channel memristor and calcium-sensitive potassium ion-channel memristor, respectively.
Small-signal circuit model of the mixed Ion Channel nonlinear resistor

  
Table 6: Explicit formulas for computing the coefficients of Y(s;Vm(Q)).

Frequency Response
A convenient way to find the total admittance Y(s; Vm(Q)) by recasting (11) into a rational function of the complex frequency variable s, is as follows: where the explicit formulas for computing the coefficients b3, b2, b1, b0, a2, a1, and a0 are summarized in Table 6.
Substituting si  = in (12a), we obtain the following small-signal admittance function at the equilibrium voltage Vm(Q): The corresponding real part Re Y (iω; Vm(Q)) and imaginary part Im Y (iω; Vm(Q)) from (12b) are given by,

Pole-zero diagram of the small-signal admittance function Y(s; Vm(Q)) and eigen values of the jacobian matrix
The location of the poles and zeros of the small signal admittance function Y(s; Vm(Q)) of ( 12a) is computed by factorizing it's denominator and numerators as

Local Activity, Edge of Chaos and Hopf-Bifurcation in Memristive Chay Model
Local Activity and edge of chaos are the powerful mathematical quantitative theories to predict whether the nonlinear system exhibits complexity or not.Local activity refers to a characteristic of nonlinear systems wherein infinitesimal fluctuations in energy are amplified, leading to the emergence of complex dynamical behavior in the system [34][35][36][37][38][39].
This section presents an extensive analysis of the memristive Chay model using the principle of local activity, edge of chaos and Hopf-bifurcation theorem to predict the mechanism of generating the complicated electrical signals in an excitable cell.

Locally active regime
The local activity theorem developed by Chua reveals that a nonlinear system must satisfy at least one of the following conditions, concerning its local transfer function about a given operating point in order to support the emergence of complexity [36].
is either a negative real number, or a complex number.
(iv) ReY(iω;Vm(Q))<0 for some ωϵ[-∞, +∞] In another words, the emergence of action potentials, oscillations, chaos, burstings or spikes in neurons are impossible unless the cells are locally active.Therefore, restricting the behavior of a nonlinear system to its local activity operating regime reduces the considerable time necessary to identify the complex phenomena, which may emerge across its physical medium as compared to a standard trial-and-error numerical μA), respectively.

Hopf-bifurcation
Hopf-bifurcation namely, super-critical and sub-critical bifurcations are local bifurcation phenomenon in which an equilibrium point changes its stability as the parameter of the nonlinear system changes under certain conditions.When an unstable equilibrium point surrounded by a stable limit cycle results to a super-critical Hopf bifurcation whereas a subcritical Hopf bifurcation refers to a qualitative change in the behavior of a system where a stable equilibrium point transitions to instability, giving rise to sustained oscillations or limit cycles as a parameter is varied.Our careful simulation at Hopf-bifurcation point 1 at Vm= =-48.763mV(resp.I=-66.671μA) 8 shows that stimulus current I should be chosen within very small edge of chaos domain 1, where the real part of the eigen values are negative, the result converges to DC equilibrium for any initial conditions.Likewise, I is selected within the bifurcation point 1, where the real part of eigen values are positive, the result converges to a stable limit cycle.Therefore, it follows from the bifurcation theory that bifurcation point 1 is a super-critical Hopf bifurcation.The super critical Hopf bifurcation point 1 and point 2 observed in this paper are just for the parameters listed in Table 2.The bifurcations phenomenon may vary for different parameters.Table 7 illustrates the computation of the potassium ion-channel activation n, calcium concentration Ca and eigen values (λ1, λ2 and λ3) as a function of the DC stimulus current I (resp.membrane potential Vm) at the DC equilibrium point Q.It is observed from Table 7 and Fig. 17

Figure 1 .
Figure 1.Typical Chay neuron model of an excitable cell [30].(a) Electrical circuit model, following conventional assumption as time varying conductances [1].(b) Equivalent memristive Chay model based on Chua's memristive theory [2-4].The potential ECa for Ca 2+ ion given in the rate of the calcium concentration in (1c) is not an external battery source and not shown in external Fig. 1(a) and Fig. 1(b), respectively.
Fig. 1(b).To verify the voltage-sensitive mixed ion-channel is a nonlinear resistor, an extensive numerical simulation for a sinusoidal input voltage source vI = 100sin(2πft) mV

Figure 2 .
Figure 2. Output waveform plotted on iI vs. vI plane when the input voltage vI = 100sin(2πft) mV is applied with three different frequencies, namely f = 100 Hz, 200 Hz, 1 KHz to the voltage-sensitive mixed ionchannel.The output nonlinear waveform observed in Fig. 2 for different frequencies confirm the mixed ion channel is a nonlinear resistor.

Fig. 4
Fig.4that, all the zero crossing pinched hysteresis loops shrink as the frequencies of the input signal increase and tend to a straight line for the frequency f=200 Hz.All of the pinched hysteresis fingerprint confirm that the calcium-sensitive potassium ion-channel is a generic memristor.

Figure 5 .
Figure 5. (a) Memristive DC Chay model at equilibrium voltage Vm.(b) DC V-I curve of mixed ionchannel nonlinear resistor at equilibrium voltage VI=Vm-EI.(c) DC V-I curve of voltage sensitive potassium ion-channel memristor at equilibrium voltage VK,V=Vm-EK.(d) DC V-I curve of calcium sensitive potassium ion-channel memristor at equilibrium voltage VK,Ca=Vm-EK.(e) DC V-I curve of leakage channel at equilibrium voltage VL=Vm-EL.(f)Plot of DC V-I curve of memristive Chay model at membrane voltage Vm.

)
the individual DC V-I curve of voltage sensitive mixed ion channel non-linear resistor, voltage sensitive potassium ion-channel memristor, calcium sensitive potassium ionchannel memristor and leakage channel at equilibrium voltage VI, VK,V, VK,Ca and VL as shown in the Fig. 5(b), Fig. 5(c), Fig. 5(d) and Fig. 5(e) respectively.Fig. 5(f) shows DC V-I curve of Fig. 5(a) over the range of DC voltage -50 mV <Vm <-24 mV.For any DC value of Vm, we calculated the corresponding value of I as the vertical axis.Our extensive calculations show that, the two Hopf bifurcation points occur at Vm=-48.763 mV (resp., I=-66.671µA) and Vm=-27.984mV (resp., I=433.594µA) respectively.Details of these two bifurcation points will be discussed in upcoming section

Figure 6 .
Figure 6.(a) Small-signal circuit model of the voltage-sensitive mixed ion-channel nonlinear resistor about the DC equilibrium point QI (VI, II).(b) Plot of the coefficient a12 and resistance R1,I as a function of the DC equilibrium voltage VI=Vm-EI where EI=100 mV.R1,I<0 over the range of local activity, edge of chaos 1 and edge of chaos 2 of the mixed ion channel nonlinear resistor is identified with respect to VI of the entire Chay circuit in Fig. 1(b) and Fig. 13 .

From
(8h), it is followed that the small-signal admittance function of the mixed ionchannel nonlinear resistor is equivalent to a linear resistor.The corresponding small-signal equivalent circuit and a plot of the coefficient a12(QI) and resistance R1,I as a function of the DC equilibrium voltage VI =Vm-EI where EI=100 mV are shown in Fig.6(a) and Fig.6(b),

Figure 7 .
Figure 7. Small-signal equivalent circuit model of the voltage-sensitive potassium ion-channel memristor about the DC equilibrium point QK,V (VK,V, IK,V).

Figure 8 .
Figure 8. Plot of coefficients (a) a11 (b) a12 (c) b11 and (d) b12 of the voltage-sensitive potassium ionchannel memristor as a function of the DC equilibrium voltage VK,V.

Figure 9 .
Figure 9. (a) Inductance LK,V (b) resistance R1K,V and (c) resistance R2K,V of the voltage-sensitive potassium ion-channel memristor as a function of DC equilibrium voltage VK,V =Vm-EK where EK=-75mV.LK,V>0, R1KV>0 and R2KV>0 shown in figures over the local activity, edge of chaos 1 and edge of chaos 2 are just corresponding range of the voltage with respect to VK,V of the entire connected Chay small signal equivalent circuit of Fig. 1(b) and Fig. 13.
=V mV where E =-75 mV = − + Local Activity Regime with respect to V K,V =V m -E K 26.455 mV<Vm<51.685mV Edge of Chaos domain 1 with respect to V K,V =V m -E K 26.455 mV<Vm<27.763mV Edge of Chaos domain 2 with respect to V K,V =V m -E K 48.984 mV <Vm<51.685mV Range of Edge of Chaos 2 of the entire Chay circuit in Fig. 1(b) and Fig. 13 with respect to V K,V Range of Edge of Chaos 1 of the entire Chay circuit in Fig. 1(b) and Fig. 13 with respect to V K,V Range of the Locally Active regime of the entire Chay circuit in Fig. 1(b) and Fig. 13 with respect to V K,

)
It follows from (10r)-(10t) that the small-signal admittance function of the first-order calcium-sensitive potassium ion-channel memristor is equivalent to the series connection of an inductor and a resistor in parallel with another resistor as shown in Fig.10.The corresponding coefficients a11, a12, b11, b12 and inductance LK,Ca, resistance R1K,Ca, and resistance R2K,Ca as a function of the DC equilibrium voltage VK,Ca are shown in Figs.11 and Figs. 12, respectively.The small-signal inductance and resistances( i.e.LK,Ca>0, R1K,Ca>0 and R2KCa>0) over the edge of chaos 1 and edge of chaos 2 with respect to the VK,Ca are shown in Fig. 12(a), Fig. 12(b) and Fig. 12(c) respectively.Please note that the local activity, edge of chaos 1 and edge of chaos 2 shown in Fig. 12(a), Fig. 12(b) and Fig. 12(c) are not the local activity, edge of chaos 1 and edge of chaos 2 of the individual calcium sensitive potassium ion channel memristor.The local activity, edge of chaos domains are just an information showing the corresponding range of voltage with respect to VK,Ca when measured across the individual calcium sensitive potassium ion channel memristor of the entire connected Chay small-signal equivalent circuit of Fig. 1(b) and Fig. 13.For the readers' convenience, the explicit formulas for computing the coefficients a11(Q,K,Ca), a12(QK,Ca), b11(Q,K,Ca), b12(Q,K,Ca)

Figure 10 .
Figure 10.Small-signal equivalent circuit model of the calcium-sensitive potassium ion-channel memristor about the DC equilibrium point QK,Ca (VK,Ca, IK,Ca).

Figure 11 .
Figure 11.Plot of coefficients (a) a11 (b) a12 (c) b11 and (d) b12 of the calcium-sensitive potassium ionchannel memristor as a function of the DC equilibrium voltage VK,Ca.

Figure 12 .
Figure 12.(a) Inductance LK,Ca (b) resistance R1K,Ca and (c) resistance R2K,Ca of the calcium-sensitive potassium ion-channel memristor as a function of DC equilibrium voltage VK,Ca.LK,Ca>0, R1K,Ca>0 and R2KCa>0 over the edge of chaos 1 and edge of chaos 2 with respect to VK,Ca of the entire connected Chay small-signal equivalent circuit of Fig.1(b) and Fig. 13. .
of Chaos 1 of the entire Chay circuit in Fig. 1(b) and Fig. 13 with respect to V K,Ca Range of Edge of Chaos 2 of the entire Chay circuit in Fig. 1(b) and Fig. 13 with respect to V K,Ca Range of the Locally Active regime of the entire Chay circuit in Fig. 1(b) and Fig. 13 with respect to V K,Ca us replace the voltage-sensitive mixed ion-channel nonlinear resistor, the voltagesensitive potassium ion-channel memristor, and the calcium-sensitive potassium ionchannel memristor in the memristive Chay neuron circuit of Fig. 1(b) with their smallsignal models about DC operating voltages VI=Vm-EI, VK,V=Vm-EK, and VK,Ca=Vm-EK, respectively.Short-circuiting all the batteries, the equivalent small-signal circuit model of the third-order neuron circuit from Fig. 1(b) about the operating point Vm(Q is found to be composed of one capacitor, two inductors, and six resistors as shown in Fig. 13.The local admittance Y(s;Vm(Q)) of this linear circuit seen from the port and formed by the capacitor terminals about Q is given by

Figure 13 .
Figure 13.Small-signal equivalent circuit model of the memristive Chay model.The DC equilibrium voltage Vm is computed at Vm=VI+EI for mixed ion channel non-linear resistor, Vm=VK,V+EK for voltage sensitive potassium ion-channel memristor and Vm=VK,Ca+EK for calcium sensitive potassium ionchannel memristor, respectively.

13 )
Fig. 15(b) that the two poles Re(p1), Re(p2) are negative while Im(p1), Im(p2) remain consistently zero for the specified DC input Vm.This observation confirms that the two poles of the admittance function possess no complex frequencies.Fig.16(a) shows the Nyquist plot, i.e. loci of the imaginary part Im(zi) versus the real part Re(zi) of the zeros as a function of the input voltage Vm over the interval −55 mV ≤ Vm ≤ 25mV.Observe that the real part of the two zeros z2 and z3 are zero at Vm= -48.763 mV(resp., I= -66.671 μA) and Vm=-27.984mV(resp., I=433.594μA), respectively.The corresponding points when Re(zi)=0 are known as Hopf bifurcation points in bifurcation theory.Fig. 16(b) and Fig. 16(c) show the zoomed version of Fig. 16(a) near to the two bifurcation points respectively.It is also observed that the Re(z2) and Re(z3) lie in open right half plane(RHP) between the bifurcation points -48.763 mV <Vm<-27.984mV(resp.-66.671 μA <I < 433.594 μA).Observe from Fig. 17 that the eigenvalues, computed from the Jacobian matrix, associated to the ODE (1a)-(1c) are identical to the zeros of the neuron local admittance Y(s; Vm(Q)), as inferable from Fig. 16, and expected from the Chua theory [3]-[4].

Figure 15 .Figure 16 .
Figure 15.Poles diagram of the small-signal admittance function Y(s; Vm(Q)) as a function of Vm over -200 mV<Vm<200 mV (a) Top and bottom figures are the plot of the real part of the pole 1 Re(p1) and Imaginary part of pole 1 Im(p1) respectively.(b) Top and bottom figures are the plot of the real part of the pole 2 Re(p2) and Imaginary part of pole 2 Im(p2) respectively.

Figure 17 .
Figure 17.Plot of the loci of the eigen values of the Jacobian Matrix (a) Nyquist plot of the eigen values λ1, λ2, λ3 in Im(λi) vs. Re(λi) plane Nyquist plot near the Hopf-bifurcation point 1, Vm=-48.763 mV(resp., I= -66.671 μA).(c) Nyquist plot near the Hopf-bifurcation point 2, Vm=-27.984mV(resp., I=433.594μA).Our numerical computations confirm the zeros of the admittance functions Y(s; Vm(Q)) obtained in Fig. 16 are identical to the eigen values of the Jacobian matrix (i) The zero of the admittance function Y(s; Vm(Q)) lie in open-right plane where Re(sz)>0 (ii) Y(s; Vm(Q)) has multiple zero on the imaginary axis (iii) Y(s; Vm(Q)) has simple zero on the imaginary axis s = iωz on the imaginary axis and

Figure 18 .
Fig. 18(f) that Re Y (iω; Vm(Q))>0, and the memrisitve Chay model is no more locally active confirming the cell is locally passive at this equilibrium.Therefore, the local activity Edge of chaos is a tiny subset of the locally-active domain where the zeros of the admittance function Y(s; Vm(Q))(equivalent to the eigen values of Jacobian matrix) lie in the open left-half plane, i.e.Re(zp)<0( eigen values λi<0) as well as ReY(iω; Vm(Q))<0.Fig.17(a) and Fig. 17 (b)  show the real part of the eigen values vanish at Vm=-48.7631 mV (resp. (c) that the real part of the eigen values vanish at λ2,3= ± 85.606i at DC equilibrium voltage Vm=-27.984mV(resp.I= 433.594 μA).It follows that the corresponding equilibrium point Vm(Q) is no longer asymptotically stable below this equilibrium point, therefore confirming the existence of a 2nd edge of chaos regime over the following interval: The nonlinear dynamical behavior of the memristive Chay model in this paper is controlled as the function the input stimulus I. The local activity, edge of chaos 1 and edge of chaos 2 regime computed in this paper under the assumption of departing the input parameter I from lower stimulus to higher stimulus (resp.low DC equilibrium voltage Vm(Q) to high equilibrium voltage Vm(Q)).

Fig. 19 (
Fig. 19(a) and Fig. 19(b) show the numerical simulations at I= -68.118 μA and I= -65.077 μA respectively.Observe, from Fig. 19(a) and Fig. 19(b) that I= -68.118 μA lying within the tiny subset of edge of chaos domain 1 converges to DC equilibrium and I= -65.077 μA lying in open right half-plane(RHP) converges to a spikes, respectively, confirming the bifurcation point 1 is a super-critical Hopf bifurcation.

Figure 19 .
Figure 19.Numerical simulations to confirm the super-critical Hopf bifurcation at bifurcation point 1.Plot of membrane potential Vm at (a) I= -68.118 μA which lies inside the tiny subset of edge of chaos domain 1 and beyond bifurcation point 1 converges to the DC equilibrium, (b) and I= -65.077 μA, chosen just to the right of bifurcation point 1, where the real parts of two zeros of the neuron local admittance lie on the open right half plane (RHP) converges to the spikes

Figure 20 .
Figure 20.Numerical simulations to confirm the super-critical Hopf bifurcation at bifurcation point 2. (a) Plot of membrane potential Vm which converges to stable DC equilibrium when I=440 µA chosen inside the tiny subset of edge of chaos domain 2 and, near and beyond the bifurcation point 2. (b) Membrane potential converging to oscillation as predicted by Hopf bifurcation theorem when I=430.884µA is chosen inside the bifurcation point (open right-half pane).
(a)  to Fig.17 (c) that the two Hopf bifurcations points 1 and 2 occur at Vm= =-48.763mV(resp.I=-66.671μA) and Vm= -27.984 mV (resp.I= 433.594 μA) respectively, where the eigen values are purely imaginary at these two equilibria.As I decreases (resp.Vm decreases) from the Hopf bifurcation point 1, the eigen values migrated to the left-hand side confirming the real parts of the eigen values are no longer positive and thereby confirming the first negative real eigen values regime exists over the following interval.
Similarly, as I increases(resp.Vm increases) from the second bifurcation points, the positive real part of the eigen values migrated from open right half to the open left half, there by confirming the second negative real eigen values regime over the following

Figs. 21 .
Figs. 21.Fig. 21(a) shows the transient waveform of membrane potential Vm converging to a stable DC equilibrium at I= -90 µA, confirming the Hopf bifurcation theorem no longer holds at this equilibrium.Similarly, when DC simulus currents I= -50 µA and -10µA are chosen inside the two bifurcation points I=-66.671μA and I= 433.594 µA, we observed different patterns of oscillations as shown in Fig. 21(b) and Fig. 21(c), confirming the bifurcation theorem holds in this regime.Likewise, when DC stimulus currents I=10 μA and I=2000μA are applied within the bifurcation points I=−66.671μA and I=433.594, respectively, oscillation patterns emerge as depicted in Fig. 22(a) and Fig. 22(b).Similarly, Fig. 22(c) illustrates the transient waveform of the membrane potential Vm, indicating its convergence to a stable DC equilibrium at I=500 μA.This observation suggests that the Hopf bifurcation theorem no longer holds at this equilibrium point.

Figure 21 .Figure 22 .
Figure 21.Patterns of oscillations when stimulus current I is chosen beyond and inside the bifurcations points.(a) DC pattern observed when I= -90 µA chosen beyond bifurcation point 1(I=-66.671μA).Different patterns of oscillatiions when I is chosen between the two bifurcation points I=-66.671μA and I= 433.594 μA, at (b ) I= -50 μA, (c) I= -10 μA.
with the predictions of the Hopf bifurcation theorem, we observed the presence of an oscillating regime between two bifurcation points within the voltage range of -48.763 mV to -27.984 mV.Our numerical simulations confirmed the super-critical Hopf bifurcation with complex conjugates of eigenvalues coincide on the purely imaginary axis at ±0.557i and ±85.606i respectively.It was also observed that a tiny change in external stimulus current I in excitable cells, far from the bifurcation points no longer holds the Hopf bifurcation theorem as it crosses the imaginary axis from right to left confirming that the real part of the eigenvalues becomes negative and converges to a DC equilibrium point.Our comprehensive comparison of the HH, FitzHugh-Nagumo, ML, and the Chay models presented in Table2along with their individual strengths and limitations reveals distinct advantages and drawbacks making them suitable for different research contexts and questions.The selection of the particular model depends on the specific objectives.We primarily focused to advance the understanding of excitable cells by modeling with the networks of memristors and predicting their responses with the concept of memristor theory, DC steady state analyses, small signal equivalent circuit, local activity principle, edge of chaos theorem and hopf bifurcations.In Conclusion, the theoretical framework outlined in this paper confirms the significance of memristors in simulating action potentials in excitable cells and also establishes a foundation for their application in neuron modeling, artificial intelligence, and brain-like machine interfaces.Our proposed model offers potential for enhancing adaptive neural networks, neuroprosthetics, neuromorphic computing architectures, and the broader scope of artificial intelligence, thereby aiding in the development of brain-like information processing systems.Appendix: Abbreviations of the Model ParametersC m =Membrane Capacitance E K =Potential across K + ion channel memristor E I =Potential across mixed ion channel memristor E L =Potential across leakage channel E Ca =Potential across Ca +2 ion channel memristor g K,V =Voltage-sensitive K+ ion-channel conductance g I =Voltage-sensitive mixed ion channel conductance g L =Leakage channel conductance g KCa =Calcium activated potassium conductance k Ca =Rate constant for the efflux of the intracellular Ca +2 ions ρ =Proportionality constant λ n =Rate constant for k + ion-channel opening m ∞ = Probability of activation of the mixed ion channel in steady state α m = The rate at which the activation of the mixed ion channel closed gates transition to an open state(s -1 ) β m = The rate at which the activation of the mixed ion channel open gates transition to the close state(s -1 ) h ∞ = Probability of inactivation of the mixed ion channel in steady state α h =The rate at which the inactivation of the mixed ion channel closed gates transition to an open state(s -1 ) β h =The rate at which the inactivation of the mixed ion channel open gates transition to the close state(s -1 ) n=Probability of n opening of the K+ ion channel memristor n ∞ =Steady state value of n α n =The rate at which K + ion channel closed gates transition to an open state(s -1 ) β n =The rate at which K + ion channel opened gates transition to an close state(s -1 )

Table 1 .
Parameters of the Chay neuron model of an exctiable cell

Table 2 .
Comparision analyses of HH, FitzHugh-Nagumo, ML and Chay models [3]nnel in the HH model are represented with generic memristors[3]It is a framework to understand the emergence of action potential propagation in neuron based on the experimental data of squid giant axon.

Table 3 .
Explicit formulas for computing the coefficients a12(QI) of the voltage-sensitive mixed ion channel nonlinear resistor
I=-66.671 μA) with pair of complex eigen values λ2,3= ± 0.557i.It follows from the edge of chaos theorem that the corresponding equilibrium point is no longer asymptotically stable, and becomes unstable thereafter confirming the 1st edge of chaos regime over the following small interval:

Table 7 .
Computation of the potassium ion-channel activation n, calcium concentration Ca and eigen values (λ1, λ2 and λ3) as a function of the stimulus current I (resp.membrane potential Vm).Rows 5 to 7 pertain to the edge of chaos 1, rows 8 to 16 pertain to the unstable local activity domain and rows 17 to 20 pertain to the edge of chaos 2. Rows 7 and 17 pertain to Hopf bifurcation point 1 and Hopf bifurcations points 2, respectively for the memristive Chay neuron model.
Observe from Table7and Fig.17(a)-Fig.17(c) that two eigenvalues of the Jacobian matrix associated to the ODE set (1a)-(1c) lie on the open RHP for each operating point Q corresponding to a DC current I value between Hopf bifurcation point 1 and Hopf bifurcation point 2. Therefore, the generation of periodic, bursting, spikes and chaos signals