# Rebuttal

## Reviewer #1

> The authors show that by adding a pump current and Kout dynamics to a simple spiking model results in a squarewave burster; that is a burster whose quiescent state terminates at a saddle node and whose firing state terminates at a homoclinic.

> There have been many previous models of bursting using Kout as a parameter so it is difficult for me to see what is new in their paper. The most salient paper is:

> Barreto, Ernest, and John R. Cressman. "Ion concentration dynamics as a mechanism for neuronal bursting." Journal of biological physics 37 2011): 361-373

> which is esentially the present paper. See, in particular figure 6 in their paper.

We thank the reviewer for pointing our this imporant article. We now substantially discuss it already in the introduction of our manuscript (lines 68 onwards). We disagree, however, with the assesment that this is "essentially the present paper". The burst mechanism discussed in our manuscript is not among the burst mechanisms discussed in Barreto and Cressman (2011), all of which are slow wave bursters. In fact, the present hysterisis loop square wave burster can not be found in their model, because there is no hyperpolarising pump current in the voltage equation of Barreto and Cressman (2011).
In order to make our point clearer, we have included simulations in the results section in Figure 5 (A and B). These simulations demonstrate that by integrating the pump current into the membrane potential equation and subsequently inducing a shear in the fast subsystem bifurcation, the model could produce the burst phenomenon introduced in this paper. In contrast, when the pump is solely attributed to its effect on potassium concentration but its current is ommited from the membrane potential equation as in Cressman (2011) paper, the square wave burst observed in our paper cannot be replicated.

> Many followup papers are closely related:

> Ma, Kaihua, Huaguang Gu, and Zhiguo Zhao. "Fast–slow variable dissection with two slow variables: a case study on bifurcations underlying bursting for seizure and spreading depression." International Journal of Bifurcation and Chaos 31.06 (2021): 2150096.
>
> here the 2 slow are Nai/Ko

This is indeed an interesting article (which analyzed the Cressman 2011 paper more throughly) along with the above to contrast to the mechanism described here. As the reviewer points out the slow dynamics is two-dimenstional as in Barreto and Cressmann (2011). The rhythm is created via an interplay between potassium and sodium dynamics. In contrast, our mechanism does not rely on two slow variables generating the rhythm, but in fact only one slow variable, which we chose to be extracellular K^+^, because it is easier to manipulate and quantify experimentally.

The articles above do not provide a bifurcation analysis of the complete system or the slow system. Only the bifurcations of the fast subsystem with fixed concentration parameters are provided. We feel that the complete analysis provides additional insight in the the system particularly how the burst to tonic spiking transition manifestes. We have provided additional material to explain the transition. Specifically, the limit cycle of tonic spiking destabilises via a cascade of period doubling bifurcations of the complete system, see Fig 7 and the discussion in lines 568 onwards.

> Hübel, Niklas, and Markus A. Dahlem. "Dynamics from seconds to hours in Hodgkin-Huxley model with time-dependent ion concentrations and buffer reservoirs." PLoS computational biology 10.12 (2014): e1003941.

> Hübel, Niklas, Eckehard Schöll, and Markus A. Dahlem. "Bistable dynamics underlying excitability of ion homeostasis in neuron models." PLoS computational biology 10.5 (2014): e1003551.

These article are indeed key to the field of concentration induced bursting. We had cite them in our maniscuript already, but have now included it in multiple other sections including the discussion of period doubling.
> [name=Mahraz Behbood] Check that

> I dont see why you need the SN at all in the bursting mechanism. This is besides the point since in that case, it is the current that determines the firing rate. You just need bistability.

This is a very good point. The saddle-node loop (SNL) is not required in all mechanisms for bursting as the examples in the literature above show. But as the reviewer pointed out it is a gate to bistability in the onset bifurcation, which the specific square wave burster we focus on make use of. The unfolding of a Bautin bifurcaftion is an other example but likely to be less common. However, as it was previously shown that the SNL point exists in a broad class of neuron models (with SNIC onset bifurcations) our proposed burst mechanism is available to a substantial amount of cortical neurons and can be controlled by the pump expression level. We have clarified this point in the mansucript (see lines 624 and following).

> Given geometric factors, by my estimates Ko will be 2-3 orders of magnitude slower (I assume currents are muA/cm^2) and K is mM. So, you have a nice slow variable built right in. In Ermentrout & Terman's book, in an exercise, they make a little oscillator with Ko and an inward rectifier K current, so this is all pretty straightforward.

Thank you very much for pointing us to the example in Ermentrout and Terman.
We think that the exersise that you are pointing to is exercise number 9 chapter 4. The exercise is not particularly abut bursting. This exercise shows that in a two dimensional system with a modification time scale separation, one can generate osciliations. Although it is very educational exercise and the oscillations are of the relaxation oscillator type, it is neither a bursting system, nor does it investigate the pump current's role as a switching into a bursting mode. We now stress in the manuscrtip that there are many previously published burst mechansims using concentration dynamics. It is the fact that adding an electrogenic Na/K-ATPases that induces a transformation of the onset bifurcation structure enabeling a burst mechanism in a generic condcutance-based model that is new here.

## Reviewer #2

> The authors present a very comprehensive and illustrative study into the dynamical structure underlying slow bursting in models endowed with ionic dynamics. The paper is very well written and easy to read. Although the basic minimalistic dynamics have been presented previously [Cressman2009, Barreto2011], the present model is sufficiently different(Wang-Buzaki and electrogenic pump), and the in depth dynamical analysis, and insights into the interplay between fast-slow dynamics are sufficient to warrant publication.

We thank the reviewer for the assesment and have now added a more comprehensive literature discussion into our manuscript.

> However, the authors have not referenced the Barreto paper and appear to be unaware of the results presented there. The authors should reference that work at a number of points, particularly regarding the bifurcation structure. A proper discussion of its relevance to your work is the only ‘major’ change required.

We thank the reviewer for the important literature advice. We have added an detailed disussion on the Cressman and Barreto articles into the introductions from lines 68 onwards. The proposed bursting mechanisms in these articles are fundamentally different as they are of the slow-wave bursing type with a two dimensional slow system and the articles do not include the electrogentic pump current into current balance equation of the voltage evolution. Our burster, on the other hand, relies on a single slow variable and an essential shear in the fast subsystem bistability. As the reviewer has correctly pointed out the Cresmann article ommits the pump current in the voltage equation, which is a key component in the mechanisms described in the present manuscript. To hightliht the importance of the pump current we have added a didactic figure in which we omit the pump current from our model in the Cressman style and how a lack of bursting in that case, see Fig. 5A,B.

> I’m not sure how familiar the readers will be with ghosts, and there are other ghosts, like Faddeev-Popov ghosts, so it may behoove you to give a brief description of them, perhaps even consider invoking the language of a Crisis and destabilization of an attractor.

We have defined what is ment by ghost in line 278 of a fixpoint more carefully and also point to further literature.

> 253. Stopped
> 331. “shear transformed”, reads better as sheared.

We have change the text accordingly.

> [Cressman 2009]Cressman, J.R., Ullah, G., Ziburkus, J. et al. The influence of sodium and potassium dynamics on excitability, seizures, and the stability of persistent states: I. Single neuron dynamics. J Comput Neurosci 26, 159–170 (2009). https://doi.org/10.1007/s10827-008-0132-4
> [Barreto2011]Barreto, E., Cressman, J.R. Ion concentration dynamics as a mechanism for neuronal bursting. J Biol Phys 37, 361–373 (2011). https://doi.org/10.1007/s10867-010-9212-6

## Reviewer #3:

> The authors construct a minimal model of a neuron with an external potassium compartment and potassium concentration dynamics. The Na/K pump is included and the potassium reversal potential is expressed in terms of the concentrations via the Nernst equation. They show that this model bursts.

> The paper is devoted to an extensive bifurcation analysis using time scale separation. The ion concentration oscillations are the slow component, and the neuronal spiking is the fast component.

> First, the extracellular potassium concentration is fixed and treated as a parameter, and bifurcations in the fast system are identified. Superimposing the full dynamics onto the resulting bifurcation diagram reveals that the bursting mechanism is mediated by a saddle-node and a homoclinic bifurcation, and the bistable region in between, where a stable equilibrium coexists with a stable limit cycle.

> Then the slow system dynamics is explored with the fast system replaced by its averaged input. The slow dynamics is explained by the interplay of potassium outflow due to spiking and potassium returned to the cell by the pump. This periodic behavior is seen as a hysteresis loop between the SN and HOM bifurcations.

> Interestingly the authors show that the region of bistability bounded by the SN and the HOM exists because of the pump. Without the pump, the unfolding of the saddle-loop bifurcation is aligned with the [K+]_out axis such that a suitable bounded bistable region for the hysteresis loop to exist in is not formed.

> Next the diagram of the pump-skewed SL unfolding is zoomed out, placing the bursting region in a larger context. Behaviors in other regions are investigated; some are mediated by additional bifurcations that are found at the top of Figure 6A. These behaviors are realized with different values of the fixed parameter I_app. The new behaviors are solutions that converge to tonic spiking or to a stable equilibrium, which is interpreted as depolarization block.

> Finally, the role of the maximal pump current (interpreted as pump density) is investigated. As was noted, the pump skews and distorts the bifurcation diagram in Fig. 5a. Fig 7 shows that bursting is restricted to a region bounded by curves corresponding to the SNL point and the homoclinic bifurcation curve inflection point (viz. Fig 6); this essentially captures the skewing as a function of the pump. That’s ok, however, I started having trouble understanding the manuscript at this point. For example, I didn’t understand the text in lines 474-479. I think they try to explain the relationship of Fig. 7 to the other parameters, and the mention of chaotic behavior is neither here nor there.

Thank you ver much for pointing us to this indeed very confusingly written part of the manuscript. We have rewritten this paragraph (see lines 504 onwards in the new manuscript) and also added new pannesl to Fig 7 explaining the period doubling nature of the boundary between stonic spiking and bursting. The period doubling cascalde is a route to chaos as depicted in the phase planes on the bottom of Fig. 7.

> But continuing, I don’t understand what defines the upper curved border of the yellow bursting region. The caption says that above this boundary, the inter-spike interval distributions are no longer bimodal, whereas the text suggests that tonic spiking occurs in region E because a fixed point in the slow system appears (viz. Fig. 7C). Are these the same phenomenon?

We agree it was very hard to understand the nature of this boundary in the old version of the manuscript. We decided to explain it not only as properties of the slow system, but additionally as a property of the complete system. Taking this perspective, the upper curve delimiting the bursting area is in fact a period doubling cascade of the complete system. We kept the old explanation, for additonal insight on why the system in the region between HOM inflection point and SNL point of the fast subsystem is not bursting, given the fast subsystem has all the necessary bursting condition. We kept the slow fast analesys for a point in this region (point E in Fig 7A). The condition is that the slow subsystem should not have any stable fixed point in the bistable region of the fast subsystem so that the hystersis loop can form.

Technically, the priod doubling line of the complete subsystem can be predicted by finding the fixed point of the slowsubsystem. However, finding this point is numerically challenging. The assumption for the slow-fast analysis are being violated near the bifurcation point of the fast subsystem due to a critical slowing down, for example near the HOM of the fast subsystem. Therefore, finding out that the intersaction of the slow subsystem branch with the zero line is really a fixed point or if it is due to a breaking of the timescale sepration is challenging. So to find the actual boundary of the tonic spiking and bursting in Fig 7A, using the bifurcation analysi of the whole system is more robust. We did in fact both to corroborate the diagram.

>[mahraz] we also checked it we simulations. so the last sentece should be refer to it also

> Finally, it is briefly noted that by adding dynamic sodium concentrations, similar behavior is seen, but it “can introduce more intricate behaviours”.

> However, the bifurcation analysis of simplified systems with ion concentration dynamics isn’t new, and it has been known for a long time that ion concentration dynamics provides a mechanism for slow passage through fast-system bifurcations that leads to bursting. Thus I think the larger context and proper novelty of the present work is not presented clearly. It would be appropriate to frame the larger context with a (brief) discussion of the contributions from the authors’ references 12, 31, 33-37 and how the current work relates to these.

We have added a more extensive literature discussion in which we highlight the differences in our article, particularly recarding the influence of the electrogenic Na^+^/K^+^-ATPase on the known bifurcation structures.

> More specifically, the slow-fast analysis presented here has been reported before in a very, very similar model of a single neuron with ion concentration dynamics. The first part of the results section (“Fast subsystem bifurcation analysis”) reports an essentially identical approach as that reported in Barreto and Cressman, Journal of Biological Physics 37, 361-373 (2011). In addition, the analysis in the next section (“Slow system analysis”) was carried out for the same model in the authors’ reference 23. Of course, this does not necessarily disqualify the current manuscript; slow-fast analysis is a widely applicable tool. However, I would ask the authors to include a discussion comparing and contrasting their work to the reports cited above. For example, the earlier model includes potassium diffusion to a bath and both potassium and sodium dynamics, but leaves out the pump in the voltage equation. In contrast, the study in the current manuscript leaves out the potassium bath and the sodium dynamics but emphasizes the important role of the pump in generating bursting. Thus, there are notable similarities and differences in the models and their bifurcation structures, and there are similar and different insights about the dynamical bursting mechanisms identified.

We thank the reviewer for pointing us to the Barreto and Cressman article. We included a discussion of it in the introduction, lines 68 onwards. There are, however, substantial differences in the discussed mechanism. The burst mechanism discussed in our manuscript is not among the burst mechanisms in Barreto and Cressman (2011), all of which are slow wave bursters. In fact, the present hysterisis loop square wave burster is notfound in their model. This is mostly due to the fact that Barreto and Cressman (2011) ommitted the electrogenic pump current from their current balance equation that governs the memebrane voltage. In order to make our point clearer, we have included new simulations with and without (as in Barreto and Cressmann, 2011) in the results section in Figure 5 (A and B). These simulations demonstrate that by integrating the pump current into the membrane potential equation and subsequently inducing a shear in the fast subsystem bifurcation, the model could produce the burst phenomenon introduced in this paper. In contrast, when the pump only effects the potassium concentration dynamics as in Barreto and Cressman (2011), the square wave burst observed in our paper cannot be replicated.

> Minor edits:
> Line 287: fix “Solow-fast method”
> Line 300: “left to right” should be “right to left”, as shown in Fig. 4.
> In the text discussing Figure 6, it’s odd that panel C is discussed before panel B. Consider swapping these panels to better align with the sequence of ideas presented in the text.
> Line 652 and 762: remove apostrophe
> Line 693: State the I_app is fixed at a constant value
> Lines 773-774: “The Rung-Kunta” should be “The Runge-Kutta”

Thank you for highlighting these mistakes. We changed the manuscript according to the suggestions.