Nonstandard numerical methods for a mathematical model for influenza disease

https://doi.org/10.1016/j.matcom.2008.04.008Get rights and content

Abstract

In this paper we construct and develop two competitive implicit finite difference scheme for a deterministic mathematical model associated with the evolution of influenza A disease in human population. Qualitative dynamics of the model is determined by the basic reproduction number, R0. Numerical schemes developed here are based on nonstandard finite difference methods. Our aim is to transfer essential properties of the continuous model to the discrete schemes and to obtain unconditional stable schemes. The proposed numerical schemes have two fixed points which are identical to the critical points of the continuous model and it is shown that they have the same stability properties. Numerical simulations with different initial conditions, parameters values and time step sizes are developed for different values of parameter R0, convergence to the disease free equilibrium point when R0<1 and to the endemic equilibrium point when R0>1 are obtained independent of the time step size. These numerical integration schemes are useful since can reproduce the dynamics of original differential equations.

Introduction

Influenza is caused by a virus that attacks mainly the upper respiratory tract, nose, throat and bronchi and rarely also the lungs. Most people recover within 1–2 weeks without requiring any medical treatment. In the very young, the elderly and people suffering from medical conditions such as lung diseases, diabetes, cancer, kidney or heart problems, influenza poses a serious risk. In these people, the infection may lead to severe complications of underlying diseases, pneumonia and death. In annual influenza epidemics 5–15% of the population are affected with upper respiratory tract infections. Hospitalization and deaths mainly occur in high risk groups (elderly, chronically ill). Annual epidemics are possibly between three and five million cases of severe illness and between 250,000 and 500,000 deaths every year around the world [26].

Many problems in mathematical epidemiology are modeled by autonomous systems of nonlinear ordinary differential equations, which implies the assumption that the parameters of the model are independent of time. In these models the variables represent subpopulations of susceptibles, infected, recovered, transmitted diseases vectors, and so forth. Thus, the solutions of ODE system describe the evolution of the different classes of subpopulations in the model for different times, therefore one important task of the mathematical modeling is to obtain accurate numerical solutions. The ideal scenario occurs when it is possible to find analytical solutions, but in most of the cases this is not possible to achieve, so it is necessary to turn to numerical methods to obtain accurate solutions.

It is well known that traditional schemes like forward Euler, Runge–Kutta and others, sometimes fail generating oscillations, bifurcations, chaos and false steady states (see [16]), moreover some methods, despite using adaptative step size, still fail (see [22]). One alternative to prevent all this class of numerical instabilities is the construction of schemes using nonstandard finite-difference method. This technique, developed by Mickens [19], [18] have brought a creation of new numerical schemes preserving the physical properties, especially the stability properties of equilibria, of the approximated system [21], [2], [1], [24], [12], [5], [6], [7], [8], [14]. Piyanwong et al. [24] and Jansen and Twizell [14] have designed positive and unconditionally stable schemes for the SIR and SEIR models, respectively.

In this paper, we develop two nonstandard finite difference schemes to obtain numerical solutions of a influenza A disease model presented by Casagrandi et al. [4]. To our knowledge these are the first nonstandard finite difference schemes applied to this model. The construction of the nonstandard difference schemes are developed using nonlocal approximation and replacing the traditional denominator of derivatives by a nonnegative function ψ(h) giving rise to fixed point schemes, which preserves the local stability properties of disease free equilibrium point and allow us to compute numerically the right endemic equilibrium point. Furthermore, the proposed numerical schemes are used with arbitrarily large time step sizes, saving computational cost when integrating over long time periods, in fact Euler method and other well-known methods produce bad approximations simulating the model for large time step sizes.

The continuous model is based on the partition of the host population into four subpopulations, susceptibles S(t), infectious I(t), recovered R(t) and cross-immune C(t). This model is represented by an ODE system with four equations, which allows to discuss how the different epidemiological parameters influence the global behavior in the evolution of the influenza disease in the human population. The subpopulation C(t) included in the model is the cross-immune individuals, i.e., those recovered after being infected by different strains of the same viral subtype in the past years. This last subpopulation is not considered in the classical SIR models for influenza disease [4]. The extension of the SIR model with the new class C(t) which allows to model the fact that, when the cross-immune hosts are exposed to antigenically similar (but different) strains, their pre-existing immune responses are boosted.

The organization of this paper is as follows. In Section 2, some well-known mathematical preliminaries on difference equations are given, which will be used in the study of systems of difference equations associated to the SIRC model for influenza disease. In Section 3 we introduce the mathematical model developed by Casagrandi et al. [4] for the evolution of influenza A disease in human population. The construction of proposed nonstandard numerical schemes is carried out in Section 4. In Section 5 analysis of the schemes is completely developed when R0<1 and, although unmanageable analytic expression of the eigenvalues of the Jacobian matrix corresponding to the linearization at the endemic point when R0>1, complicates the theoretical analysis, numerical simulations developed in Section 6 for different parameters values and step size values of both numeric schemes show the convergence properties and performance versus other well-known schemes. Discussion and conclusions are presented in Section 7.

Section snippets

Preliminaries on the difference equations

Steady state equilibria provide essential reference points for the characterization of nonlinear dynamical systems. Generally, a nonlinear system may be characterized by the existence of a unique equilibrium, the existence of multiplicity of (distinct) steady state equilibria, the existence of chaotic behavior, or the non-existence of a steady state equilibria. Stability analysis of the nonlinear system in general requires its linear approximation in the vicinity of its steady state

Mathematical model

In this section stability properties of the mathematical model of the influenza disease proposed in [4] are studied in order to know the correct behavior of the numerical solutions coming out from the discretization of the model using an appropriate nonstandard scheme.

The model presented in [4] for the spread influenza disease in the human population classifies the population in four groups or classes: S(t) is the proportion of susceptibles at time t (individuals that do not have the virus), I(t

Discretizations

In this section we construct the proposed schemes for the systems (1), (3). The main idea of these schemes choice is to obtain unconditionally stability and positivity in the variables representing the subpopulations S(t),I(t),R(t) and C(t). The first motivation is important since large time step sizes can be used, saving computational cost when integrating over long time periods. Second motivation is important due to the fact that variables representing subpopulations must never take negative

Analysis of the schemes

In this section we show that the proposed schemes (6), (8) have the properties of asymptotic stability using the Linearized Stability Theorem in the DFE. First at all, we present a Lemma which gives some information about roots of polynomials of degree two.

Lemma 5.1

[3], p.82

For the quadratic equation λ2aλ+b=0 both roots satisfy |λi|<1,i=1,2, if and only if the following conditions are satisfied:

  • (1)

    1a+b>0,

  • (2)

    1+a+b>0,

  • (3)

    b<1.

Let us consider the linearization of system (6) at the DFE point F=(1,0,0). The Jacobian

Numerical results

In this section numerical results illustrate the advantages of the proposed nonstandard finite difference schemes. In order to test convergence and stability properties of the schemes, we perform several numerical simulations varying time step size and parameters of the influenza model. As contact rate parameter β is very important from the dynamics point of view, different values were taken in the numerical simulations. Numerical results are presented in different subsections, one for the DFE

Discussion and conclusions

In this paper we use nonstandard finite difference numerical method to develop two nonstandard schemes for the solution of a SIRC model for influenza disease which is represented by a nonlinear system of ordinary differential equations. This epidemic model with temporary partial immunity or variable susceptibility is interesting to understand the evolution of this disease as well as other diseases of similar characteristics. The aforementioned model have two biological equilibrium points, one

Acknowledgements

The authors are grateful to the anonymous reviewers for their valuable comments and suggestions which improved the quality and the clarity of the paper.

References (26)

  • W. Piyawong et al.

    An unconditionally convergent finite-difference scheme for the SIR model

    Appl. Math. Comput.

    (2003)
  • R. Anguelov et al.

    Contributions to the mathematics of the nonstandard finite difference method and applications

    Numer. Methods Partial Differ. Eq.

    (2001)
  • F. Brauer et al.

    Mathematical Models in Population Biology and Epidemiology

    (2001)
  • Cited by (0)

    View full text