Analysis of Interspike-Intervals for the General Class of Integrate-and-Fire Models with Periodic Drive

We study one-dimensional integrate-and-fire models of the general type ẋ = F (t, x) and analyze properties of the firing map which iterations recover consecutive spike timings. We impose very week constraints for the regularity of the function F (t, x), e.g. often it suffices to assume that F is continuous. If additionally F is periodic in t, using mathematical study of the displacement sequence of an orientation preserving circle homeomorphism, we provide a detailed description of the regularity properties of the sequence of interspike-intervals and behaviour of the interspike-interval distribution.


Introduction
We consider one-dimensional integrate-and-fire (IF) models of the general type: where for the equation (1.1) we assume that for any initial condition x(t 0 ) = x 0 , t 0 , x 0 ∈ R, it has the property of existence and uniqueness of the solution x(t) (not necessarily defined on the whole R) and (1.2) is the resetting mechanism added to this continuous dynamics given by (1.1), explained below. In these systems, most commonly used as simplified models of neuron's activity (cf. [4,6,7,10,16]), the continuous dynamics of the variable x(t), standing for the membrane voltage of the cell, is interrupted by the so-called thresholdreset behaviour (1.2), which is supposed to mimic spiking (generation of action potential). However, IF systems (and circle mappings induced by them in case of periodic forcing) can also be used in modeling of cardiac rhythms and arrhythmias [2], in some engineering applications (e.g. electrical circuits of certain type, see [5]) or as models of many other phenomena, which involve accumulation and discharge processes that occur on significantly different time scales.
Integrate-and-fire models belong to the so-called hybrid systems, coupling continuous-time nonlinear dynamics with discrete events, which are ubiquitous in applied science. Beyond neuroscience, they have a broad range of applications (impact physical systems, control theory, economics. . . ) and the theory of hybrid dynamical systems has been the subject of recent advances in dynamical systems (see, for instance, [9]).
In integrate-and-fire models the dynamical variable x(t) evolves according to the differential equation (1.1) until it reaches the threshold-value x = x T , say at some time t 1 . Next it is immediately reset to the resting value x = x r and the system continues again from the new initial condition (t 1 , x r ) until possibly next time t 2 when the threshold is reached again, etc. Thus the forward-solutions x(t), t ≥ t 0 , of the system (1.1)-(1.2) with the initial condition (t 0 , x 0 ), x 0 < x T , are piece-wise differentiable functions x(t), which in the intervals (t n−1 , t n ) satisfy (1.1) and at the discontinuity points t n fulfill x(t n ) = x T and lim t→t + n x(t) = x r . Of course, it might happen that the solution starting from (t 0 , x 0 ) never spikes, i.e. always stays below x T and in this case the solution of the spiking system (1.1)-(1.2) is the continuous solution of (1.1).
We set for simplicity x r = 0 and x T = 1. It is also possible to consider varying (i.e. time-dependant) threshold and reset, which allows to introduce some other biologically realistic phenomena (such as refractory periods and threshold modulation [10]). However, from the mathematical point of view, often analysis of models with varying threshold and reset can be reduced through the appropriate nonlinear change of variables to studying the case of constant x r and x T (see e.g. [4]).
Below we introduce basic concepts necessary for our further considerations.
Note that the firing map Φ(t) does not need to be well defined for every t ∈ R since for some t it might happen that the solution x(·; t, 0) never reaches the value x = 1. Thus the natural domain of Φ is the set (compare with [6]): D Φ = t ∈ R : there exists s > t such that x(s; t, 0) = 1 .
Later on we will give sufficient conditions for the firing map Φ : R → R of the models considered to be well-defined (Lemma 1, Lemma 2 and Proposition 1).
The basic quantity associated with the integrate-and-fire systems is the firing rate: .
Its multiplicative inverse is in turn the average interspike-interval: Obviously, in general the limits above might not exist or depend on the initial condition (t 0 , 0).
In [16] the following observation for periodically driven models was made (the remark was not directly formulated in this way but it is a well-know fact): for all x and t), then the firing map Φ has periodic displacement Ψ := Φ − Id. In particular for T = 1 we have Φ(t + 1) = Φ(t) + 1 and thus Φ is a lift of a degree one circle map under the standard projection p : t → exp(2πıt).
In case of periodic forcing, the underlying circle map ϕ : S 1 → S 1 such that Φ is a lift of ϕ, is referred to as the firing phase map.
One of the most popular linear IF models is Leaky Integrate-and-Fire model which for σ = 0 reduces to the Perfect Integrator (PI): There are few purely analytical works on one-dimensional IF models and they include e.g. [4,6,10,20]. Combination of analytical and numerical approach to the firing map was taken, for example, in [7] (phase-locking and Arnold tongues), [16] (LIF model with sinusoidal input), [21] (LIF with periodic input) and [25] (LIF with periodic input and noise).
Analytical results concerning the firing map Φ were obtained assuming that F (t, x) is regular enough (always at least continuous) and often periodic in t. However, in [20] (and in the thesis [24]) it was proved that for the linear models, LIF and PI, many required properties (e.g. continuous dependence of the firing map on the input function f in appropriate topology), hold even if f is just locally-integrable (and thus not-necessary continuous).
In this paper we consider the whole class of IF models (1.1)-(1.2), thus including also non-linear systems. We show that even for such a general case, rigorous results on the behaviour of the firing map and interspike-intervals can be proven. We impose not so restrictive constraints for the regularity of the function F (t, x) (as for instance the authors of [6], who considered analytical functions F ) since usually, as we will see, it suffices to assume that F is continuous and such that (1.1) has the property of existence and uniqueness of solutions (thus e.g. that F is Lipschitz in x).
Firstly, we prove general properties of the firing map arising from the systeṁ x = F (t, x), where the function F is continuous and positive. Then assuming that F is also periodic in t, we give a detailed description of the sequence of interspike-intervals and interspike-interval distribution. Here we make use of mathematical result concerning displacement sequence of an orientation preserving circle homeomorphism proved in [19] and [24].
Although for the firing map of systems with continuous (and sometimes also periodic) drive some rigorous results have been proved (e.g. in [4,6,10]), the sequence of interspike-intervals even in such a case has not been investigated in detail yet. It is worth mentioning that sometimes the sequence of interspike-intervals might be of greater importance than the exact spiking times themselves [23] since interspike-intervals are said to take part in information encoding by neurons (see e.g. [11] and references therein).
This work together with our previous paper [20] provides detailed and complete description of interspike-intevals for periodically driven one-dimensional integrate-and-fire models. This is the natural extension of [20] where we studied only linear models.
The main results, under assumptions (1.)-(3.) x) is periodic in t (with period T = 1 without the loss of generality), are Theorems A and B on the interspike-interval distribution, proved in part 4: where the functions F and F n satisfy (1.)-(3.), have continuous and uniformly bounded partial derivatives with respect to x on U = R × [0, 1] (i.e. there exists K such that ∂F ∂x (t, x) < K and also ∂Fn ∂x (t, x) < K for (t, x) ∈ U, at least for n sufficiently large) and moreover F (t, x) > m > 0 for some m. Suppose that all the induced firing maps Φ n and Φ have irrational rotation numbers, n and , respectively. By µ (n) ISI and µ ISI denote the interspike-interval distributions, correspondingly for Φ n and Φ, with respect to the corresponding invariant measures µ (n) and µ of their firing phase maps. If Theorem B. Consider the integrate-and-fire systemsẋ = F Θ1 (t, x) andẋ = F Θ2 (t, x), where F Θi again satisfy (1.)-(3.) and F Θi and ∂x are continuous and bounded on U = R×[0, 1], i = 1, 2. By Φ Θ1 and Φ Θ2 denote the firing maps emerging from the corresponding systems. Suppose that the rotation number associated with Φ Θ1 is irrational and that F Θ1 (t, x) > m > 0.
For any ε > 0 there exists a neighborhood N of F Θ1 in C 0 (R 2 ) (in sup R 2topology) such that if F Θ2 ∈ N , then for every initial condition (t, 0) we have where ω (Θ2) n,t is the empirical interspike-interval distribution for the run of the systemẋ = F Θ2 (t, x) starting from (t, 0) and µ (Θ1) ISI is the interspike-interval distribution forẋ = F Θ1 (t, x) with respect to its invariant measure µ (Θ1) .
We remark that for the other choice of values x r and x T , the above theorems remain true under replacement of U = R × [0, 1] by U := R × [x r , x T ]. It is also worth pointing out that although we assume that the partial derivatives of F and F n (at least with respect to x) are continuous, we only require the convergence F n → F in C 0 (R × R). Similarly, it is not necessary that the functions F n and F are bounded on the whole plane R 2 .
These results are at the end illustrated by the numerical example involving quadratic integrate-and-fire model.

General Properties of the Firing Map
For the linear integrate-and-fire models such asẋ = −σx + f (t) we are able to easily prove the required properties of the firing map Φ. However, when we deal with non-linear ones, different methods shall be used since usually we cannot write down the solution x(t) of the differential equation in the explicit form.
One of the most popular non-linear models is the quadratic integrate-andfire (QIF) model:ẋ = x 2 + f (t), (2.1) in which the threshold value x T should be ∞ (since x(t) blows to ∞ in finite time) but in practise we set x T big but finite value. Since the equation (2.1) is a Riccati type equation, sometimes we can solve it, get the implicit formula for the firing map and use similar technics as in case of the LIF-model in [20] to prove desired properties of the firing map and the sequence of interspikeintervals. However, we will not discuss QIF in details. Instead, we consider the general case ofẋ = F (t, x).
To this purpose we can use, for example, a number of theorems on the continuous dependence of solutions on parameters. In this work we extent results for the LIF and PI models obtained for locally integrable input functions f (t) in [20] to the general modelẋ = F (t, x) at the cost of more restrictive constraints on F . Here the technical problem remains often to assure that the firing map Φ(t) is defined for all t ∈ R.
The following theorem was proved in [6] ("Injectivity theorem" and "Continuity theorem"): Theorem 2. [6] Assume that Ω is a region that contains the strip {−∞ < t < ∞, 0 ≤ x ≤ 1} and that the function F : Ω ⊂ R 2 → R is analytic. Then This theorem allows to determine whether the firing map Φ is continuous and injective or not. However, the assumption that F is analytic seems too restrictive and we will show that if we assume instead that F is continuous and positive, then Φ is continuous and injective as well.
Another important thing is to assert that the firing map Φ : R → R is well defined, which is necessary if we want later on Φ to be a lift of the circle homeomorphism. In [4] the integrate-and-fire modelsẋ = F (t, x) were investigated under the hypothesis I or II: Calling a run a spiking trajectory x(·; t 0 , 0) of the systemẋ = F (t, x) starting from (t 0 , 0), the author proves that when at least one of the conditions I or II is satisfied, then Notice that the fact that for every t 0 ∈ R the run x t0 := x(·, t 0 , 0) has infinitely many spikes, means in particular that Φ(t 0 ) ∈ R is defined for all t 0 ∈ R. However, in our approach in order to assure that Φ : R → R is well defined, in general we do not assume that I. or II. are satisfied. Instead we provide the following sufficient condition for sustained firing: x) is continuous on R 2 and such that the eq. (1.1) for every initial condition (t 0 , x 0 ) ∈ R 2 has the unique solution x(t; t 0 , x 0 ) defined on the whole R. Further, assume that F (t, x) satisfies Then the firing map Φ : R → R is well-defined.
Proof. Choose t 0 ∈ R. The solution of (1.1) with the initial condition (t 0 , 0) satisfies ) is continuous and thus also locally integrable. Now if lim sup t→∞ Remark 1. We added additional constraint that the solution x(t; t 0 , x 0 ) of the initial value problem is defined on the whole R in order to guarantee that it does not blow up in finite time (to +∞ or −∞) and hence that x(·; t 0 , x 0 ) : R → R is well-defined and continuous on the whole real line, which enabled us to make use of the assumption (2.2).
Note that, for instance, by Picard-Lindelöf Theorem (see e.g. Theorem 1.1 in [14]), every solution of (1.1) with given initial condition can be extended onto We remark that the condition (2.2) is not equivalent to Indeed, is an example of a continuous function which satisfies (2.3) but does not satisfy (2.2) (take for instance x(t) = t). These two conditions are not equivalent even for functions F periodic in t. The function The condition (2.2) might not be so easy to verify but we formulate another condition which is more effective: In particular, ω(u) might be any locally integrable function on [0, ∞), which is positive (at least for u > 0) and goes to 0 monotonically with u → ∞ (at least for u > 0), but slowly enough so that ∞ 0 ω(u) du = ∞. We do not demand that ω is continuous. Examples include is non-decreasing and, consequently, either x(t) blows to +∞ in finite time or x(t) is defined on the whole R. The first possibility immediately asserts that where the last equality is explained by the fact that for ω locally integrable, lim t→∞ t 0 ω(u) du = ∞ means that lim t→∞ t a ω(u) du = ∞ for any a > 0. Thus lim t→∞ x(t) = ∞ and consequently there exists t * such that x(t * ) = 1. This contradicts the assumption that Φ(t 0 ) is not defined and the proof is completed.
However, if additionally F is periodic in t the sufficient condition for the firing map is even simpler: Proof. Let x(t) be the solution of the differential equationẋ = F (t, x) with the initial condition (t 0 , 0). Suppose that Φ(t 0 ) is not defined. Again, this means that 0 ≤ x(t) < 1 for every t ≥ t 0 . Thus due to the periodicity of F in t and the fact that x(t) is bounded we have: ). Consequently, the solution x(·) increases at least by ς on the interval [t 0 , t 0 + 1]. Similarly, we justify that x(·) increases at least by nς on [t 0 , t 0 + n], n ∈ N, which contradicts the assumption that x(t) is bounded for t > t 0 . It follows that x(t) is not bounded (in particular not bounded from above) and Φ(t 0 ) is defined.
by definition of the firing map and the fact that the increasing solutions x(·; t 1 , 0) and x(·; t 2 , 0) cannot cross. If Φ(t 1 ) ≤ t 2 then also Φ(t 2 ) > Φ(t 1 ) since always Φ(t) > t. This shows that the firing map is increasing in D Φ and ends the proof.
Proposition 3. Let F (t, x) be continuous and positive on R 2 . Suppose that the firing map Φ : R → R is well-defined. Then Φ is continuous.
In proving the continuity of the firing map we will make use of the continuous dependence of the solutions on initial conditions which we formulate as the following lemma: and that for every t 0 ∈ J, the unique solution x(·; t 0 , 0) with the initial condition (t 0 , 0), exists on some open interval Ω t0 ⊃ I. Then the map t → x(·; t, 0) is continuous from J to C(I), where C(I) denotes the space of all continuous functions on I with the uniform norm.
Proof. The statement is a consequence of the continuous dependence of the solutions of differential equations on initial conditions (which is uniform continuity on closed intervals, see e.g. Theorem 2.1 in [14] or, alternatively, the proof of Lemma 8.10 in [6]).
By definition of the firing map Now fix ε > 0. From Lemma 3 it follows that x(Φ(t 0 + δ); t 0 , 0) − x(Φ(t 0 + δ); t 0 + δ, 0) < ζε for δ small enough. Finally, Φ(t 0 + δ) − Φ(t 0 ) < ε and Φ is right continuous at t 0 . In the same way we show that Φ is left continuous and thus Φ is continuous at t 0 . Since t 0 was arbitrary, the statement follows. Proof. Under these assumptions the firing map Φ : R → R is well-defined. It is continuous, increasing and satisfies Φ(t+T ) = Φ(t)+T on the account of Fact 1. Consequently, it is a lift of an orientation preserving circle homeomorphism (cf. [15,Chapter 11]).
In what follows we assume that is periodic in t (with period T = 1 without the loss of generality).

Regularity properties of the sequence of interspike-intervals
By ϕ : S 1 → S 1 we denote the underlying circle homeomorphism of the firing map Φ, i.e. ϕ ∼ Φ mod 1 and ϕ is called the firing phase map.
Following the classical Poincaré rotation theory, the rotation number n of Φ at t exists and does not depend on t (cf. Proposition 11.1.1 in [15]). Moreover, we can set := Φ (t) mod 1 to be the unique rotation number associated with the firing phase map ϕ. It is then uniquely defined, irrespectively of the choice of t and the lift Φ. If the rotation number of the firing map Φ (which is exactly the average interspike-interval) is rational, then ϕ necessarily has periodic orbits (cf. Proposition 11.1.4 in [15]). However, unless it is not strictly conjugated to the rational rotation, there are also non-periodic orbits. Nevertheless, these non-periodic orbits are attracted to the periodic ones (Proposition 11.2.2 in [15]) which leads to the following theorem: Proposition 5. If the rotation number = p/q of the firing map is rational, then the sequence of interspike intervals is asymptotically periodic with period q, i.e. for every initial condition (t 0 , 0) we have We only remark that when ϕ is conjugated to the rational rotation, then the sequence ISI n (t 0 ) is purely periodic for every t 0 . However, determining for which parameter values λ a given modelẋ = F λ (t, x) is conjugated to the rotation is often not so easy to settle (see discussion in [4] for the LIF model).
The asymptotic q-periodicity with the winding number p which is stable under small change of parameters corresponds to the so-called phase-locking [7,10,21].
The case of irrational rotation number is more interesting from the mathematical point of view but before we will discuss it, we make the simple observation: Remark 2. If ϕ is a rigid rotation r by = (Φ), then the sequence ISI n (t) is constant: ISI n (t) = .
In this special case, which holds when F (t, x) =:F (x) does not depend on t, the rotation number does not influence periodicity of the sequence of interspikeintervals, i.e. the sequence is periodic, no matter whether is rational or not. The above remark follows from the fact that the rotation is an isometry and in this case Φ(t) = t + (up to mod 1) for every t.
In order to make the forthcoming considerations more clear we remind that when ϕ : S 1 → S 1 has irrational rotation number and is transitive (i.e. there exists an orbit {ϕ n (z 0 )} n∈N dense in S 1 ), then it is conjugated to the rotation by , meaning that for some orientation preserving circle homeomorphism γ we have γ • ϕ = r • ϕ (in fact, γ is unique up to an additive constant in the lift). When it is not transitive, it is still semi-conjugated to the rotation, which means that γ • ϕ = r • ϕ holds but the lift Γ of γ is a continuous, surjective, non-decreasing and non-invertible function (thus γ is not a homeomorphism; see Theorem 11.2.7 in [15]). Independently of its transitivity, ϕ is metrically isomorphic to the irrational rotation. The unique (up to normalization) invariant ergodic measure µ is the Lebesgue measure transported by γ where γ is the (semi-)conjugacy: µ(A) = Λ(γ(A)). In non-transitive case µ is concentrated on the minimal set ∆ ⊂ S 1 and the semi-conjugacy γ becomes a strict conjugacy on the set ∆ ⊂ ∆ which is the set ∆ \ E for some countable subset E (we might think of ∆ as of ∆ with identified endpoints of the intervals complementary to ∆; cf. [15, pp. 399]).
Theorem 5. Consider the modelẋ = F (t, x) where F is as above and (Φ) ∈ R \ Q. Then 5.1 if ϕ is transitive, then the sequence of interspike-intervals {ISI n (t 0 )} for every t 0 ∈ R is almost strongly recurrent, i.e.
Proof. Since the firing phase map ϕ is a circle homeomorphism with the irrational rotation number, the statement follows from the fact that every point z ∈ S 1 for ϕ transitive (similarly, every point z ∈ ∆ for non-transitive case) is almost periodic under the dynamical system (ϕ, S 1 ). This was shown in [19] basing on the fact (the classical result proved in [13]) that every point x ∈ X is almost-periodic under ϕ, if X is a compact metric space and X is minimal for ϕ (the set is minimal if it is non-empty, closed, invariant and such that no proper subset of it shares all these properties).
The property of almost strong recurrence of the sequence ISI n (t 0 ) means that for every ε and every n the sequence of indexes r k = k + i such that the element ISI n+r k (t 0 ) equals ISI n (t) with ε-accuracy, has gaps bounded by some N and that this constant N can be chosen uniformly for all n.
Notice also that ϕ with irrational rotation number is transitive, whenever F ∈ C 2 (R 2 ) (positive, periodic) since then also ϕ ∈ C 2 (S 1 ) and on the account of Denjoy'a Theorem ϕ is conjugated to the rotation r (compare with the proof of Proposition 8).

Interspike-interval distribution
In the remaining part of the article we will work with the case of irrational rotation number. Let then µ be the unique invariant probability measure of the firing phase map ϕ. Note that µ gives the distribution of phases of firing times, i.e. it is the limit distribution of points Φ n (t 0 ) with n → ∞ for a.e. initial condition t 0 (by distribution we mean a normalized measure).
Definition 2. The interspike-interval distribution µ ISI with respect to the invariant measure µ is defined as follows: Note that since Φ mod 1 is periodic with period 1, we consider only t ∈ [0, 1]. Moreover, although the measure µ has support contained in [0, 1], as it is the invariant measure for Φ mod 1 : [0, 1] → [0, 1], the measure µ ISI has support equal to Ψ ([0, 1]), which in general might not be contained in [0, 1] but it is always contained is some interval of length not greater than 1 (Φ maps intervals of length 1 onto intervals of length 1 due to the fact that Φ(t + 1) = Φ(t) + 1 for every t and the resulting interval, containing supp(µ ISI ), is shifted by a from its projection mod 1 into [0, 1], where a > 0 is such that Φ(0) = a). We can consider µ Ψ (A), where A is an arbitrary subset of R, if we adopt the convention that µ ISI is defined on the whole R but it simply vanishes everywhere in R outside its support.
In [19] we proved many properties of the displacement distribution µ Ψ of an orientation preserving circle homeomorphism ϕ with an irrational rotation number which is defined exactly as in Definition 2, up to taking mod 1. Thus interspike-interval distribution µ ISI of the firing map Φ has exactly the same properties as the displacement distribution µ Ψ .
Moreover, when one takes t 0 ∈ R \ ∆, then for every t ∈ ∆ there exists an increasing sequence n k , k ∈ N such that for every l ∈ N we have Proof. The theorem is an immediate consequence of Proposition 2.1 in [19].
With the use of Birkhoff Ergodic Theorem and unique ergodicity of the firing phase map ϕ = Φ mod 1 we justify Proposition 7. Under the assumptions of Proposition 6 (regardless the transitivity of ϕ), for A ⊂ R we have where denotes the number of elements of the set, and the above convergence is uniform (with respect to t ∈ R). The average interspike interval aISI (which equals the rotation number (Φ)) is the mean of the distribution µ ISI : Under stronger assumption on the regularity of F we prove: Proof. By the Implicit Function Theorem, if F is C k , then Φ is C k in the neighbourhood of any t such that F (Φ(t), 1) > 0 (compare with Theorem 2 in [4]). Thus under stated assumptions Φ is locally C 2 but since Φ is also invertible and F is periodic in t, Φ is in fact a lift of C 2 -diffeomorphism ϕ : S 1 → S 1 . Due to Denjoy Theorem [8] ϕ is transitive and on the grounds of Proposition 6 the statement follows.
Proof. Since x(t) and y(t) are the solutions of the corresponding differential equations, for α < t < β we have where the integral t t0 F (u, x(u))du is to be understood as a vector (correspondingly for t t0 G(u, x(u))du). Hence Using Gronwall's inequality we obtain Now we can prove the following Theorem 7. Let some constant K > 0 be fixed. Consider the integrate-and-fire systemẋ = F (t, x), where the function F : R 2 → R satisfies the assumptions (1.)-(3.), the hypothesis of Theorem 6 with K and U = R × [0, 1]. The mapping F → Φ, from a subset of such functions F in C 0 (R × R), is continuous from sup R 2 -topology into sup R -topology at every point F , where F (t, x) > m > 0 on R 2 for some m > 0.
By sup R 2 -topology (and similarly for sup R ) we mean the topology induced by the metric d sup R 2 (F, G) := sup (t,x) |F (t, x) − G(t, x)|. The above theorem says that the firing maps Φ F and Φ G arising from the systemsẋ = F (t, x) andẋ = G(t, x), respectively, are arbitrarily close in uniform (sup R ) topology whenever F and G are sufficiently close, provided that F and G fulfill certain conditions.
Proof of Theorem 7. Suppose that F and G satisfy (1.)-(3.) and the assumptions of Theorem 6 with the constant K and additionally that F (t, x) > m > 0. Let sup (t,x)∈R 2 |F (t, x) − G(t, x)| < δ ≤ min{1/m, m/2}. Choose t 0 ∈ R and assume firstly that Φ F (t 0 ) > Φ G (t 0 ). Since x(t) = t t0 F (u, x(u)) du and y(t) = t t0 G(u, y(u)) du, where x(t) and y(t) are the solutions of the equationṡ x = F (t, x) andẋ = G(t, x) satisfying x(t 0 ) = y(t 0 ) = 0, by the definition of the firing map we get as follows from (3.2) since it is enough to consider the solutions x(t) and y(t) on the common interval of length less than 2/m, which includes all the points t 0 , Φ G (t 0 ) and Φ F (t 0 ) (we might estimate that Φ F (t 0 ) − t 0 < 1/m for every t 0 ∈ R using the definition of Φ and the fact that F (t, x) > m for all t and x; by the same argument we have Φ G (t 0 ) − t 0 < 2/m). Combining this with (3.3) we obtain that Thus δ ≤ min{ 1 m , m 2 , m 2 ε 4e K/m } yields |Φ F (t 0 )−Φ G (t 0 )| < ε and the statement follows from arbitrary choice of t 0 .

Main Results: Approximation of the Interspike-Interval Distribution
In this part we are concern with interspike-interval distribution for parameterdependant IF systems. We will see what kind of assumptions guarantee that the distribution µ ISI changes continuously with parameters. For this purpose we recall the notion of the weak convergence of measures: Definition 3 [see e.g. [3]]. Let X be a complete separable metric space and M(X) the space of all finite measures defined on the Borel σ-field B(X) of subsets of X.
A sequence µ n of elements of M(X) is called weakly convergent to µ ∈ M(X) if for every bounded and continuous function f on X We denote the weak convergence as µ n =⇒ µ.
One can show (cf. [22]) that µ n =⇒ µ if and only if for each continuity set A of µ, lim µ n (A) = µ(A) (a Borel set A is said to be a continuity set for µ if A has µ-null boundary, i.e. µ(∂A) = 0).
Next theorems can be proved by the arguments involving the continuity of the mapping F → Φ within some class of functions F , asserted by Theorem 7.
Theorem A. Consider the systemsẋ = F (t, x) andẋ = F n (t, x), n ∈ N, where the functions F and F n satisfy (1.)-(3.), have continuous and uniformly bounded partial derivatives with respect to x on U = R × [0, 1] (i.e. there exists K such that ∂F ∂x (t, x) < K and also ∂Fn ∂x (t, x) < K for (t, x) ∈ U, at least for n sufficiently large) and moreover F (t, x) > m > 0 for some m. Suppose that all the induced firing maps Φ n and Φ have irrational rotation numbers, n and , respectively. By µ (n) ISI and µ ISI denote the interspike-interval distributions, correspondingly for Φ n and Φ, with respect to the corresponding invariant measures µ (n) and µ of their firing phase maps. If Proof. By Theorem 7, if F n → F in C 0 (R × R), then also Φ n → Φ, where Φ n and Φ are the corresponding firing maps, being the lifts of corresponding orientation preserving circle homeomorphisms (firing phase maps) ϕ n → ϕ with irrational rotation numbers. But since ϕ n → ϕ in C 0 (S 1 ), also the displacement distributions converge weakly: µ (n) Ψn → µ Ψ (this is directly asserted by Proposition 2.8 in [19]). Since distributions µ (n) Ψn and µ Ψ correspond to interspike-interval distributions, the statement follows.

Empirical approximation of the interspike-interval distribution
In next, we would like to tackle the case when the IF system with irrational firing rate is approximated by another system, for which we do not demand that the rotation number is irrational too. To formulate the result in a mathematically rigorous way, we introduce the concepts of the Fortet-Mourier metric and of the empirical interspike-interval distribution: Definition 4. Let µ and ν be the two Borel probability measures on a measurable space (Ω, F), where Ω is a compact metric space. Then the distance between the measures µ and ν is defined as Definition 5. Let Φ be the firing map arising from the IF systemẋ = F (t, x). Choose the initial condition (t, 0) (x r = 0). Then the empirical interspikeinterval distribution for the run of length n (i.e. having n-spikes) starting at (t, 0) equals where δ ISIi(t) is the Dirac delta centered at the point ISI Thus Note that if ϕ with rotation number is close in C 0 (S 1 )-metric to ϕ with rotation number , then the rotation numbers and are also close due to the continuity of the rotation number in C 0 (S 1 ) (cf. Proposition 11.1.6 in [15]) but, of course, in general might be rational, even if we take irrational .
For any ε > 0 there exists a neighborhood N of F Θ1 in C 0 (R 2 ) (in sup R 2topology) such that if F Θ2 ∈ N , then for every initial condition (t, 0) we have where ω (Θ2) n,t is the empirical interspike-interval distribution for the run of the systemẋ = F Θ2 (t, x) starting from (t, 0) and µ (Θ1) ISI is the interspike-interval distribution forẋ = F Θ1 (t, x) with respect to its invariant measure µ (Θ1) .
Proof. Now this result is a direct consequence of Theorem 2.17 in [19] about approximation in the Fortet-Mourier metric by the sample displacement distribution of a homeomorphism which is just close enough to the given homeomorphism.
As the convergence under Fortet-Mourier metric implies weak convergence of measures [12], we conclude: Corollary 1. Under the notation as in Theorem B, for every t ∈ R we have Theorem A says that the interspike-interval distribution with respect to the invariant measure behaves continuously (in terms of the weak convergence of measures) with the small change of the input in uniform topology (in fact the general form of the equationẋ = F (t, x) in our considerations allows to change the input and the intrinsic dynamics of the cell as well, i.e. simply we can make an acceptable change of the right-hand side of the differential equation). Note that we can speak formally of the interspike-interval distribution µ ISI only when the unique invariant probability measure µ exists, which means only in case of the irrational rotation number. However, basically we are often not able to compute directly the rotation number and due to numerical approximations we work in fact with systems with rational rotation number (rational firing rate), which well approximates the entire irrational rotation number. Nevertheless, Theorem B assures that empirically obtained interspike-interval distribution by tracking a trajectory of the "close" system (for which we admit irrational and rational rotation number, as well), well approximates (in Fortet-Mourier metric) the theoretical interspike-interval distribution for the given system.
In order to illustrate this phenomenon we encourage the reader to look, for example, at the histograms of firing phases and interspike-intervals in section 5 of [16] for the model du dτ = −σu + S(1 + B cos τ ) with fixed σ and S and varying B, such that the induced firing maps are homeomorphisms. Indeed, for B = 0 we have irrational rotation (implying, in particular, that the rotation number for B = 0 is irrational) and for B = 0.1 and B = 0.2 we see that the histograms evolve continuously.
Below we provide another simple example when the mathematical scheme developed in this work allows to anticipate the phenomena observed in the model.
Example. Let us consider the quadratic integrate-and-fire with periodic forcing: where a, λ and κ are real parameters. We fix a = 0.2 and κ = 0.5 and consider λ ∈ [0, 0.5]. Notice that for any λ < 0.5 it holds that F λ (t, x) := ax 2 + λ sin(2πt) + κ ≥ κ − λ > 0, thus the function F λ is separated from 0. Simultaneously, for any choice of reset x r and threshold We take x r = 0 and x T = 1. However, one can also consider taking bigger value of the threshold potential since the spiking solution of (4.1) has many subthreshold oscillations of the potential x(t) and thus taking higher value of x T allows to better distinguish between spikes and oscillations of relatively smaller amplitude. For any λ < 0.5 the function F λ fulfills conditions (1.)-(3.). Moreover, F λ and its partial derivatives are bounded on the sets R × J , where J is an arbitrary bounded interval. We notify that the parameter value λ = 0.5 is a border value for these assumptions since at λ = 0.5 F λ (t, x) can vanish for x = 0. Nevertheless, we decided to include this case in the below simulations, where we plotted histograms of interspike-intervals for different values of λ. In particular, for λ = 0 by solving the equation directly we are able to compute exactly the rotation number = √ 10 arctan( √ 2/ √ 5) ≈ 1.7833, which is equal to the average interspike-interval. For this parameter value the induced firing map Φ is simply a lift of the irrational rotation by mod 1 and thus phases of firing times are uniformly distributed in [0, 1], whereas the distribution of interspike-intervals is a Dirac delta centred at (see Figure 1). According to Theorems A and B as λ → 0 we should have ω λ t =⇒ µ 0 ISI , where ω λ t := lim n→∞ ω λ n,t are the empirical interspike-interval distributions obtained for the spiking solution ofẋ = F λ (t, x) starting at (t, x r ). Note that for λ = 0 we can define the formal distribution µ 0 ISI as we are able to assure that the rotation number is irrational so that the unique invariant measure µ 0 , giving distribution of points of the orbits, exists.
In Figure 1 we include numerically computed histograms of interspikeintervals for different values of λ which indeed evolve "continuously" when varying λ.
Remark 3. Note that due to the properties of an orientation preserving circle homeomorphism, corresponding statements of Proposition 5, Theorem 5 and Theorems A and B are also true for the sequence of firing phases {Φ n (t 0 ) mod 1}.
Another question which one can ask is concern with the absolute continuity of the interspike-interval distribution with respect to the Lebesgue measure, i.e. with the existence of the density of measure µ ISI . We recall that even if the invariant measure µ, which gives the distribution of firing phases, is absolutely continuous, it might not be so for the interspike-interval distribution µ ISI as happens for example for the systemẋ =F (x) + I with constant input I, where the firing map Φ is a lift of the rotation and the interspike-interval distribution is Dirac delta centered at , being the rotation number of Φ (and thus µ ISI is singular). We formulate the following theorem which is just a paraphrase of its correspondence for the displacement sequence of a circle homeomorphism (diffeomorphism), Theorem 2.14 in [19], and thus will be stated without the proof: Theorem 8. Suppose that the firing map Φ arising from the systemẋ = F (t, x) is a C 1 -diffeomorphism with irrational rotation number , which is conjugated with the translation by via a C 1 -diffeomorphism Γ and that the set where the latter is well-defined almost everywhere in supp(µ ISI ), i.e. in supp(µ ISI ) \ V (C Ψ ), with V (C Ψ ) denoting the set of critical values of Ψ [0, 1].
The above formula is most effective for the simplest linear model, i.e. the Perfect Integratorẋ = f (t), where we have the exact formula for the conjugating homeomorphism Γ (see [4]). However, the real content of this theorem, which is the existence of the density for µ ISI , might be applied for every type of IF model for which we are able to verify the condition on the zero measure of the set of critical points of the displacement function of the firing map.

Discussion
We have shown many specific properties of the interspike-interval sequence arising from general class of periodically driven integrate-and-fire models for which the emerging firing phase map is a circle homeomorphism. Although the obtained results might not seem surprising to the specialists in computational neuroscience, their implications, especially Theorems A and B on stability of the interspike-interval distribution, extend some numerical results on this distribution presented by some other authors (see e.g. ISI histograms in [16]) that so far have lacked analytical explanation.
The natural extension of this research includes first of all obtaining such rigorous results on interspike-intervals for periodically driven integrate-and-fire models for which the firing phase map ϕ is not necessary a homeomorphism, but just a circle map (even for continuous circle mappings one might expect significantly different properties due to the fact that in this case we have rotation intervals instead of the unique rotation number). Further development, mathematically more challenging, would be a detailed description of the interspikeinterval sequence for IF systems with an almost periodic input, for which a formal framework was prepared in [18] and [24]: Namely, we proved that in linear systemsẋ = −σx + f (t) the almost periodic input function f (t) (in the sense of Stepanov, to include also not continuous input functions) induces the firing map with uniformly almost periodic displacement. This is generalization of the Fact 1 which enables, perhaps, analysis of the firing map and interspikeintervals arising from almost-periodically driven integrate-and-fire models with the use of locally compact topological groups (with the special case of S 1 for purely periodic input) due to the Bohr compactification of the reals [17]. This is not only a theoretical task, since almost periodic inputs include sums of periodic inputs with rationally incommensurable periods and thus could correspond to different periodic signals which a neuron is receiving from different sources (from presynaptic neurons or injected electrodes).
Systems which recently attract a lot of attention, are bidimensional IF models (see, for instance, [26]). They could be another field of investigation in terms of exact results on the interspike-interval sequence.