Abstract
We propose a theoretical model consisting of coupled differential equation of membrane potential phase and temperature for describing the neuronal signal in mammals cold receptor. Based on the results from previous work by Roper et al., we modified a nonstochastic phase model for cold receptor neuronal signaling dynamics in mammals. We introduce a new set of temperature adjusted functional parameters which allow saturation characteristic at high and low steady temperatures. The modified model also accommodates the transient neuronal signaling process from high to low temperature by introducing a nonlinear differential equation for the “effective temperature” changes which is coupled to the phase differential equation. This simple model can be considered as a candidate for describing qualitatively the physical mechanism of the corresponding transient process.
1. Introduction
Mammals complex thermoreceptor systems consisting of free nerve ending fibers are located in the dermis, muscle, skeleton, liver, and hypothalamus [1]. It is a phasic receptor which is active when there is a change in environmental temperature and rapidly becomes steady when reaching the stable temperature. Based on its characteristics with respect to the temperature level, it can be classified into warm or cold receptor [2, 3], which is, respectively, sensitive to high or low temperature relative to the normal body temperature, characterized by its way in delivering the neuronal signals. The corresponding neuronal signals are delivered in the form of bursting, that is, rhythmic of action potential consisting of spikes and punctuated by periods of inactivity [4, 5]. Their characteristics depend strongly on the associated temperature levels.
In this report, we focus our discussion on the dynamics of mammals cold receptor. In a low temperature condition, the corresponding neuronal signals produce periodic bursts with uniform duration and slow oscillation characteristic, but with nonuniform spike frequencies for each burst. When the temperature is raised up by a quasistatic process, the amount of spikes per burst tends to decrease forming a periodic single spike or beating. At a relatively higher temperature, the spike pattern becomes aperiodic; namely, it can also exhibit either double spike or stochastically phase-locked spike (skipping) phenomenon [3]. An experimental study on the static and dynamic discharge of a specific mammals cold receptor, that is, cat’s lingual nerve, has been comprehensively conducted by Braun et al. [4]. In particular, they showed that the dynamical response of the associated cold receptor is different for various temperature transitions between 10°C and 40°C.
Nowadays, many models have been proposed to explain the dynamical characteristics of mammals cold receptor. One of the most profound models is the conductance-based model which relies on the conductance voltage-dependent phenomenon due to the existence of Na+ and K- ions. For example, Braun et al. [6] in their report had discussed a Hodgkin-Huxley voltage-conductance type equation in their attempt to understand the role of nonlinearity and noise on the dynamics of nerve cell membrane through mammals cold receptor data. Another conductance-based cold receptor model was also discussed in different reports [7, 8].
In the meantime, there is a certain type of ion channel called transient receptor potential melastatin 8 (TRPM8) that plays an important role in delivering the cold receptor neuronal signal (see [9] for review). The role of TRPM8 has been shown experimentally in thermosensation mechanism in mice as discussed in [10–12]. Very recently, a conductance-based model which includes the role of TRPM8 ion channel has been proposed by Olivares et al. and showed a good agreement with the experimental data found from the cold receptor of mice [13]. The corresponding Olivares model successfully resembled the experimental data of increasing the firing rate for quasistatically increasing exposed temperature protocol.
Apart from those conductance-based models, a fully ionic model has been proposed by Longtin and Hinzer [5], which discussed the stochastic action potential phase model specifically for a cat’s lingual cold receptor. This model was further simplified by Roper et al. [14], namely, by introducing a simplified phase differential equation. It was demonstrated that the corresponding model was able to approximate the Longtin-Hinzer model for temperature interval 17.8°C to 40°C. Compared to the conductance-based model, Roper’s model [14] offered a relatively simple mathematical description. However, we discovered that this model did not lead to a realistic description on phenomena that occurred in higher or lower temperature conditions.
Based on this fact, in the present report we discuss a possible modification on the corresponding Roper’s model for the nonstochastic limit by introducing a new functional form of parameters that appeared in the corresponding model. Furthermore, we also discuss an extension of the corresponding modified model to accommodate the dynamical response of neuronal signals during a transition process from high to low temperature condition. This dynamical model is able to explain the phenomenon of sudden increasing amount of spikes per burst due to decreasing temperature, which is followed by a gradual decreasing of the corresponding amount of spikes per burst until the receptor reaches a steady condition at the lower temperature [4, 15]. We explain this phenomenon by considering an additional differential equation to describe the temperature dynamics, which is coupled to the associated phase differential equation.
We organize the report as follows: Section 2 discusses the phase model for the case of the steady temperature condition. The modified models for static temperature and dynamic transient process from high to low temperature are given in Section 3, namely, by defining a new set of functional parameters in the corresponding phase differential equation and introducing a new differential equation of temperature coupled to the phase differential equation and we focus our discussion on the characteristics of spike per burst, burst period, and interspike interval. We end this report with a conclusion in Section 4. Comprehensive discussions regarding the biological and chemical related properties of the corresponding cold receptor have been given in detail previously [4, 5, 14], such that in this report we only focus on the modified mathematical model.
2. Model and Method
The corresponding nonstochastic phase differential equation for steady condition of neuronal signaling at a specific temperature developed previously by Roper et al. [14] is given as follows:withwhereHere, the symbol represents the phase of membrane potential in the trigger region, in which its full rotation describes the generation of an action potential [14]. The parameter is related to the modulation of the mean potential of the cell, while the term is a zero mean periodic term that oscillates with the frequency , with as the corresponding magnitude. The function , with an inverse time unit, describes the dynamics of the corresponding neuronal signal bursting [14].
It is seen that there are two important terms in (2), namely, and functions, as given by (3) and (4), respectively. The burst occurs when , where the average amount of spikes in each burst is proportional to the maximum width of the overlap area of both curves, denoted by , as exemplified in Figure 1(a) along with the corresponding phase of membrane potential , which is found by solving (1), and neuronal signal bursting functions as shown in Figures 1(b) and 1(c), respectively. It is obvious that the amount of spikes per burst can be controlled by changing the value of and as well as which also determine the period between two consecutive bursts. It was assumed previously [14] that these parameters are of linear functional forms of temperature as follows:with , , , , , and being constants to be determined. This assumption was aimed at yielding a decreasing period between two consecutive bursts when the temperature increases through a quasistatic process. In their work, Roper et al. defined the value of each parameter as follows: , , , , , and . We used all these parameters at = 35°C to depict the example shown in Figure 1.
(a)
(b)
(c)
Based on the above formulation, it is clearly seen that these linear assumptions will lead to an unrealistic scenario at the high and low temperature conditions, since all those parameters are not saturated at these limits. Therefore, it is reasonable to assume that phenomenologically at those temperatures the neuronal signals become saturated since in that range the receptor becomes less sensitive [16]. It is interesting to note that this model can also be further developed to describe the transient transition from high to low temperatures as previously reported by Ring and de Dear [15]. For this, we propose assuming that the corresponding temperature should be considered as a function of time with Morse-like characteristic described by a differential equation which is coupled to the corresponding phase differential equation given by (5). To study the dynamical characteristics of this model, we numerically solve the related coupled differential equations by means of standard Runge-Kutta method.
3. Results and Discussion
3.1. Modified Phase Model for Steady Temperature Condition
To develop a more realistic model, we consider the modification of , and parameters by introducing the following tanh functional forms: which exhibit sigmoidal saturation characteristic at relatively high and low temperatures. Here, and are parameters to be adjusted. We denote the parameter as an “effective temperature” for describing the dynamics of neuronal signal during the transient process from high to low environmental temperature. The meaning of this parameter will become clear in later discussion (see Section 3). Obviously, the functional parameter forms given by (8) will lead to saturated behavior of neuronal signal at high and low temperatures. By considering the same values with that used previously [14] for a temperature interval of 40°C to 15°C and at → −∞ we found /°C and = 33.75°C, while the other parameters are set to , , , , , and .
Demonstrated in Figure 2 is the comparison between the previous set of functional parameter forms and the new ones. It should be realized that, for the actual cases, those functional forms should be adjusted using experimental data of the associated spiking and bursting neuronal signaling phenomena in low and high temperature conditions. Indeed, it is also important to note that one can choose different type of functional forms. A nonsigmoidal functional form was previously proposed in [17] which was aimed at mimicing the bifurcation characteristics of the conductance-based Huber-Braun model [7].
(a)
(b)
(c)
The corresponding neuronal signals at steady temperatures for = 40°C to 15°C, which are compared to the previous reported results [14], are depicted in Figure 3. Here the function of (2) is calculated by solving first function in (1) and then the solution is inserted into the corresponding equation. The average amount of spikes per burst (SB) for both functional forms is given in Figure 4. It is observed that, at low temperature condition, the tanh functional forms exhibit a more reasonable amount of spike per burst than the linear functional forms of Roper’s model [14] that exhibits higher amount as demonstrated in Figure 4. Indeed, there are some discrepancies between both linear and tanh functional forms since both assumptions do not perfectly coincide as clearly shown in Figure 2. But indeed, both forms share qualitatively similar bursting and spiking characteristics. As shown by the figures, it is important to realize that, for decreasing , the amount of spikes per burst is increasing.
The other important characteristic, namely, the interspike interval histogram (ISIH) of the neuronal signal at the corresponding different temperatures for the tanh functional forms, is given in Figure 5(a) for both Roper and the present modified model, along with its plot as a multivalued function of in Figure 5(b). The ISIH shows the existence of beating and skipping at high temperature condition which are indicated by the presence of large intervals. Clearly, the modified model exhibits a bit different characteristics than the original Roper’s model.
(a)
(b)
3.2. Model for Transient Transition Process
During a transient transition from high to low temperature, the existence of a peak response with relatively large amount of spikes per burst at a certain time was shown experimentally as the transition process begins as reported in [15, 18]. Based on the burst characteristic at steady temperatures as exemplified by Figure 3, we suspect that this condition might be perceived by the brain to occur at a temperature lower than the final temperature. In the meantime, another experimental result demonstrated that the period between two consecutive bursts during the transient transition process is higher than the period at final temperature [4]. To model the corresponding dynamical response, we consider first a conjecture that the effective temperature in our formulation is a function of time in the following Morse-like function [19]:where and are effective temperature at time and its lowest value, respectively, whereas denotes the time when Graphically, the parameter determines the depth of Morse-like curve as illustrated in Figure 6, while denotes its effective width. The corresponding Morse-like function is chosen because it is a mathematically well-defined function with no singularity. Phenomenologically, it is likely to be the best geometrical shape to describe the corresponding transient response characteristics among other similar functional forms such as the Lennard-Jones [20], the Buckingham Exponential-6 [21], and the Mie potential functions [22]. All these functions are commonly used to describe the molecular interactions [22]. In contrast to the Morse function, the other three functions contain a singularity.
It should be emphasized that the existence of the abovementioned peak response with large amount of spike per burst during the transient transition is the reason to define the term “effective temperature” as a tuning factor in our formulation based on the following argument: as shown in Figure 3, the amount of SB for low temperature is larger than the higher one. At the same time, a sudden increasing amount of SB occurs due to decreasing temperature, which is followed by a gradual decrease of SB until the receptor reaches a steady condition at the lower temperature [4, 15]. From all these facts, we therefore propose that the functional parameter, with its curve given in Figure 6, should be considered as a dynamical tuning factor that is needed to describe the dynamics of the related neuronal signal propagation.
It is easy to prove that the function given by (9) satisfies the following differential equation:where In the ensuing discussion, we choose to solve (10) numerically rather than using (9) in order to explain the corresponding peak response phenomenon. It should be noted that, to ensure the corresponding numerical solution of (10) is the Morse-like function as given by (9), one should consider a negative initial condition for ; that is, .
We expect the parameters and can be determined experimentally. However, in our calculation, we assume that the parameter is fixed to where and denote the final and initial effective temperature, respectively, such that in (9) and (11) is fixed toIndeed, one can assume different values for this parameter and it is clear that different leads to a different in (9). On the other hand, it is reasonable to assume that the parameter in (10) should be expressed as a function of , that is, , because it is natural to think that the shape of the associated Morse-like function is different in different transition process. To formulate the corresponding expression, first we plot functions and given by (3) and (4) and adjust the value of in (10) to meet a matching condition, which is indicated by the coincidence between the first overlaps width of both functions (denoted by in Figure 1) and the lowest effective temperature . Exemplified in Figure 7 is the associated matching condition for the dynamical response from = 40°C to = 15°C. For this transition, we found that the matching condition occurs at The calculation result for this parameter from = 40°C to various is given in Figure 8. Using a standard fitting procedure, it is found that all those values can be approximated by the following function:with and , while is given by (12). Therefore, it is clear that the differential equation (10) should be rewritten as follows:which is coupled to (1) through parameters given by (8). It is important to note that the parameter and function phenomenologically correspond to the characteristics of the burst period and amount of SB around the matching condition, respectively. Therefore, as mentioned previously, it is obvious that these two parameters should be determined experimentally by observing the spiking and bursting characteristics similar to what was done in [4].
Using this new model, the simulation results of transition processes from = 40°C to = 35°C, 30°C, 25°C, 20°C, and 15°C are depicted in Figure 9. Given in Figure 10 are the SB for transition to = 35°C and 15°C, along with the corresponding burst period (BP) as defined in Figure 1(c), which is qualitatively similar to that given in [4]. The parametric plot in phase-plane of and to give a more clear description is also given in Figure 10. It is shown that the SB and BP characteristics exhibit pronounced transient feature for transition to = 15°C. In the meantime, the occurrence of dense patterns of bursting process when is demonstrated, justifying the existence of previously discussed peak response phenomenon as shown experimentally in [4] during a transient time. In contrast, it is interesting to note that a monotonous response characteristic is exhibited during
(a)
(b)
An example of the approximate Morse-like function found from (15) is given in Figure 7, which shows the same position of the related matching condition compared to Morse-like function found by solving (10) numerically. It is clearly shown that, right after the transition process begins, namely, at the first burst, the amount of the corresponding spikes is larger than the next bursts. As a consequence of choosing the saturated tanh functional forms in the model parameters as given by (8), we found that the change of the SB at matching condition is at a reasonable level, especially in the case of = 15°C, where
To validate this modified model with experimental data, we focus on comparing qualitatively the SB and BP characteristics with the results reported by Braun et al. [4]. Given in Figure 11 are the spiking and bursting phenomena for the transition similar to what was discussed in [4] along with the associated function. We calculate the SB of four different transition intervals, as well as the BP parameter. The results are shown in Figure 12. It is interesting to note that qualitatively the corresponding SB characteristic exhibits fairly similar trend with the experimentally found SB figure in Figure in [4], while it is seen that graphically that BP exhibits similar characteristics with SB. Note that the transient characteristics are demonstrated significantly in the cases of low temperature transitions, that is, = 25°C to = 20°C and = 20°C to = 15°C. This feature can be explained as a consequence of function with wider profile due to larger value as shown in Figure 11.
(a)
(b)
(c)
(d)
(a)
(b)
On the other hand, in comparison with Olivares’s model [13], taking into account the role of TRPM8 ion channel transient current to explain the phenomenon of increasing firing rate (which indicates the increasing of SB) as exposed temperature decreases, this modified model offers a different perspective to the transition mechanism from high to low temperature, namely, by introducing the “effective temperature” functional parameter as a dynamic tuning factor that coupled to the phase of membrane potential. As discussed previously, although the physical meaning of this dynamical parameter is not clearly understood at this moment, we proposed that this parameter might be interpreted as an actual temperature being perceived by the mammal brain and it is likely reasonable to assume that the corresponding transient function is related to the complex role of TRPM8 ion transient channel. Indeed, this hypothesis should be separately investigated.
Furthermore, although our model is able to describe the existence of peak response at matching condition, however, it should be noted that the model leads to the increasing pause duration or the time distance between two consecutive bursts at the corresponding condition, while in reality this is not the case as reported previously [4]. This problem is a bit complicated to be solved and we suggest it can be overcome by defining new functional forms for parameters given by (8). The other problem with our model is related to the value of in (11) where in our calculation it was considered to be fixed. We expect this parameter can be determined experimentally, and this is beyond the scope of our study.
To this end, apart from the above mentioned problems, it is also realized that this modified model should be improved further, since the related effective temperature differential equation given by (15) does not take into account the influence of phase of membrane potential. We suggest that fully coupled differential equations that accommodate this feature will likely be able to give a good quantitative explanation of the dynamics and characteristics of neuronal signals of the corresponding cold receptors. This issue could be a challenging topic for future investigation.
4. Conclusion
We have discussed a modified Roper’s model for describing the characteristics of neuronal signaling in mammals cold receptor, especially for the temperature transition processes. The model consists of coupled phase-temperature nonlinear differential equations equipped with a set of functional parameters that saturate at low and high temperature. It was shown that our modified model is able to describe the experimental fact that the characteristics of neuronal signal in a transient transition process from high to low temperature exhibit the existence of large amounts of spikes per burst right after the process initiated, namely, by introducing the new functional parameter “effective temperature,” which plays a role as a dynamical tuning factor to explain the corresponding phenomenon. We propose that this dynamical tuning factor might be interpreted as a perceived temperature by the mammal brain in which its perception of temperature at has the lowest value, while and are coincides with the environmental temperatures. Certainly, it is intriguing to further examine experimentally whether this interpretation is correct or not. For instance, by observing the related mammal brain activity that corresponds to the temperature perception. Further studies should be conducted in order to overcome a few problems that still exist. However, this modified model can be considered as a dynamic simple alternative candidate to complex ionic models to describe qualitatively the transient transition from high to low temperature of the mammals cold receptor.
Disclosure
Firman Ahmad Kirana and Husin Alatas are co-first authors.
Competing Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
The authors would like to thank Yessie W. Sari from the Department of Physics, Bogor Agricultural University, for a useful discussion.