When I_app = −0.05 μA/cm², the theoretical neuron model exhibits resting state when g_h = 0.05mS/cm², 0.04 mS/cm² (in the presence of I_h), and 0.0 mS/cm² (blockage of I_h). When a hyperpolarization square current with duration 100 ms is injected to the theoretical model, the membrane potential exhibits different characteristics for different g_h levels, as shown in Figure 2.

For gh=0.05 mS/cm2, when the strength of the hyperpolarization square current is small (−0.4μA/cm², red), both voltage sag, which appears at the initial phase of the hyperpolarization current, and post-inhibitory rebound phenomenon, which appears after the termination of the hyperpolarization current, are induced, as shown in the 1st panel of Figure 2A. When the strength increases to a middle level (−0.8μA/cm², blue), the sag and rebound become strong and the rebound spike appears. When the strength increases to a high level (−1.2μA/cm², black), both the voltage sag and rebound spike become strong, and the rebound spike appears earlier than that of the middle strength (−0.8μA/cm²).

For gh=0.04 mS/cm2, both voltage sag and rebound (spike) appear, and the rebound (spike) phenomenon becomes strong with increasing the current strength, which is similar to those of gh=0.05 mS/cm2. However, compared with gh=0.05 mS/cm2, both voltage sag and rebound (spike) become weak for gh=0.04 mS/cm2, as shown in the 1st panel of Figure 2B), which shows that the voltage sag and rebound (spike) become weak with decreasing g_h. For example, when the current strength is middle (−0.8μA/cm², blue), rebound spike appears for gh=0.05 mS/cm2 while not spike but rebound appears for gh=0.04 mS/cm2.

For gh=0.0 mS/cm2, which corresponds to the blockage of I_h, both voltage sag and post-inhibitory rebound spike cannot be induced by the hyperpolarization square currents with the same strengths as Figures 2A,B, as depicted in the 1st panel of Figure 2C.

The results are consistent with the previous experimental investigations of I_h, which were performed on the hippocampal pyramidal neurons and stratum oriens interneurons (Manseau et al., 2008; Ascoli et al., 2010; Zemankovics et al., 2010; Gastrein et al., 2011; Borel et al., 2013).

The changes of the membrane potential for the resting state and of the firing frequency for the firing behavior in (I_app, g_h) plane are shown in Figures 3A,B, respectively. The resting state locates at down-left half and the firing at up-right half of (I_app, g_h) plane. When I_app is fixed, with increasing g_h, the membrane potential of the resting state increases, as shown in Figure 3A. The blank area corresponds to the firing behavior. When I_app is fixed, with increasing g_h, the firing frequency increases, which means that the ISIs of firing decreases, as shown in Figure 3B. The blank area corresponds to the resting state. Both results are consistent with the biological experiments (Elgueta et al., 2015). In addition, the firing threshold for I_app decreases with increasing g_h, which shows that g_h induces the decrease of the firing threshold for I_app. Both enhancement of I_h and increase of I_app play the same positive roles in evoking firing, therefore, the increase of g_h leads to the decrease of firing threshold for I_app.

When I_app = 0.17 μA/cm², the firing for both gh=0mS/cm2 (blockage of I_h) and gh=0.02mS/cm2 (in presence of I_h) in the deterministic model is period-1 and the CV of ISIs is 0.

When noise is introduced, the ISIs of the firing exhibit variability, which leads to a nonzero CV, for example, as shown in Figure 4A. With increasing noise intensity D, although CV for both gh=0mS/cm2 (red) and gh=0.02mS/cm2 (black) increases, CV for gh=0.02mS/cm2 is always higher than that for gh=0mS/cm2 in a large range of noise intensity (0 < D < 1 μA/cm²), which shows that spike-timing precision in the presence of I_h is higher than that of the blockage of I_h.

The more spike-timing precision induced by I_h can also be found from Figure 4B. When noise intensity D is fixed, for example D = 0. 6 μA/cm²,CV decreases with increasing g_h. The dependence of CV on both I_app and g_h are shown by the color scale in (I_app, g_h) plane in Figure 4C (D = 0.2μA/cm²) and Figure 4D (D = 0.6μA/cm²). For both noise intensities, CV decrease with increasing g_h, which shows that I_h includes more spike-timing precision. No firing locates in down-left half of (I_app, g_h) plane (the blank area) and no CV is calculated.

The detailed difference of the spike-timing precision between the presence of I_h and the absence of I_h can be found from the spike trains of the firing and the ISI histogram (ISIH). For example, when I_app = 0.17 μA/cm² and D = 0.2μA/cm², the spike trains for gh=0mS/cm2 and gh=0.02mS/cm2 are shown in Figures 5A,B), respectively, and the ISIH for gh=0mS/cm2 and gh=0.02mS/cm2 is illustrated in Figures 5C,D), respectively, which closely matches the experimental observations on the hippocampal pyramidal neurons and stratum oriens interneurons (Gastrein et al., 2011). The mean value of ISIs is 208.99 ms for gh=0mS/cm2 and 76.98 ms for gh=0.02mS/cm2, respectively. The standard deviation (STD) of ISIs is 102.24 ms for gh=0mS/cm2 and 14.09 ms for gh=0.02mS/cm2, respectively. The CV of the ISIs is 0.494 for gh=0mS/cm2 and 0.183 for gh=0.02mS/cm2, respectively. The result indicates that I_h reduces the variability of ISIs, i.e., enhances the precision of the spike-timing, which can also be found from the ISIH depicted in Figures 5C,D). Compared with Figure 5D (blockage of I_h), ISIH becomes broader and lower for Figure 5C (in the presence of I_h), which also shows that I_h can enhance the spike-timing precision.

Two cases of the subthreshold resonance for different I_app values are simulated.

When Iapp=-0.05μA/cm2, the theoretical model exhibits resting state when g_h = 0.0 mS/cm², 0.03 mS/cm², and 0.05 mS/cm². When a subthreshold “ZAP” current (Figure 6D) is applied, the membrane potential exhibits different responses at different g_h values, as shown in Figure 6A (g_h = 0.05 mS/cm²), Figure 6B (g_h = 0.03 mS/cm²), and Figure 6C (g_h = 0.0 mS/cm²). When g_h = 0.05 mS/cm², the amplitude of the oscillation increases firstly and then decreases with respect to the increase of time/frequency, correspondingly, the impedance profile exhibits an inverse “U” shape, as shown by the first line from the top (red) in Figure 6E, showing that a “typical”-resonance appears. The resonance frequency is around 3.1 Hz. When g_h decreases to 0.03 mS/cm², the amplitude of the oscillation exhibits a slight increase at first and then a decrease with increasing time/frequency, and the height of the inverse “U” shape of the impedance profile becomes short, which shows that the resonance becomes weak, as shown by the third line (blue) from the top in Figure 6E. In the present paper, such a resonance is called as “weak”-resonance. When g_h decreases to 0.0 mS/cm², the amplitude of the oscillation decreases with increasing time/frequency and the phenomenon of resonance disappears, called as non-resonance, as shown by the last line (black) from the top in Figure 6E. With decreasing g_h from 0.05 mS/cm², to 0.04mS/cm², to 0.03 mS/cm², to 0.02 mS/cm², and to 0.0 mS/cm², the maximal impedance decreases, as shown by the black square in Figure 6E. The result is consistent with the experimental observations in many types of neurons such as the stratum oriens interneuron (Zemankovics et al., 2010; Gastrein et al., 2011; Pavlov et al., 2011; Borel et al., 2013; Gonzalez et al., 2015).

Another representative similar to I_app = −0.05μA/cm² is shown in Figure 6F (I_app = −0.3 μA/cm²). The upper 2 lines (g_h = 0.12 and 0.1 mS/cm²), middle 2 lines (g_h = 0.08 and 0.06 mS/cm²), and lowest line (g_h = 0.0 mS/cm²) for impedance vs. frequency exhibit “typical”-, “weak”-, and non-resonances.

When Iapp=0.08μA/cm2, the response of the membrane potential to “ZAP” current and resonance exhibits characteristic different from Iapp=-0.05μA/cm2 and I_app = −0.3 μA/cm², as shown in Figure 7. From g_h = 0.02 mS/cm² to 0.01 mS/cm², and to 0.0 mS/cm², the responses of the membrane potential are shown in Figures 7A–C). The amplitude of the oscillation of the membrane potential decreases with increasing time/frequency, which shows that not “typical”-resonance but non- or “weak”-resonance appears, as shown by the red (g_h = 0.02 mS/cm²), blue (g_h = 0.01 mS/cm²), and black lines (g_h = 0.0 mS/cm²) in Figure 7D. With decreasing g_h, the maximal impedance decreases, as shown the black square in Figure 7D.

Considering that the subthreshold resonance for Iapp=0.08μA/cm2 exhibits characteristics different from that of Iapp=-0.05μA/cm2, the bifurcations with respect to g_h for Iapp=0.08μA/cm2 and I_app = −0.05 are calculated and manifest different dynamics, as shown in Figure 8.

For Iapp=0.08μA/cm2, the bifurcation of the theoretical model with respect to g_h are shown in Figures 8A,B (Figures 8B is the partial enlargement of Figures 8A). The theoretical neuron model exhibits firing/stable limit cycle (blue bold lines) when gh≥0.0229919mS/cm2. The upper and lower blue lines represent the maximal and minimal values of the membrane potential of firing. The equilibrium contains three branches (dotted, dashed, or thin solid lines). The upper branch (dotted line) and middle branch (dashed line) are unstable focus and saddle (unstable), respectively, and are not associated with the resting state. The lower branch (thin solid line) is related the resting state and composed of three parts, a stable node (gh<0.0169329mS/cm2, black solid line left to the red solid line), a stable focus (0.0169329<gh<0.0229915mS/cm2, red solid line), and a stable node (0.0229915<gh<0.0229919mS/cm2, black solid line right to the red solid line). Therefore, the resting state contains two dynamical behaviors, stable node and stable focus. There is a SNIC bifurcation at gh≈0.0229919mS/cm2, wherein a stable node (right part of the lower branch, black solid line right to the red line) contacts with a saddle (middle branch, dotted line) to disappear and a limit cycle (blue) appears. As g_h changes across the SNIC bifurcation point, the transition between the resting state corresponding to a stable node and the firing corresponding to a stable limit cycle happens.

Therefore, with increasing g_h, there are four stable behaviors and 3 bifurcation/transition points. The four behaviors are stable node, stable focus, stable node, and firing (stable limit cycle), and the 3 bifurcation/transition points are the transition point from the stable node to stable focus, transition point from the stable focus to stable node, and the SNIC bifurcation point.

For Iapp=-0.05μA/cm2, the bifurcations are different from those of Iapp=0.08μA/cm2, as shown in Figures 8C–E) (Figure 8D is the partial enlargement of Figure 8C, and Figure 8E is the partial enlargement of Figure 8D). Except for the firing/stable limit cycle (blue bold lines) and 3 branches of the equilibrium points (dashed or thin solid lines), a small unstable limit cycle (green) appears. The lower branch of the equilibrium point associated with the resting state is different from that of Iapp=0.08μA/cm2. With increasing g_h, the equilibrium (lower branch) is stable node (gh<0.0454454mS/cm2, black solid line), stable focus (0.0454454<gh<0.0620557mS/cm2, red solid line), unstable focus (0.0620557<gh<0.0623584mS/cm2, red dashed line), and unstable node (0.0623584<gh<0.0623686mS/cm2, black dashed line). Furthermore, there are more types of bifurcations for Iapp=-0.05μA/cm2. For example, the firing/stable limit cycle related to a big saddle homoclinic orbit (BHom) appears at gh≈0.0610595mS/cm2, wherein a stable limit cycle with large amplitude (blue) contacts with a saddle (middle branch of the equilibrium). Other three types of bifurcations are subcritical Hopf (SubH) bifurcation at gh≈0.0620557mS/cm2, saddle-node (SN) bifurcation at gh≈0.0623686mS/cm2, and small saddle homoclinic orbit (SHom) at g_h ≈ 0.061832mS/cm². The stable focus (red solid line) changes to unstable focus (red dashed line) via the SubH bifurcation and an unstable limit cycle with small amplitude (green) appears left to the SubH. The unstable limit cycle with small amplitude (green) contacts with the middle branch of the equilibrium (saddle) to form a SHom and disappears. The unstable node at the lower branch of the equilibrium contacts with the saddle (middle branch) to form a SN bifurcation and both equilibriums disappear.

As g_h increases across the SubH bifurcation point (gh≈0.0620557mS/cm2), the resting state corresponding to a focus (red solid line) changes to the firing corresponding to a stable limit cycle (blue line). As g_h decreases across the BHom bifurcation point (gh≈0.0610595mS/cm2), the firing corresponding to a stable limit cycle (blue line) changes to the resting state corresponding to a focus (red solid line). Therefore, bistability or coexistence of the resting state and firing lies between the BHom and SubH bifurcation points, which closely matches the experimental observations in the Purkinje cells (Engbers et al., 2013).

With increasing g_h, there are 6 bifurcation/transition points, the transition point from the stable node to stable focus, the BHom point, the SHom point, the SubH point, the transition point from the unstable focus to unstable node, and the SN bifurcation point related to the unstable node. There are 4 stable behaviors, the stable node, the stable focus, the coexistence of the resting state and firing, and the firing. The 3 bifurcation/transition points related to the stable dynamical behaviors, the transition point from the stable node to focus, the BHom point, and the SubH point, form the border of the parameter region between the four stable dynamical behaviors.

To investigate the comprehensive roles of I_h on the dynamical behaviors, bifurcations, and transitions between different kinds of bifurcations, the bifurcation diagram in (I_app, g_h) plane is acquired, as shown in Figure 9. In (I_app, g_h) plane, the codimension-1 bifurcation or transition points connect to form codimension-1 bifurcation or transition curves, which is the border of the parameter region between different dynamical behaviors. There are multiple codimension-1 bifurcation or transition curves, such as the SubH curve, BHom curve, SN curve, SNIC curve, and transition curve between node and focus which is labeled as NF curve. NF_U (dashed red line) and NF_S (solid red line) curves represent the part of the NF curve related to stable focus/node and unstable focus/node, respectively. SN_S curve is the part of the SN curve related to the stable node. The intersection point between different kinds of codimension-1 bifurcation curves is the codimension-2 bifurcation point, which is the key point to characterize the transition between different kinds of bifurcations. There are 2 codimmension-2 bifurcation points. One is Bogdanov-Takens (BT) point, which is labeled with a star and locates at Iapp≈0.0432μA/cm2 and gh≈0.03413mS/cm2. The other is the saddle-node homoclinic orbit (SNHO) bifurcation point, which is labeled with a black circle point and locates at Iapp≈0.07118μA/cm2 and gh≈0.02549mS/cm2. In addition, the NF_S curve (red line) is composed of the upper and lower branches and the connection point between two branches is called as T point (I_app ≈ 0.15751μA/cm² and g_h ≈ 0.00001mS/cm²), which is labeled with a black triangle (Figure 9D).

The enlargement of 3 part of the bifurcations from resting state to firing, which are labeled with arrows in Figure 9A, are shown in Figure 9B (Left to BT point), Figure 9C (between BT and SNHO points), and Figure 9D (right to SNHO point). The (I_app, g_h) plane is divided into 4 parts by the codimension-1 bifurcation or transition curves. From top to bottom, the stable behaviors within the 4 sub-regions are firing (gray), coexistence of firing and resting state (blank area in Figures 9B,C, too small to be seen in Figure 9A), stable focus (red horizontal lines), and stable node (vertical lines), as shown in Figure 9. The border of the parameter region of the 4 stable behaviors composed of codimension-1 bifurcation or transition curve and codimension-2 bifurcation points is shown in Table 1. From left to right, the bottom border of the firing is SubH curve (black line, Figures 9A,B), BT point, SN curve (black dashed line in Figures 9A–C, very close to the upper branch of the NF_S curve in Figure 9C), SNHO point, SNIC curve (black dashed line in Figures 9A–D). This border left to the SNHO point is the top border of the coexistence behavior, and the bottom border of the coexistence is the BHom curve (blue line in Figures 9A–C) ended at SNHO point. The top border of the stable focus region is composed of the BHom curve (blue line in Figures 9A–C) and the upper branch of the NF_S curve (red line in Figures 9C,D) between the BT point and T point, and the bottom border is the lower branch of the NF_S curve (red line in Figures 9A–D). The stable node locates in the remained area (Vertical line in Figures 9A,B,D).

The upper branch of the transition curve between the stable node and focus (NF_S curve) intersects with the codimension-1 SN bifurcation curve, SubH curve, and SHom curve at the BT point, which shows that BT point is related to both Hopf and SN/SNIC bifurcations, i.e., the transition point between Hopf and SN/SNIC bifurcations. The SNHO point is the intersection point between the bifurcation curves of SN, BHom, and SNIC. Other bifurcation (SHom and NF_U) curves not related to the stable behaviors, which are close to some of the bifurcation curves mentioned above, are not addressed in the present paper.

In summary, two obvious characteristics can be found from the bifurcation diagram. One is that the codimension-1 bifurcation curve or the border between the firing and resting state show a negative slope in the (I_app, g_h) plane, which correspond to the firing threshold shown in Figure 3. The firing behavior and the resting state locate up-right half and down-left half of (I_app, g_h) plane, respectively. The result shows that increasing g_h means far away from the bifurcation point for the firing behavior and approaching bifurcation point for the resting state corresponding to a stable focus or stable node. The other is that the distribution of the dynamical behaviors in the parameter region left and right to the BT point is different. Left to the BT point, the stable focus locates upper to the stable node. With increasing I_app, the parameter region of the stable focus becomes narrow. Down-right to the BT point, the parameter region of the resting state and stable focus becomes very narrow and at last the stable focus disappears with increasing I_app. The behavior upper, lower, and right to the stable focus is stable node (Figure 9D).

Firstly, the parameter regions for the two cases of the subthreshold resonance can be identified. Two different cases of subthreshold resonance shown in subsection 3.1.4 locate at left (I_app = −0.05 μA/cm² and I_app = −0.3 μA/cm², Figure 6) and right (I_app = 0.08 μA/cm², Figure 7) to the BT point, respectively, which implies that the BT point plays very important roles in the identification of dynamical mechanism of the subthreshold resonance.

Secondly, with help of the distribution of the stable focus and node left to the BT point, the “typical” -, “weak” -, and non-resonance can be well understood. The “typical”-resonance appears at large g_h value (g_h = 0.05 mS/cm² for I_app = −0.05 μA/cm², g_h = 0.1 mS/cm² and 0.12 mS/cm² for I_app = −0.3 μA/cm²), i.e., within the parameter region of the stable focus. The non-resonance appears at g_h = 0 mS/cm² (I_app = −0.05 μA/cm² and I_app = −0.3 μA/cm²), i.e., within the parameter region of the stable node. The “weak” -resonance appears near the transition curve between the stable focus and node (g_h = 0.02, 0.03 and 0.04 mS/cm²for I_app = −0.05 μA/cm², g_h = 0.06 and 0.08 mS/cm² for I_app = −0.3 μA/cm²). The results show that the “typical”-resonance in the presence of I_h is related to the stable focus, and the non-resonance in the absence of I_h (blockage of I_h) associated with the stable node. Therefore, I_h induces a transition from the stable node to the stable focus. In nonlinear dynamics, Hopf bifurcation related to the stable focus/firing and SN/SNIC bifurcation associated with the stable node/firing. Therefore, the transition point between the Hopf bifurcation and SN/SNIC bifurcation, i.e., the BT point, exists in the neuronal system with I_h.

Thirdly, the “typical”-resonance for the stable focus or non-resonance for the stable node and the difference between them are origin from the distinction of the frequency response due to eigenvalue. The eigenvalues for the stable focus are complex numbers and for the stable node are real numbers. When disturbed, the behavior of a stable focus is a damped subthreshold oscillation with a fixed frequency which is related to the imaginary part of the complex eigenvalues, as shown in Figure 10A (gh=0.05mS/cm2,Iapp=-0.05μA/cm2). However, the behavior of a stable node after perturbed exhibits no oscillations, as shown in Figure 10B (gh=0.01mS/cm2,Iapp=-0.05μA/cm2). Such two different behaviors exhibit different frequency responses to the periodic stimulus with single frequency. For a stable focus, the oscillation amplitude of the membrane potential reaches maximal when the stimulus frequency is close to the fixed frequency of the focus. For example, when stimulus frequency is 3.1 Hz, the oscillation amplitude is higher than those of other frequencies, as shown in Figure 10C. For a stable node, the oscillation amplitude of the membrane potential remains nearly unchanged at low frequency and then decreases with increasing the stimulus frequency, as shown in Figure 10D. The periodic stimulus current is asin2πf₁t (a = 0.01μA/cm², f₁ is the frequency). The different frequency responses between the stable focus and stable node are the fundamental base of the “typical”-resonance and non-resonance phenomenon induced by the “ZAP” current.

Fourthly, another important characteristic of the subthreshold resonance (Figures 6E,F,D) is that the maximal impedance increases with increasing g_h when I_app is fixed, which can be well understood with the bifurcation theory and the changes of the membrane potential shown in Figure 3A. With increasing g_h, the membrane potential increases. When I_app is fixed, increasing g_h means that approaching the bifurcation point from resting state to firing, and the oscillation amplitude of the membrane potential to a “ZAP” current stimulation increases, as shown by the difference between the maximal and minimal values of the membrane potential in Figure 11A (Iapp=-0.05μA/cm2) and Figure 11B (Iapp=0.08μA/cm2). The increase of the oscillation amplitude of the membrane potential leads to the increase of the maximal impedance with increasing g_h. Such a result is consistent with bifurcation theory wherein as a stable dynamical behavior approaches the bifurcation point, its capability of resisting disturbance such as regular stimulus becomes weaker. For the stable node or stable focus which should keep unchanged and remain at a fixed level, the weakness of the capability of resisting disturbance with increasing g_h, i.e., approaching the bifurcation point, is characterized by the increase of the oscillation amplitude of the membrane potential to the “ZAP” current stimulus.

Last, complex or “strange” phenomenon related to the subthreshold resonance appears near the bifurcation point. (1) In the small parameter region of the resting state or focus down-right to the BT point or SNHO point, not “typical”-resonance, but non-resonance or “weak”-resonance (Figure 7), appears for the stable focus, which is “strange” because a “typical”-resonance should appear. The dynamical mechanism is very complex and may be from the narrow region of focus and the close distance from the codimension-2 bifurcation point BT with double zero eigenvalues. Such a complex dynamical behavior is similar to those simulated in the parameter region near a special equilibrium with double zero eigenvalues in a linear model (Rotstein, 2013; Rotstein and Nadim, 2014). The dynamical mechanism underlying the nonlinear system used in the present paper and the difference to the linear system should be studied in future. (2) The “weak”-resonance appears for the stable node within the neighborhood of the transition curve between the stable node and focus, which is also “strange” due to a non-resonance should appear for a stable node. Such a result also shows that the complex behaviors appear near the bifurcation/transition point. Except for the bifurcation/transition point, the external stimulus signal, the “ZAP” current, may be an important factor to induce the “weak”-resonance for the stable node. (3) There is complex phenomenon in the very narrow region of the coexistence behavior. We do not address the phenomenon because the parameter region is very narrow. The complex dynamics within these three narrow parameter regions (Figure 9D) near the bifurcation/transition point should be further studied in future. Except for the three narrow regions, the “typical”-resonance and non-resonance appear in the parameter region of the stable focus and the stable node far away from the bifurcation/transition point, respectively.

The more spike-timing precision induced by I_h can also well be understood with the bifurcation theory. For the firing behavior, the increase of g_h means that the firing becomes farther away from the bifurcation point. According to the bifurcation theory, as a stable dynamical behavior becomes farther away from the bifurcation point, its capability of resisting disturbance such as noise becomes stronger (Izhikevich, 2000, 2007; Tateno et al., 2004). For a neuronal firing, the enhancement of the capability of resisting disturbance is characterized by the further decrease of the standard deviation (STD) of ISIs than the mean value of ISIs, with increasing g_h, as shown in Figures 12, 13, which leads to the decrease of CV of ISIs, i.e., the enhancement of the spike-timing precise.

The faster decrease of STD than the mean value can be found from Figure 12. For D = 0.2μA/cm², the mean value decreases from 290 ms to 20 ms, as shown in Figure 12A, and the STD decreases from 128 ms to 0.5 ms, as shown in Figure 12B. 290/20 (14.5) is much smaller than 128/0.5(256). For D = 0.6μA/cm², the mean value decreases from 140ms to 20ms, as depicted in Figure 12C, and the STD decreases from 64ms to 4ms, as shown in Figure 12D. 140/20 (7) is smaller than 64/4(16). The results are similar in a large noise intensity (0 < D < 1.0μA/cm²).

The faster decrease of the STD than the mean value can also be found from two representatives shown in Figure 13. If the mean value for different g_h values are normalized by the largest value corresponding to the smallest g_h, the decrease of the normalized values with increasing g_h is shown by the red line in Figure 13A (I_app = 0.17 μA/cm² and D = 0.6μA/cm²) and Figure 13B (I_app = −0.08 μA/cm² and D = 0.2μA/cm²). Similarly, the decrease of the normalized STD with increasing g_h is depicted by the black line in Figures 13A,B). For both representatives, the STD (black) decreases faster than the mean value (red).

The results show that the enhancement of the spike-timing precision is mainly caused by the fast decrease of standard deviation of ISIs with increasing g_h, which obeys the bifurcation theory wherein a stable dynamical behavior far away from the bifurcation point exhibits stronger capability of resisting disturbance such as noise.

The voltage sag is mainly induced by the hyperpolarization activation of I_h, which can be found from the 2nd−4th panels of Figures 2A–C. At the initial phase of the hyperpolarization current, I_h increases with increasing time, which leads to the increase of the membrane potential, i.e., the voltage sag. If the strength of the hyperpolarization current is stronger, I_h becomes higher. Correspondingly, the voltage sag becomes stronger. If g_h is larger, I_h is higher, which also leads to the stronger voltage sag. If there is no I_h, no voltage sag appears.

The rebound (spike) is mainly induced by both I_h and Na⁺ current, as shown by the 2nd to 4th panels of Figure 2. At the termination time of the hyperpolarization current, I_h is much higher than that of the resting state, which leads to the depolarization of the membrane potential. The depolarization of the membrane potential induces the decrease of I_h and the increase of Na⁺ current which plays the dominant role and leads to further depolarization of the membrane potential. The larger the g_h, the higher the I_h; the larger the strength of the hyperpolarization square current, the higher the I_h. When g_h and/or the strength of the hyperpolarization current are lower, I_h is lower and the depolarization is not strong enough to induce a spike, which corresponds to the rebound phenomenon. When g_h and/or the strength of the hyperpolarization current are larger, I_h is higher and the depolarization is strong enough to induce a spike, which corresponds to the rebound spike phenomenon.

The rebound (spike) becomes strong with increasing g_h, which can also be well understood with the bifurcation theory. According to the bifurcation theory (Izhikevich, 2000, 2007), as a stable dynamical behavior approaches the bifurcation point, its capability of resisting disturbance such as regular stimulus becomes weaker. For the resting state, the increase of g_h means that it approaches the bifurcation point, and the weakness of the capability of resisting disturbance such as the hyperpolarization square current is characterized by the stronger rebound (spike).

For the resting state for other I_app values, the results are similar to those obtained with I_app = −0.05 μA/cm². If there is no I_h, no rebound (spike) is induced.