Analysis of a SIR model with pulse vaccination and temporary immunity: Stability, bifurcation and a cylindrical attractor

Abstract

A time-delayed SIR model with general nonlinear incidence rate, pulse vaccination and temporary immunity is developed. The basic reproduction number is derived and it is shown that the disease-free periodic solution generically undergoes a transcritical bifurcation to an endemic periodic solution as the vaccination coverage drops below a critical level. Using numerical continuation and a monodromy operator discretization scheme, we track the bifurcating endemic periodic solution as the vaccination coverage is decreased and a Hopf point is detected. This leads to a bifurcation to an attracting, invariant cylinder. As the vaccination coverage is further decreased, the geometry of the cylinder contracts along its length until it finally collapses to a periodic orbit when the vaccination coverage goes to zero. In the intermediate regime, phase locking on the cylinder is observed.

Introduction

Pulse vaccination is a disease control policy under which at certain times, a portion of the population is vaccinated en-masse. It has been argued empirically and verified analytically that pulse vaccination might be more effective than continuous vaccination in preventing epidemics that exhibit seasonality, such as measles [1], [2]. Since then, the impact on pulse vaccination has been studied in ever more complex models of disease transmission. For instance, finite infectious periods [3], saturation incidence with latent period and immune period [4], incubation period [5], force of infection by distributed delay [6], nonlinear vaccination [7], quarantine measures [8] and stochastic effects [9] have been considered.

Dynamical analysis of these pulsed vaccination models often include stability criteria for the disease-free equilibrium or periodic orbit, effectively providing a proxy for the basic reproduction number. However, due to the presence of the impulse effect, establishing the existence of an endemic periodic orbit is much more difficult. When there are no delayed terms, methods of bifurcation theory for discrete time systems have been used to prove the existence of endemic periodic orbits from bifurcations at disease-free states; see [7], [10], [11], [12], [13] for some recent examples. In contrast, when delays are present, most analytical studies prove only permanence when $R_0>1$, which means that the disease persists. Numerical simulations are needed to obtain further detail, and this provides only a heuristic description of the orbit structure at a possible bifurcation point. We refer the reader to [3], [4], [5], [8] for examples. It would therefore appear that it is not for lack of interest that no authors have studied bifurcations in pulsed vaccination models involving delays, but rather that bifurcation theory techniques such as centre manifold reduction [14] have only recently been developed for impulsive functional differential equations.

Restricting to SIR models without pulse vaccination specifically, there are many papers that consider various forms of population dynamics and their interplay with delayed effects. Since endemic equilibrium points are often analytically available, Hopf bifurcations can often be studied analytically without the aid of numerical methods. One may consult [15], [16], [17], [18], [19] for some recent examples of this.

It is our goal to use centre manifold theory for impulsive delay differential equations [14], [20] to obtain more precise information about the orbit structure in a particular pulsed SIR vaccination model involving delay. Our starting point is the model of Kyrychko and Blyuss [21]:

$$\dot{S}=\mu -\mu S-\eta f\left(I\left(t\right)\right)S\left(t\right)+\gamma I(t-\tau )e^{-\mu \tau }$$$$\dot{I}=\eta f\left(I\left(t\right)\right)S\left(t\right)-(\mu +\gamma )I\left(t\right)$$$$\dot{R}=\gamma I\left(t\right)-\gamma I(t-\tau )e^{-\mu \tau }-\mu R\left(t\right).$$

Here, $f\left(I\right)$ is a general nonlinear incidence rate, infected individuals clear their infection at rate $\gamma$ and acquire temporary immunity of length $\tau$, $\eta$ is a recruitment rate and $\mu$ is a natural death rate, with birth rate scaled accordingly so that $N\left(t\right)=S\left(t\right)+I\left(t\right)+R\left(t\right)$ approaches unity as $t\to \infty$. The incidence rate is assumed to satisfy the properties: $f\left(0\right)=0$, $f^{\prime }\left(0\right)>0$, $f^{\prime \prime }\left(0\right)<0$ and ${lim}_{I\to \infty }f\left(I\right)=c<\infty$. Kyrychko and Blyuss [21] proved global stability of the disease-free equilibrium when $R_0<1$ for arbitrary nonlinear incidence satisfying the previous conditions, and considered the existence and stability of an endemic equilibrium for the particular incidence $f\left(I\right)=I∕(1+I)$. They numerically observed Hopf bifurcations at this equilibrium upon varying the immunity period $\tau$. Soon after, Jiang and Wei [17] proved that the endemic equilibrium may indeed undergo a Hopf bifurcation, by taking $\eta$ as a bifurcation parameter.

We here extend the model of Kyrychko and Blyuss to include pulse vaccination. To do this, we make the following assumptions.

• At specific instants of time $t_k$ for $k\in \mathbb{Z}$, any individuals that received their vaccine at time $t_k-\tau$ and are still alive lose their immunity and re-enter the susceptible cohort, at which point a fraction $v\in [0,1)$ of the total susceptible cohort is vaccinated.

• Vaccinated individuals are immune to infection for a period $\tau$ (the same immunity period as having recovered from infection) and are subject to the same natural death rate $\mu$.

• The sequence of vaccination times is periodic with shift of $\tau$: there exists $q>0$ such that $t_{k+q}=t_k+\tau$ for all $k\in \mathbb{Z}$.

The interpretation of (3) is that the period of the pulse vaccination schedule is synchronized with the immunity period. This seems reasonable for seasonal flu epidemics, for example, where vaccination effort should typically be focused near the beginning of flu season [22]. We therefore make the approximation that the vaccination period is modelled as a single pulse. With these assumptions in place, the pulse vaccination model takes the following form, where we will ignore the recovered ($R$) component since it is decoupled from the remaining equations:

$$\dot{S}=\mu -\mu S-\eta f\left(I\left(t\right)\right)S\left(t\right)+\gamma I(t-\tau )e^{-\mu \tau },t\ne t_k$$$$\dot{I}=\eta f\left(I\left(t\right)\right)S\left(t\right)-(\mu +\gamma )I\left(t\right),t\ne t_k$$$$\Delta S=-vS\left(t^{-}\right)+vS(t-\tau )e^{-\mu \tau },t=t_k.$$

A derivation of the impulse condition $\Delta S$ using assumptions (1)–(3) is available in Appendix A.

It is known [17] that the model of Kyrychko and Blyuss can exhibit Hopf bifurcation. Numerically, it appears as though the bifurcating periodic orbit may be globally (excluding the other two equilibria) asymptotically stable. In the presence of impulse effects, Hopf points are known to generically lead to bifurcations to invariant cylinders [14]. The ramifications of this result to the present model are that, in the presence of pulse vaccination, we expect a bifurcation from an endemic periodic solution to an invariant cylinder. Verifying this hypothesis is our primary goal.

The structure of the paper is as follows. Section 2 contains some necessary background material and notation that will be used throughout, as well as a reformulation of the model that will be needed for some of the bifurcation analysis. We study the disease-free periodic solution in Section 3. A numerical analysis of the cylinder bifurcation is completed in Section 4. The biological interpretation of the results is provided in Section 5. We end with the concluding Section 6.