Next Article in Journal
Digital PCR: What Relevance to Plant Studies?
Next Article in Special Issue
Multipath: An R Package to Generate Integrated Reproducible Pathway Models
Previous Article in Journal
Biochemical, Physiological, and Productive Response of Greenhouse Vegetables to Suboptimal Growth Environment Induced by Insect Nets
Previous Article in Special Issue
Fluid-Structure Interaction Simulation of an Intra-Atrial Fontan Connection
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Computational Model Informs Effective Control Interventions against Y. enterocolitica Co-Infection

by
Reihaneh Mostolizadeh
1,2,3,4,* and
Andreas Dräger
1,2,3,4
1
Computational Systems Biology of Infections and Antimicrobial-Resistant Pathogens, Institute for Bioinformatics and Medical Informatics (IBMI), University of Tübingen, 72076 Tübingen, Germany
2
Department of Computer Science, University of Tübingen, 72076 Tübingen, Germany
3
German Center for Infection Research (DZIF), Partner Site Tübingen, 72076 Tübingen, Germany
4
Cluster of Excellence ‘Controlling Microbes to Fight Infections’, University of Tübingen, 72076 Tübingen, Germany
*
Author to whom correspondence should be addressed.
Biology 2020, 9(12), 431; https://doi.org/10.3390/biology9120431
Submission received: 16 October 2020 / Revised: 19 November 2020 / Accepted: 23 November 2020 / Published: 30 November 2020
(This article belongs to the Special Issue Computational Biology)

Abstract

:

Simple Summary

Medical control strategies for infectious diseases remain enormously important. One germ that can cause gastrointestinal infections is Yersinia enterocolitica. This study investigates and analyzes a computational model to identify the occurrence of disease-free and co-infection states. Thereby, the reproduction number R 0 informs us about the germ’s ability to spread disease. Suppose this fundamental quantity takes a value between zero and one. In that case, every infectious strain will cause less than one secondary infection, so the strain will disappear. In contrast, if R 0 exceeds one, every infectious strain causes more than one secondary infection, and Yersinia infection strains will persist. A disease-free state occurs when the commensal bacteria’s growth rate exceeds the maximum immune action and the rate at which the intestines release the bacteria. With a large enough commensal bacteria growth rate, this state can be stable. Co-infection occurs when the maximum growth rates of the wild-type and mutant strains become unequal. Studying the immune system’s behavior can result in an infection’s disappearance from hosts with a healthy microbiota immune system. In this case, Yersinia strains do not spread in the lumen when the commensal bacteria’s growth rate exceeds the growth rate of wild-type and mutant Yersinia.

Abstract

The complex interplay between pathogens, host factors, and the integrity and composition of the endogenous microbiome determine the course and outcome of gastrointestinal infections. The model organism Yersinia entercolitica (Ye) is one of the five top frequent causes of bacterial gastroenteritis based on the Epidemiological Bulletin of the Robert Koch Institute (RKI), 10 September 2020. A fundamental challenge in predicting the course of an infection is to understand whether co-infection with two Yersinia strains, differing only in their capacity to resist killing by the host immune system, may decrease the overall virulence by competitive exclusion or increase it by acting cooperatively. Herein, we study the primary interactions among Ye, the host immune system and the microbiota, and their influence on Yersinia population dynamics. The employed model considers commensal bacterial in two host compartments (the intestinal mucosa the and lumen), the co-existence of wt and mut Yersinia strains, and the host immune responses. We determine four possible equilibria: disease-free, wt-free, mut-free, and co-existence of wt and mut in equilibrium. We also calculate the reproduction number for each strain as a threshold parameter to determine if the population may be eradicated or persist within the host. We conclude that the infection should disappear if the reproduction numbers for each strain fall below one, and the commensal bacteria growth rate exceeds the pathogen’s growth rate. These findings will help inform medical control strategies. The supplement includes the MATLAB source script, Maple workbook, and figures.

1. Introduction

Forecasting the evolution of infectious diseases is the primary motivation behind the use of mathematical models in biology. Determining those threshold values that predict whether the disease will spread within the host or can be contained is crucial to inform medical control strategies.
Yersinia entercolitica (Ye) is a Gram-negative enteropathogen causing foodborne gastrointestinal infections. Within the small intestine (SI), Ye can adhere to and invade the intestinal epithelial lining mainly via the adhesins Yersinia adhesin A (YadA) [1] and Invasin [2,3,4]. YadA is among the essential virulence factors as a YadA-deficient strain is impaired 97 colonizations, and systemic spread in a mouse model of infection has been shown Upon attachment, Ye can engage its Type III secretion system (T3SS) to inject effector proteins, the Yersinia outer proteins (Yops), into host cells, thus evading the host immune response and establishing a productive infection [5,6,7]. Indeed, the first line of host defense against invading Ye is a massive infiltration of phagocytotic cells. However, Ye can counteract phagocytic killing via its T3SS [6,7]. Together, both the T3SS and YadA contribute to the efficient colonization of the intestinal tract, where Ye induces an inflammatory response that likely accounts for a reduction in density and complexity of the commensal microbiome [8,9]. Several Ye serotypes have been isolated from animal reservoirs and the human gastrointestinal tract, but only a few of them cause disease in humans. Although causing severe distress, gastrointestinal infections are generally self-contained. Typically, healthy persons will only receive symptomatic treatment aimed at avoiding dehydration. However, in individuals with underlying disease, elderly persons, or newborn children, gastrointestinal infections can cause high morbidity and even fatal outcomes. Especially in such patients, identifying those at high risk of developing fatal diseases could help tailor personalized therapeutic interventions to improve the outcome of infection. To this aim, we recently developed a computational model that can calculate Ye population dynamics during gastrointestinal infection and predict pathogen expansion, gut colonization, and infection course [10].
Previously published models have focused on the virulence evolution during co-infection and superinfection by more than one pathogen [11,12,13,14,15,16,17]. Furthermore, Nurtay et al. [18] investigated theoretical conditions in co-existence strains in their research focused on viral populations of wt and mut [18]. Herein, we apply our computational model to predict the behavior of Ye population dynamics during the co-infection of mice by different Ye mutant (mut) strains, lacking distinct virulence factors, and a Ye wt strain in bacterial gastroenteritis. We perform a bifurcation analysis to show how the dynamics of the system change as multiple parameters are varied [19,20].
The existence of a backward bifurcation, i.e., the co-existence of stable disease-free equilibrium with one or more stable endemic equilibria, has significant consequences on the process of producing medical policies aimed at eradicating or controlling an infectious disease within the host [21]. A fundamental parameter in models displaying a backward bifurcation is the primary reproduction number R 0 . R 0 represents the expected number of new infections caused within the host or between hosts when a pathogen enters into the host. If R 0 > 1 , after an initial introduction, the infection spreads within the host/or between hosts, thus creating an infectious disease. If R 0 < 1 , small initial introductions are not sufficiently transmissible to spread the infection within the host (or between hosts), and the cells get infected within a host. This is called an endemic disease that will fade out. Thus, many control policies like medication/vaccination have focused on reaching coverage levels sufficient to reduce R 0 below 1 [15,16,22,23]. Therefore, in this analysis, we compute the system’s equilibrium points, study their stability and the center manifold, and investigate how the different points would affect the infection rate, virulence, and mutation rate. Next, we calculate the primary reproductive number for scenarios in which multiple strains are introduced, and an R 0 number is calculated for each strain. Finally, we discuss the biological implications of our results.

2. Materials and Methods

2.1. Model Description

Our model for Ye, as described in Geißert et al. [10], considers three different sites with their individual population dynamics. Figure 1 depicts an overview of the model in the form of a Systems Biology Graphical Notation (SBGN) [24] PD diagram [25] based on an illustration by Geißert et al. [10]. The lumen and the mucosal site of the small intestine are illustrated as two separate compartments by abbreviated notations L and M in Table 1. The mucosa and lumen sites include three different population dynamics: commensal bacteria, wt Yersinia, and mut Yersinia.
Besides the bacterial populations, we have the strength of the host immune response. For the sake of clarity, the complex immune cell population dynamics, which is made up of innate and adaptive immunity, are implemented in our model as a single abstract immune response in Table 1. An advantage of this rather abstract immune response is that it is easily adjustable. Only the Ye populations in the mucosa activate the immune response, whereas tolerance to commensal bacteria exists. However, this immune response affects all bacterial populations at the mucosal site. Upon oral co-infection, a Ye population enters the small intestinal lumen. Due to their particular virulence traits, some Ye cells are also able to enter the mucosal site.
As shown in Models (1)–(7), the population dynamics of bacterial species consist of both growth, i.e., an increase in populations due to a distinct growth rate α , and a reduction through peristalsis, which moves the bacterial populations towards the colon where they finally end up in feces. Peristalsis will only influence the populations in the lumen. Furthermore, bacterial populations in the mucosa can additionally be reduced through killing by the host immune attack.
Both compartments (the lumen and the mucosal) are colonized by commensal bacterial species but to different extents. The colonization capacity of the lumen considerably exceeds that of the mucosal site. This is due to the host’s natural mechanisms to limit bacterial contact to the epithelium through physical barriers, such as the adherent mucus layer, and a high concentration of anti-microbial peptides (AMPs). Hence, we assume a much lower carrying capacity of the mucosal site compared to the luminal compartment. Another assumption is that as soon as bacterial numbers exceed the mucosa’s carrying capacity, they spill-over to the luminal site, feeding the luminal populations [10].
The model’s immune system represented by I (see Table 1) is activated as soon as at least one Ye cell enters the small intestinal mucosal compartment. This immune activation increases with the number of Ye cells in the mucosa. By triggering the killing of all bacterial cells at the mucosal site, this leads to a reduction in bacterial populations spilling over to the lumen [10].
In contrast to commensal bacteria, Ye exhibits a number of virulence traits to evade the host immune response. This is implemented in our model by different immunity adjustment factors for either the wt strain or the different mut strains. Consequently, the number of commensal species is more affected by the host immune attack than the number of Ye mut strains, which are, of course, more affected than the wt strain. Our model’s final output is the number of bacteria, or colony-forming units (CFUs), finally ending up in feces. These population dynamics are represented as the following system of differential equations [10]:
d d t B M = α ( B ) σ M L ( B ) γ · I · B M
d d t Y M ( w t ) = α ( w t ) σ M L ( w t ) γ · f I ( w t ) · I · Y M ( w t )
d d t Y M ( m u t ) = α ( m u t ) σ M L ( m u t ) γ · f I ( m u t ) · I · Y M ( m u t )
d d t B L = α L ( B ) β · B L + σ M L ( B ) · B M
d d t Y L ( w t ) = α L ( w t ) β · Y L ( w t ) + σ M L ( w t ) · Y M ( w t )
d d t Y L ( m u t ) = α L ( m u t ) β · Y L ( m u t ) + σ M L ( m u t ) · Y M ( m u t )
d d t I = Y M ( w t ) + Y M ( m u t ) · κ · C I I C I .
For a detailed description of the parameters, see Table 2. Besides the seven differential equations describing the population dynamics of the Ye strains, the commensal bacteria and the immune system at the mucosal site, and all three populations (commensal bacteria, wt Yersinia, and mut Yersinia) at the luminal site, the model consists of different mass-action terms describing the spill-over rates from the mucosa into the lumen and describes the growth rates of the bacterial populations, which are limited through a given capacity “C” and are defined as follows:
σ M L ( B ) = α ( B ) B M + Y M ( w t ) + Y M ( m u t ) C M Spill-over rate commensal
σ M L ( w t ) = α ( w t ) B M + Y M ( w t ) + Y M ( m u t ) C M Spill-over rate wt
σ M L ( m u t ) = α ( m u t ) B M + Y M ( w t ) + Y M ( m u t ) C M Spill-over rate mut
α L ( B ) = α ( B ) C L B L + Y L ( w t ) + Y L ( m u t ) C L Bacteria growth rate in the lumen
α L ( w t ) = α ( w t ) C L B L + Y L ( w t ) + Y L ( m u t ) C L wt growth rate in the lumen
α L ( m u t ) = α ( m u t ) C L B L + Y L ( w t ) + Y L ( m u t ) C L mut growth rate in the lumen
Model parameters are shown in Table 2. Since the system models the Yersinia population in the mucosa and lumen, we assume that all state variables and parameters of the model are non-negative t 0 .

2.2. The Basic Reproduction Number R 0

The basic reproduction ratio R 0 is calculated by the fraction of the transmission rate for each strain (spill-over) and the average infectious period for each strain in compartments, as defined by Diekmann et al. [28].
In addition, van den Driessche and Watmough [29,30] defined R 0 as a general compartmental disease transmission model suited to heterogeneous populations that can be modeled by a system of ordinary differential equations. Suppose there are n disease compartments and m non-disease compartments, and let x R n and y R m be the subpopulations in each of these compartments. In this splitting, the authors thus represented the diseased compartment with F i ( x ) as the rate of the appearance of new infections in compartment i, and V i as the rate disease progression, death, and recovery reduced in compartment i, which can be written as x ˙ i = f i ( x ) = F i ( x ) V i , i = 1 , , n .
To understand compartment V i better, this can be written as V i = V i + V i , where V i + is defined as the rate of transfer of individuals into compartment i by all other means, and V i is the rate of transfer of individuals out of the compartment i. Therefore, the disease transmission model consisting of non-negative initial conditions will be as follows:
d d t x i = F i ( x , y ) V i ( x , y ) i { 1 , , n }
d d t y j = G j ( x , y ) j { 1 , , m } .
It is assumed that functions F i ( x ) and V i = V i + V i are continuously differentiable at least twice in each variable, and they satisfy the assumptions A(1)–A(5) described below [29,31]:
A(1)
F i ( 0 , y ) = 0 , V i ( 0 , y ) = 0 : y > 0 and i = 1 , , n (no immigration of individuals into the disease compartments)
A(2)
F i ( x , y ) 0 : x i 0 y i 0 and i = 1 , , n (the new infections will be represented by F , so it cannot be negative)
A(3)
V i ( x , y ) 0 : whenever x i = 0 , i = 1 , , n (if the compartment is empty, it can only have inflow, and the net outflow from the compartment must be negative)
A(4)
i V i ( x , y ) 0 : x i 0 y i 0 (sum is net outflow)
A(5)
The system y ˙ = G ( 0 , y ) has a unique asymptotically stable equilibrium, y * (all solutions with initial conditions of the form ( 0 , y ) approach a point ( 0 , y * ) as t )
Theorem 1.
Consider the disease transmission model given by X ˙ = f ( X ) with f ( X ) satisfying conditions A(1)–A(5). If X 0 is a DFE of the model, then X 0 is locally asymptotically stable if R 0 < 1 , but unstable if R 0 > 1 , where R 0 acts as a threshold parameter.
Proof. 
See [29]. □
By defining F = F i x j ( 0 , y * ) and V = V i x j ( 0 , y * ) as n × n matrices, the basic reproduction number R 0 is computed by R 0 = ρ F V 1 , where ρ ( A ) denotes the spectral radius of a matrix A.
Several models found in the literature have been used to show that when R 0 crosses the threshold, R 0 = 1 , it can act as a bifurcation parameter and transcritical bifurcation takes place. That is, local asymptotical stability is transferred from the disease-free state to the new (positive) equilibria. In ordinary differential equations, we will encounter bifurcations of equilibrium and periodic orbits, which are typical in the sense that they occur when a small smooth change is made to the threshold parameter values (the bifurcation parameters μ ) of a system. It causes a sudden qualitative or topological change in the behavior of the system. Consider the continuous dynamical system described by the Ordinary Differential Equation (ODE)
X ˙ = f ( X , μ ) with f : R n × R R n .
A local bifurcation occurs at ( X 0 , μ 0 ) if an eigenvalue with zero real part is included in the Jacobian matrix of the system. If the eigenvalue is equal to zero, the bifurcation is a steady-state bifurcation. Therefore, we now recall the analysis of the center manifold near the critical ( X = X 0 , R 0 = 1 ) , which allows clarifying the direction of the bifurcation near the bifurcation point using the Center Manifold Theorem 2. This theory describes the local stability at the non-hyperbolic equilibrium (linearization matrix has at least one eigenvalue with zero real parts) and the existence of another equilibrium (bifurcated from the non-hyperbolic equilibrium).

2.3. Center Manifold

The center manifold theorem provides a means for systematically reducing the dimensions of the state spaces, which need to be considered when analyzing bifurcations of a given type.
Theorem 2 (Center Manifold Theorem for Flows).
Let f be a C r vector field on R n vanishing at the origin f ( 0 ) = 0 , and let A = D f ( 0 ) . Divide the spectrum of A into three parts, σ s , σ c , a n d σ u , with
R e λ < 0 if λ σ s = 0 if λ σ c > 0 if λ σ u .
Let the (generalized) eigenspaces of σ s , σ c , and σ u be E s , E c , and E u , respectively. Then, there exists C r stable and unstable invariant manifolds W u and W s tangent to E u and E s at 0 and a C r 1 center manifold W c tangent to E c at 0. The manifolds W u , W s , and W c are all invariant for the flow of f. Stable and unstable manifolds are unique, but W c need not be.
Proof. 
See [32,33]. □
To achieve this, consider the general system
x ˙ = J x + F ( x ) , x R n
where J x is the linear part of the system. We must first find the differential equations on its center manifold and then reduce the system to its normal form. Without loss of generality, we assume that x = 0 is the fixed point of interest for the system.
Suppose J has n c eigenvalues with zero real-parts and n s eigenvalues with negative and positive real parts, and n = n c + n s . Using the eigenvectors of J to form a transformation matrix, the system can be rewritten in block matrix form as
x ˙ c = A x c + f ( x c , x s ) x ˙ s = B x s + g ( x c , x s )
where ( x c , x s ) R n c × R n s , A R n c × R n c , and B R n s × R n s . With the eigenvalues of zero real parts, the Center Manifold Theorem 2 guarantees that there exists a smooth manifold W c = ( x c , x s ) x s = q ( x c ) , the equilibrium point such that the local behavior in the center direction of the system is qualitatively the same as that on the manifold. By differentiating x s = q ( x c ) , we get x ˙ s = D q ( x c ) x ˙ c . Substituting (19) into the previous identity and rearranging the equation, we get
D q ( x c ) A x c + f x c , q ( x c ) B q ( x c ) g x c , q ( x c ) = 0 .
By solving for q ( x c ) , we get a function describing the center manifold. In general, q ( x c ) cannot be solved explicitly. Instead, substituting a Taylor expansion q ( x c ) = a x c 2 + O x c 3 into (20), we can find the coefficients for the expansion by balancing the lower-order terms. Based on q ( x c ) , we now have a system in the reduced form:
x ˙ c = A x c + f x c , q ( x c ) .
The following proposition helps us understand if the type of transcritical bifurcation is forward or backward.
Proposition 1.
Assume that conditions A(1)–A(5) are satisfied. Furthermore, assume that the following hypotheses are satisfied by the system in (14) and (15):
H(1) 
In the balance equations for the infected compartments, nonlinear terms are present only in the rate of the appearance of new infections;
H(2) 
Nonlinear terms are bilinear;
H(3) 
There is no linear transfer from infected to uninfected compartments.
Then, the transcritical bifurcation of the system in (14) and (15) at R 0 = 1 is forward, and if at least one of these features is not present in the model structure, then a backward bifurcation may occur [31].

3. Results

We used the model in (1)–(7) to understand how the dynamics change following variations of different parameters. We calculated the equilibria of (1)–(7), conducted a linear stability analysis, and identified the analytical conditions that lead to a transcritical bifurcation.

3.1. Existence of Equilibria

For mathematical convenience, we divide the model in (1)–(7) such that the first four equations correspond to infected individuals. The distinction between infected and uninfected populations must be determined from the model’s epidemiological interpretation and cannot be deduced from the structure of the equations alone. Thus, we have
d d t Y M ( w t ) = α ( w t ) σ M L ( w t ) γ · f I ( w t ) · I · Y M ( w t )
d d t Y M ( m u t ) = α ( m u t ) σ M L ( m u t ) γ · f I ( m u t ) · I · Y M ( m u t )
d d t Y L ( w t ) = α L ( w t ) β · Y L ( w t ) + σ M L ( w t ) · Y M ( w t )
d d t Y L ( m u t ) = α L ( m u t ) β · Y L ( m u t ) + σ M L ( m u t ) · Y M ( m u t )
d d t B M = α ( B ) σ M L ( B ) γ · I · B M
d d t B L = α L ( B ) β · B L + σ M L ( B ) · B M
d d t I = Y M ( w t ) + Y M ( m u t ) · κ · C I I C I ,
while we assume that the growth rate α ( B ) of the endogenous commensal bacteria is higher than the Ye growth rates α ( w t ) and α ( m u t ) , respectively.
On the other side, the Yersinia model has three compartments (mucosa, lumen, immune system), which are analyzed separately. The model system is analyzed in a suitable, feasible region. All forward solutions of the system are feasible t 0 if they enter the invariant region Ω = Ω L × Ω M × Ω I where
Ω L = B L , Y L ( w t ) , Y L ( m u t ) R + 3 : B L + Y L ( w t ) + Y L ( m u t ) C L
Ω M = B M , Y M ( w t ) , Y M ( m u t ) R + 3 : B M + Y M ( w t ) + Y M ( m u t ) C M
Ω I = ( I ) R + 1 : I C I .
The existence of equilibrium points for system (1)–(7) would be as follows:
  • The trivial equilibrium point is as an origin equilibrium 0 , 0 , 0 , 0 , 0 , 0 , 0 . This solution appears when all populations are extinct. For all parameters, this point never becomes stable due to the positivity of eigenvalues in (A2).
  • The first equilibrium point appears in the absence of Yersinia Y M ( w t ) = Y M ( m u t ) = Y L ( w t ) = Y L ( m u t ) = 0 . System (1)–(7) has a disease-free equilibrium, which is given by
    P 0 = B M 0 , Y M 0 ( w t ) , Y M 0 ( m u t ) , B L 0 , Y L 0 ( w t ) , Y L 0 ( m u t ) , I 0 = B M 0 , 0 , 0 , B L 0 , 0 , 0 , I 0
    and
    B M 0 = C M 1 γ α ( B )
    B L 0 = C L 1 β α ( B )
    I 0 = 1 .
    It describes a disease-free state whereby only the commensal bacteria persist. In order for the disease-free state P 0 to be biologically meaningful, the conditions γ < α ( B ) and β < α ( B ) must hold. These conditions correspond to the maximal growth rate of intestinal bacteria exceeding the rate at which intestines are charged and the maximal immunity action, which is not that strong in the absence of Yersinia strains. However, the population of the immune system is at its maximum carrying capacity (in health, not in fighting with any infection).
  • A second equilibrium corresponds to the commensal bacteria’s persistence and the Yersinia mut strain in the absence of the wt strain. Without loss of generality, the commensal bacteria are supposed to be zero because they are not infective. This point is obtained by setting Y M ( w t ) = Y L ( w t ) = 0 :
    P ( m u t ) = P 1 = B M 1 , Y M 1 ( w t ) , Y M 1 ( m u t ) , B L 1 , Y L 1 ( w t ) , Y L 1 ( m u t ) , I 1 = 0 , 0 , Y M 1 ( m u t ) , 0 , 0 , Y L 1 ( m u t ) , C I
    with
    Y M 1 ( m u t ) = C M 1 C I γ f I ( m u t ) α ( m u t )
    Y L 1 ( m u t ) = 1 2 C L 1 β α ( m u t ) + 4 C M C L C I γ f I ( m u t ) α ( m u t ) 1 2 + C L 2 1 β α ( m u t ) 2 .
  • The other equilibrium corresponds to the persistence of commensal bacteria and the Yersinia wt strain in the absence of the mut strain. Without loss of generality, the commensal bacteria are supposed to be zero because they are not infective. This point is obtained by setting Y M ( m u t ) = Y L ( m u t ) = 0 :
    P ( w t ) = P 2 = B M 2 , Y M 2 ( w t ) , Y M 2 ( m u t ) , B L 2 , Y L 2 ( w t ) , Y L 2 ( m u t ) , I 2 = 0 , Y M 2 ( w t ) , 0 , 0 , Y L 2 ( w t ) , 0 , C I
    with
    Y M 2 ( w t ) = C M 1 C I γ f I ( w t ) α ( w t )
    Y L 2 ( w t ) = 1 2 C L 1 β α ( w t ) + 4 C M C L C I γ f I ( w t ) α ( w t ) 1 2 + C L 2 1 β α ( w t ) 2 .
  • Finally, the last equilibrium point corresponds to a state of the co-existence of wt and mut Yersinia strains. This point is achieved by supposing B M = B L = 0 :
    P ( w t ) ( m u t ) = P 3 = B M 3 , Y M 3 ( w t ) , Y M 3 ( m u t ) , B L 3 , Y L 3 ( w t ) , Y L 3 ( m u t ) , I 3 = 0 , Y M 3 ( w t ) , Z , 0 , Y L 3 ( w t ) , Y L 3 ( m u t ) , C I
    with
    Y M 3 ( w t ) = C M 1 C I γ f I ( w t ) α ( w t ) Z
    Y L 3 ( w t ) = α ( m u t ) Z 1 C I γ f I ( w t ) α ( w t ) γ f I ( w t ) 1 C I α ( w t ) 2 β 1 α ( m u t ) α ( w t )
    Y L 3 ( m u t ) = α ( m u t ) 1 C I γ f I ( w t ) α ( w t ) Z 1 C I γ f I ( w t ) α ( w t ) β 1 α ( m u t ) α ( w t )
    where Z is defined as 1 α ( B ) f I ( w t ) α ( w t ) α ( B ) f I ( m u t ) α ( m u t ) to make the equilibrium point biologically meaningful.

3.2. Analysis of the Disease-Free Equilibrium Point

The system’s behavior close to the equilibrium points was determined through linear stability analysis and bifurcations. It was carried out by calculating the Jacobian matrix of the model in equilibrium points. The Jacobian matrix for (1)–(7) is given by Equation (A1) in the Appendix A. The disease-free equilibrium of the model is the steady-state solution of the model in the absence of the disease. The eigenvalues of the Jacobian matrix (A1) at this point were calculated as follows:
λ 1 = 0
λ 2 = α ( B ) γ
λ 3 = α ( B ) β
λ 4 = γ f I ( w t ) α ( w t ) α ( B )
λ 5 = γ f I ( m u t ) α ( m u t ) α ( B )
λ 6 = β 1 α ( w t ) α ( B )
λ 7 = β 1 α ( m u t ) α ( B )
in which it is assumed that α ( w t ) < α ( B ) and α ( m u t ) < α ( B ) . Based on the general compartmental model describing an infectious disease transmission within a heterogeneous population, the host population is grouped into two general classes: the infected and uninfected compartments. Therefore, the system (1)–(7) is divided into two infected and uninfected compartments.
In classical epidemic models, it is common to observe that a disease-free equilibrium loses its stability for R 0 = 1 , and a transcritical bifurcation occurs. A transcritical bifurcation is when a fixed point exists for all parameter values and is never destroyed. However, as the parameter values vary, this fixed point transitions from a stability region to an instability region. Biologically speaking, as R 0 < 1 , the disease-free equilibrium would stay stable. Thus, only the commensal bacteria persist, and wt and mut Yersinia strains cannot pass the invasion boundary. Therefore, as R 0 increases, wt and mut Yersinia strains, fed by commensal bacteria, can appear. We mathematically analyze this aspect within the structure of the model to assess which parts of the structure might be responsible for the occurrence of the transcritical bifurcation.
Let us consider system (1)–(7) with the above-calculated eigenvalues. Using the eigenvectors (A9) shown in the Appendix A as a basis for a new coordinate system, we set the transformation matrix whose columns are the eigenvectors; T = [ v 1 , v 2 , v 3 , v 4 , v 5 , v 6 , v 7 ] t . We thus have
B M Y M ( w t ) Y M ( m u t ) B L Y L ( w t ) Y L ( m u t ) = T u 1 u 2 u 3 u 4 u 5 u 6 u 7 and u 1 u 2 u 3 u 4 u 5 u 6 u 7 = T 1 B M Y M ( w t ) Y M ( m u t ) B L Y L ( w t ) Y L ( m u t ) .
Under this transformation, the system (1)–(7) changes to
u ˙ 1 u ˙ 2 u ˙ 3 u ˙ 4 u ˙ 5 u ˙ 6 u ˙ 7 = A u 1 u 2 u 3 u 4 u 5 u 6 u 7 + f 1 ( u 1 , , u 7 ) f 2 ( u 1 , , u 7 ) f 3 ( u 1 , , u 7 ) f 4 ( u 1 , , u 7 ) f 5 ( u 1 , , u 7 ) f 6 ( u 1 , , u 7 ) f 7 ( u 1 , , u 7 )
where
A = 0 0 0 0 0 0 0 0 γ α ( B ) f I ( w t ) α ( w t ) α ( B ) 0 0 0 0 0 0 0 γ α ( B ) f I ( m u t ) α ( m u t ) α ( B ) 0 0 0 0 0 0 0 β α ( B ) α ( w t ) α ( B ) 0 0 0 0 0 0 0 β α ( B ) α ( m u t ) α ( B ) 0 0 0 0 0 0 0 α ( B ) + β 0 0 0 0 0 0 0 α ( B ) + γ
so that the linear part is now in a standard (diagonal) form. In the ( u 1 , u 2 , u 3 , u 4 , u 5 , u 6 , u 7 ) coordinates, the center manifold is a curve tangent to the u 1 axis . The projection of the system onto the u 1 axis is obtained by setting u 2 = u 3 = u 4 = u 5 = u 6 = u 7 = 0 in the equation for u ˙ 1 . It yields u ˙ 1 = 0 . Since E c is one dimension, we can approximate h i ( u 1 ) by Taylor expansion such that u i = h i ( u 1 ) , i { 2 , , 7 } , satisfying the following equations:
u ˙ 2 = D h 2 ( u 1 ) u ˙ 1 , u ˙ 3 = D h 3 ( u 1 ) u ˙ 1 , u ˙ 4 = D h 4 ( u 1 ) u ˙ 1 , u ˙ 5 = D h 5 ( u 1 ) u ˙ 1 , u ˙ 6 = D h 6 ( u 1 ) u ˙ 1 , u ˙ 7 = D h 7 ( u 1 ) u ˙ 1 , h i ( 0 ) = D h i ( 0 ) = 0
Thus, the reduced system is
u ˙ 1 = α ( B ) C M u 1 2 + O u 1 3 .
The advantages of a center manifold are clear from this calculation. We may study a one-dimensional system instead of a seven-dimensional system. That is, the method of center manifolds enables one to reduce the dimensions of the system by studying the flow restricted to the center manifold, in which the transients associated with the nonzero eigenvalues have decayed. As long as α ( w t ) < α ( B ) and α ( m u t ) < α ( B ) , the disease-free equilibrium point will persist. When the infection rate of wt and mut strains increases, the system becomes unstable at the disease-free equilibrium point, and the trajectory of the system approaches asymptotically to the bacterial population of wt or mut and, finally, to the co-existence of the strains.

3.3. Computing and Analysis of the System through the Basic Reproduction Number

For scenarios where multiple strains (subtypes) of an infectious disease exist, an R 0 number is calculated for each strain. Hence, in our model, we must define two different R 0 , one for the wt and the other for the mut strain, rather than an R 0 for the whole model. This is a significantly more difficult task. Herein, we introduce the following variations to the basic reproduction number, which are calculated by the average number of secondary infections:
  • R 0 w t : the basic reproductive numbers for wt strain = α ( w t ) α ( B ) f I ( w t )
  • R 0 m u t : the basic reproductive numbers for mut strain = α ( m u t ) α ( B ) f I ( m u t ) .
To understand the role of the basic reproduction number, we also define a single reproduction number for commensal bacteria R 0 c o m , computed as an expected number of secondary cases of commensal bacteria produced by a single commensal bacteria:
R 0 c o m = max α ( B ) β , α ( B ) γ .
Due to biologically meaningful disease-free state P 0 and holding the conditions β < α ( B ) and γ < α ( B ) , the disease-free state would be meaningful if R 0 c o m > 1 . Therefore, the commensal bacteria appear, R 0 c o m slightly reaches above one, and a stable disease-free state happens. Due to the presence of commensal bacteria and instability of the trivial solution, the solution corresponds to the extinction of all populations that never appear because if R 0 c o m < 1 , all components of the disease-free equilibrium will lie out of the biologically meaningful region.
For R 0 w t < 1 and R 0 m u t < 1 , thus, λ 1 = 0 is a simple zero eigenvalue, and the other eigenvalues are real and negative. This implies that transcritical bifurcation occurs in the disease-free equilibrium for R 0 w t < 1 and R 0 m u t < 1 , and the uninfected state is stable (not asymptotically stable) as none of the strains persist. As soon as the two equilibria collide non-destructively, exchanging their stability and resulting in the Jacobian matrix having a single eigenvalue that is equal to zero, then either the mut strain or the wt appears and persists for R 0 m u t or R 0 w t , respectively.
As shown in Figure 2, when the disease-free equilibrium loses its stability, different scenarios can occur. Herein, we discuss the different cases.
(I)
If R 0 w t = 1 and R 0 m u t = 1 , then α ( w t ) = α ( B ) f I ( w t ) and α ( m u t ) = α ( B ) f I ( m u t ) . Thus, the intersection of the transcritical curves R 0 w t and R 0 m u t results in a triple transcritical bifurcation. As shown in (A9), the Jacobian has a triple zero eigenvalue at this point ( λ 2 = 0 , λ 3 = 0 ). Kuznetsov [34] has proved that such a point would be an indicator of the onset of a non-degenerate or degenerate Bogdanov–Takens bifurcation [34,35]. The disease-free equilibrium P 0 loses its stability, and the wt-free and mut-free include one simple zero eigenvalue ( λ 3 = 0 ), meaning that the dynamics of the model change as the target parameter is within the threshold value.
(II)
If R 0 w t > 1 , the wt strain equilibrium in region II will persist when R 0 m u t < 1 . The wt strain will spread and possibly persist within the host population. In general, for a strain to persist, its basic reproduction number has to be strictly greater than one. Therefore, in this region, the disease-free, mut strain, and co-existence state exchange stability: P 0 becomes unstable, P 2 becomes locally asymptotically stable, and P 1 and P 3 remain unstable. This means that the immune system could kill one of the strains more efficiently.
(III)
If R 0 m u t > 1 , the mut strain equilibrium in region III will persist when R 0 w t < 1 . The mut strain will spread and possibly persist within the host population since its basic reproduction number is greater than one. Therefore, in this region, the disease-free, wt strain, and co-existence state exchange stability: P 0 becomes unstable, P 1 becomes locally asymptotically stable, and P 2 and P 3 remain unstable. This means that the immune system could defeat the wt strains. However, the risk of this situation to happen is low because the mut strains are influenced more efficiently than wt strains by immune action.
(IV)
If R 0 w t > 1 and R 0 m u t > 1 , the co-existence population spreads, and both strains persist. The overall R 0 can be defined as R 0 c = d R 0 m u t + ( 1 d ) R 0 w t . A mut with d = 1 is thoroughly dominant, while one with d = 0 is completely recessive; scenarios of incomplete dominance ( 0 < d < 1 ) , under-dominance ( d < 0 ) , and over-dominance ( d > 1 ) are possible as well. For instance, a mut could achieve a higher R 0 than the wt via a higher growth rate that increases transmission. In a co-infection, the faster-growing mut strain would outcompete the wt and reach its maximum capacity. This situation would change the co-infection to the conditions where a single infection happens. Thus, the overall R 0 of the co-infection would be similar to that of the mut by itself, making the mut a dominant one. Furthermore, the effort for having a co-existence equilibrium and analysis of the co-infection model will fail. By contrast, let us assume a mut strain achieves a higher R 0 . Nevertheless, the virulence of the wt strain neutralizes the higher R 0 value of the mut. This would make the mut a recessive one. In summary, virtually any two-strain co-infection model can be mapped to a set of values for d, allowing scenarios of particular interest to be explored in a context broader than the one possible with typical models.
The mathematical model (1)–(7) calculates the number of pathogen-specific characteristics in different layouts, e.g., when colonization resistance mediated by the endogenous microbiome is lacking or when the immune response is partially impaired. In this paper, we use the experimental data obtained in the lab upon infection of a host either lacking a microbiome (mimicking antibiotic treatment of patients) or a fully functional immune system [10].
To challenge the numerical simulation with experimental data, we tested whether the numerical simulation from the mathematical analysis could fit the experimental data.
Correspondingly, we simulated the process of Yersinia enterocolitica co-infection in specific-pathogen-free (SPF) (i.e., wt), germ-free (GF), and MyD88-deficient (MyD88−/−) mut mice. Using the model values in Table 3, which was estimated through experimental data [10], we first examined the influence of the infection rates α ( B ) , α ( w t ) , and α ( m u t ) on the dynamics beginning by constructing one-dimensional bifurcation diagrams using α ( B ) as the bifurcation parameter and fixing α ( w t ) and α ( m u t ) . The model values in Table 3 are the output of estimation parameters in Table 2 of the (1)–(7) by experimental data [10]. In theory, some parameter values were defined based on the biological aspect, and the rest were estimated by the optimization problem with the maximum likelihood estimation [10].
As shown in Figure 3, R 0 c o m should be larger than one to achieve the disease-free equilibrium. As long as R 0 w t and R 0 m u t are smaller than one, the disease-free equilibrium is stable.
Following the analysis of the model with respect to α ( B ) , two cases are considered. The first case is when α ( w t ) and α ( m u t ) are equal. In this case, we do not have co-existence of wt and mut (the denominators of Y L 3 ( w t ) and Y L 3 ( m u t ) will be equal to zero). The second one is when α ( w t ) and α ( m u t ) have different values as the co-existence equilibrium is achieved and is biologically meaningful.
Thus, we start with values of the wt and mut-strain infections rate given by α ( w t ) = 4.44 E 01 and α ( w t ) = 4.44 E 01 in Table 3. Corresponding to those infections rates, the basic reproduction numbers R 0 ( w t ) and R 0 ( m u t ) in terms of α ( B ) are computed.
Through the definitions of R 0 c o m , R 0 w t , and R 0 m u t , we display an x-scale as α ( B ) in Figure 4. This shows the appearance of a different equilibrium of Model (14) in terms of the basic reproduction numbers.
As shown in Figure 4, we expect no co-existence equilibrium when α ( w t ) = α ( m u t ) , and only wt equilibrium exists (mut-free equilibrium) when 1.76 < α ( B ) < 11.20 . However, the disease-free equilibrium can happen when α ( B ) is large enough to surpass α ( w t ) and α ( m u t ) . These findings are displayed in Figure 5.
Secondly, to show the appearance of the co-existence, we consider the parameter values of in Ye wt/T3S0 from Table 3 with the assumption of 1.86 = α w t α m u t = 1.83 . Thus, we compute R 0 w t and R 0 m u t with respect to parameter α ( B ) . This results in the appearance co-existence when 1 < α ( B ) < 4.90 as Figure 6 shows.
Furthermore, we analyze the immune system’s influence. In our model, at least one Ye within the mucosal compartment activates the immune system. This activation increases proportionally to the number of Ye. The immune system is assumed to influence the bacterial populations primarily at the mucosal site compared to bacteria within the lumen. Herein, we simulate the process of Yersinia enterocolitica with an immune response that was experienced in SPF, GF, and MyD88−/− mice. As the responses of the immune system in these mice are different, we conclude different behaviors. Let us denote the numerical simulation:
  • The effect of the maximum rate of immune growth κ on wt Yersinia strain in the mucosa, Figure 7;
  • The effect of the maximum rate of immune growth κ on mut Yersinia strain in the mucosa, Figure 8;
  • The effect of the maximum rate of immune growth κ on wt Yersinia strain in the lumen, Figure 9;
  • The effect of the maximum rate of immune growth κ on mut Yersinia strain in the lumen, Figure 10,
The experimental data obtained in the lab showed that the Ye mut strain was eliminated more efficiently than Ye wt in the mucosal compartment. However, since GF mice lack a microbiome, both the Ye wt and the mut strains can colonize the GIT at high numbers. In the SPF-colonized MyD88 mice, the immune response is weaker [10]. Therefore, as shown in Figure 7, Figure 8, Figure 9 and Figure 10, the following results conclude:
  • Figure 7d and Figure 8d show a slow reduction in the wt strain in comparison to mut strain as κ changes.
  • Figure 7b, Figure 8b, Figure 9b and Figure 10b in comparison with similar situations in SPF and MyD88−/− show elimination of wt and mut strains is less efficient.
  • Figure 7c, Figure 8c, Figure 9c and Figure 10c show that wt and the mut strains increased faster. Besides, the influence of κ on wt and the mut strains cannot project properly (Figure 9f and Figure 10f) at the same speed of producing strains.
Our results allow us to state that if mice possess a healthy microbiota and immune system, as long as the growth rate of commensal bacteria is larger than the growth rate of wt and mut Yersinia, then the infection will not spread, as Ye strains cannot enter the lumen compartment. Therefore, the disease-free equilibrium, P 0 exists and is stable when R 0 w t and R 0 m u t are less than one while α ( B ) is large enough. To challenge this situation, we tested the data values Ye SPF wt/A0 from Table 3 with the assumed boundary α ( B ) > 11.20 to face a disease-free state, Figure 4. Our computer simulations ran for 366 h and are shown in Figure 11. As predicted for the disease-free state, the commensal fractions B M and B L approached their carrying capacities C M and C L , respectively.
Therefore, we note that the disease-free equilibrium always exists and is (asymptotically) stable whenever the reproductive numbers for both strains of the disease are below unity. Simulations support these findings. The two single-strain equilibria, where one strain persists while the other strain dies out, are symmetrical and exist when at least one of the strains has a reproductive number above unity. However, in our model, (14) with parameter values from Table 3, we could not observe a state that corresponds to Ye mut strain equilibrium. As shown in Figure 4, there is no region where R 0 m u t > 1 and R 0 w t < 1 . Whenever R 0 m u t > 1 , R 0 w t is already above one; this facilitates the situation for a co-existence equilibrium. This is where both strains persist, existing when both reproductive numbers are above a certain threshold. However, we could not analytically solve the stability criteria for the equilibrium due to complexity, but simulations show that these stability criteria exist.

4. Discussion and Conclusions

In this paper, we proved that our computational model of the Yersinia enterocolitica [10] population in co-infections of mice with Ye wt and mut strains always showed stable results in terms of R 0 , such that if we can keep or reduce R 0 s to less than one, then wt and mut strains cannot spread, and the infection resolves. The medical control strategies for this infectious disease, like other infectious diseases, are frequently predicted through the basic reproduction number. To know what the difference is, we should be aware that the decomposition in dynamics and the designation of compartments as infected or uninfected are unique, and the basic reproduction number is achieved in the same way for all different compartmental models. This also raises the question of how the infection becomes persistent or is resolved based on the calculation of R 0 as a significant measure in control policy.
The answer to this question depends on the magnitude of the basic reproductive number, R 0 , since calculating the stability of an equilibrium is very complicated. Therefore, we computed R 0 for our computational model based on different strains. This resulted in three R 0 s : one R 0 c o m , R 0 w t , and R 0 m u t .
As it is obvious in any compartmental model, the disease-free equilibrium is locally (asymptotically) stable when 0 < R 0 < 1 and unstable if R 0 > 1 . In other words, when 0 < R 0 < 1 , every infectious strain will cause less than one secondary infection; hence, the strain will eliminate. When R 0 > 1 , every infectious strain will cause more than one secondary infection; hence, Yersinia infection strains will persist. However, all public health control measures can usually be based on methods that, if practical, would lower R 0 to below unity. On the other hand, the co-infection equilibrium is locally stable when R 0 > 1 and unstable when 0 < R 0 < 1 .
This trivial result is essential, but we had to adjust our multi-strain model with different R 0 s . To more smoothly analyze R 0 s , we looked for a common parameter in all R 0 s . This resulted in the parameter α ( B ) . We thus analyzed at which level of sensibility the system’s parameter α ( B ) changes the model’s invasion behavior. This came to an end when a disease-free state occurs, when α ( B ) satisfies the conditions β < α ( B ) and γ < α ( B ) . In addition, if α ( B ) is large enough, then the disease-free state can be stable, and the infection with other strains does not appear or spread. A mut strain equilibrium does not occur for parameter values estimated through experimental data [10] since there is no region with only R 0 m u t > 1 . Although there is a region with two threshold parameters in the model, R 0 w t and R 0 m u t stay above one for a co-infection setup. This region includes a small range for the parameter α ( B ) . In addition, it is important that the maximum growth rates associated with the wt and mut strains should be unequal. Otherwise, the co-existence scenario cannot take place.
Furthermore, we looked for conditions in the immune system to see different scenarios of a weaker and more robust immune system in fighting against infection. This resulted in the analysis of systems with respect to parameter κ . We found that if mice possess a healthy microbiota and immune system, as long as the growth rate of commensal bacteria is larger than the growth rate of wt and mut Yersinia, then the infection will not spread, as Ye strains cannot enter the lumen compartment.
Our results indicate that the within-host dynamics of immunity can, in principle, have significant consequences on population-level dynamics. However, immunity alone never creates a backward bifurcation of the disease-free steady state under biologically realistic hypotheses, and this would require some other complementary conditions.
Altogether, Yersinia is a great model system and can be used to predict the spread of infections for more clinically relevant bugs, i.e., enteropathogenic or enterohemorrhagic Escherichia coli. To control the spread of infection, efforts must be made to ensure that the co-infection equilibrium is unstable, i.e., 0 < R 0 < 1 , and in our model R 0 w t < 1 and R 0 m u t < 1 .
Furthermore, an additional factor to take into account is whether antibiotics are used from the beginning of the infection, in which case the spread of Yersinia may or may not be extinguished. We also believe that further research into the elaborate makeup of the strains and how the immune system differentiates between strains could be incorporated to more accurately represent control policy.

Author Contributions

Conceptualization, R.M. and A.D.; methodology, R.M.; software, R.M.; validation, R.M. and A.D.; formal analysis, R.M.; investigation, A.D.; resources, A.D.; data curation, R.M. and A.D.; writing—original draft preparation, R.M.; writing—review and editing, A.D.; visualization, R.M.; supervision, A.D.; project administration, A.D.; funding acquisition, R.M. and A.D. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the German Center for Infection Research (DZIF, doi:10.13039/100009139) within the Deutsche Zentren der Gesundheitsforschung (BMBF-DZG, German Centers for Health Research of the Federal Ministery of Education and Research), grant № 8020708703 and supported by infrastructural funding from the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation), Cluster of Excellence EXC 2124 Controlling Microbes to Fight Infections. Parts of this work were funded by the Federal Ministry of Education and Research (BMBF, Germany) and the Baden-Württemberg Ministry of Science as part of the Excellence Strategy of the German Federal and State Governments. The authors acknowledge support by the Open Access Publishing Fund of the University of Tübingen (https://uni-tuebingen.de/de/58988).

Acknowledgments

The authors thank Libera Lo Presti for feedback on this manuscript.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

AMPantimicrobial peptide
BMBFFederal Ministry of Education and Research
BMBF-DZGDeutsche Zentren der Gesundheitsforschung
CFUcolony-forming unit
DFEdisease-free equilibrium
DFGDeutsche Forschungsgemeinschaft
DZIFGerman Center for Infection Research
GFgerm-free
GITgastrointestinal tract
mutmutant
MyD88−/−MyD88-deficient mice
ODEOrdinary Differential Equation
PDProcess Description
RKIRobert Koch Institute
SBGNSystems Biology Graphical Notation
SIsmall intestine
SPFspecific-pathogen-free
T3SSType iii@ secretion system
wtwild-type
YadAYersinia adhesin A
YeYersinia entercolitica
YopYersinia outer protein

Appendix A. Mathematical Calculation

Appendix A.1. Jacobian

The Jacobian of the system (1)–(7) is calculated as follows:
J = α ( B ) 1 2 B M + Y M ( w t ) + Y M ( m u t ) C M γ I α ( B ) B M C M α ( w t ) Y M ( w t ) C M α ( w t ) 1 B M + 2 Y M ( w t ) + Y M ( m u t ) C M γ I f I ( w t ) α ( m u t ) Y M ( m u t ) C M α ( m u t ) Y M ( m u t ) C M α ( B ) 2 B M + Y M ( w t ) + Y M ( m u t ) C M α ( B ) B M C M α ( w t ) Y M ( w t ) C M α ( w t ) B M + 2 Y M ( w t ) + Y M ( m u t ) C M α ( m u t ) Y M ( m u t ) C M α ( m u t ) Y M ( m u t ) C M 0 κ ( C I I ) C I α ( B ) B M C M 0 α ( w t ) Y M ( w t ) C M 0 α ( m u t ) 1 B M + 2 Y M ( w t ) + Y M ( m u t ) C M γ I f I ( m u t ) 0 α ( B ) B M C M α ( B ) 1 2 B L + Y L ( w t ) + Y L ( m u t ) C L β α ( w t ) Y M w t C M α ( w t ) Y L ( w t ) C L α ( m u t ) B M + Y M ( w t ) + 2 Y M ( m u t ) C M α ( m u t ) Y L ( m u t ) C L κ ( C I I ) C I 0 0 0 γ B M 0 0 γ f I ( w t ) Y M ( w t ) 0 0 γ f I ( m u t ) Y M ( m u t ) α ( B ) B L C L α ( B ) B L C L 0 α ( w t ) 1 B L + 2 Y L ( w t ) + Y L ( m u t ) C L β α ( w t ) Y L w t C L 0 α ( m u t ) Y L ( m u t ) C L α ( m u t ) 1 B L + Y L ( w t ) + 2 Y L ( m u t ) C L β 0 0 0 κ ( Y M ( w t ) + Y M ( m u t ) ) C I

Appendix A.2. The Eigenvalues of Trivial Equilibrium Point

The eigenvalues for trivial equilibrium points are as follows:
λ 1 = 0
λ 2 = α ( B )
λ 3 = α ( w t )
λ 4 = α ( m u t )
λ 5 = α ( B ) β
λ 6 = α ( w t ) β
λ 7 = α ( m u t ) β

Appendix A.3. The Eigenvalues and Eigenvectors of the Disease-Free Equilibrium point

The eigenvectors corresponding to the eigenvalues for disease-free equilibrium points are as follows:
λ 1 = 0 , v 1 = γ C M α ( B ) , 0 , 0 , 0 , 0 , 0 , 1
λ 2 = γ f I ( w t ) α ( B ) α ( w t ) α ( B ) , v 2 = v 21 , v 22 , 0 , 1 , v 25 , 0 , v 27
λ 3 = γ f I ( m u t ) α ( B ) α ( m u t ) α ( B ) , v 3 = v 31 , 0 , v 33 , 1 , 0 , v 36 , v 37
λ 4 = β α ( B ) α ( w t ) α ( B ) , v 4 = 0 , 0 , 0 , 1 , α ( B ) 2 2 α ( B ) β + β α ( w t ) α ( B ) α ( B ) β , 0 , 0
λ 5 = β α ( B ) α ( mut ) α ( B ) , v 5 = 0 , 0 , 0 , 1 , 0 , α ( B ) 2 2 α ( B ) β + β α ( m u t ) α ( B ) α ( B ) β , 0
λ 6 = α ( B ) + β , v 6 = 0 , 0 , 0 , 1 , 0 , 0 , 0
λ 7 = α ( B ) + γ , v 7 = 1 , 0 , 0 , 0 , 0 , 0 , 0

Appendix A.4. The Eigenvalues wt Strain Equilibrium

The eigenvalues of the equilibrium of wt strain are calculated as
λ 1 = C I γ f I ( m u t ) α ( m u t ) f I ( w t ) α ( w t )
λ 2 = α ( w t ) 1 C I γ f I ( w t ) α ( w t )
λ 3 = C I γ 1 α ( B ) f I ( w t ) α ( w t )
λ 4 = κ C M C I 1 C I γ f I ( w t ) α ( w t )
λ 5 = α ( w t ) 4 C M C L 1 C I γ f I ( w t ) α ( w t ) 2 + 1 β α ( w t ) 2
λ 6 = 1 2 2 β α ( B ) 1 + β α ( w t ) + α ( B ) 4 C M C L 1 C I γ f I ( w t ) α ( w t ) 2 + 1 β α ( w t ) 2
λ 7 = 1 2 2 β α ( m u t ) 1 + β α ( w t ) + α ( m u t ) 4 C M C L 1 C I γ f I ( w t ) α ( w t ) 2 + 1 β α ( w t ) 2

Appendix A.5. The Eigenvalues mut Strain Equilibrium

The eigenvalues of equilibrium of mut strain are calculated as
λ 1 = C I γ f I ( w t ) α ( w t ) f I ( m u t ) α ( m u t )
λ 2 = α ( m u t ) 1 C I γ f I ( m u t ) α ( m u t )
λ 3 = C I γ 1 α ( B ) f I ( m u t ) α ( m u t )
λ 4 = κ C M C I 1 C I γ f I ( m u t ) α ( m u t )
λ 5 = α ( m u t ) 4 C M C L 1 C I γ f I ( m u t ) α ( m u t ) 2 + 1 β α ( m u t ) 2
λ 6 = 1 2 2 β α ( B ) 1 + β α ( m u t ) + α ( B ) 4 C M C L 1 C I γ f I ( m u t ) α ( m u t ) 2 + 1 β α ( m u t ) 2
λ 7 = 1 2 2 β α ( w t ) 1 + β α ( m u t ) + α ( w t ) 4 C M C L 1 C I γ f I ( m u t ) α ( m u t ) 2 + 1 β α ( m u t ) 2

Appendix A.6. The Basic Reproduction Numbers

The basic reproduction number is calculated as follows:
  • Based on Equations (14) and (15), we split (22)–(28) into two compartments: one disease compartment with populations Y M ( w t ) , Y M ( m u t ) , Y L ( w t ) , Y L ( m u t ) and one non-disease compartment with populations ( B M , B L , I ) .
  • We consider the first four Equations (22)–(28) correspond to the appearance of a new infection as well as the rate of other transitions between other infected compartments.
  • F i ( x ) corresponds to the rate of the appearance of new infections, and V i corresponds to the rate of transfer assumed as
    F 1 = α ( w t ) σ M L ( w t ) · Y M ( w t )
    F 2 = α ( m u t ) σ M L ( m u t ) · Y M ( m u t )
    F 3 = α L ( w t ) · Y L ( w t )
    F 4 = α L ( m u t ) · Y L ( m u t )
    and
    V 1 = γ · f I ( w t ) · I · Y M ( w t )
    V 2 = γ · f I ( m u t ) · I · Y M ( m u t )
    V 3 = β · Y L ( w t ) σ M L ( w t ) · Y M ( w t )
    V 4 = β · Y L ( m u t ) σ M L ( m u t ) · Y M ( m u t )
  • The basic reproduction number is achieved as the spectral radius of the matrix F V 1 ;
    as R = ρ F V 1 where F = F i ( P 0 ) x j , V = V i ( P 0 ) x j , and x j are as follows:
    x j = Y M ( w t ) , Y M ( m u t ) , Y L ( w t ) , Y L ( m u t )
    F = α ( w t ) ( γ α ( B ) ) 0 0 0 0 α ( m u t ) ( γ α ( B ) ) 0 0 0 0 α ( w t ) ( β α ( B ) ) 0 0 0 0 α ( m u t ) ( β α ( B ) )
    and
    V = γ f I ( w t ) 0 0 0 0 γ f I ( m u t ) 0 0 α ( w t ) γ + α ( B ) α ( B ) 0 β 0 0 α ( m u t ) γ + α ( B ) α ( B ) 0 β .
    Then
    R = ρ F V 1
    = max α ( w t ) α ( B ) , α ( m u t ) α ( B ) , α ( w t ) α ( B ) f I ( w t ) , α ( m u t ) α ( B ) f I ( m u t )

Appendix B. Data Availability

  • The Matlab script
  • The Maple workbook
  • The described model is available in SBML format [26] (Level 3 Version 2 [27]) from BioModels database [36] under model identifier MODEL2002070001.

References

  1. El Tahir, Y.; Skurnik, M. YadA, the multifaceted Yersinia Adhes. Int. J. Med. Microbiol. 2001, 291, 209–218. [Google Scholar] [CrossRef] [PubMed]
  2. Handley, A.J.; Koster, R.; Monsieurs, K.; Perkins, G.D.; Davies, S.; Bossaert, L. European Resuscitation Council Guidelines for Resuscitation 2005: Section 2. Adult basic life support and use of automated external defibrillators. Resuscitation 2005, 67, S7–S23. [Google Scholar] [CrossRef] [PubMed]
  3. Mühlenkamp, M.; Oberhettinger, P.; Leo, J.C.; Linke, D.; Schütz, M.S. Yersinia Adhes. A (YadA)- Beast. Int. J. Med. Microbiol. 2015, 305, 252–258. [Google Scholar] [CrossRef] [PubMed]
  4. Young, G.A.; LaVon, G.D.; Taylor, G.W. High Efficiency Absorbent Articles for Incontinence Management. US Patent 5,147,345, 15 September 1992. [Google Scholar]
  5. Pepe, G.J.; Albrecht, E.D. Actions of placental and fetal adrenal steroid hormones in primate pregnancy. Endocr. Rev. 1995, 16, 608–648. [Google Scholar] [PubMed]
  6. Cornelis, G.R. Yersinia Type III Secretion: Send Eff. J. Cell Biol. 2002, 158, 401. [Google Scholar] [CrossRef] [Green Version]
  7. Ruckdeschel, K.; Roggenkamp, A.; Schubert, S.; Heesemann, J. Differential contribution of Yersinia Enterocolitica Virulence Factors Evas. Microbicidal Action Neutrophils. Infect. Immun. 1996, 64, 724–733. [Google Scholar] [CrossRef] [Green Version]
  8. Lupp, C.; Robertson, M.L.; Wickham, M.E.; Sekirov, I.; Champion, O.L.; Gaynor, E.C.; Finlay, B.B. Host-mediated inflammation disrupts the intestinal microbiota and promotes the overgrowth of Enterobacteriaceae. Cell Host Microbe 2007, 2, 119–129. [Google Scholar] [CrossRef] [Green Version]
  9. Stecher, B.; Robbiani, R.; Walker, A.W.; Westendorf, A.M.; Barthel, M.; Kremer, M.; Chaffron, S.; Macpherson, A.J.; Buer, J.; Parkhill, J.; et al. Salmonella enterica serovar typhimurium exploits inflammation to compete with the intestinal microbiota. PLoS Biol. 2007, 5, e244. [Google Scholar] [CrossRef] [Green Version]
  10. Geißert, J.; Bohn, E.; Mostolizadeh, R.; Dräger, A.; Autenrieth, I.; Beier, S.; Deusch, O.; Eichner, M.; Schütz, M. Model-based prediction of bacterial population dynamics in gastrointestinal infection. bioRxiv 2020. [Google Scholar] [CrossRef]
  11. Nowak, M.; Robert, M. May. Superinfection and the evolution of parasite virulence. Proc. R. Soc. B 1994, 255, 81n89. [Google Scholar]
  12. Castillo-Chavez, C.; Velasco-Hernández, J.X. On the relationship between evolution of virulence and host demography. J. Theor. Biol. 1998, 192, 437–444. [Google Scholar] [CrossRef] [PubMed]
  13. Dye, C.; Williams, B.G.; Espinal, M.A.; Raviglione, M.C. Erasing the world’s slow stain: Strategies to beat multidrug-resistant tuberculosis. Science 2002, 295, 2042–2046. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  14. Boldin, B.; Geritz, S.A.; Kisdi, É. Superinfections and adaptive dynamics of pathogen virulence revisited: A critical function analysis. Evol. Ecol. Res. 2009, 11, 153–175. [Google Scholar]
  15. Liu, W.m.; Hethcote, H.W.; Levin, S.A. Dynamical behavior of epidemiological models with nonlinear incidence rates. J. Math. Biol. 1987, 25, 359–380. [Google Scholar] [CrossRef] [PubMed]
  16. Madden, L.V.; Hughes, G.; Van Den Bosch, F. The Study of Plant Disease Epidemics; Am Phytopath Society: St. Paul, MN, USA, 2007. [Google Scholar]
  17. Gross, K.L.; Porco, T.C.; Grant, R.M. HIV-1 superinfection and viral diversity. Aids 2004, 18, 1513–1520. [Google Scholar] [CrossRef] [PubMed]
  18. Nurtay, A.; Hennessy, M.G.; Sardanyés, J.; Alsedà, L.; Elena, S.F. Theoretical conditions for the coexistence of viral strains with differences in phenotypic traits: A bifurcation analysis. R. Soc. Open Sci. 2019, 6, 181179. [Google Scholar] [CrossRef] [Green Version]
  19. Mostolizadeh, R.; Shamsara, E.; Afsharnezhad, Z. Hopf Bifurcation and Chaos in a model for HTLV-I infection of CD4+ T- cells. Ann. Tiberiu Popoviciu Semin. Funct. Equ. Approx. Convexity Ser. B Math. Interdiscip. Res. 2017, 15, 59–79. [Google Scholar]
  20. Mostolizadeh, R.; Afsharnezhad, Z. Hopf bifurcation in a model for adult T-cell leukemia. Math. Methods Appl. Sci. 2018, 41, 6210–6225. [Google Scholar] [CrossRef]
  21. Shamsara, E.; Mostolizadeh, R.; Afsharnezhad, Z. Transcritical bifurcation of an immunosuppressive infection model. Iran. J. Numer. Anal. Optim. 2016, 6, 1–16. [Google Scholar]
  22. Hethcote, H.W. The mathematics of infectious diseases. SIAM Rev. 2000, 42, 599–653. [Google Scholar] [CrossRef] [Green Version]
  23. Anderson, R.M.; May, R. Infectious Diseases of Humans; Oxford Science Publication: New York, NY, USA, 1991. [Google Scholar]
  24. Touré, V.; Dräger, A.; Luna, A.; Dogrusoz, U.; Rougny, A. The Systems Biology Graphical Notation: Current Status and Applications in Systems Medicine. In Systems Medicine; Wolkenhauer, O., Ed.; Academic Press: Oxford, UK, 2020; Volume 3, pp. 372–381. [Google Scholar]
  25. Rougny, A.; Touré, V.; Moodie, S.; Balaur, I.; Czauderna, T.; Borlinghaus, H.; Dogrusoz, U.; Mazein, A.; Dräger, A.; Blinov, M.L.; et al. Systems Biology Graphical Notation: Process Description language Level 1 Version 2.0. J. Integr. Bioinform. 2019, 16. [Google Scholar] [CrossRef] [PubMed]
  26. Keating, S.M.; Waltemath, D.; König, M.; Zhang, F.; Dräger, A.; Chaouiya, C.; Bergmann, F.T.; Finney, A.; Gillespie, C.S.; Helikar, T.; et al. SBML Level 3: An extensible format for the exchange and reuse of biological models. Mol. Syst. Biol. 2020, 16, e9110. [Google Scholar] [CrossRef] [PubMed]
  27. Hucka, M.; Bergmann, F.T.; Chaouiya, C.; Dräger, A.; Hoops, S.; Keating, S.M.; König, M.; Le Novère, N.; Myers, C.J.; Olivier, B.G.; et al. Systems Biology Markup Language (SBML) Level 3 Version 2 Core Release 2. J. Integr. Bioinform. 2019, 16, 1. [Google Scholar] [CrossRef] [PubMed]
  28. Diekmann, O.; Heesterbeek, J.A.P.; Metz, J.A. On the definition and the computation of the basic reproduction ratio R 0 in models for infectious diseases in heterogeneous populations. J. Math. Biol. 1990, 28, 365–382. [Google Scholar] [CrossRef] [Green Version]
  29. Van den Driessche, P.; Watmough, J. Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission. Math. Biosci. 2002, 180, 29–48. [Google Scholar] [CrossRef]
  30. Van den Driessche, P.; Watmough, J. Further notes on the basic reproduction number. In Mathematical Epidemiology; Springer: Berlin/Heidelberg, Germany, 2008; pp. 159–178. [Google Scholar]
  31. Buonomo, B. A note on the direction of the transcritical bifurcation in epidemic models. Nonlinear Anal. Model Control 2015, 20, 38–55. [Google Scholar] [CrossRef] [Green Version]
  32. Guckenheimer, J.; Holmes, P. Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2013; Volume 42. [Google Scholar]
  33. Wiggins, S. Introduction to Applied Nonlinear Dynamical System and Chaos; Spring: New York, NY, USA, 1990. [Google Scholar]
  34. Kuznetsov, Y.A. Elements of Applied Bifurcation Theory; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2013; Volume 112. [Google Scholar]
  35. Kuznetsov, Y.A. Practical computation of normal forms on center manifolds at degenerate Bogdanov–Takens bifurcations. Int. J. Bifurc. Chaos 2005, 15, 3535–3546. [Google Scholar] [CrossRef]
  36. Malik-Sheriff, R.S.; Glont, M.; Nguyen, T.V.N.; Tiwari, K.; Roberts, M.G.; Xavier, A.; Vu, M.T.; Men, J.; Maire, M.; Kananathan, S.; et al. BioModels—15 years of sharing computational models in life science. Nucleic Acids Res. 2020, 48, D407–D415. [Google Scholar] [CrossRef] [Green Version]
Figure 1. An overview of the Yersinia enterocolitica population dynamics model. Filled circles represent the entity pool nodes for the populations of Ye in their respective compartments as well as the strength of the immune reaction I. Table 1 explains the notations used in the model and this figure. The black arrows represent processes with an impact on the population dynamics of the entity pools. Arrows pointing from empty set symbols to pool nodes denote an increase in the population size or a decrease if the process arcs point from entity pools to empty sets. Migration across compartments of the respective populations appears as vertical process arcs. Some of these processes receive stimulating effects from the immune reaction or from the size of the Ye populations within the mucosa, as colored arcs indicate. Reference [10] provides a more detailed description of the model’s structure.
Figure 1. An overview of the Yersinia enterocolitica population dynamics model. Filled circles represent the entity pool nodes for the populations of Ye in their respective compartments as well as the strength of the immune reaction I. Table 1 explains the notations used in the model and this figure. The black arrows represent processes with an impact on the population dynamics of the entity pools. Arrows pointing from empty set symbols to pool nodes denote an increase in the population size or a decrease if the process arcs point from entity pools to empty sets. Migration across compartments of the respective populations appears as vertical process arcs. Some of these processes receive stimulating effects from the immune reaction or from the size of the Ye populations within the mucosa, as colored arcs indicate. Reference [10] provides a more detailed description of the model’s structure.
Biology 09 00431 g001
Figure 2. The diagram displays three regions with different qualitative behaviors in terms of the basic reproduction numbers. Region I: infection-free state; Region II: mut-free state; Region III: wt-free state; Region IV: co-existence of all strains.
Figure 2. The diagram displays three regions with different qualitative behaviors in terms of the basic reproduction numbers. Region I: infection-free state; Region II: mut-free state; Region III: wt-free state; Region IV: co-existence of all strains.
Biology 09 00431 g002
Figure 3. The diagrams display the role of R 0 c o m in appearance/non-appearance of the trivial solution and disease-free equilibrium. (a) By fixing the parameter values Table 3 in model (14), the commensal bacteria in the mucosa appear when R 0 c o m > 1 . Therefore, as long as R 0 c o m < 1 , only the trivial solution for the model exists. Since the trivial solution is always unstable, the extinction of all populations is never achieved. (b) By fixing the parameter values Table 3 in model (14), the commensal bacteria in the lumen appear when R 0 c o m > 1 . Therefore, as long as R 0 c o m < 1 , only the trivial solution for the model exists. Since the trivial solution is always unstable, the extinction of all populations is never achieved.
Figure 3. The diagrams display the role of R 0 c o m in appearance/non-appearance of the trivial solution and disease-free equilibrium. (a) By fixing the parameter values Table 3 in model (14), the commensal bacteria in the mucosa appear when R 0 c o m > 1 . Therefore, as long as R 0 c o m < 1 , only the trivial solution for the model exists. Since the trivial solution is always unstable, the extinction of all populations is never achieved. (b) By fixing the parameter values Table 3 in model (14), the commensal bacteria in the lumen appear when R 0 c o m > 1 . Therefore, as long as R 0 c o m < 1 , only the trivial solution for the model exists. Since the trivial solution is always unstable, the extinction of all populations is never achieved.
Biology 09 00431 g003
Figure 4. The diagram displays α ( B ) regions with the appearance of a different equilibrium in terms of the basic reproduction numbers. This display is plotted by fixing the parameter values from Table 3 in Model (14).
Figure 4. The diagram displays α ( B ) regions with the appearance of a different equilibrium in terms of the basic reproduction numbers. This display is plotted by fixing the parameter values from Table 3 in Model (14).
Biology 09 00431 g004
Figure 5. The sensitivity of parameter α ( B ) with respect to commensal bacteria, wild-type, and mutant Yersinia in the mucosa and lumen for 336 h. 0 < α ( B ) < 1 , no biologically meaningful region (out of our interest). α ( B ) > 1 , all populations appear. When 1 < α ( B ) < 2.28 , the region is a region for the appearance of the co-existence equilibrium, but the hypothesis of the co-existence equilibrium is not satisfied. Therefore, wild-type strain does not grow, and the mutant strain is going down slowly. When α ( B ) > 2.28 , this is a region of the appearance of a wt equilibrium R 0 w t > 1 . Thus, (a,e) are increasing fast and stay at the maximum level as (a,d) are going back to zero. Additionally, (c,d) do not grow anymore R 0 m u t < 1 .
Figure 5. The sensitivity of parameter α ( B ) with respect to commensal bacteria, wild-type, and mutant Yersinia in the mucosa and lumen for 336 h. 0 < α ( B ) < 1 , no biologically meaningful region (out of our interest). α ( B ) > 1 , all populations appear. When 1 < α ( B ) < 2.28 , the region is a region for the appearance of the co-existence equilibrium, but the hypothesis of the co-existence equilibrium is not satisfied. Therefore, wild-type strain does not grow, and the mutant strain is going down slowly. When α ( B ) > 2.28 , this is a region of the appearance of a wt equilibrium R 0 w t > 1 . Thus, (a,e) are increasing fast and stay at the maximum level as (a,d) are going back to zero. Additionally, (c,d) do not grow anymore R 0 m u t < 1 .
Biology 09 00431 g005
Figure 6. The diagram displays the appearance of the co-existence of wt and mut Yersinia strains as α ( B ) is changing in a region where R 0 w t > 1 and R 0 m u t > 1 . (a) An immune reaction influences the wt and mut Yersinia strains in the mucosa. Additionally, another part is spilled over and moves to the lumen compartment. Therefore, they simultaneously increase or decrease to reach their maximum values. (b) The wt and mut Yersinia strains in the lumen simultaneously increase to reach the carrying capacity of the lumen compartment. In both (a,b), when α ( B ) reaches 4.90, the co-existence of wt and mut Yersinia strains disappears.
Figure 6. The diagram displays the appearance of the co-existence of wt and mut Yersinia strains as α ( B ) is changing in a region where R 0 w t > 1 and R 0 m u t > 1 . (a) An immune reaction influences the wt and mut Yersinia strains in the mucosa. Additionally, another part is spilled over and moves to the lumen compartment. Therefore, they simultaneously increase or decrease to reach their maximum values. (b) The wt and mut Yersinia strains in the lumen simultaneously increase to reach the carrying capacity of the lumen compartment. In both (a,b), when α ( B ) reaches 4.90, the co-existence of wt and mut Yersinia strains disappears.
Biology 09 00431 g006
Figure 7. The sensitivity of parameter κ in the elimination of wt Yersinia strain in the mucosa. (ac) The effect of the parameter κ on the different types of mice. (df) The projection of κ for a different types of mice.
Figure 7. The sensitivity of parameter κ in the elimination of wt Yersinia strain in the mucosa. (ac) The effect of the parameter κ on the different types of mice. (df) The projection of κ for a different types of mice.
Biology 09 00431 g007
Figure 8. The sensitivity of parameter κ in the elimination of mut Yersinia strain in the mucosa. (ac) The effect of the parameter κ on the different types of mice. (df) The projection of κ for different types of mice.
Figure 8. The sensitivity of parameter κ in the elimination of mut Yersinia strain in the mucosa. (ac) The effect of the parameter κ on the different types of mice. (df) The projection of κ for different types of mice.
Biology 09 00431 g008
Figure 9. The sensitivity of parameter κ in the elimination of wt Yersinia strain in the lumen. (ac) The effect of the parameter κ on the different types of mice. (df) The projection of κ for different types of mice.
Figure 9. The sensitivity of parameter κ in the elimination of wt Yersinia strain in the lumen. (ac) The effect of the parameter κ on the different types of mice. (df) The projection of κ for different types of mice.
Biology 09 00431 g009
Figure 10. The sensitivity of the parameter κ in the elimination of the mut Yersinia strain in the lumen. (ac) The effect of the parameter κ on the different types of mice. (df) The projection of κ for different types of mice.
Figure 10. The sensitivity of the parameter κ in the elimination of the mut Yersinia strain in the lumen. (ac) The effect of the parameter κ on the different types of mice. (df) The projection of κ for different types of mice.
Biology 09 00431 g010
Figure 11. Diagram displaying the disease-free state by fixing the parameter values Ye SPF wt/A0 from Table 3 in Model (14). As long as α ( B ) is large enough, the disease-free state persists. However, by reducing α ( B ) , the basic reproduction number R 0 w t reaches its threshold value ( R 0 w t = 1 ). This causes changes in the dynamic behavior of the model (14), as shown in Figure 4.
Figure 11. Diagram displaying the disease-free state by fixing the parameter values Ye SPF wt/A0 from Table 3 in Model (14). As long as α ( B ) is large enough, the disease-free state persists. However, by reducing α ( B ) , the basic reproduction number R 0 w t reaches its threshold value ( R 0 w t = 1 ). This causes changes in the dynamic behavior of the model (14), as shown in Figure 4.
Biology 09 00431 g011
Table 1. The variable symbols and their meaning. The lumen is abbreviated with L, the mucosa as M. We also indicate each variable with mut or wt to denote to which population they refer. An upper-case I refers to the immune system. SBML [26,27] defines the units item and dimensionless to indicate that a quantity occurs in a piece number or that its unit originates from the cancellation of other units.
Table 1. The variable symbols and their meaning. The lumen is abbreviated with L, the mucosa as M. We also indicate each variable with mut or wt to denote to which population they refer. An upper-case I refers to the immune system. SBML [26,27] defines the units item and dimensionless to indicate that a quantity occurs in a piece number or that its unit originates from the cancellation of other units.
Variable SymbolMeaningUnits
B M Commensal bacteria in the mucosaitem
Y M ( w t ) wt Yersinia in the mucosaitem
Y M ( m u t ) mut Yersinia in the mucosaitem
B L Commensal bacteria in the lumenitem
Y L ( w t ) wt Yersinia in the lumenitem
Y L ( m u t ) mut Yersinia in the lumenitem
IStrength of immune reactiondimensionless
Table 2. Definition of the parameters.
Table 2. Definition of the parameters.
ParameterDefinitionUnit
α ( B ) Maximal growth rate of intestinal bacteria1/ h
α ( w t ) Maximal growth rate of wt Yersinia1/ h
α ( m u t ) Maximal growth rate of mut Yersinia1/ h
f I ( w t ) Immunity adjustment factor for wt Yersiniadimensionless
f I ( m u t ) Immunity adjustment factor for mut Yersiniadimensionless
C M Carrying capacity of the mucosaitem
C L Carrying capacity of the lumenitem
C I Carrying capacity of the immune systemitem
γ Maximal immunity action1/ h
κ Maximal rate of immune growth1/ h
β Rate at which intestines are discharged1/ h
Table 3. Parameter values.
Table 3. Parameter values.
ParameterValues in Ye SPF wt/A0Values in Ye SPF wt/T3S0Values in Ye GF wt/A0Values in Ye MyD88−/− wt/A0
α ( B ) 4.89 × 10−12.001.995.40 × 10−1
α ( w t ) 4.44 × 10−11.861.605.78 × 10−1
α ( m u t ) 4.44 × 10−11.861.605.78 × 10−1
f I ( w t ) 3.96 × 10−19.48 × 10−31.10 × 10−16.23 × 10−2
f I ( m u t ) 1.95 × 10−13.73 × 10−11.19 × 10−11.28 × 10−1
C M 1.76 × 1056.27 × 1031.3 × 1061.28 × 105
C L 2.14 × 1076.13 × 1064.99 × 1099.98 × 109
C I 1.001.001.001.00
γ 1.001.009.97 × 10−11.00 × 10−1
κ 7.83 × 10−14.28 × 10−16.50 × 10−14.37 × 10−1
β 2.50 × 10−12.50 × 10−18.33 × 10−21.82 × 10−1
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Mostolizadeh, R.; Dräger, A. Computational Model Informs Effective Control Interventions against Y. enterocolitica Co-Infection. Biology 2020, 9, 431. https://doi.org/10.3390/biology9120431

AMA Style

Mostolizadeh R, Dräger A. Computational Model Informs Effective Control Interventions against Y. enterocolitica Co-Infection. Biology. 2020; 9(12):431. https://doi.org/10.3390/biology9120431

Chicago/Turabian Style

Mostolizadeh, Reihaneh, and Andreas Dräger. 2020. "Computational Model Informs Effective Control Interventions against Y. enterocolitica Co-Infection" Biology 9, no. 12: 431. https://doi.org/10.3390/biology9120431

APA Style

Mostolizadeh, R., & Dräger, A. (2020). Computational Model Informs Effective Control Interventions against Y. enterocolitica Co-Infection. Biology, 9(12), 431. https://doi.org/10.3390/biology9120431

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop