Population Dynamics of Patients with Bacterial Resistance in Hospital Environment.

During the past decades, the increase of antibiotic resistance has become a major concern worldwide. The researchers found that superbugs with new type of resistance genes (NDM-1) have two aspects of transmission characteristics; the first is that the antibiotic resistance genes can horizontally transfer among bacteria, and the other is that the superbugs can spread between humans through direct contact. Based on these two transmission mechanisms, we study the dynamics of population in hospital environment where superbugs exist. In this paper, we build three mathematic models to illustrate the dynamics of patients with bacterial resistance in hospital environment. The models are analyzed using stability theory of differential equations. Positive equilibrium points of the system are investigated and their stability analysis is carried out. Moreover, the numerical simulation of the proposed model is also performed which supports the theoretical findings.

1. Introduction

During the past decades, the increase of antibiotic resistance has become a major concern worldwide. Antibiotic resistance is a type of drug resistance where a microorganism is able to survive exposure to an antibiotic. Serious infections caused by bacteria that have become resistant to commonly used antibiotics have become a principle global healthcare problem in the 21st century. Superbugs with the drug-resistant genes are special kinds of bacteria which can be spread by human contact. Drug-resistant genes can be transferred between bacteria in horizontal fashion by conjugation, transduction, or transformation [1, 2]. Therefore a gene for antibiotic resistance which had evolved via natural selection may be shared. In 2010 the August issue of the journal The Lancet: Infectious Diseases, a multinational team reported the emergence and spread of 180 cases of patients infected by bacteria carrying the NDM-1, thus suggesting a widespread dissemination [3]. The horrible problem of NDM-1 is not its ability to attack a person, but the fact that the resistant gene can horizontally transfer to any other bacteria, which can change common bacteria to superbugs.

The World Health Organization announced that the infections caused by superbugs did not respond to conventional treatments, which often induced a long-term illness and greater risk of death. Patients have to spend more time and money because of the bacterial resistance. It will bring huge economic burden and higher risk of death; meanwhile infectious patients will continuously spread superbugs to other people by contact [4]. Therefore, the superbugs bring not only a medical challenge, but also a serious social problem.

The NDM-1 encoding gene is located on different large plasmids that are easily transferable to susceptible bacteria at a high frequency. These plasmids also harbor genes conferring resistance to almost all antibiotics, thus making their rapid dissemination in clinically relevant bacteria a serious threat for therapy [5]. Scientists are afraid that, once resistance genes are combined with certain dangerous bacteria in a multispecies environment, the consequences will be very serious. Due to the superbugs with new drug-resistant genes having just been found, there are mainly experimental studies and testing, and the relevant documents with mathematical models have not been given.

There are some corresponding researches for the dynamic of other resistant bacteria, such as MRSA [6-11] and VRE [10, 12–14]. The literatures about drug resistance under different circumstances, for instance, school, home [15, 16], hospital [8, 9, 11, 13, 17–21], and the community environment, have also attempted to elucidate the mechanisms. Mathematical modeling and simulations are essential approaches to understand biological phenomena. A deterministic mathematical model was developed in papers by D'Agata et al. [8, 17] to characterize the factors contributing to the replacement of hospital-acquired MRSA with CA-MRSA and to quantify the effectiveness of interventions aimed at limiting the spread of CA-MRSA in health care settings. A tremendous variety of models have been formulated recently, mathematically analyzed, and applied to bacterial resistance. Stochastic simulation [9, 13] and agent-based simulation [14, 18–20], combined with hospital environment, are also established. For social environment, researchers are using differential equation models [22-24], stochastic simulation model [25], and neural network simulation model [26] to discuss the dissemination of resistant bacteria.

Opatowski et al. gave a literature review about contribution of mathematical modeling to fight against antibiotic resistance [27]. Temime et al. illustrated the influence of mathematical models in the corresponding research work with example of resistant bacteria and indicated that the importance of mathematical modeling was gradually upgrading [28]. zur Wiesch et al. conducted a review on the generation, evolution, control of resistant bacteria, and so on [29].

The paper is organized as follows. A basic mathematical model is introduced in Section 2; the model depicts the population dynamics of patients with one disease in hospital environment; we analysis the stability properties of the model, and the stability of the disease-free equilibrium and the endemic equilibrium are also discussed. A complex model is given in Section 3, which describes the dynamics of population when two diseases exist in the system. When there are m diseases in the hospital environment, the population dynamics is studied in Section 4. Numerical result and parameter analysis are derived in Section 5. The model is analyzed using stability theory of differential equations. Positive equilibrium points of the system are investigated and their stability analysis is carried out. Moreover, the numerical simulation of the proposed model is also performed by using fourth-order Runge-Kutta method which supports the theoretical findings.

2. The Basic Model

Patients suffering from a certain disease who enter the hospital are divided into two categories: infected patients and susceptible patients. Here, I stands for infected patients, who suffer from the certain disease and are also with bacterial resistance; S stands for susceptible patients, who suffer from the certain disease but are without bacterial resistance. The rate of hospital admissions per day is α, p  (0 < p < 1) represents the transfer rate from susceptible patients to infected patients, the cure rate of susceptible patients is β, and 1/β is the average length of stay. k  (0 < k < 1) is the resistant strength coefficient; the cure rate of infected patients is reduced to kβ because of bacterial resistance. Here, the smaller the value of k is, the lower the cure rate of infected patients is, which also means longer average length of stay in hospital. The death rate of patients is μ; we assume that bacterial resistance has no effect on mortality but the cure rate and these bacteria are transmitted between patients in hospital via direct contact between patients, through contamination of the institutional environment, or with the inadvertent help of human vectors.

The model shown in Figure 1 considers the dynamics of patients with bacterial resistance in hospital environment where a single disease exists. The basic system is described by the following set of ordinary differential equations:

Figure 1

A compartment model of antibiotic resistance in a hospital setting. See text for description and equations.

2.1. Steady State Analysis

One of the most important concerns about any infectious disease is its ability to invade a population. An equilibrium point is a point at which variables of a system remain unchanged over time. System (1) possesses the following equilibria.

The model has two types of equilibrium points:

The disease-free equilibrium (DFE) is given by E 0(S 0, 0) = (α/(μ + β), 0).

The endemic equilibrium is given by E (S , I ) = ((αp − (μ + β)(μ + kβ))/(μ + kβ)p, (μ + kβ)/p).

It is obvious that the endemic equilibrium exists if and only if αp > (μ + β)(μ + kβ).

2.2. Stability of DFE

2.2.1. Local Stability of DFE

To discuss the local stability of equilibrium points we compute the variational matrix of system (1). The signs of the real parts of the eigenvalues of the variational matrix evaluated at a given equilibrium determine its stability. The entries of general variational matrix are given by differentiating the right-hand side of system (1). The matrix is given by

Define

We get Theorem 1 about the stability of disease-free equilibrium as follows.

Theorem 1 .

When R < 1, E 0 will be locally asymptotically stable and the endemic equilibrium does not exist.

Otherwise, E 0 is not stable and the endemic equilibrium exists.

Proof

It is easy to calculate the characteristic equation of E 0 which isSo λ 1 = −μ − β < 0, and λ 2 = pS 0 − μ − kβ = pα/(μ + β) − μ − kβ.

Considering S 0 = α/(μ + β) and the defining R as that in (3), then when R < 1, the eigenvalues has negative real part, which follows that the disease-free equilibrium E 0 is locally asymptotically stable under this condition. Otherwise, the DFE is unstable and the endemic equilibrium exists.

Biological Meanings. Here, R is the basic reproduction number, which is “the expected number of secondary cases produced, in a completely susceptible population, by a typical infective individual” [30]. The condition R < 1 has obviously biological meanings. The drug resistance has the possibility to be extinct in these populations only if the infection rate of drug resistance is small enough. At this moment, the resistant strength coefficient plays a minor role on the dissemination of drug resistance.

2.2.2. Global Stability of Disease-Free Steady State

Theorem 2 .

Define R as that in (3); then when R ≤ 1, the disease-free equilibrium E 0 will be globally asymptotically stable.

To prove the global stability of the disease-free equilibrium, we use the method of Castillo-Chávez et al. [31].

Set X = S, Z = I, and rewrite system (1) in the following form:in which G(X , 0) = 0.

Represent E 0 by U 0 = (X , 0). According to the theorem in Castillo-Chávez et al. [31], in order to get the global stability of E 0, system (1) should satisfy three conditions as follows:

E 0 is locally asymptotically stable.

For , X is globally asymptotically stable.

, where is an M-matrix.

For system (1), the first condition has been proven in Theorem 1. Since F(X , 0) is a limiting function of , that is, lim X = X , so the second condition is easy to get.

Now we compute A and as follows:For system (1) we get that S = α/(μ + β + pI), so S < S 0. So is always established. Then we conclude that the disease-free equilibrium E 0 is globally stable if R < 1. The proof is completed.

2.3. Global Stability of the Endemic Steady State

Now let us discuss the stability of E . About the local stability for the disease steady state, we have the following Theorem 3.

Theorem 3 .

Define R as that in (3). Then when R > 1, the endemic equilibrium E is globally asymptotically stable.

The two equations of system (1) constitute a planar system as follows:in which (S, I) ∈ D = {(S, I)∣S ≥ 0, I ≥ 0}. The stability of E is determined by the value of p 1, q in characteristic equation

When R > 1, q(E 0) = (μ + β)(μ + kβ) − αp < 0, the disease-free equilibrium E 0 is unstable. It is easy to verify q(E ) = αp − (μ + β)(μ + kβ) > 0 and p 1(E ) = αp/(μ + kβ) > 0, so the endemic equilibrium E is locally asymptotically stable.

Because D is a positive invariant set in system (1), to prove E is global steady in set D is equivalent to proving that no periodic orbit of system (1) exists in D. We choose the Dulac function B(S, I) = 1/I, to evaluate the following expression ∂(BP)/∂S + ∂(BQ)/∂I = −p − μ/I − β/I < 0. So no periodic orbit of system (1) exists and E is globally asymptotically stable in set D.

Biological Meanings. When R > 1, it means once the patients with drug resistance enter the hospital, the drug resistance will be epidemic. Ultimately the number of patients without drug resistance and patients with drug resistance will be stable in (μ + kβ)/p and

3. Bacterial Resistance Spread between Patients with Two Diseases

In this model, we assumed that the patients in the hospital are divided into four compartments:

is the population size of patients infected in disease 1 without bacterial resistance at time t.

is the population size of patients infected in disease 2 without bacterial resistance at time t.

is the population size of patients infected in disease 1 with bacterial resistance at time t.

is the population size of patients infected in disease 2 with bacterial resistance at time t.

There are two diseases that exist in the hospital environment; disease 1 and disease 2 are not infectious, so we assume that the hospital did not take the isolation precautions. Each kind of patients is divided into infected patients with bacterial resistance I 1, I 2 and susceptible patients without bacterial resistance S 1, S 2. Patients infected with disease i  (i = 1,2) enter the hospital with the rate α . The cure rate of susceptible patients is β , and 1/β is the average length of stay. p   (0 < p < 1) represents the transfer rate from susceptible patients to infected patients, k is the resistant strength coefficient, and the cure rate of patients I is reduced to k β because of bacterial resistance. The death rate of patients with disease i  (i = 1,2) is μ ; we also assume that bacterial resistance has no effect on mortality but affects the cure rate. Individuals enter the hospital in one of these states and exit via death or hospital discharge. Figure 2 shows an expanded form of system (1), in which two diseases exist, called disease 1 and disease 2.

Figure 2

The extended model, in which two diseases exist in hospital setting.

Based on the horizontal of the bacterial resistance between patients, the model is built as follows:

3.1. Steady State Analysis

System (10) has two possible steady states.

The disease-free equilibrium is given by E 0(S 1 0, S 2 0, 0,0) = (α 1/(μ 1 + β 1), α 2/(μ 2 + β 2), 0,0).

The endemic equilibrium is given by E (S 1 , S 2 , I 1 , I 2 ).

In this section, we analyze the steady states of the model. Define

3.2. Stability of DFE

3.2.1. Local Stability of DFE

Theorem 4 .

If S 1 0, S 2 0, I 1 0, I 2 0 ≥ 0, then the solutions are nonnegative and remain bounded in the positive cone of R 4. If R 0 1 + R 0 2 < 1, then the disease-free steady state E 0 is locally asymptotically stable. If R 0 1 + R 0 2 > 1, then E 0 is unstable.

It is easy to see that the solutions remain in the positive cone if the initial conditions are in the positive cone. Let T = I 1 + I 2 + S 1 + S 2. ThenThus, the solutions remain bounded in the positive cone of R 4 and the system induces a global semiflow in the positive cone of R 4.

To determine the stability of the disease-free steady state E 0, we use the results in van den Driessche and Watmough [32]. Reorder the components of E 0 as I 1 0 = 0, I 2 0 = 0, S 1 0 = α 1/(μ 1 + β 1), and S 2 0 = α 2/(μ 2 + β 2). SetThenSimilarly,Therefore,which implies that the spectral radius of the matrix FV −1 isIf R 0 1 + R 0 2 < 1, then ρ(FV −1) < 1. By Theorem 2 in van den Driessche and Watmough [32], we know that the disease-free steady state E 0 is locally asymptotically stable. E 0 is unstable if R 0 1 + R 0 2 > 1.

Remark 5 .

The case when R 0 1 + R 0 2 < 1 corresponds to the situation that there are no drug resistance strains prevailing in the hospital. Define R = R 0 1 + R 0 2, where R is the reproduction number of system (10).

3.2.2. Global Stability of Disease-Free Steady State

Let

Theorem 6 .

Define R as that in (18); then when R ≤ 1, the disease-free equilibrium E 0 will be globally asymptotically stable.

To prove the global stability of the disease-free equilibrium, we use the method of Castillo-Chávez et al.

Set X = (S 1, S 2), Z = (R 1, R 2), and rewrite system (10) in the following form:in which G(X , 0) = 0. Represent E 0 by U 0 = (X , 0). According to the theorem in Castillo-Chávez et al. [31], in order to get the global stability of E 0, system (10) should satisfy three conditions as follows:

E 0 is locally asymptotically stable.

For , X is globally asymptotically stable.

, , where is an M-matrix.

For system (10), the first condition has been proven in Theorem 1. Since F(X , 0) is a limiting function of , that is, lim X = X , so the second condition is easy to get.

Now we compute A and as follows:From system (10) we get that S 1 = α 1/(p 1 I 1 + p 2 I 2 + μ 1 + β 1), so S 1 < S 1 0. In the same way we get that S 2 < S 2 0. So is always established. Then we conclude that the disease-free equilibrium E 0 is globally stable if R ≤ 1. The proof is completed.

3.3. Stability of the Endemic Equilibrium

Proposition 7 .

Define R as that in (18); then the endemic equilibrium E (S 1 , S 2 , I 1 , I 2 ) is locally asymptotically stable.

We only numerically investigated the system's behavior around the interior feasible equilibrium point E and provide the necessary numerical proof in the next section.

Theorem 8 .

No periodic orbit of system (10) exists in

If a periodic orbit of model (10) in Ω exists, its projection onto some two-dimensional subspace of R 4 should also be periodic. Therefore, we have to investigate if any periodic solution exists or not in an all two-dimensional subspace. There are six different two-dimensional subsystems of (10).

For the subsystem,

We choose the Dulac function B 1 = 1/S 1 S 2, to evaluate the following expression:

For the other five subsystems, we choose the Dulac functionsSimilarly

Now, using the Bendixson-Dulac negative criterion, no periodic solution in these two dimensions can exist. Therefore, the solution of (10) in R 4 also cannot oscillate persistently.

4. Bacterial Resistance Spreads between Patients with m Diseases

There are m diseases that exist in the hospital environment; because these diseases are not infectious, we assume that the hospital did not take the isolation precautions. Each kind of patients is divided into infected patients with bacterial resistance I and susceptible patients without bacterial resistance S   (i = 1,…, m). Patients infected with disease i  (i = 1,…, m) enter the hospital with the rate α . The cure rate of patients without bacterial resistance is β , and 1/β is the average length of stay. k is the resistant strength coefficient; the cure rate of patients I is reduced to k β because of bacterial resistant. The death rate of patients with disease i  (i = 1,…, m) is μ ; we also assume that bacterial resistance has no effect on mortality but the cure rate. Each person may come into contact with an infected type; individuals enter the hospital in one of S   (i = 1,…, m) states and exit via death or hospital discharge.

The patients in the hospital are divided into 2m compartments, based on the horizontal of the bacterial resistance between patients; the model is built as follows:for i = 1,…, m, where x = (I 1,…, I , S 1,…, S ). The incidence, p , depends on individual behavior, which determines the amount of mixing between the different groups.

The DFE for this model is x 0 = (0,…, 0, S 1 0,…, S 0), where S 0 = α /(μ + β ).

Linearizing (26) about x = x 0 giveswhere δ is one if i = j, but zero otherwise. Thus, F has rank one, and the basic reproduction number isThat is, the basic reproduction number of the disease is the sum of the “reproduction numbers” for each group.

5. Numerical Result and Parameter Analysis

5.1. Numerical Results

The stability and instability of the equilibrium points of the system are studied using the linear stability approach. For further analysis of the steady state of equilibrium points and the parameter effects of p 1 and p 2, we illustrate some key numerical solutions in system (10).

When the parameter values are fixed at α 1 = 4, α 2 = 6, μ 1 = 0.06, μ 2 = 0.06, p 1 = 0.004, p 2 = 0.002, β 1 = 0.2, β 2 = 0.15, k 1 = 0.7, and k 2 = 0.5, we can calculate the equilibrium points as The basic reproduction number is R = 1.1436, and as seen in Proposition 7 the equilibrium point E = (16.9873,32.3503,2.6222,7.4905) is locally asymptotically stable.

The numerical proof of Proposition 7 is as follows.

For the interior equilibrium point E = (16.9873,32.3503,2.6222,7.4905), the Jacobi matrix is The eigenvalues of J are λ 1 = −0.0175, λ 2 = −0.2261, λ 3 = −0.1397, and λ 4 = −0.1801, and all the real parts of the eigenvalues are negative. Hence, the Routh-Hurwitz criteria are satisfied.

Therefore, E = (16.9873,32.3503,2.6222,7.4905) is locally asymptotically stable.

5.2. Role of Parameter p 1 in Model (10)

p 1  (0 < p 1 < 1) represents the conversion rate from susceptible patients to infected patients. Therefore, as p 1 plays a major role in the outcome of the model, we will discuss its effect on the system. The following initial values are used: S 1(0) = 60, S 2(0) = 80, I 1(0) = 4, I 2(0) = 0, and the parameters except p 1 are α 1 = 4, α 2 = 6, μ 1 = 0.06, μ 2 = 0.06, p 2 = 0.001, β 1 = 0.2, β 2 = 0.15, k 1 = 0.7, and k 2 = 0.5.

In order to clearly show population dynamics for each case, two numerical results are given in the following for different values of p 1.

(a) p 1 = 0.001. In this case, the basic reproduction number is R = 0.7973 < 1 and so the equilibrium E 0 = (19.0476,37.5000,0, 0) is locally asymptotically stable (Figure 3). Both of the infected patients are extinct, and there are only susceptible patients. In the wide range 0 < p 1 < 0.0028, the behavior of the system is qualitatively the same.

Figure 3

Evolution of the system populations for p 1 = 0.001.

(b) p 1 = 0.005. In this case, the basic reproduction number is R = 1.2590 > 1 and so the interior equilibrium E = (15.5828,29.0286,4.4097,12.3220) is locally asymptotically stable (the proof is provided above). The system tends toward the coexistence equilibrium (Figure 4). When increasing p 1 is from 0.0028 to 1, there is only quantitative change, and the system's qualitative behavior remains the same. Four populations can coexist.

Figure 4

Evolution of the system populations for p 1 = 0.005.

5.3. Role of Parameter p 2 in Model (10)

Besides the parameter p 1, we are also concerned about the role of parameter p 2. The initial values are S 1(0) = 60, S 2(0) = 80, I 1(0) = 4, and I 2(0) = 0, and the parameters except p 2 are α 1 = 4, α 2 = 6, μ 1 = 0.06, μ 2 = 0.06, p 1 = 0.004, β 1 = 0.2, β 2 = 0.15, k 1 = 0.7, and k 2 = 0.5.

We describe the influence of p 2 on the system under two different situations.

(a) p 2 = 0.001. In this case, the basic reproduction number is R = 0.8027 < 1 and so the equilibrium E 0 = (19.0476,37.5000,0, 0) is locally asymptotically stable (Figure 5). Both of the infected patients are extinct, and there are only susceptible patients. In the wide range 0 < p 2 < 0.0016, the behavior of the system is qualitatively the same.

Figure 5

Evolution of the system populations for p 2 = 0.001.

(b) p 2 = 0.002. In this case, the basic reproduction number is R = 1.1436 > 1 and so the interior equilibrium E = (16.9873,32.3503,2.6222,7.4905) is locally asymptotically stable (the proof is provide above). The system tends toward the coexistence equilibrium (Figure 6). When increasing p 2 from 0.0016 to 1, there is only quantitative change, and the system's qualitative behavior remains the same. Four populations can coexist.

Figure 6

Evolution of the system populations for p 2 = 0.002.

5.4. Role of Parameter p 1 and p 2 in Model (10)

In the following section, we discuss the effect of p 1 and p 2 together, while maintaining the other parameters fixed at α 1 = 6, α 2 = 8, μ 1 = 0.06, μ 2 = 0.06, β 1 = 0.2, β 2 = 0.15, k 1 = 0.7, and k 2 = 0.5. The resulting image is shown in Figure 7, in which two different regions were obtained after 500 time steps, which are indicated by different colors.

Figure 7

Changing p 1 and p 2 while maintaining the other parameters and the initial values fixed.

In region A, four populations coexist. In region B, both of the infected patients are extinct, and there are only susceptible patients. As shown in Figure 7, for a fixed p 1, as the value of p 2 gets smaller, the possibility of coexistence will also be smaller, even becoming impossible.

6. Discussion

Two steady states were obtained from the model described in system (10): in the first state, there are only susceptible patients and both of the infected patients are extinct; in the second state, four populations coexist. Among these two cases, only the first is a favorable outcome for humans. As the result shown above, fix one transfer rate, as the value of another transfer rate gets smaller; the possibility of coexistence will also be smaller, even becoming impossible. Meanwhile, the basic production number is determined by the value of transfer coefficients and resistant strength coefficients. Because the resistant strength coefficients are difficult to control, horizontal transfer coefficients are the key parameters that not only can affect the basic reproduction number but also can be controlled by human. There are some papers about how to reduce acquisition of antimicrobial-resistant bacteria [10, 33, 34], such as improving hand hygiene, designing of an efficient sentinel hospital surveillance system, and controlling the connections and number of connections that a given hospital has with other hospitals. In the following work, we will try to add some control measures in the model.