- Research
- Open Access
- Published:

# Temperature-dependent bursting pattern analysis by modified Plant model

*Molecular Brain*
**volume 7**, Article number: 50 (2014)

## Abstract

Many electrophysiological properties of neuron including firing rates and rhythmicaloscillation change in response to a temperature variation, but the mechanismunderlying these correlations remains unverified. In this study, we analyzed variousaction potential (AP) parameters of bursting pacemaker neurons in the abdominalganglion of *Aplysia juliana* to examine whether or not bursting patterns arealtered in response to temperature change. Here we found that the inter-burstinterval, burst duration, and number of spike during burst decreased as temperatureincreased. On the other hand, the numbers of bursts per minute and numbers of spikesper minute increased and then decreased, but interspike interval during burst firstlydecreased and then increased. We also tested the reproducibility oftemperature-dependent changes in bursting patterns and AP parameters. Finally weperformed computational simulations of these phenomena by using a modified Plantmodel composed of equations with temperature-dependent scaling factors tomathematically clarify the temperature-dependent changes of bursting patterns inburst-firing neurons. Taken together, we found that the modified Plant model couldtrace the ionic mechanism underlying the temperature-dependent change in burstingpattern from experiments with bursting pacemaker neurons in the abdominal ganglia of*Aplysia juliana*.

## Introduction

To date, researchers have investigated the effect of temperature on the electricalactivity and firing patterns in neurons from many animals, including *Aplysiajuliana*, crabs, crayfish, frogs, lobsters, locusts, snails, and squids [[1]-[9]]. Especially the functions and properties of bursting pacemaker neuron R15 inthe abdominal ganglion of *Aplysia* have been extensively studied [[10]-[12]] and mathematical simulation of bursting activity has been successfullyperformed [[13]-[19]]. However, few studies on the temperature dependence of action potential (AP)parameters in the R15 bursting pacemaker neuron have been reported [[20],[21]]. They analyzed typical changes of AP parameters in burst-firing neurons byinvestigating the effect of heat on R15 bursting pacemaker neuron activity as thetemperature increased and reported temperature-dependent changes in inter-burstinterval, burst duration, number of spike per burst, intra-burst spike broadening andspike height [[21]]. However, the reproducible properties and mechanism of temperature-dependentchanges of AP parameters remain unknown yet.

Changes in temperature can produce numerous effects on the neural tissue of mostorganisms. Indeed it has been reported that hyperthermic temperature may inducedepolarization and spontaneous firing of pyramidal neurons leading to enhancedexcitability of hippocampus [[22],[23]]. On the other hand, a body of evidence indicating a therapeutic effect ofhypothermia has accumulated in several conditions. Orlowski and colleagues publishedthat in refractory status epilepticus unresponsive to conventional treatment, systemichypothermia (30-31°C) was an effective therapeutic method leading to burstsuppression on electroencephalography (EEG) [[24]]. Despite a lot of works aiming at demonstrating the neural effect oftemperature, the mechanism of temperature-dependent electrophysiological changecontaining alteration of bursting and firing pattern in neuron remains unclear. Here, weexamined the effects of temperature changes on the neuronal activity and burstingpatterns during several consecutive heating − cooling cycles by usingbursting pacemaker neurons which are a proper specimen with capability of long-lastingrecording for mathematical modeling. Next, we sought to identify the mechanismunderlying temperature-dependent bursting patterns of these neurons by analyzing andcomparing the experimental data to computational simulation data calculated by modifiedPlant equations with temperature-dependent scaling factors, ρ(T) and ϕ(T).

## Materials and methods

### Animals and dissection

Animals (*A. juliana*) were collected locally in Seogwipo, South Korea and allexperimental procedures were approved by the Jeju National University Animal Care andUse Committee. These animals were dissected to observe the temperature-dependency ofbursting patterns generated by pacemaker neurons R15 in the abdominal ganglia of*A. juliana* as described in previously published paper [[1]]. Briefly, the animals were anesthetized with an injection of 0.38 MMgCl_{2} amounting to half of each specimen’s weight before theabdominal ganglia were removed. Each abdominal ganglion was incubated at 34°Cfor 40 minute in a solution containing equal volumes of artificial sea water(ASW; in mM: 460 NaCl, 10 KCl, 11 CaCl_{2} , 55 MgCl_{2} , and 10HEPES; pH 7.6) and isotonic Leibovitz’s L-15 media (Cat. No. L-4386;Sigma) containing 1% protease (type IX; Sigma) (hereafter, ASW:L-15(1:1)). They werethen washed several times with ASW and were placed in a low-temperature incubator(VS-1203PIN; Hanback Co., Daejeon, Korea) at 18°C. At last they were removedonto Petri dishes (50x9 mm), pinned down on Sylgard plates (Dow Corning, USA) filledwith ASW:L-15(1:1), and desheathed.

### Data acquisition and analysis

A PT100 platinum resistance temperature sensor connected to a digital thermometer(TRM-006, Toho, Japan) was placed near the abdominal ganglia of *Aplysia*soaked in ASW:L-15 (1:1, v/v) to measure the medium temperature, which could beincreased or decreased by activating the temperature controller system (HMN 3940,Acetec Co., Korea). We used a low flow rate peristaltic pump (BJ100-2 J; BaodingLonger Precision Pump Co., China) to maintain a viable solution with ASW:L-15 (1:1)media. The flow rate and the speed were 0.14 mL/min and 1.0 rpm,respectively. Intracellular recordings of bursting pacemaker neurons were performedto measure the membrane potential. The glass intracellular electrodes were filledwith 3 M KCl and the membrane potentials recorded were simultaneously savedusing a DAQ card (NI PCI-6221, National Instruments) and the Labview program(National Instruments). The electrical signals of the APs were identified using adigital oscilloscope (54622A; Agilent, Colorado Springs, Colo., USA). These data wererecorded at a rate of 3 kHz to reproduce the APs on the computer because thebursting cells generated APs at a rate of 0–4 Hz; each component of thisdata is composed of 180,000 pairs of temperature and membrane potential values withunique file names, and the dataset was composed of 1,713 files as shown inTable 1. With this selected dataset, the average valuesof each of the AP parameters of the burst-firing neuron could be calculated usingOrigin 6.0 (Microcal Software, Inc.) and a computer program we designed. For eachexperiment, the bursting neuron was held at room temperature for one hour as shown inTable 2 and then changed with continuous ramp fromabout 16°C to 30°C. Then, the data recorded by increments (or decrements)of 2°C from about 16°C to 30°C (16-18-20-22-24-26-28-30°C) wereaveraged over. To reduce experimental variability, in each experiment the burstingneuron was incubated at least 5 min to adapt to a new steady-state temperaturebefore recording rhythmic change.

## Results

### Definition of AP parameters in burst-firing neuron

The graph in Figure 1A represents the typical burstingpatterns of electrical signals shown in *Aplysia* bursting pacemaker neuron.The intraburst interspike interval, interburst interval, and burst duration arerepresented by symbols (1), (2), and (3), respectively. The membrane potentials atthe positive peak, V_{pp} (mV), and negative peak, V_{np} (mV), aredefined as the values of membrane potentials at points T and B shown inFigure 1B, respectively. The definition of the firsthalf of the rising phase of AP, Δt_{r1} (ms), and the following APparameters are very similar to those defined in our previous study [[1]]. The last half of the rising phase of AP, Δt_{r2} (ms) andthe first half of the falling phase of AP, Δt_{f1} (ms) are timeintervals shown in Figure 1B. The last half of the fallingphase Δt_{f2} (ms) is defined as the values of the time interval betweentwo points H2 and B shown in Figure 1B. The AP half-widthduration, Δt_{AP, 1/2} (ms), is defined asΔt_{r2} + Δt_{f1} (ms). The interspikeinterval, ISI (ms), is defined as Δt_{r1} + Δt_{AP,1/2} + Δt_{f2}. The spontaneous firing frequency,simply referred to as Frequency (s^{−1}), is defined asISI^{−1} (s^{−1}). For convenience, the average anglesθ_{1} and θ_{2} are defined to compare the slope of thesecond half of the rising phase of the AP with that of the first half of the fallingphase. The angle θ_{1} is defined as the inverse tangent of the ratio ofhalf A_{AP} to Δt_{r2}, and the angle θ_{2} is alsodefined as the inverse tangent of the ratio of half A_{AP} toΔt_{f1}: ${\mathrm{\theta}}_{1}\equiv \mathrm{Arctan}\phantom{\rule{0.12em}{0ex}}\left[\frac{\mathrm{V}\left({\mathrm{t}}_{\mathrm{T}}\right)\u2010\mathrm{V}\left({\mathrm{t}}_{\mathrm{H}1}\right)}{{\mathrm{t}}_{\mathrm{T}}\u2010{\mathrm{t}}_{\mathrm{H}1}}\right]$ and ${\mathrm{\theta}}_{2}\equiv \mathrm{Arctan}\phantom{\rule{0.12em}{0ex}}\left[\frac{\mathrm{V}\left({\mathrm{t}}_{\mathrm{T}}\right)\u2010\mathrm{V}\left({\mathrm{t}}_{\mathrm{H}2}\right)}{{\mathrm{t}}_{\mathrm{H}2}\u2010{\mathrm{t}}_{\mathrm{T}}}\right]$. Although angles are usually defined as dimensionless quantities,these angles are defined differently for convenience. The AP amplitude,A_{AP} (mV), is defined as the absolute value of difference between thetwo values of membrane potentials at points T and B, respectively. More detailedanalytical techniques regarding the definition and computational analysis of APparameters and others were described previously [[1],[25]]. Data were processed using the scientific data analysis tools,C^{++}, Origin 6.0, and Mathematica 5.1, and the resulting data werepresented as mean ± SE (standard error).

### Structure of the selected dataset for AP parameter analysis

Eight experiments (A − H) were conducted using *A. juliana*specimens, and Table 1 shows the animal weights, selecteddata for analysis, total number of spike, and total recording time.

Animal weights were between 130 and 310 g, and the average weight was226 g. The average total number of spike and average total recording time ofeach experiment were 15,039 spikes and 607 min, respectively. It was necessaryto select the middle portion of each recorded dataset for analysis. The averagevalues of time intervals for analysis, numbers of cycles of temperature change, andnumber of spike selected for analysis of these selected data were 214 min,2 cycles, and 5,433 spikes, respectively. The dataset gathered for analysis wascomposed of these selected data from each experiment. The total time interval foranalysis and the total number of spike of this dataset were 1,713 minutes and43,468 spikes, respectively. The numbers written in parentheses in the middle columnof Table 1 represent the total number of cycles oftemperature change for each experiment.

In each experiment, the temperature was maintained at room temperature for the firstone hour. As shown in Table 2, average values oftemperature and four AP parameters of burst-firing neurons were calculated using theselected dataset saved by intracellular recordings during maintaining at roomtemperature for 30 minutes; all data were selected after 5 to 10 minutesfrom temperature change onset. Because the standard errors of the values of theseparameters were low, the bursting patterns of these experiments might be consideredas regular while the temperatures were held constant at room temperature.

### Temperature dependence of AP parameters in burst-firing neurons

We selected burst trains in a continuous time series from the intracellular recordingdata of experiment A during the fourth falling phase (from 440 to 487 min)(Figure 2) out of the eight falling phases intemperature recorded during eight consecutive heating − coolingcycles: it was composed of 48 panels. Temperature and membrane potentials arerepresented as the upper and lower traces in Figure 2,respectively. These figures demonstrated that interburst interval, burst duration,and number of spike during burst decreased as temperature increased. Next, weextracted six burst signals (6 panels) among many burst signals to examine thetemporal change of instantaneous ISI during burst at specific temperature;30.0°C (at 440 min), 26.9°C (at 454 min), 23°C (at464 min), 21.3°C (at 468 min), 19.1°C (at 474 min), and16.4°C (487 min). As shown in Figure 3, ISIduring burst versus time displayed parabolic pattern, suggesting thattemperature-dependent bursting patterns are parabolic bursts characterized by lowerspike frequency at the beginning and end of the burst. The analyzing dataset wascomposed of 1,713 files corresponding to the sum of each selected time interval foranalysis as shown in the 4th column of the Table 1.

With the selected data shown in Table 1, thetemperature-dependent properties of AP parameters were analyzed by using techniquessimilar to those described in Hyun et al. [[1]]. As shown in Figure 4A and C, all values ofISI and Δt_{r1} decreased and then increased as the temperature wasraised, but frequency shown in Figure 4B rose and thenfell as the temperature increased. Figure 4D, E and G showthat the parameters Δt_{f2}, Δt_{1/2},Δt_{r2}, Δt_{f1}, and V_{np,ave} decreased as thetemperature was raised with small standard errors for each value. On the other hand,the values of angles shown in Figure 4F increased as thetemperature was raised, and the values of angle θ_{1} were larger thanthose of θ_{2.} In addition, both A_{AP, max} and A_{AP,ave} shown in Figure 4H increased and thendecreased as the temperature was raised. These temperature dependencies of APparameters (except A_{AP, max}, A_{AP, ave}, and angles) were similarto the case of those analyzed using a dataset obtained from beating cells: in thatcase, A_{AP,ave} decreased as the temperature was raised between 16°Cand 28°C, but A_{AP, ave}, and angles were not shown [[1]].

### Temperature dependence of bursting patterns in burst-firing neurons

The values of six temperature-dependent burst parameters were averaged by using theseselected dataset shown in Table 1. As shown inFigure 5A, B, and C respectively, interburst interval,burst duration, and number of spike during burst decreased as temperature increased.ISI during burst shown in Figure 5D firstly decreased andthen increased as temperature increased. In contrast, the numbers of burst per minuteand of spike per minute shown in Figure 5E and Frespectively increased and then decreased as the temperature was raised.

However, since these six burst parameters were not independent parameters, someequations representing an interrelationship between these parameters of burst-firingneurons could be driven; number of spike during burst = (burstduration)/(inter spike interval during burst), number of spike perminute = (number of spike during burst) × (number of burst perminute), and number of burst per minute = 60/(interburstinterval + burst duration).

The values of these six temperature-dependent burst parameters during lowering thetemperature (●) and raising the temperature (○) are shown inFigure 6. In the cases of interburst interval, burstduration, number of spike during burst, and ISI during burst shown inFigure 6A, B, C, and D respectively, there were nosignificant differences in these parameter values between the cooling-off and heatperiods, which could be considered that there were reproducible properties oftemperature dependencies in these burst parameters. Moreover, the number of burst perminute and spike per minute during heating shown in Figure 6E and F were slightly larger than those during cooling between 19°Cand 25°C. However, it might be thought that the reproducible property oftemperature-dependent changes in bursting patterns could not be undermined by these.It is because the values of the experimental percentage error of interspike intervalduring burst were between 0.1% and 14.1%; experimental percentageerror ≡ 100% × | value of bursting parameter during heating(cooling) - average value of bursting parameter |/ average value of burstingparameter.

Thus, each mean value of functional forms of reproducible temperature-dependentbursting parameters in Figures 5 and 6 could be calculated by using these dataset and the computer program wedesigned. It was shown that these six bursting parameters changed in response totemperature variation, but the mechanism underlying these correlations remainsunverified. So, it was suggested that performing computational simulations of thesephenomena by using a modified Plant model which were composed of equations withtemperature-dependent scaling factors was necessary to mathematically clarify thetemperature-dependent changes of bursting patterns in burst-firing neurons.

### Simulation of temperature-dependent bursting patterns

In order to simulate various temperature-dependent spiking patterns of burstingpacemaker neuron, we set up nonlinear differential equations by modifying the Plantmodel [[15]] with temperature-dependent scaling factors [[25]-[27]]. The Plant model studied by Rinzel and Lee was designed to assessparabolic bursters by analyzing the fast and slow processes to show how spikeactivities were generated by their mutual interaction [[16]].

The cell membrane of modified Plant model contains sodium channels carrying a fastsodium current, I_{Na}, and potassium channels underlying a fast potassiumcurrent, I_{K}. The slow processes include the slow inward calcium current,I_{Ca}, and the slow changes in intracellular free calcium concentration,Ca. The accumulation of calcium turns on outward calcium-activated potassium current,I_{K(Ca)}, and undermines the slow inward calcium current, and finallybrings to an end of bursting activity. This system also contains a leak current,I_{L}.

Where

and

and

Where ${\overline{\mathrm{g}}}_{\mathrm{Na}}$, ${\overline{\mathrm{g}}}_{\mathrm{Ca}}$, ${\overline{\mathrm{g}}}_{\mathrm{K}}$, and ${\overline{\mathrm{g}}}_{\mathrm{L}}$ are maximal conductances for the Na^{+}, Ca^{2+},K^{+}, and Cl^{−} currents, respectively, andV_{Na}, V_{Ca}, V_{K}, and V_{L} are the reversalpotentials for the respective currents: ${\overline{\mathrm{g}}}_{\mathrm{K}\left(\mathrm{Ca}\right)}$ is maximal conductance for calcium-activated potassium current.Here, voltage, V indicates membrane potential (mV).

The voltage-dependent activation and inactivation variables for sodium channels are*m* and *h* respectively, and the activation variable for potassiumchannels is *n*. A slowly activating conductance for calcium current isdenoted by χ, and a slow change in intracellular free calcium concentration Ca,is treated as parameters. The maximal relaxation time constants of h and n aredefined as $\frac{1}{\mathrm{\lambda}}$. It is taken that ρ^{−1} is an estimate for thetime-constant of the Ca-equation. The temperature-dependent scaling factors,ρ(T) and ϕ(T), are defined as $\mathrm{\rho}\left(\mathrm{T}\right)\equiv {1.3}^{\frac{{\mathrm{T}\u2010\mathrm{T}}_{0}}{10\phantom{\rule{0.24em}{0ex}}\xb0\mathrm{C}}}$ and $\mathrm{\varphi}\left(\mathrm{T}\right)\equiv {3}^{\frac{\mathrm{T}-{\mathrm{T}}_{0}}{10\xb0\mathrm{C}}}$, respectively [[25]]. The steady-state values of activation or inactivation variablesm_{∞}, h_{∞}, χ_{∞}, andn_{∞} are functions of voltage. The relaxation time constants arerepresented by τ_{h}, τ_{χ}, andτ_{n}.

### Comparison of experimental and simulation results

In order to obtain good computer simulation of temperature-dependent burstingpatterns generated respectively by these eight bursting pacemaker neurons, it wasnecessary to select three good data files of each experiment. The figures drawn byusing data files selected during the same rising (or falling) phase of temperaturechange from dataset of experiment A at temperatures below 20°C, between20°C and 25°C, above 25°C, were expressed by symbols L1, M1, H1,respectively, shown on Figure 7. Similarly, the figuresshown by using data files selected from datasets of experiments B – H at thesame temperature ranges were symbolized by L2- L8, M2 - M8, and H2 - H8,respectively, shown on Figure 7.

Then, it was necessary to describe how to draw three figures shown on the panels L1,M1, and H1 in Figure 7. Firstly, temperature values shownon these figures should be substituted for T in the temperature-dependent scalingfactors, ρ(T) and ϕ(T) involved in modified Plant model. Secondly, computersimulation should be carried out until percentage errors of all parameters had to becalculated below 50%: the results of these calculations were shown at first threerows in Table 3. However, percentage error could bedefined with the following formula shown below:

Finally, numerical values of each of the 24 parameters involved in equations (1) to(4) would be fixed: C_{m} = 1 μF/cm^{2}, ${\overline{\mathrm{g}}}_{\mathrm{KCa}}=0.018\phantom{\rule{0.24em}{0ex}}\mathrm{mmho}/{\mathrm{cm}}^{2}$, ${\overline{\mathrm{g}}}_{\mathrm{Na}}=4.0\phantom{\rule{0.24em}{0ex}}\mathrm{mmho}/{\mathrm{cm}}^{2}$, ${\overline{\mathrm{g}}}_{\mathrm{Ca}}=0.007\phantom{\rule{0.24em}{0ex}}\mathrm{mmho}/{\mathrm{cm}}^{2}$, ${\overline{\mathrm{g}}}_{\mathrm{K}}=0.60\phantom{\rule{0.12em}{0ex}}\mathrm{mmho}/{\mathrm{cm}}^{2}$, ${\overline{\mathrm{g}}}_{\mathrm{L}}=0.017\phantom{\rule{0.24em}{0ex}}\mathrm{mmho}/{\mathrm{cm}}^{2}$, V_{Na} = 40 mV,V_{Ca} = 140 mV,V_{K} = −75 mV,V_{L} = −40 mV, λ = 0.18,ρ = 0.000074 ms^{−1},τ_{χ} = 1, 500,K_{c} = 0.0275 mV^{−1},α = 127/105, β = 8265/105,γ = 0.3, δ = −18,μ_{m} = 0.1, μ_{h} = 0.08,μ_{n} = 0.016, ν_{n} = 0.1, ${\overline{\mathrm{\tau}}}_{\mathrm{n}}=1$, T_{0} = 23.0°C. Here, we took a value of23°C as the reference temperature T_{0} because bursts are usuallyactivated from 22°C to 25°C, and this represented the middle value oftemperature range of the experiment A; from 16.0°C to 30.0°C.

To simulate the other figures shown in Figure 7, it wasnecessary to fix the same numerical values of each of the 22 parameters as those usedin experiment A, except the values of ρ and τ_{x}. And then,computer simulation should be carried out by changing parameter values ρ andτ_{x} until percentage errors of all parameters had to be calculatedbelow 50%: interburst interval, burst duration, and number of spike per burst wereincreased (number of burst per minute decreased) as ρ, and τ_{x}increased. The results of these calculations were shown from 4th row to the last rowin Table 3: the values of all calculated percentage errorswere between 0 and 47.7%. Here, the results about the number of spike per minuteamong burst parameters was excluded from Table 3, becauseall percentage errors in the number of burst per minute were 0.0 and then the valuesof number of spike during burst would be the same as the values of number of spikeper minute.

The results of computer simulation with the values of temperature shown on each panelin Figure 7 are represented as the time series of burstingactivity. Although the results of simulation (red solid line) of AP amplitude werenot completely optimal, simulation data for number of burst per minute, number ofspike during burst, inter-burst interval, burst duration, and interspike intervalduring burst were well reflected by modified Plant model when these simulationresults were compared to the experimental results (blue dotted line).

Therefore, analysis and simulation of the experimental data by our equation modelmight be helpful for understanding the mechanisms underlying the change oftemperature-dependent bursting activity in neurons*.*

## Discussion

In the previous study we analyzed the temperature-dependent change of 14 AP parameters(such as the AP amplitude, membrane potential at the positive peak, ISI, first half andlast half of the temperature rising phase, first half and last half of the temperaturefalling phase, absolute value of the membrane potential at negative peak, absolute valueof the maximum slope of the AP during the temperature rising and falling phases, andspiking frequency) in *Aplysia* neurons [[1]]. With these findings, in this study we tried to examine thetemperature-dependent change of bursting patterns in pacemaker neuron at abdominalganglion of *Aplysia juliana.* Furthermore, we attempted to identify themechanism underlying the temperature dependence of busting by developing an equationmodel predicting neuronal electrophysiological activity with a function oftemperature.

### Temperature-dependent change of bursting patterns

Firstly we conducted real experiments with bursting pacemaker neurons in order toexplore the functional relation between temperature and bursting patterns as well asAP parameters. We examined the temperature-dependent activity of bursting pacemakerneurons during several consecutive heating − cooling cycles. Wetested the temperature dependencies of AP parameters in bursting pacemaker neurons inthe temperature range of 16–30°C. The temperature-dependent properties ofAP parameters were analyzed with the selected data shown in Table 1 by using techniques similar to those described in Hyun et al. [[1]].

The values of inter-burst interval, burst duration, and number of spike during burstdecreased as the temperature increased, and it was possible to confirm functionalproperties of temperature dependencies in six burst parameters. When compared to datareported by other research group, Fletcher and Ram analyzed changes of AP parametersof burst-firing neurons during heating by investigating *Aplysia* R15pacemaker neuron activity in the temperature range of 23 – 37°C; theyfound that interburst interval, burst duration, and number of spike per burst alldecreased within the temperature range of 23 to 32°C [[21]]. We found that our data are consistent with results reported by Fletcherand Ram [[21]] within the temperature range of 23 to 30°C. Furthermore, the presentstudy revealed that the values of ISI and Δt_{r1} decreased and thenincreased, frequency increased and then decreased, and the parametersΔt_{f2}, Δt_{1/2}, Δt_{r2},Δtf_{1}, and V_{np} all decreased as the temperature wasraised. These results were similar to those shown in previously published paper [[1]]. Fletcher and Ram also showed that intra-burst spike broadening and spikeheight decreased as temperature increased [[21]]. These results were also in good accordance with our experimental resultswithin the temperature range of 23 to 30°C. Because Δt _{1/2} andA_{AP,ave} decreased as temperature increased in the temperature range of23 to 30°C, as shown in Figure 4D and H: twoparameters, spike-broadening and spike height, could be compared with two APparameters, Δt_{AP,1/2} and A_{AP,ave}, respectively. Moreoverthe present study generated quantitative data regarding the temperature-dependentproperties of the other three AP parameters of burst-firing neurons (ISI duringburst, number of burst per minute, and number of spike per minute) and nine other APparameters (ISI, Δt_{r1}, Δt_{r2}, Δt_{f1},Δt_{f2}, θ_{1}, θ2, V_{np}, and A_{AP,max}).

### Mechanism underlying the temperature dependence of busting: choice of themodel

Numerous revised model of bursting in R15 [[17]-[19],[28]] have been developed since the mid-1980s. Canavier et al. [[17]] developed the model of bursting in R15 with 10 currents and 11 staticvariables; 10 currents (the fast Na^{+}-current I_{Na}, the fastCa^{2+}-current I_{Ca,} the delayed rectifierK^{+}-current I_{K}, the slow inward Ca^{2+}-currentI_{SI}, the nonspecific cation current I_{NS}, the anomalousrectifier current I_{R}, the leak current I_{L},K^{+}-current I_{K}, the Na^{+}-Ca^{2+} exchangercurrent I_{NaCa}, the Na^{+}-K^{+} pump currentI_{NaK}, and the Ca^{2+} pump current I_{CaP}), and 11static variables (membrane potential V, intracellular concentration ofCa^{2+} [Ca]_{i}, the occupancy of the intracellularCa^{2+}-buffer O_{c}, and eight voltage-dependent state-variables(m, h, d, f, n, l, s, b)). The parabolic bursting pattern results from two competingprocesses, the slow voltage-dependent activation of the slow inwardCa^{2+}-current (I_{SI}), and the still more slower calcium-dependentinactivation of this same current. Butera et al. [[19]] also developed the models of bursting in R15 with 10 currents and 12static variables; 10 currents (the fast Na^{+}-current I_{Na}, thefast Ca^{2+}-current I_{Ca,} the delayed rectifierK^{+}-current I_{K}, the slow inward C^{2+}-currentI_{SI}, the nonspecific cation current I_{NS}, the anomalousrectifier current I_{R}, the leak current I_{L},K^{+}-current I_{K}, the Na^{+}-Ca^{2+} exchangercurrent I_{NaCa}, the Na^{+}-K^{+} pump currentI_{NaK}, and the Ca^{2+} pump current I_{CaP}), and 12static variables (membrane potential V, intracellular concentration ofCa^{2+} [Ca^{2+}]_{i}, intracellular concentration of cAMP[cAMP], the occupancy of the intracellular Ca^{2+}-buffer O_{c}, andeight voltage-dependent state-variables (m, h, d, f, n, l, s, b)). The slow inwardCa^{2+} current (I_{SI}) in this model is the key currentresponsible for leading to bursting phenomena.

In addition to these models, Bertram [[18],[29]] developed a mathematical model of bursting neuron R15 including 8currents and 11 static variables; 8 currents (excitatory sodium currentI_{Na}, excitatory calcium current I_{Ca}, inhibitory potassiumcurrent I_{K1} and I_{K2} , small leakage current I_{L},calcium current which initiate the burst I_{NSR}, and cation-nonspecificcurrent I_{D}, and a potassium current which is activated only near thepotassium equilibrium potential I_{R'}) and 11 static variables (membranepotential V, intracellular concentration of Ca^{2+} c, and ninevoltage-dependent state-variables (m, h, n, j, q, y, μ, r, x)). Two subthresholdburst currents, the calcium current which initiate the burst I_{NSR} and thecation-nonspecific current I_{D}, induce the busting oscillation.

To mathematically clarify the temperature-dependent changes of bursting patterns inburst-firing neurons, it was suggested to perform computational simulations of thesephenomena by using a proper model composed of equations with temperature-dependentscaling factors.

We have tried to conduct computational simulations by using modified models: inmodifying the models that were developed by Canavier et al. [[17]] and Butera et al. [[19]], the maximal conductances of currents ${\overline{\mathrm{g}}}_{\mathrm{i}}$ were multiplied by temperature-like scaling factor for currents $\mathrm{\rho}\left(\mathrm{T}\right)\phantom{\rule{0.12em}{0ex}}({\overline{\mathrm{g}}}_{\mathrm{i}}\times \mathrm{\rho}\left(\mathrm{T}\right)={\overline{\mathrm{g}}}_{\mathrm{i}}\times 1{.3}^{\frac{{\mathrm{T}\u2010\mathrm{T}}_{0}}{10\xb0\mathrm{C}}}$, i = Na, Ca, K, SI, NS), and the relaxation timeconstant T_{j} were multiplied by temperature-like scaling factor for ionickinetics ${\mathrm{T}}_{\mathrm{j}}\left(\left(1/{\mathrm{T}}_{\mathrm{j}}\right)\times \mathrm{\varphi}\left(\mathrm{T}\right)\right.=\left(1/{\mathrm{\tau}}_{\mathrm{j}}\right)\times {3}^{\frac{{\mathrm{T}\u2010\mathrm{T}}_{0}}{10\xb0\mathrm{C}}}$, j = m, h, d, f, n, l, s, b). In order to obtaincomputer simulation of temperature-dependent bursting patterns, it was necessary toselect three good data files obtained from experiment A at different temperatures;18.1°C, 22.1°C, and 27.9°C. The figures drawn by using these datafiles were expressed by symbols L1, M1 and H1, respectively, shown on the panels atFigure 8A and B. After temperature values shown onthese panels should be substituted for T in the temperature-dependent scalingfactors, ρ(T) and ϕ(T) involved in these modified models(T_{0} = 23°C), we could get these figures drawn in redsolid lines and these were shown on these panels. In the results of simulations withthese modified models, it could be found that burst duration and number of spikeduring burst increased as temperature increased. But these results were notconsistent with our experimental results. So, although the model presented byCanavier et al. [[17]] and Butera et al. [[19]] were based on experimental data to a great extent, it might be thoughtthat these models are not appropriate for investigation about the ionic mechanismunderlying the temperature-dependent change in bursting pattern from experiments withbursting pacemaker neurons in the abdominal ganglia of *Aplysia juliana*.

In a similar way, we could get three figures shown on Figure 8C by simulation using the modified models that were developed by Bertram [[18]]. Although the fact that the values of burst duration and number of spikeduring burst that were obtained from simulation by using this modified modeldecreased as temperature increased was consistent with our experimental results in asmall temperature range, it was very difficult to simulate the temperature-dependentchanges of bursting patterns in a large temperature range; from 18.1°C to27.9°C. So we tried to simulate temperature-dependent bursting patterns in asmall temperature range from 22.2°C to 25.8°C and the results were shown onFigure 8D. But these results were also not consistentwith our experimental results. Namely, the values of interspike interval during burstobtained from simulations were larger than those got from experiments, and thesevalues could not be reduced: the values of the other bursting parameters might begiven rise to inconsistencies with our experimental results. Thus, it might bethought that this model was not also appropriate for investigation about that.

Any way, it was shown that that I_{K(Ca)} is not involved in bursttermination in *Aplysia* bursting neurons Kramer and Zucker [[30]]. Even though calcium-activated potassium channels are located in the somaof R15 neuron, it was not necessarily proved that this channel was important forbursting pattern. Next, we have to discuss on the following question: are there anyother experimental evidences suggesting that I_{K(Ca)} is involved inbursting pattern of *Aplysia* neurons, particularly R15 neuron? Although wedidn’t yet find any one of research papers that were published after year 1985,it might be necessary to discuss whether the figures that we have suggested in thispaper could be given us any chance to run back over the experimental evidence onit.

Now, it might be necessary to compare the temperature-dependent changes in burstfiring patterns of subthreshold currents (or burst currents) involved in Plant modelwith those included in the other models. Time series of temperature-dependentsubthreshold current activities (ρ(T)∙I_{SI}) simulated by usingthe equation of slow inward calcium current based on the model developed by Canavieret al. [[17]] at temperature values 18°C, 23°C, and 28°C are showntogether in blue lines (L), black lines (M), and red lines (H), respectively on panelA in Figure 9. Those (ρ(T)∙I_{SI}) and(ρ(T)∙(I_{Ca} + I_{K(Ca)})) obtained fromsimulation data by using the model made by Butera et al. [[19]] and Plant [[15]] at the same temperature values are shown together on panel B and D inFigure 9, respectively. Those(I_{D} + ρ(T)∙I_{NSR}) taken fromsimulation data calculated by using the model made by Bertram [[18]] at temperature values 22.3°C, 23°C, and 23.7°C are showntogether on panel C in Figure 9. Next, it might be wantedto calculate the maximum values of temperature-dependent burst currents underlyingthe hyperpolarization of the inter-burst intervals, and then compare the values thatwere calculated by using Plant developed model with those by using the other models.The maximum values of these calculated by using Plant model at temperature values18°C, 23°C, and 28°C were 0.28 μA, 0.12 μA, and0.081 μA, respectively. The maximum values calculated by using Canavier etal. (Butera et al.) developed model at the same temperature values were−1.62 nA, −0.9 nA, and −0.75 nA (−1.5 nA,−1.1 nA, and −0.77 nA), respectively. The maximum valuescalculated by using the model that Bertram developed at the temperature values22.3°C, 23°C, and 23.7°C were −0.18 nA,−0.38 nA, and −0.58 nA, respectively. These absolute maximumvalues calculated by using the models that Plant, Canavier, and Butera developeddecreased as temperature increased, but those calculated by Bertram developed modelincreased as temperature increased. But more detailed analysis of these data remainedas a challenge for future study.

Anyway, these facts might show that it was hard to clarify mathematically thetemperature-dependent changes of bursting patterns at burst-firing neurons of*Aplysia juliana* with these three models involving many currents and alarge number of static variables. So, it was necessary to investigate much simplermodel with small number of currents and static variables such as Plant model.

At this stage, it might be useful to remind that temperature-dependent impulsepatterns of mammalian cold receptors could be well simulated by using nonlineardifferential equations involving I_{K(Ca)}. Different types of impulsepatterns of mammalian cold receptors can be observed as a function of skintemperature: irregular and less frequent burst discharges, regular and frequentbursting patterns, and irregular single spike patterns are observed from low to hightemperatures. These patterns could be simulated by Braun et al., and the Huber-Brauncold receptor model has been described in detail with 5 currents and 5 staticvariables (membrane potential and four voltage-dependent state-variables) [[26],[27],[31]-[36]]. This model consisted of two minimal sets of ionic conductances operatingat different voltage levels with different delays and the leakage current. The firsttwo voltage-dependent currents that generate the action potentials mean depolarizingcurrent I_{Na} and repolarizing current I_{K}. The next twovoltage-dependent slow currents that generate subthreshold potential oscillationswere slow depolarizing noninactivating Na^{+}-current I_{Nap} andslow repolarizing Ca^{2+}-dependent K^{+}-current I_{K(Ca)}with voltage dependent activation of Ca^{2+}-current. The temperaturedependences were given by temperature-like scaling factors for the maximumconductances and the time constants, with reference temperatureT_{0} = 25°C. Q_{10} of 3.0 for the activationvariable and Q_{10} of 1.3 for temperature dependences of maximumconductances were chosen.

Now, it would be necessary to look into another simple model to investigate themechanisms of temperature-dependent changes of bursting patterns which share a fewsimilarities with mammalian cold receptors. Plant model consisted of 5 currents and 6static variables: 5 currents (fast sodium current I_{Na}, fast potassiumcurrent I_{K}, inward calcium current I_{Ca}, calcium-activatedpotassium current I_{K(Ca)}, and leakage current, I_{L}) and 6 staticvariables (membrane potential, intracellular concentration of Ca^{2+}, andfour voltage-dependent state-variables). There are two fast and slow processes inthis model. The fast process had three components: the activation and inactivationvariables for Na^{+} channels “m” and “h”,respectively, and the activation of K^{+} channels n. The slow process hadtwo components: a slow conductance for Ca^{2+}-current “X”, andintracellular free calcium concentration “Ca”. Not only Plant model had asimple structure, but also theoretical analyses of Plant’s model were performedalready [[16],[37],[38]].

However, it is not easy to justify the reason of our choice of Plant model, becausesome papers excluded I_{K(Ca)} as a key bursting current in R15 and focusedon the inward current with slow activating and inactivating components. Plant modelwas based on the works of Gorman and Thomas [[39]], who demonstrated that I_{K(Ca)} was linearly dependent onincreasing concentration of intracellular calcium ions (Ca^{2+}) injectedinto the cytoplasm and expected to be activated during a burst. Chay [[40]] constructed mathematical model applicable to the *Aplysia*bursting neurons with involving the properties of Ca^{2+}-activatedK^{+} channel: intracellular Ca^{2+} concentration increasesgradually during activity to levels where it activates the I_{K(Ca)} [[10]]. Although Kramer and Zucker [[30]] concluded that I_{K(Ca)} did not play a role in burst terminationin bursting neurons of *Aplysia*, Adams and Levitan [[41]] did not exclude a role for I_{K(Ca)} in the repolarization ofaction potential. But they asked the question: why a long-lasting I_{K(Ca)}did not be activated by action potential? The plausible answer was the amount ofcalcium ions that entered into the cell during the activity of action potential wasinsufficient to generate it. Even if Canavier et al. [[17]] developed the model without I_{K(Ca),} it was suggested as alimitations of this model that more experimental data of I_{K(Ca)} wererequired for further investigation. Bertram [[18]] revised the model with I_{K(Ca)} conductance depended only on theconstant calcium concentration of the soma. Besides this, it has been known that thecalcium activated potassium channels are located in R_{15} soma [[42]] and electrode was inserted into the soma of the cell in the abdominalganglion of *A. juliana* to measure the membrane potentials in ourexperiments. Thus, we cannot underestimate the impact of I_{K(Ca),} but wedid not here want to claim the quantitative application of Plant model, rather, weused it fr understanding the mechanism underlying the temperature dependence ofbusting patterns by computer simulations: it might be challenging to work later witha Rinzel and Lee’s model based on the hypothesis that it is necessary for burstgeneration to consider an inward current which slowly inactivates as Ca^{2+}accumulates during the burst to the exclusion of I_{K(Ca)} [[16]]. By making a comparison between experimental data and simulation resultscalculated with modified Plant equations, it was suggested that the mechanismunderlying temperature-dependent bursting patterns of bursting pacemaker neurons atabdominal ganglion of *Aplysia juliana* might be derived from temperature-likescaling factors ρ(T) and ϕ(T) for the maximum conductances and the timeconstants, respectively; together with reference temperatureT_{0} = 23°C, Q_{10} of 3.0 for the activationvariable, and Q_{10} of 1.3 for temperature dependences of maximumconductances.

Taken together, it was suggested that a modified Plant model could be used tosimulate the temperature-dependent bursting activity of bursting pacemaker neurons inthe abdominal ganglia of *Aplysia juliana* and to unravel the mechanism oftemperature dependence in bursting patterns.

### Nomenclature

A_{AP}, action potential amplitude

A_{AP, ave}, averaged action potential amplitude

A_{AP, max}, maximum action potential amplitude

Ca, slow change in intracellular free calcium concentration

Frequency, spontaneous firing frequency (=1/ISI)

ISI, interspike interval (=Δt_{r1} + Δt_{AP,1/2} + Δt_{f2})

I_{Ca}, slow inward calcium current

I_{K}, fast potassium current

I_{K(Ca)}, calcium-activated potassium current

I_{L}, leak current

I_{Na}, fast sodium current

T_{0}, reference temperature (=23°C)

V, membrane potential

V_{Ca}, reversal potential for calcium current

V_{K}, reversal potential for potassium current

V_{L}, reversal potential for leak current

V_{Na}, reversal potential for sodium current

V_{np}, membrane potential at the negative peak

V_{pp}, membrane potential at the positive peak

Δt_{r1}, first half of the rising phase of AP

Δt_{r2}, last half of the rising phase of AP

Δt_{f1}, first half of the falling phase of AP

Δt_{f2}, last half of the falling phase of AP

Δt_{AP, 1/2}, AP half-width duration(=Δt_{r2} + Δt_{f1})

θ_{1,} inverse tangent of the ratio of half A_{AP} toΔt_{r2}

θ_{2}, inverse tangent of the ratio of half A_{AP} toΔt_{f1}

${\overline{\mathrm{g}}}_{\mathrm{Ca}}$, maximal conductance for calcium ion current

${\overline{\mathrm{g}}}_{\mathrm{K}}$, maximal conductance for potassium ion current

${\overline{\mathrm{g}}}_{\mathrm{L}}$, maximal conductance for leak current

${\overline{\mathrm{g}}}_{\mathrm{Na}}$, maximal conductance for sodium ion current

*h,* voltage-dependent inactivation variable for sodium channel;h_{∞} = steady state value of inactivation h

*m,* voltage-dependent activation variable for sodium channel;m_{∞} = steady state value of activation m

*n,* activation variable for potassium channel;n_{∞} = steady state value of activation n

μ_{m}, constant value involved in the defining equation form_{∞}

μ_{h}, constant value involved in the defining equation forh_{∞}

μ_{n}, *v*_{n}, constant values involved in the defining equation forn_{∞}

α, β, γ, δ, constant values involved in the defining equation forsteady state values of activation (inactivation) variables and relaxationconstants

1/λ, maximal relaxation time constant of h and n

${\overline{\mathrm{\tau}}}_{\mathrm{n}}$, maximal relaxation time constant of n

τ_{h,} relaxation constant of h

τ_{χ,} relaxation constant of χ

τ_{n,} relaxation constant of n

χ, slowly activating conductance for calcium current;χ_{∞} = steady state value of activation χ

ρ^{−1}, estimation for the time-constant of the Ca-equation

ρ(T), temperature-dependent scaling factor; $\rho \left(\mathrm{T}\right)\equiv {1.3}^{\frac{{\mathrm{T}\u2010\mathrm{T}}_{0}}{{10}^{\xb0\mathrm{C}}}}$

ϕ(T), temperature-dependent scaling factor; $\varphi \left(\mathrm{T}\right)\equiv {3}^{\frac{{\mathrm{T}\u2010\mathrm{T}}_{0}}{{10}^{\xb0\mathrm{C}}}}$

I_{SI,} slow inward Ca^{2+}-current

I_{NS,} nonspecific cation current

I_{R,} anomalous rectifier current

I_{NaCa,} Na^{+}-Ca^{2+} exchanger current

I_{NaK,} Na^{+}-K^{+} pump current

I_{CaP,} Ca^{2+} pump current

I_{K1,} inhibitory potassium current

I_{NSR,} calcium current which initiate the burst

I_{D,} cation-nonspecific current

I_{R,} potassium current which is activated only near the potassiumequilibrium potential

I_{Nap,} slow depolarizing noninactivating Na^{+}-current

## References

- 1.
Hyun NG, Hyun KH, Lee K, Kaang BK: Temperature dependence of action potential parameters in Aplysia neurons. Neurosignals. 2012, 20: 252-264. 10.1159/000334960.

- 2.
Murray RW: The effect of temperature on the membrane properties of neurons in the visceralganglion of Aplysia. Comp Biochem Physiol. 1966, 18: 291-303. 10.1016/0010-406X(66)90188-5.

- 3.
Stephens PJ, Frascella PA, Mindrebo N: Effects of ethanol and temperature on a crab axon action potentials: a possiblemechanism for peripheral spike generation. J Exp Biol. 1983, 103: 289-301.

- 4.
Heitler WJ, Edwards DH: Effect of temperature on a voltage-sensitive electrical synapse in crayfish. J Exp Biol. 1998, 201: 503-513.

- 5.
Frankenhaeuser B, Moore LE: The effect of temperature on the sodium and potassium permeability changes inmyelinated nerve fibres of

*Xenopus laevis*. J Physiol. 1963, 169: 431-437. - 6.
Dalton JC, Hendrix DE: Effects of temperature on membrane potentials of lobster giant axon. Am J Physiol. 1962, 202: 491-494.

- 7.
Burrows M: Effects of temperature on a central synapse between identified motor neurons inlocust. J Comp Physiol A. 1989, 165: 687-695. 10.1007/BF00611000.

- 8.
Kerkut GA, Ridge RM: The effect of temperature changes on the activity of the neurons of the snailHelix aspersa. Comp Biochem Physiol. 1962, 5: 283-295. 10.1016/0010-406X(62)90057-9.

- 9.
Hodgkin AL, Katz B: The effect of temperature on the electrical activity of the giant axon of thesquid. J Physiol. 1949, 109: 240-249.

- 10.
Gorman ALF, Hermann A, Thomas MV: Ionic requirements for membrane oscillations and their dependence on the calciumconcentration in a molluscan pace-maker neurone. J Physiol. 1982, 327: 185-217.

- 11.
Adams WB: Slow depolarizing and hyperpolarizing currents which mediate bursting in Aplysianeuron R15. J Physiol. 1985, 360: 51-68.

- 12.
Adams WB, Benson JA: The generation and modulation of endogenous rhythmicity in the Aplysia burstingpacemaker neuron R15. Prog Biophys molec Biol. 1985, 46: 1-49. 10.1016/0079-6107(85)90011-2.

- 13.
Plant RE, Kim M: On the mechanism underlying bursting in the Aplysia abdominal ganglion R15cell. Math Biosci. 1975, 26: 357-375. 10.1016/0025-5564(75)90022-X.

- 14.
Plant RE, Kim M: Mathematical description of a bursting pacemaker neuron by a modification of theHodgkin-Huxley equations. Biophys J. 1976, 16: 227-224. 10.1016/S0006-3495(76)85683-4.

- 15.
Plant RE: Bifurcation and resonance in a model for bursting nerve cells. J Math Biology. 1981, 11: 15-32. 10.1007/BF00275821.

- 16.
Rinzel J, Lee YS: Dissection of a model for neuronal parabolic bursting. J Math Biol. 1987, 25: 653-675. 10.1007/BF00275501.

- 17.
Canavier CC, Clark JW, Byrne JH: Simulation of the bursting activity of Neuron R15 in Aplysia: role of ioniccurrents, calcium balance, and modulatory transmitters. J Neurophysiol. 1991, 66: 2107-2124.

- 18.
Bertram R: A computational study of the effects of serotonin on a molluscan neuron. Biol Cybern. 1993, 69: 257-269. 10.1007/BF00198966.

- 19.
Butera RJ, Clark JW, Canavier CC, Baxter DA, Byrne JH: Analysis of the effects of modulatory agents on a modeled bursting neuron: dynamicinteractions between voltage and calcium dependent systems. J Comput Neurosci. 1995, 2: 19-44. 10.1007/BF00962706.

- 20.
Carpenter DO: Temperature effects on pacemaker generation, membrane potential, and criticalfiring threshold in Aplysia neurons. J Gen Physiol. 1967, 50: 1469-1484. 10.1085/jgp.50.6.1469.

- 21.
Fletcher SD, Ram JL: High temperature induces reversible silence in Aplysia R15 bursting pacemakerneuron. Comp Biochem Physiol A. 1991, 98: 399-405. 10.1016/0300-9629(91)90422-9.

- 22.
Wu J, Fisher RS: Hyperthermic spreading depressions in the immature rat hippocampal slice. J Neurophysiol. 2000, 84: 1355-1360.

- 23.
Kim JA, Connors BW: High temperatures alter physiological proper-ties of pyramidal cells andinhibitory interneurons in hippocampus. Front Cell Neurosci. 2012, 6: 27-

- 24.
Orlowski JP, Erenberg G, Lueders H, Cruse RP: Hypothermia and barbiturate coma for refractory status epilepticus. Crit Care Med. 1984, 12: 367-372. 10.1097/00003246-198404000-00006.

- 25.
Hyun NG, Hyun KH, Hyun KB, Han JH, Lee K, Kaang BK: A computational model of the temperature-dependent changes in firing patterns inAplysia neurons. Korean J Physiol Pharmacol. 2011, 15: 371-382. 10.4196/kjpp.2011.15.6.371.

- 26.
Braun HA, Huber MT, Dewald M, Schäfer K, Voigt K: Computer simulation of neural signal transduction: the role of nonlinear dynamicsand noise. Int J Bifur Chaos. 1998, 8: 881-889. 10.1142/S0218127498000681.

- 27.
Finke C, Freund JA, Rosa E, Braun HA, Feudel U: On the role of subthreshold currents in the Huber-Braun cold receptor model. Chaos. 2010, 20: 045107–045107–10-10.1063/1.3527989.

- 28.
Chay TR, Cook DL: Endogenous bursting patterns in excitable cells. Math Biosci. 1988, 90: 139-1530. 10.1016/0025-5564(88)90062-4.

- 29.
Bertram R: Reduced-system analysis of the effects of serotonin on a molluscan bursterneuron. Biol Cybern. 1994, 70: 359-368. 10.1007/BF00200333.

- 30.
Kramer RH, Zucker RS: Calcium-induced inactivation of calcium current causes the inter-bursthyperpolarization of Aplysia bursting neurons. J Physiol. 1985, 362: 131-160.

- 31.
Braun W, Eckhardt B, Braun HA, Huber M: Phase-space structure of thermoreceptor. Phys Rev E. 2000, 62: 6352-6360. 10.1103/PhysRevE.62.6352.

- 32.
Braun HA, Huber MT, Anthes N, Voigt K, Neiman A, Pei X, Moss F: Interaction between slow and fast conductances in the Huber/Braun model ofcold-receptor discharges. Neurocomputing. 2000, 32–33: 51-59. 10.1016/S0925-2312(00)00143-0.

- 33.
Braun HA, Huber MT, Anthes N, Voigt K, Neiman A, Pei X, Moss F: Noise-induced impulse pattern modifications at different dynamical period-onesituations in a computer model of temperature encoding. BioSystems. 2001, 62: 99-112. 10.1016/S0303-2647(01)00140-X.

- 34.
Braun H, Vogit K, Huber M: Oscillations, resonances and noise: basis of flexible neuronal patterngeneration. BioSystems. 2003, 71: 39-50. 10.1016/S0303-2647(03)00108-4.

- 35.
Feudel U, Neiman A, Pei X, Wojtenek W, Braun H, Huber M, Moss F: Homoclinic bifurcation in a Hodgkin-Huxley model of thermally sensitiveneurons. Chaos. 2000, 10: 231-239. 10.1063/1.166488.

- 36.
Mosekilde E, Sosnovtseva OV, Postnov D, Braun H, Huber MT: Noise-activated and noise-induced rhythms in neural systems. Nonlin Stud. 2004, 11: 449-467.

- 37.
Honerkamp J, Mutschler G, Seitz R: Coupling of a slow and a fast oscillator can generate bursting. Bull Math Biol. 1985, 47: 1-21. 10.1007/BF02459643.

- 38.
Soto-Trevino C, Kopell N, Watson D: Parabolic bursting revisited. J Math Biol. 1996, 35: 114-128. 10.1007/s002850050046.

- 39.
Gorman ALF, Thomas MV: Potassium conductance and internal calcium accumulation in a molluscan neuron. J Physiol. 1980, 308: 287-313.

- 40.
Chay TR: Electrical bursting and intracellular Ca

^{2+}oscillations in excitablecell models. Biol Cybern. 1990, 63: 15-23. 10.1007/BF00202449. - 41.
Adams WB, Levitan IB: Voltage and ion dependences of the slow currents which mediate bursting in Aplysianeuron R15. J Physiol. 1985, 360: 69-93.

- 42.
Lewis DV: Calcium-activated inward spike after-current in bursting neuron R15 of Aplysia. J Physiol Lond. 1988, 395: 285-302.

## Acknowledgements

This work was supported by a research grant from the Sochun Academic Research Fund ofJeju National University in 2012. K.L. was supported by Basic Science ResearchProgram NRF-2013R1A1A3010216 and NRF-2013R1A3A1072702 funded by the Ministry ofEducation, Science and Technology.

## Author information

## Additional information

### Competing interests

The authors declare that they have no competing interests.

### Authors’ contributions

N.G.H. carried out electrophysiological experiments, performed computer simulation ofexperimental data, and drafted the manuscript. K.H.H. and K.B.H. participated incomputer programming for simulation of experimental data. K.L. drafted and revised themanuscript. All authors read and approved the final manuscript.

## Authors’ original submitted files for images

Below are the links to the authors’ original submitted files for images.

## Rights and permissions

## About this article

#### Received

#### Accepted

#### Published

#### DOI

### Keywords

- Temperature dependence
- Bursting patterns
- Modified Plant model
- R15 pacemaker neuron
- Aplysia juliana