# Effects of Coinfection on the Dynamics of Two Pathogens in a Tick-Host Infection Model

Authors: Bei Sun; Xue Zhang; Marco Tosato
Journal: Complexity
 (2020)

---

## Abstract

As both ticks and hosts may carry one or more pathogens, the phenomenon of coinfection of multiple tick-borne diseases becomes highly relevant and plays a key role in tick-borne disease transmission. In this paper, we propose a coinfection model involving two tick-borne diseases in a tick-host population and calculate the basic reproduction numbers at the disease-free equilibrium and two boundary equilibria. To explore the impact of coinfection, we also derive the invasion reproduction numbers which indicate the potential of a pathogen to persist when another pathogen already exists in tick and host populations. Then, we obtain the global stability of the system at the disease-free equilibrium and the boundary equilibrium, respectively, and further demonstrate the existence conditions for uniform persistence of the two diseases. The final numerical simulations mainly verify the theoretical results of coinfection.


---

## Body

## 1. Introduction

Tick-borne diseases, which mainly comprise Lyme disease, human babesiosis, tick-borne encephalitis, and human granulocytic, are becoming an increasingly significant danger to people living in countrysides or near woodlands all over the world. Specifically, Lyme disease, which is one of the most widespread tick-born diseases in the world, is caused by the pathogenBorrelia burgdorferi, and medical reports prove that there were more than 20,000 cases already as early as 2013 in America [1] and that the number has continued to increase year after year. As for human babesiosis, the parasite Babesia microti is responsible for the occurrence of this disease, and the number of human cases is in steady increase in northeastern America [2]. Tick-borne encephalitis, instead, is an arthropod-transmitted viral infection caused by TBE virus and frequently occurs in some European countries. Human granulocytic anaplasmosis is caused by Anaplasma phagocytophilum that is a bacterium, which has done great harm to human health over the last decade [3].The main species that transmit all the tick-borne diseases mentioned above areIxodes ticks, so their dynamics plays a significant role in studying the transmission of tick-borne diseases. For instance, a black-legged tick named Ixodes scapularis can be infected by a large number of different kinds of pathogens including B. burgdorferi, B. microti, and A. phagocytophilum or any combination of them concurrently [4]. An early study reported occurrences of coinfection of different pathogens in tick and host populations [5], and several case studies described the coinfection among tick-borne diseases from biotic experiments and ecological perspectives.Diuk-Wasser et al. [6] analyzed the promoting effect of coinfection of B. burgdorferi and B. microti on disease transmission in conformity with epidemiological, ecological, and clinical effects. Horowitz et al. [7] applied blood culture and serology to detect how human granulocytic anaplasmosis in Lyme disease infections played a role in the apparent rate of coinfection and in the severity to illness by studying concrete cases to explore their interaction. Welc-Faleciak et al. [8] performed a retrospective study concerning tick-exposed hosts living in southeastern Poland, aimed at probing into the risk of coinfection between Borrelia species and Anaplasma phagocytophilum or Babesia spp.. These articles have clearly shown that coinfection among tick-borne pathogens or diseases is a real phenomenon through actual case studies or experimental analysis and demonstrated that coinfection can indeed have a certain impact on the disease transmission.Mathematical models have been applied to the coinfection of infectious diseases. Gao et al. [9] established a SIS epidemic model that includes coinfection between two infectious diseases and proved a sufficient condition for coexistence of both diseases by introducing the invasion reproduction number. Tang et al. [10] integrated the relevant transmission dynamics describing coinfection of dengue and Zika into a mathematical model and focused on how dengue vaccine affected the outbreak of Zika. Wang et al. [11] proposed a Zika-dengue model emphasizing the joint dynamics between dengue and Zika and discussed the impact of vaccination against dengue and antibody-development enhancement. These studies are all related to mosquito-borne diseases; however, to our best knowledge, reports of mathematical models regarding tick-borne coinfection are relatively scarce. Exclusively, Lou et al. [2] built a tick-borne pathogen model with coinfection and discussed the promotion effect of coinfection on two diseases transmission. Unfortunately, they just considered coinfection with two determinated diseases, Lyme disease and human babesiosis, and ignored the effects of invasion on coexistence of the two pathogens.Based on the previous studies mentioned above, our objective here is principally to formulate a tick-borne coinfection model comprising tick dynamics and host dynamics with two pathogens and obtain the coexistence conditions of two tick-borne diseases. This paper is organized as follows. We establish a coinfection model involving tick and host populations in Section2. Then, we calculate basic reproduction numbers and invasion reproduction numbers and present their explicit expressions in Section 3. We further discuss the threshold dynamics at the disease-free equilibrium and two boundary equilibria, as well as the coexistence conditions of two tick-borne diseases in Section 4. Numerical simulations validate the coinfection theories and reveal the coexistence conditions of two diseases, and finally, we provide some discussions in Section 5.


## 2. The Model

We formulate a SIS-type model consisting of two tick-borne pathogens to describe the disease dynamics of coinfection as well as the effect on the spread of two tick-borne diseases. Note that we only consider one life stage of tick, nymph, which can account for the infection on the host to a great extent [12]. Both ticks and hosts can be infected with two or more pathogens through tick bites and blood meals.LetTt and Ht be the total tick population and host population, respectively. The tick population is divided into four subclasses: susceptible ticks to two diseases, infected ticks with only disease 1, infected ticks with only disease 2, and infected ticks with both diseases, which can be described as Tst, T1it, T2it, and T3it, respectively. Here, the total tick population is(1)Tt=Tst+T1it+T2it+T3it.The host population is also partitioned into four categories: susceptible hosts to two diseases, infected hosts with only disease 1, infected hosts with only disease 2, and infected hosts with both diseases, which are expressed asHst, H1it, H2it, and H3it, respectively. In the same way,(2)Ht=Hst+H1it+H2it+H3it.In combination with the above statements, a transmission diagram describing coinfection of two diseases transmission among ticks and hosts is explicitly depicted in Figure1. Then, we construct the following model:(3)dTstdt=ΛT−β1H1it+β2H2it+β3H3it+β31H3it+β32H3itTstHt−μTTst,dT1itdt=β1H1itTstHt+β31H3itTstHt−β12H2itT1itHt−β13H3itT1itHt−μTT1it,dT2itdt=β2H2itTstHt+β32H3itTstHt−β21H1itT2itHt−β23H3itT2itHt−μTT2it,dT3itdt=β3H3itTstHt+β12H2itT1itHt+β13H3itT1itHt+β21H1itT2itHt+β23H3itT2itHt−μTT3it,dHstdt=ΛH−α1T1it+α2T2it+α3T3it+α31T3it+α32T3itHstHt−μHHst,dH1itdt=α1T1itHstHt+α31T3itHstHt−α12T2itH1itHt−α13T3itH1itHt−μHH1it,dH2itdt=α2T2itHstHt+α32T3itHstHt−α21T1itH2itHt−α23T3itH2itHt−μHH2it,dH3itdt=α3T3itHstHt+α12T2itH1itHt+α13T3itH1itHt+α21T1itH2itHt+α23T3itH2itHt−μHH3it,where all parameters involved in the model are listed in Table 1. We assume the summations of transmission probabilities β1+β31 and α1+α31 are greater than original β1 and α1 without coinfection, respectively, and which also hold for the transmission probabilities of disease 2. Susceptible ticks Tst recruited at a constant rate ΛT are infected and can move into T1itT2it by blood feeding on H1it or H3it (H2it or H3it) with rates of disease transfer β1, β31≥0 (β2, β32≥0). T1itT2it can move to T3it by feeding on H2it or H3it (H1it or H3it) with infection probability β12, β13≥0 (β21, β23≥0). In addition, susceptible ticks Tst may also be infected by taking a blood meal from H3it directly with infection probability β3 and move into compartment T3it. We also consider μT which represents the exit rate from the current compartment of the tick population. Similarly, the infection process with two diseases in the host population can be derived.Figure 1
A transmission diagram of coinfection of two tick-borne pathogens among tick population (a) and host population (b), respectively.
(a)(b)Table 1
Model parameters and descriptions.
ParameterDescriptionΛTRecruitment rate of ticksΛHRecruitment rate of hostsμTExit rate of ticksμHMortality rate of hostsβ1,β2Probability of infection for susceptible ticks that feed on infectious hostsH1i or H2i, respectivelyβ3Probability of infection for susceptible ticks that feed on infectious hostsH3i and become infected by both pathogensβ31,β32Probability of infection for susceptible ticks that feed on infectious hostsH3i and become infected by only disease 1 or only disease 2, respectivelyβ12,β13Probability of infection for infectious ticks with only disease 1 that feed onH2i or H3i, respectively, and also become infected by disease 2β21,β23Probability of infection for infectious ticks with only disease 2 that feed onH1i or H3i, respectively, and also become infected by disease 1α1,α2Probability of infection for susceptible hosts that are fed by infectious ticksT1i or T2i, respectively.α3Probability of infection for susceptible hosts that are fed by infectious ticksT3i and become infected by both pathogensα31,α32Probability of infection for susceptible hosts that are fed by infectious ticksT3i and become infected by only disease 1 or only disease 2, respectivelyα12,α13Probability of infection for infectious hosts with only disease 1 that are fed byT2i or T3i, respectively, and also become infected by disease 2α21,α23Probability of infection for infectious hosts with only disease 2 that are fed byT1i or T3i, respectively, and also become infected by disease 1For the total populations of ticks and hosts, it can be shown that(4)dTtdt=ΛT−μTTt,dHtdt=ΛH−μHHt.It follows that the asymptotic equilibria of the two populations satisfylimt⟶+∞Tt=ΛT/μT and limt⟶+∞Ht=ΛH/μH. Here, Lemma 1 illustrates the basic properties of solutions in system (3), and the proof corresponding to this lemma can be derived in the Appendix of Wang [11].Lemma 1.
Letχt=Tst,T1it,T2it,T3it,Hst,H1it,H2it,H3it. All solutions χt of system (3) with nonnegative initial values χ0 remain nonnegative for t≥0. Specifically, the solutions χt satisfy χt>0 under the condition that the initial value χ0>0 for t>0. Furthermore, the solution set Γ=Tst,T1it,T2it,T3it,Hst,H1it,H2it,H3itTst≥0,T1it≥0,T2it≥0,T3it≥0,Hst≥0,H1it≥0,H2it≥0,H3it≥0,0≤Tst+T1it+T2it+T3it≤ΛT/μT,0≤Hst+H1it+H2it+H3it≤ΛH/μH is positively invariant.


## 3. The Reproduction Numbers

The reproduction numberℛ0, which is an essential threshold concept in epidemiological studies, denotes the average number of secondary cases caused by an infectious individual when it is introduced in an entirely susceptible population [11]. Mathematically, ℛ0 can measure the maximum reproductive capacity and determine whether the disease will become endemic or die out [13]. In this section, we calculate basic reproduction numbers and invasion reproduction numbers.


### 3.1. Basic Reproduction Number

In system (3), there always exists a disease-free equilibrium E0=T0s,0,0,0,H0s,0,0,0 in which T0s=ΛT/μT and H0s=ΛH/μH. Based on the calculation of the relevant next generation matrices [14], we can easily obtain the specific form of the basic reproduction number, that is,(5)ℛ0=maxℛ1,ℛ2,ℛ3,where ℛi=αiβiΛT/ΛHμT2 for i=1,2,3. Biologically, ℛ1 indicates the average number of secondary cases infected by an infectious tick or host with disease 1 when this infection is drawn into a susceptible population. Concretely, it is the product of the number of new ticks infected with disease 1, which are those produced by susceptible ticks feeding on an infectious host with disease 1 at a rate β1T0s/H during this host’s average lifetime 1/μH and the number of new hosts infected with disease 1 generated by an infectious tick with disease 1 when taking blood meals of susceptible hosts at a rate α1H0s/H during the period 1/μT. The biological interpretations for ℛ2 and ℛ3 are analogous.It is clear that there are two boundary equilibriaE1 and E2, which represent the cases in which just disease 1 and just disease 2 are endemic in the population, respectively.Therefore, atE1, we have that T2i, T3i, H2i, and H3i are all zero, and(6)E1=Ts∗,T1i∗,0,0,Hs∗,H1i∗,0,0,where(7)Ts∗=ΛTμT+λT∗,T1i∗=λT∗ΛTμTμT+λT∗,Hs∗=ΛHμH+λH∗,H1i∗=λH∗ΛHμHμH+λH∗,in which(8)λT∗=β1H1i∗H∗,λH∗=α1T1i∗H∗,with H∗=Hs∗+H1i∗.We substitute (7) into (8) and get(9)λH∗=α1λT∗ΛTμTμT+λT∗H∗.Moreover, there also exists(10)λH∗=λT∗μH2H∗β1ΛH−λT∗μHH∗.Combining (9) and (10), we obtain(11)λT∗=μHμT2H∗ℛ12−1α1ΛT+μTμHH∗.Therefore, (11) is positive if ℛ1>1, but when ℛ1<1, no positive root can be found. Thus, the equilibrium E1 exists provided ℛ1>1.Similarly, the boundary equilibriumE2 takes the following form:(12)E2=Ts∗∗,0,T2i∗∗,0,Hs∗∗,0,H2i∗∗,0,where(13)Ts∗∗=ΛTμT+λT∗∗,T2i∗∗=λT∗∗ΛTμTμT+λT∗∗,Hs∗∗=ΛHμH+λH∗∗,H2i∗∗=λH∗∗ΛHμHμH+λH∗∗,in which(14)λT∗∗=β2H2i∗∗H∗∗,λH∗∗=α2T2i∗∗H∗∗,with H∗∗=Hs∗∗+H2i∗∗. Also, the equilibrium E2 exists if ℛ2>1.


### 3.2. Invasion Reproduction Number

The basic reproduction number works as a threshold that only includes information about diseases transmission in one population containing a single pathogen and focuses on the static state, but it cannot represent the number of secondary cases infected by an infectious individual with one disease when this infection is led into a population where another disease already exists. Thus, it is necessary to introduce a new threshold quantity: the invasion reproduction number [15].In the following, we calculate two reproduction numbers which are the invasion reproduction number for disease 1 when disease 2 is already endemic and the invasion reproduction number for disease 2 when disease 1 is already endemic, respectively.Considering that the derivation of the invasion reproduction number for disease 2, namely,ℛ21, is based on a population in which disease 1 is already endemic, we concentrate on the boundary equilibrium E1 to reckon the matrices ℱ2 and V2 as follows:(15)ℱ2β2H2iTs∗H∗+H2i+H3i+β32H3iTs∗H∗+H2i+H3iβ3H3iTs∗H∗+H2i+H3i+β12H2iT1i∗H∗+H2i+H3i+β13H3iT1i∗H∗+H2i+H3iα2T2iHs∗H∗+H2i+H3i+α32T3iHs∗H∗+H2i+H3iα3T3iHs∗H∗+H2i+H3i+α12T2iH1i∗H∗+H2i+H3i+α13T3iH1i∗H∗+H2i+H3i,V2μTT2i+β21T2iH1i∗H∗+H2i+H3i+β23H3iT2iH∗+H2i+H3iμTT3i−β21T2iH1i∗H∗+H2i+H3i−β23H3iT2iH∗+H2i+H3iμHH2i+α21H2iT1i∗H∗+H2i+H3i+α23T3iH2iH∗+H2i+H3iμHH3i−α21H2iT1i∗H∗+H2i+H3i−α23T3iH2iH∗+H2i+H3i.Next, taking the derivative ofℱ2 and V2 with respect to T2i,T3i,H2i,H3iT and substituting the boundary equilibrium E1 into the variables, we have(16)F2=00β2Ts∗H∗β32Ts∗H∗00β12T1i∗H∗β3Ts∗+β13T1i∗H∗α2Hs∗H∗α32Hs∗H∗00α12H1i∗H∗α3Hs∗+α13H1i∗H∗00,V2=μT+β21H1i∗H∗000−β21H1i∗H∗μT0000μH+α21T1i∗H∗000−α21T1i∗H∗μH.Then, the next generation matrix is given by(17)F2V2−1=00ab00cdef00kg00,where(18)a=μHH∗β2Ts∗+α21T1i∗β32Ts∗μHH∗μHH∗+α21T1i∗,b=β32Ts∗μHH∗,c=β12T1i∗μHH∗+α21T1i∗β3Ts∗+β13T1i∗α21T1i∗μHH∗μHH∗+α21T1i∗,d=β3Ts∗+β13T1i∗μHH∗,e=α2Hs∗μTH∗+α32Hs∗β21H1i∗μTH∗μTH∗+β21H1i∗,f=α32Hs∗μTH∗,k=α12H1i∗μTH∗+α3Hs∗β21H1i∗+α13H1i∗β21H1i∗μTH∗μTH∗+β21H1i∗,g=α3Hs∗+α13H1i∗μTH∗,in which Ts∗, T1i∗, Hs∗, H1i∗, and H∗ are expressed in (7). The characteristic equation is(19)fλ=λ4−cf+ae+dg+bkλ2+eg−kfad−bc=0.Let(20)A=cf+ae+dg+bk,B=eg−kfad−bc.It follows from A2−4B>0 that the dominant eigenvalue of equation (19) satisfies λ2=A+A2−4B/2, and then we can acquire the invasion reproduction number for disease 2, that is, ℛ21=ρF2V2−1=A+A2−4B/2.Similarly, the invasion reproduction number for disease 1, namely,ℛ12, can also be derived following the above process, which is given by(21)ℛ12=ρF1V1−1=A˜+A˜2−4B˜2,where(22)A˜=c˜f˜+a˜e˜+d˜g˜+b˜k˜,B˜=e˜g˜−k˜f˜a˜d˜−b˜c˜,(23)a˜=μHH∗∗β1Ts∗∗+α12T2i∗∗β31Ts∗∗μHH∗∗μHH∗∗+α12T2i∗∗,b˜=β31Ts∗∗μHH∗∗,c˜=β21T2i∗∗μHH∗∗+α12T2i∗∗β3Ts∗∗+β12T2i∗∗α23T2i∗∗μHH∗∗μHH∗∗+α12T2i∗∗,d˜=β3Ts∗∗+β23T2i∗∗μHH∗∗,e˜=α1Hs∗∗μTH∗∗+α31Hs∗∗β12H2i∗∗μTH∗∗μTH∗∗+β12H2i∗∗,f˜=α31Hs∗∗μTH∗∗,k˜=α21H2i∗∗μTH∗∗+α3Hs∗∗β12H2i∗∗+α23H2i∗∗β12H2i∗∗μTH∗∗μTH∗∗+β12H2i∗∗,g˜=α3Hs∗∗+α23H2i∗∗μTH∗∗.Note thatTs∗∗, T2i∗∗, Hs∗∗, H2i∗∗, and H∗∗ are defined in (13).


## 3.1. Basic Reproduction Number

In system (3), there always exists a disease-free equilibrium E0=T0s,0,0,0,H0s,0,0,0 in which T0s=ΛT/μT and H0s=ΛH/μH. Based on the calculation of the relevant next generation matrices [14], we can easily obtain the specific form of the basic reproduction number, that is,(5)ℛ0=maxℛ1,ℛ2,ℛ3,where ℛi=αiβiΛT/ΛHμT2 for i=1,2,3. Biologically, ℛ1 indicates the average number of secondary cases infected by an infectious tick or host with disease 1 when this infection is drawn into a susceptible population. Concretely, it is the product of the number of new ticks infected with disease 1, which are those produced by susceptible ticks feeding on an infectious host with disease 1 at a rate β1T0s/H during this host’s average lifetime 1/μH and the number of new hosts infected with disease 1 generated by an infectious tick with disease 1 when taking blood meals of susceptible hosts at a rate α1H0s/H during the period 1/μT. The biological interpretations for ℛ2 and ℛ3 are analogous.It is clear that there are two boundary equilibriaE1 and E2, which represent the cases in which just disease 1 and just disease 2 are endemic in the population, respectively.Therefore, atE1, we have that T2i, T3i, H2i, and H3i are all zero, and(6)E1=Ts∗,T1i∗,0,0,Hs∗,H1i∗,0,0,where(7)Ts∗=ΛTμT+λT∗,T1i∗=λT∗ΛTμTμT+λT∗,Hs∗=ΛHμH+λH∗,H1i∗=λH∗ΛHμHμH+λH∗,in which(8)λT∗=β1H1i∗H∗,λH∗=α1T1i∗H∗,with H∗=Hs∗+H1i∗.We substitute (7) into (8) and get(9)λH∗=α1λT∗ΛTμTμT+λT∗H∗.Moreover, there also exists(10)λH∗=λT∗μH2H∗β1ΛH−λT∗μHH∗.Combining (9) and (10), we obtain(11)λT∗=μHμT2H∗ℛ12−1α1ΛT+μTμHH∗.Therefore, (11) is positive if ℛ1>1, but when ℛ1<1, no positive root can be found. Thus, the equilibrium E1 exists provided ℛ1>1.Similarly, the boundary equilibriumE2 takes the following form:(12)E2=Ts∗∗,0,T2i∗∗,0,Hs∗∗,0,H2i∗∗,0,where(13)Ts∗∗=ΛTμT+λT∗∗,T2i∗∗=λT∗∗ΛTμTμT+λT∗∗,Hs∗∗=ΛHμH+λH∗∗,H2i∗∗=λH∗∗ΛHμHμH+λH∗∗,in which(14)λT∗∗=β2H2i∗∗H∗∗,λH∗∗=α2T2i∗∗H∗∗,with H∗∗=Hs∗∗+H2i∗∗. Also, the equilibrium E2 exists if ℛ2>1.


## 3.2. Invasion Reproduction Number

The basic reproduction number works as a threshold that only includes information about diseases transmission in one population containing a single pathogen and focuses on the static state, but it cannot represent the number of secondary cases infected by an infectious individual with one disease when this infection is led into a population where another disease already exists. Thus, it is necessary to introduce a new threshold quantity: the invasion reproduction number [15].In the following, we calculate two reproduction numbers which are the invasion reproduction number for disease 1 when disease 2 is already endemic and the invasion reproduction number for disease 2 when disease 1 is already endemic, respectively.Considering that the derivation of the invasion reproduction number for disease 2, namely,ℛ21, is based on a population in which disease 1 is already endemic, we concentrate on the boundary equilibrium E1 to reckon the matrices ℱ2 and V2 as follows:(15)ℱ2β2H2iTs∗H∗+H2i+H3i+β32H3iTs∗H∗+H2i+H3iβ3H3iTs∗H∗+H2i+H3i+β12H2iT1i∗H∗+H2i+H3i+β13H3iT1i∗H∗+H2i+H3iα2T2iHs∗H∗+H2i+H3i+α32T3iHs∗H∗+H2i+H3iα3T3iHs∗H∗+H2i+H3i+α12T2iH1i∗H∗+H2i+H3i+α13T3iH1i∗H∗+H2i+H3i,V2μTT2i+β21T2iH1i∗H∗+H2i+H3i+β23H3iT2iH∗+H2i+H3iμTT3i−β21T2iH1i∗H∗+H2i+H3i−β23H3iT2iH∗+H2i+H3iμHH2i+α21H2iT1i∗H∗+H2i+H3i+α23T3iH2iH∗+H2i+H3iμHH3i−α21H2iT1i∗H∗+H2i+H3i−α23T3iH2iH∗+H2i+H3i.Next, taking the derivative ofℱ2 and V2 with respect to T2i,T3i,H2i,H3iT and substituting the boundary equilibrium E1 into the variables, we have(16)F2=00β2Ts∗H∗β32Ts∗H∗00β12T1i∗H∗β3Ts∗+β13T1i∗H∗α2Hs∗H∗α32Hs∗H∗00α12H1i∗H∗α3Hs∗+α13H1i∗H∗00,V2=μT+β21H1i∗H∗000−β21H1i∗H∗μT0000μH+α21T1i∗H∗000−α21T1i∗H∗μH.Then, the next generation matrix is given by(17)F2V2−1=00ab00cdef00kg00,where(18)a=μHH∗β2Ts∗+α21T1i∗β32Ts∗μHH∗μHH∗+α21T1i∗,b=β32Ts∗μHH∗,c=β12T1i∗μHH∗+α21T1i∗β3Ts∗+β13T1i∗α21T1i∗μHH∗μHH∗+α21T1i∗,d=β3Ts∗+β13T1i∗μHH∗,e=α2Hs∗μTH∗+α32Hs∗β21H1i∗μTH∗μTH∗+β21H1i∗,f=α32Hs∗μTH∗,k=α12H1i∗μTH∗+α3Hs∗β21H1i∗+α13H1i∗β21H1i∗μTH∗μTH∗+β21H1i∗,g=α3Hs∗+α13H1i∗μTH∗,in which Ts∗, T1i∗, Hs∗, H1i∗, and H∗ are expressed in (7). The characteristic equation is(19)fλ=λ4−cf+ae+dg+bkλ2+eg−kfad−bc=0.Let(20)A=cf+ae+dg+bk,B=eg−kfad−bc.It follows from A2−4B>0 that the dominant eigenvalue of equation (19) satisfies λ2=A+A2−4B/2, and then we can acquire the invasion reproduction number for disease 2, that is, ℛ21=ρF2V2−1=A+A2−4B/2.Similarly, the invasion reproduction number for disease 1, namely,ℛ12, can also be derived following the above process, which is given by(21)ℛ12=ρF1V1−1=A˜+A˜2−4B˜2,where(22)A˜=c˜f˜+a˜e˜+d˜g˜+b˜k˜,B˜=e˜g˜−k˜f˜a˜d˜−b˜c˜,(23)a˜=μHH∗∗β1Ts∗∗+α12T2i∗∗β31Ts∗∗μHH∗∗μHH∗∗+α12T2i∗∗,b˜=β31Ts∗∗μHH∗∗,c˜=β21T2i∗∗μHH∗∗+α12T2i∗∗β3Ts∗∗+β12T2i∗∗α23T2i∗∗μHH∗∗μHH∗∗+α12T2i∗∗,d˜=β3Ts∗∗+β23T2i∗∗μHH∗∗,e˜=α1Hs∗∗μTH∗∗+α31Hs∗∗β12H2i∗∗μTH∗∗μTH∗∗+β12H2i∗∗,f˜=α31Hs∗∗μTH∗∗,k˜=α21H2i∗∗μTH∗∗+α3Hs∗∗β12H2i∗∗+α23H2i∗∗β12H2i∗∗μTH∗∗μTH∗∗+β12H2i∗∗,g˜=α3Hs∗∗+α23H2i∗∗μTH∗∗.Note thatTs∗∗, T2i∗∗, Hs∗∗, H2i∗∗, and H∗∗ are defined in (13).


## 4. Stability and Coinfection Dynamics

To study the dynamics of two diseases transmission, it is essential to explore the stability of system (3) by studying that of the disease-free equilibrium and the one-disease equilibria for two subsystems and yield dynamic conditions for a second disease invasion and coinfection.


### 4.1. Stability of Disease-Free Equilibrium

Firstly, we discuss the stability of the disease-free equilibriumE0=T0s,0,0,0,H0s,0,0,0, where T0s=ΛT/μT and H0s=ΛH/μH.Theorem 1.
The disease-free equilibriumE0 for system (3) is globally asymptotically stable if ℛ0<1 but unstable if ℛ0>1.Proof.
As mentioned above, the tick and host populations satisfy the following two equations:(24)dTdt=ΛT−μTT,⟹T0s=ΛTμT,dHdt=ΛH−μHH,⟹H0s=ΛHμH,which means that system (3) is globally attractive for the two populations [16]. Then, there exists small enough ϵ>0 and t1=tϵ>0 such that T0s−ϵ≤Ts≤T0s+ϵ, H0s−ϵ≤Hs≤H0s+ϵ, 0≤T1i≤ϵ, 0≤T2i≤ϵ, 0≤T3i≤ϵ, 0≤H1i≤ϵ, 0≤H2i≤ϵ, and 0≤H3i≤ϵ for any t≥t1. So, for t≥t1, we can obtain(25)dT1itdt≤β1T0s+ϵH0s−ϵH1it+β31T0s+ϵH0s−ϵH3it−μTT1it,dT2itdt≤β2T0s+ϵH0s−ϵH2it+β32T0s+ϵH0s−ϵH3it−μTT2it,dT3itdt≤β3T0s+ϵH0s−ϵH3it+β12+β13+β21+β23ϵ2H0s−ϵ−μTT3it,dH1itdt≤α1H0s+ϵH0s−ϵT1it+α31H0s+ϵH0s−ϵT3it−μHH1it,dH2itdt≤α2H0s+ϵH0s−ϵT2it+α32H0s+ϵH0s−ϵT3it−μHH2it,dH3itdt≤α3H0s+ϵH0s−ϵT3it+α12+α13+α21+α23ϵ2H0s−ϵ−μHH3it.
We consider the auxiliary linear system(26)dxtdt=M0ϵxt,where vector xt=x1t,x2t,x3t,x4t,x5t,x6tT and(27)M0ϵ=−μT00β1T0s+ϵH0s−ϵ0β31T0s+ϵH0s−ϵ0−μT00β2T0s+ϵH0s−ϵβ32T0s+ϵH0s−ϵ00−μT00β3T0s+ϵH0s−ϵα1H0s+ϵH0s−ϵ0α31H0s+ϵH0s−ϵ−μH000α2H0s+ϵH0s−ϵα32H0s+ϵH0s−ϵ0−μH000α3H0s+ϵH0s−ϵ00−μH.
Whenℛ0<1, it follows from Theorem 2 in [14] that sM00<0, which means that M00 is stable. Thus, we have limt⟶+∞xit=0 for i=1,2,…,6. In accordance with the comparison principle, we have(28)limt⟶+∞T1it,T2it,T3it,H1it,H2it,H3it=0,0,0,0,0,0,and by the theory of asymptotically semiflow in [17], we get(29)limt⟶+∞Tst=T0sandlimt⟶+∞Hst=H0s.
This proves thatE0 is globally asymptotically stable when ℛ0<1. Similarly, there is small enough τ>0 and t2=tτ>0 such that T0s−τ≤Ts≤T0s+τ, H0s−τ≤Hs≤H0s+τ, 0≤T1i≤τ, 0≤T2i≤τ, 0≤T3i≤τ, 0≤H1i≤τ, 0≤H2i≤τ, and 0≤H3i≤τ for any t≥t2. Thus, for t≥t2, we also have(30)dT1itdt≥β1T0s−τH0s+τH1it+β31T0s−τH0s+τH3it−β12+β13τ2H0s+τ−μTT1it,dT2itdt≥β2T0s−τH0s+τH2it+β32T0s−τH0s+τH3it−β21+β23τ2H0s+τ−μTT2it,dT3itdt≥β3T0s−τH0s+τH3it−μTT3it,dH1itdt≥α1H0s−τH0s+τT1it+α31H0s−τH0s+τT3it−α12+α13τ2H0s+τ−μHH1it,dH2itdt≥α2H0s−τH0s+τT2it+α32H0s−τH0s+τT3it−α21+α23τ2H0s+τ−μHH2it,dH3itdt≥α3H0s−τH0s+τT3it−μHH3it.
And we can construct the auxiliary linear system(31)dxtdt=M0′τxt,where(32)M0′τ=−μT00β1T0s−τH0s+τ0β31T0s−τH0s+τ0−μT00β2T0s−τH0s+τβ32T0s−τH0s+τ00−μT00β3T0s−τH0s+τα1H0s−τH0s+τ0α31H0s−τH0s+τ−μH000α2H0s−τH0s+τα32H0s−τH0s+τ0−μH000α3H0s−τH0s+τ00−μH.
It is shown obviously thatsM0′0>0 if ℛ0>1, which implies that yit⟶+∞ when t⟶+∞ for i=1,2,…,6. Then, the disease-free equilibrium E0 is unstable.


### 4.2. Stability of Boundary Equilibria

We first analyze the two subsystems containing only one disease prior to exploring the stability of the boundary equilibriaE1 and E2.The subsystem with only disease 1 is given by(33)dTstdt=ΛT−β1H1itTstHt−μTTst,dT1itdt=β1H1itTstHt−μTT1it,dHstdt=ΛH−α1T1itHstHt−μHHst,dH1itdt=α1T1itHstHt−μHH1it.The feasible region of this subsystem (33) is Γ1=Tst,T1it,Hst,H1itTst≥0,T1it≥0,Hst≥0,H1it≥0,0≤Tst+T1it≤ΛT/μT,0≤Hst+H1it≤ΛH/μH and is positively invariant.There is a disease-free equilibriumE01=T0s,0,H0s,0 in subsystem (33), in which T0s=ΛT/μT and H0s=ΛH/μH. If ℛ1>1, then a unique endemic equilibrium E11=Ts∗,T1i∗,Hs∗,H1i∗, where Ts∗, T1i∗, Hs∗, and H1i∗ are expressed as (7), exists in subsystem (33). We aim to study the global dynamics of system (33) as follows.Theorem 2.
The following statements are valid:(i)
Ifℛ1<1, then subsystem (33) at E01 is globally asymptotically stable in Γ1(ii)
Ifℛ1>1, then there exists a unique positive solution E11 that is globally asymptotically stable for subsystem (33) in Γ1\E01Proof.
In subsystem (33), the tick and host populations satisfy(34)dTdt=ΛT−μTT,dHdt=ΛH−μHH,which imply that the subsystem is globally attractive, and there is E01=T0s,0,H0s,0, where T0s=ΛT/μT and H0s=ΛH/μH.
Ifℛ1<1, there exists small enough ζ>0 and t3=tζ>0, and when t≥t3, we take the limiting system for T1i and H1i. Considering the auxiliary system(35)dwtdt=M1ζwt,with vector wt=w1t,w2t and matrix(36)M1ζ=−μTβ1T0s+ζH0s−ζα1H0s+ζH0s−ζ−μH,we have sM10<0 if ℛ1<1 and finally obtain that E01 is globally asymptotically stable when ℛ1<1 according to the comparison principle.
Similarly, statement (ii) can also be proved via taking the limiting system forT1i and H1i, a detailed proof process similar to Lemma 2.3 in Lou and Zhao [18] and Theorem 3.2 in Gao et al. [19].Analogously, the other subsystem containing only disease 2 can be written as follows:(37)dTstdt=ΛT−β2H2itTstHt−μTTst,dT2itdt=β2H2itTstHt−μTT2it,dHstdt=ΛH−α2T2itHstHt−μHHst,dH2itdt=α2T2itHstHt−μHH2it.In view of the feasible regionΓ2=Tst,T2it,Hst,H2itTst≥0,T2it≥0,Hst≥0,H2it≥0,0≤Tst+T2it≤ΛT/μT,0≤Hst+H2it≤ΛH/μH, we can find that Γ2 is positively invariant. In this subsystem (37), there exists a disease-free equilibrium E02=T0s,0,H0s,0 in which T0s=ΛT/μT, H0s=ΛH/μH, and a unique nontrivial equilibrium E22=Ts∗∗,T2i∗∗,Hs∗∗,H2i∗∗ in which Ts∗∗, T2i∗∗, Hs∗∗, and H2i∗∗ are expressed by (13). Then, we obtain the global behaviors at E02 and E22 for subsystem (37).Theorem 3.
The following statements are valid:(i)
Ifℛ2<1, then subsystem (37) at E02 is globally asymptotically stable in Γ2(ii)
Ifℛ2>1, then there exists a unique positive solution E22 that is globally asymptotically stable for subsystem (37) in Γ2\E02Next, we discuss the relevant dynamics of general system (3) in Γ, which comprise not only these basic reproduction numbers describing the stability of the boundary equilibria E1 and E2 but also the invasion reproduction numbers explaining how a disease can invade a tick-host population where another disease is already prevalent.Theorem 4.
Whenℛ1>1, the boundary equilibrium E1 is locally asymptotically stable for system (3) if ℛ21<1 but unstable if ℛ21>1.Proof.
DefineWSE1 as the stable manifold of E1, and(38)Γ3=Ts,T1i,T2i,T3i,Hs,H1i,H2i,H3i∈ΓT2i=T3i=H2i=H3i=0,Γ4=Ts,T1i,T2i,T3i,Hs,H1i,H2i,H3i∈ΓTs=Ts∗,T1i=T1i∗,Hs=Hs∗,H1i=H1i∗.
Note thatΓ3 and Γ4 are all invariant sets for system (3) with respect to the boundary equilibrium E1. According to the center manifold theorem, we consider two invariant manifolds W1 and W2 tangent to Γ3 and Γ4, respectively, which satisfy dimW1=dimΓ3 and dimW2=dimΓ4.
Likewise,Γ3\E0 is also an invariant set for system (3), and there is an invariant manifold expressed as W3 such that dimW3=dimΓ3\E0=dimΓ3. Since E1 is globally stable in Γ3\E0 by Theorem 2, we have W3⊂WSE1. Then, we focus on the linearization of system (3) with respect to the subspace Γ4. Calculating the respective Jacobian matrix J2, we obtain(39)J2=F2−V2,where F2 and V2 are shown in (16). In accordance with Theorem 2 in [14],(40)sJ2<0>0⟺ℛ21<1>1,in which s· represents the maximum real part of the eigenvalues of J2.
Ifℛ21>1, that is, sJ2>0, the boundary equilibrium E1 is unstable. If ℛ21<1, the real parts of the eigenvalues of J2 are all negative; therefore, Γ4 is stable. At the same time, the manifold W2 is also stable. Then, we can prove that dimWSE1≥dimW3+dimW2=8; therefore, dimWSE1=8, and the equilibrium E1 is stable. Finally, when ℛ1>1, the equilibrium E1 is locally asymptotically stable for system (3) if ℛ21<1 but unstable if ℛ21>1.The proof method of Theorem4 can also be applied to the following theorem, which shows how the stability of E2 depends on the invasion reproduction number ℛ12.Theorem 5.
Whenℛ2>1, the boundary equilibrium E2 is locally asymptotically stable for system (3) if ℛ12<1 but unstable if ℛ12>1.


### 4.3. Coexistence of Two Diseases

We have shown that the boundary equilibriumE1 (or E2) is unstable in the tick-host population where disease 1 (or 2) is already prevalent once the invasion reproduction number ℛ21 (or ℛ12) is greater than one, which means that disease 2 (or 1) is able to invade in the population. From this perspective, the instability of E1 and E2 may be an advantage to study the coexistence of these two diseases. Here, we present Theorem 6 involving disease persistence.Theorem 6.
Ifℛ1>1, ℛ2>1, ℛ21>1, and ℛ12>1, then system (3) remains uniformly persistent, namely, there exists a constant ε>0 such that every solution χt≡Tst,T1it,T2it,T3it,Hst,H1it,H2it,H3it for system (3) with positive initial value χ0≡Ts0,T1i0,T2i0,T3i0,Hs0,H1i0,H2i0,H3i0 satisfies(41)lim inft⟶+∞Ts,T1i,T2i,T3i,Hs,H1i,H2i,H3i≥ε,ε,ε,ε,ε,ε,ε,ε.Proof.
Define(42)X0=Ts,T1i,T2i,T3i,Hs,H1i,H2i,H3i∈ΓTs>0,T1i>0,T2i>0,T3i>0,Hs>0,H1i>0,H2i>0,H3i>0,(43)∂X0=Γ\X0.
It can be shown that system (3) remains uniformly persistent regarding X0,∂X0.
By Lemma1, first of all, we know that Γ and X0 are positively invariant in system (3). Also, X0 is relatively closed in Γ, and system (3) is point dissipative. Let(44)M∂=χtχtsatisfies system3as well asχt∈∂X0,∀t≥0,and set(45)M¯∂0=Ts,0,0,0,Hs,0,0,0Ts≥0,Hs≥0,M¯∂1=Ts,T1i,0,0,Hs,H1i,0,0Ts≥0,Hs≥0,T1i≥0,H1i≥0,M¯∂2=Ts,0,T2i,0,Hs,0,H2i,0Ts≥0,Hs≥0,T2i≥0,H2i≥0.
AssumeM¯∂=M¯∂0∪M¯∂1∪M¯∂2. Then, we prove that(46)M∂=M¯∂.
Firstly,M¯∂⊆M∂ clearly holds, so we concentrate on the proof of M∂⊆M¯∂.
Proving by contradiction, we suppose thatM∂⊆M¯∂ does not hold. As χt denotes the solution of system (3) with positive initial value χ0, then χt∈M∂ for t≥0, and χt∉M¯∂. Thus, χt≡Tst,T1it,T2it,T3it,Hst,H1it,H2it,H3it is in accordance with one of these conditions:(1)
There existst0≥0 such that at least one of T3it0 and H3it0 is not zero(2)
IfT3it0=H3it0=0 as t0≥0, then the solution satisfies the following: (i) there exists t0≥0 such that at least one of T1it0 and H1it0 is not zero; (ii) there exists t0≥0 such that at least one of T2it0 and H2it0 is not zero
We supposeχt satisfies condition (1). Without loss of generality, we may as well presume that T3it0>0 and H3it0=0. According to Lemma 1, the solution T3it>0 can be given for any t>t0. For system (3) and for t>t0, we have(47)H1it=H1it0exp−∫t0tα12T2isHds−∫t0tα13T3isHds−μHt−t0+exp−∫t0tα12T2isHds−∫t0tα13T3isHds−μHt−t0∫t0tα1T1isHssH+α31T3isHssHexp∫stα12T2iξH+α13T3iξHdξ+μHt−sds>0,H2it=H2it0exp−∫t0tα21T1isHds−∫t0tα23T3isHds−μHt−t0+exp−∫t0tα21T1isHds−∫t0tα23T3isHds−μHt−t0∫t0tα2T2isHssH+α32T3isHssHexp∫stα21T1iξH+α23T3iξHdξ+μHt−sds>0,T1it=T1it0exp−∫t0tβ12H2isHds−∫t0tβ13H3isHds−μTt−t0+exp−∫t0tβ12H2isHds−∫t0tβ13H3isHds−μHt−t0∫t0tβ1H1isTssH+β31H3isTssHexp∫stβ12H2iξH+β13H3iξHdξ+μTt−sds>0,T2it=T2it0exp−∫t0tβ21H1isHds−∫t0tβ23H3isHds−μTt−t0+exp−∫t0tβ21H1isHds−∫t0tβ23H3isHds−μTt−t0∫t0tβ2H2isTssH+β32H3isTssHexp∫stβ21H1iξH+β23H3iξHdξ+μTt−sds>0,H3it=H3it0exp−μHt−t0+exp−μHt−t0∫t0tα3T3isHssH+α12T2isH1isH+α13T3isH1isH+α21T1isH2isH+α23T3isH2isHexpμHt−sds>0.
This implies thatχt∉∂X0 when t>t0. Likewise, we have χt∉∂X0 for t>t0 using condition (2), which contradicts the hypothesis χt∈∂X0. Therefore, it follows that M∂⊆M¯∂, and then M∂=M¯∂.
We know that there exist three equilibriaE0, E1, and E2 in which all are isolated invariant sets in Γ. In addition, we can show that WSE0∩X0=∅, WSE1∩X0=∅, and WSE2∩X0=∅, in which WSEi represents the stable manifold of Ei for i=0,1,2, namely, there is κ>0 such that, for every solution Φtχ0 of system (3) with χ0∈X0, we obtain(48)lim supt⟶+∞dΦtχ0,E0≥κ,lim supt⟶+∞dΦtχ0,E1≥κ,lim supt⟶+∞dΦtχ0,E2≥κ,in which d is a distance function for X0.
Without loss of generality, we only prove thatlimsupt⟶+∞dΦtχ0,E1≥κ. Suppose it does not hold; then, limsupt⟶+∞dΦtχ0,E1<ε for any ε>0, that is, there is t4=tε>0 such that Ts∗−ε≤Ts≤Ts∗+ε, T1i∗−ε≤T1i≤T1i∗+ε, Hs∗−ε≤Hs≤Hs∗+ε, H1i∗−ε≤H1i≤H1i∗+ε, 0≤T2i≤ε, 0≤T3i≤ε, 0≤H2i≤ε, and 0≤H3i≤ε for any t>t4. In system (3) for t>t4, we have(49)dT2itdt≥β2H2iTs∗−εH∗+β32H3iTs∗−εH∗−β21T2iH1i∗+εH∗−β23ε2H∗−μTT2i,dT3itdt≥β3H3iTs∗−εH∗+β12H2i+β13H3iT1i∗−εH∗+β21T2iH1i∗−εH∗−μTT3i,dH2itdt≥α2T2iHs∗−εH∗+α32T3iHs∗−εH∗−α21H2iT1i∗+εH∗−α23ε2H∗−μHH2i,dH3itdt≥α3T3iHs∗−εH∗+α12T2i+α13T3iH1i∗−εH∗+α21H2iT1i∗−εH∗−μHH3i.
Then, we consider the following auxiliary system:(50)dudt=M2εu,in which the vector u=u1,u2,u3,u4T and(51)M2ε=−β21H1i∗+εH∗−μT0β2Ts∗−εH∗β32Ts∗−εH∗β21H1i∗−εH∗−μTβ12T1i∗−εH∗β3Ts∗−ε+β13T1i∗−εH∗α2Hs∗−εH∗α32Hs∗−εH∗−μH−α21T1i∗+εH∗0α12H1i∗−εH∗α3Hs∗−ε+α13H1i∗−εH∗α21T1i∗−εH∗−μH.
It is clear thatM20=F2−V2, where F2 and V2 are given by (16). By referring to Theorem 2 in [14], we get that if ℛ21>1, then sM20>0. As M2ε is continuous for small ε, there exists sufficiently small ε such that sM2ε>0. Consequently, we obtain a positive eigenvalue of M2ε with a positive eigenvector. Let ut=u1t,u2t,u3t,u4tT denote a positive solution for auxiliary system (50), which is strictly increasing, namely, uit⟶+∞ as t⟶+∞ for i=1,2,3,4. Using the comparison principle, we have(52)limt⟶+∞T2it=+∞,limt⟶+∞T3it=+∞,limt⟶+∞H2it=+∞,limt⟶+∞H3it=+∞,and this conflicts with the above assumption.
Thus,E1 is not only an isolated invariant set for Γ, but it also satisfies WSE1∩X0=∅. Analogously, we can also verify that E2 is not only an isolated invariant set for Γ, but it also satisfies WSE2∩X0=∅ when ℛ12>1. In addition, when ℛ1>1 and ℛ2>1, E0 is not only an isolated invariant set for Γ, but it also satisfies WSE0∩X0=∅. According to Thieme [20], system (3) is uniformly persistent if ℛ1>1, ℛ2>1, ℛ21>1, and ℛ12>1.


## 4.1. Stability of Disease-Free Equilibrium

Firstly, we discuss the stability of the disease-free equilibriumE0=T0s,0,0,0,H0s,0,0,0, where T0s=ΛT/μT and H0s=ΛH/μH.Theorem 1.
The disease-free equilibriumE0 for system (3) is globally asymptotically stable if ℛ0<1 but unstable if ℛ0>1.Proof.
As mentioned above, the tick and host populations satisfy the following two equations:(24)dTdt=ΛT−μTT,⟹T0s=ΛTμT,dHdt=ΛH−μHH,⟹H0s=ΛHμH,which means that system (3) is globally attractive for the two populations [16]. Then, there exists small enough ϵ>0 and t1=tϵ>0 such that T0s−ϵ≤Ts≤T0s+ϵ, H0s−ϵ≤Hs≤H0s+ϵ, 0≤T1i≤ϵ, 0≤T2i≤ϵ, 0≤T3i≤ϵ, 0≤H1i≤ϵ, 0≤H2i≤ϵ, and 0≤H3i≤ϵ for any t≥t1. So, for t≥t1, we can obtain(25)dT1itdt≤β1T0s+ϵH0s−ϵH1it+β31T0s+ϵH0s−ϵH3it−μTT1it,dT2itdt≤β2T0s+ϵH0s−ϵH2it+β32T0s+ϵH0s−ϵH3it−μTT2it,dT3itdt≤β3T0s+ϵH0s−ϵH3it+β12+β13+β21+β23ϵ2H0s−ϵ−μTT3it,dH1itdt≤α1H0s+ϵH0s−ϵT1it+α31H0s+ϵH0s−ϵT3it−μHH1it,dH2itdt≤α2H0s+ϵH0s−ϵT2it+α32H0s+ϵH0s−ϵT3it−μHH2it,dH3itdt≤α3H0s+ϵH0s−ϵT3it+α12+α13+α21+α23ϵ2H0s−ϵ−μHH3it.
We consider the auxiliary linear system(26)dxtdt=M0ϵxt,where vector xt=x1t,x2t,x3t,x4t,x5t,x6tT and(27)M0ϵ=−μT00β1T0s+ϵH0s−ϵ0β31T0s+ϵH0s−ϵ0−μT00β2T0s+ϵH0s−ϵβ32T0s+ϵH0s−ϵ00−μT00β3T0s+ϵH0s−ϵα1H0s+ϵH0s−ϵ0α31H0s+ϵH0s−ϵ−μH000α2H0s+ϵH0s−ϵα32H0s+ϵH0s−ϵ0−μH000α3H0s+ϵH0s−ϵ00−μH.
Whenℛ0<1, it follows from Theorem 2 in [14] that sM00<0, which means that M00 is stable. Thus, we have limt⟶+∞xit=0 for i=1,2,…,6. In accordance with the comparison principle, we have(28)limt⟶+∞T1it,T2it,T3it,H1it,H2it,H3it=0,0,0,0,0,0,and by the theory of asymptotically semiflow in [17], we get(29)limt⟶+∞Tst=T0sandlimt⟶+∞Hst=H0s.
This proves thatE0 is globally asymptotically stable when ℛ0<1. Similarly, there is small enough τ>0 and t2=tτ>0 such that T0s−τ≤Ts≤T0s+τ, H0s−τ≤Hs≤H0s+τ, 0≤T1i≤τ, 0≤T2i≤τ, 0≤T3i≤τ, 0≤H1i≤τ, 0≤H2i≤τ, and 0≤H3i≤τ for any t≥t2. Thus, for t≥t2, we also have(30)dT1itdt≥β1T0s−τH0s+τH1it+β31T0s−τH0s+τH3it−β12+β13τ2H0s+τ−μTT1it,dT2itdt≥β2T0s−τH0s+τH2it+β32T0s−τH0s+τH3it−β21+β23τ2H0s+τ−μTT2it,dT3itdt≥β3T0s−τH0s+τH3it−μTT3it,dH1itdt≥α1H0s−τH0s+τT1it+α31H0s−τH0s+τT3it−α12+α13τ2H0s+τ−μHH1it,dH2itdt≥α2H0s−τH0s+τT2it+α32H0s−τH0s+τT3it−α21+α23τ2H0s+τ−μHH2it,dH3itdt≥α3H0s−τH0s+τT3it−μHH3it.
And we can construct the auxiliary linear system(31)dxtdt=M0′τxt,where(32)M0′τ=−μT00β1T0s−τH0s+τ0β31T0s−τH0s+τ0−μT00β2T0s−τH0s+τβ32T0s−τH0s+τ00−μT00β3T0s−τH0s+τα1H0s−τH0s+τ0α31H0s−τH0s+τ−μH000α2H0s−τH0s+τα32H0s−τH0s+τ0−μH000α3H0s−τH0s+τ00−μH.
It is shown obviously thatsM0′0>0 if ℛ0>1, which implies that yit⟶+∞ when t⟶+∞ for i=1,2,…,6. Then, the disease-free equilibrium E0 is unstable.


## 4.2. Stability of Boundary Equilibria

We first analyze the two subsystems containing only one disease prior to exploring the stability of the boundary equilibriaE1 and E2.The subsystem with only disease 1 is given by(33)dTstdt=ΛT−β1H1itTstHt−μTTst,dT1itdt=β1H1itTstHt−μTT1it,dHstdt=ΛH−α1T1itHstHt−μHHst,dH1itdt=α1T1itHstHt−μHH1it.The feasible region of this subsystem (33) is Γ1=Tst,T1it,Hst,H1itTst≥0,T1it≥0,Hst≥0,H1it≥0,0≤Tst+T1it≤ΛT/μT,0≤Hst+H1it≤ΛH/μH and is positively invariant.There is a disease-free equilibriumE01=T0s,0,H0s,0 in subsystem (33), in which T0s=ΛT/μT and H0s=ΛH/μH. If ℛ1>1, then a unique endemic equilibrium E11=Ts∗,T1i∗,Hs∗,H1i∗, where Ts∗, T1i∗, Hs∗, and H1i∗ are expressed as (7), exists in subsystem (33). We aim to study the global dynamics of system (33) as follows.Theorem 2.
The following statements are valid:(i)
Ifℛ1<1, then subsystem (33) at E01 is globally asymptotically stable in Γ1(ii)
Ifℛ1>1, then there exists a unique positive solution E11 that is globally asymptotically stable for subsystem (33) in Γ1\E01Proof.
In subsystem (33), the tick and host populations satisfy(34)dTdt=ΛT−μTT,dHdt=ΛH−μHH,which imply that the subsystem is globally attractive, and there is E01=T0s,0,H0s,0, where T0s=ΛT/μT and H0s=ΛH/μH.
Ifℛ1<1, there exists small enough ζ>0 and t3=tζ>0, and when t≥t3, we take the limiting system for T1i and H1i. Considering the auxiliary system(35)dwtdt=M1ζwt,with vector wt=w1t,w2t and matrix(36)M1ζ=−μTβ1T0s+ζH0s−ζα1H0s+ζH0s−ζ−μH,we have sM10<0 if ℛ1<1 and finally obtain that E01 is globally asymptotically stable when ℛ1<1 according to the comparison principle.
Similarly, statement (ii) can also be proved via taking the limiting system forT1i and H1i, a detailed proof process similar to Lemma 2.3 in Lou and Zhao [18] and Theorem 3.2 in Gao et al. [19].Analogously, the other subsystem containing only disease 2 can be written as follows:(37)dTstdt=ΛT−β2H2itTstHt−μTTst,dT2itdt=β2H2itTstHt−μTT2it,dHstdt=ΛH−α2T2itHstHt−μHHst,dH2itdt=α2T2itHstHt−μHH2it.In view of the feasible regionΓ2=Tst,T2it,Hst,H2itTst≥0,T2it≥0,Hst≥0,H2it≥0,0≤Tst+T2it≤ΛT/μT,0≤Hst+H2it≤ΛH/μH, we can find that Γ2 is positively invariant. In this subsystem (37), there exists a disease-free equilibrium E02=T0s,0,H0s,0 in which T0s=ΛT/μT, H0s=ΛH/μH, and a unique nontrivial equilibrium E22=Ts∗∗,T2i∗∗,Hs∗∗,H2i∗∗ in which Ts∗∗, T2i∗∗, Hs∗∗, and H2i∗∗ are expressed by (13). Then, we obtain the global behaviors at E02 and E22 for subsystem (37).Theorem 3.
The following statements are valid:(i)
Ifℛ2<1, then subsystem (37) at E02 is globally asymptotically stable in Γ2(ii)
Ifℛ2>1, then there exists a unique positive solution E22 that is globally asymptotically stable for subsystem (37) in Γ2\E02Next, we discuss the relevant dynamics of general system (3) in Γ, which comprise not only these basic reproduction numbers describing the stability of the boundary equilibria E1 and E2 but also the invasion reproduction numbers explaining how a disease can invade a tick-host population where another disease is already prevalent.Theorem 4.
Whenℛ1>1, the boundary equilibrium E1 is locally asymptotically stable for system (3) if ℛ21<1 but unstable if ℛ21>1.Proof.
DefineWSE1 as the stable manifold of E1, and(38)Γ3=Ts,T1i,T2i,T3i,Hs,H1i,H2i,H3i∈ΓT2i=T3i=H2i=H3i=0,Γ4=Ts,T1i,T2i,T3i,Hs,H1i,H2i,H3i∈ΓTs=Ts∗,T1i=T1i∗,Hs=Hs∗,H1i=H1i∗.
Note thatΓ3 and Γ4 are all invariant sets for system (3) with respect to the boundary equilibrium E1. According to the center manifold theorem, we consider two invariant manifolds W1 and W2 tangent to Γ3 and Γ4, respectively, which satisfy dimW1=dimΓ3 and dimW2=dimΓ4.
Likewise,Γ3\E0 is also an invariant set for system (3), and there is an invariant manifold expressed as W3 such that dimW3=dimΓ3\E0=dimΓ3. Since E1 is globally stable in Γ3\E0 by Theorem 2, we have W3⊂WSE1. Then, we focus on the linearization of system (3) with respect to the subspace Γ4. Calculating the respective Jacobian matrix J2, we obtain(39)J2=F2−V2,where F2 and V2 are shown in (16). In accordance with Theorem 2 in [14],(40)sJ2<0>0⟺ℛ21<1>1,in which s· represents the maximum real part of the eigenvalues of J2.
Ifℛ21>1, that is, sJ2>0, the boundary equilibrium E1 is unstable. If ℛ21<1, the real parts of the eigenvalues of J2 are all negative; therefore, Γ4 is stable. At the same time, the manifold W2 is also stable. Then, we can prove that dimWSE1≥dimW3+dimW2=8; therefore, dimWSE1=8, and the equilibrium E1 is stable. Finally, when ℛ1>1, the equilibrium E1 is locally asymptotically stable for system (3) if ℛ21<1 but unstable if ℛ21>1.The proof method of Theorem4 can also be applied to the following theorem, which shows how the stability of E2 depends on the invasion reproduction number ℛ12.Theorem 5.
Whenℛ2>1, the boundary equilibrium E2 is locally asymptotically stable for system (3) if ℛ12<1 but unstable if ℛ12>1.


## 4.3. Coexistence of Two Diseases

We have shown that the boundary equilibriumE1 (or E2) is unstable in the tick-host population where disease 1 (or 2) is already prevalent once the invasion reproduction number ℛ21 (or ℛ12) is greater than one, which means that disease 2 (or 1) is able to invade in the population. From this perspective, the instability of E1 and E2 may be an advantage to study the coexistence of these two diseases. Here, we present Theorem 6 involving disease persistence.Theorem 6.
Ifℛ1>1, ℛ2>1, ℛ21>1, and ℛ12>1, then system (3) remains uniformly persistent, namely, there exists a constant ε>0 such that every solution χt≡Tst,T1it,T2it,T3it,Hst,H1it,H2it,H3it for system (3) with positive initial value χ0≡Ts0,T1i0,T2i0,T3i0,Hs0,H1i0,H2i0,H3i0 satisfies(41)lim inft⟶+∞Ts,T1i,T2i,T3i,Hs,H1i,H2i,H3i≥ε,ε,ε,ε,ε,ε,ε,ε.Proof.
Define(42)X0=Ts,T1i,T2i,T3i,Hs,H1i,H2i,H3i∈ΓTs>0,T1i>0,T2i>0,T3i>0,Hs>0,H1i>0,H2i>0,H3i>0,(43)∂X0=Γ\X0.
It can be shown that system (3) remains uniformly persistent regarding X0,∂X0.
By Lemma1, first of all, we know that Γ and X0 are positively invariant in system (3). Also, X0 is relatively closed in Γ, and system (3) is point dissipative. Let(44)M∂=χtχtsatisfies system3as well asχt∈∂X0,∀t≥0,and set(45)M¯∂0=Ts,0,0,0,Hs,0,0,0Ts≥0,Hs≥0,M¯∂1=Ts,T1i,0,0,Hs,H1i,0,0Ts≥0,Hs≥0,T1i≥0,H1i≥0,M¯∂2=Ts,0,T2i,0,Hs,0,H2i,0Ts≥0,Hs≥0,T2i≥0,H2i≥0.
AssumeM¯∂=M¯∂0∪M¯∂1∪M¯∂2. Then, we prove that(46)M∂=M¯∂.
Firstly,M¯∂⊆M∂ clearly holds, so we concentrate on the proof of M∂⊆M¯∂.
Proving by contradiction, we suppose thatM∂⊆M¯∂ does not hold. As χt denotes the solution of system (3) with positive initial value χ0, then χt∈M∂ for t≥0, and χt∉M¯∂. Thus, χt≡Tst,T1it,T2it,T3it,Hst,H1it,H2it,H3it is in accordance with one of these conditions:(1)
There existst0≥0 such that at least one of T3it0 and H3it0 is not zero(2)
IfT3it0=H3it0=0 as t0≥0, then the solution satisfies the following: (i) there exists t0≥0 such that at least one of T1it0 and H1it0 is not zero; (ii) there exists t0≥0 such that at least one of T2it0 and H2it0 is not zero
We supposeχt satisfies condition (1). Without loss of generality, we may as well presume that T3it0>0 and H3it0=0. According to Lemma 1, the solution T3it>0 can be given for any t>t0. For system (3) and for t>t0, we have(47)H1it=H1it0exp−∫t0tα12T2isHds−∫t0tα13T3isHds−μHt−t0+exp−∫t0tα12T2isHds−∫t0tα13T3isHds−μHt−t0∫t0tα1T1isHssH+α31T3isHssHexp∫stα12T2iξH+α13T3iξHdξ+μHt−sds>0,H2it=H2it0exp−∫t0tα21T1isHds−∫t0tα23T3isHds−μHt−t0+exp−∫t0tα21T1isHds−∫t0tα23T3isHds−μHt−t0∫t0tα2T2isHssH+α32T3isHssHexp∫stα21T1iξH+α23T3iξHdξ+μHt−sds>0,T1it=T1it0exp−∫t0tβ12H2isHds−∫t0tβ13H3isHds−μTt−t0+exp−∫t0tβ12H2isHds−∫t0tβ13H3isHds−μHt−t0∫t0tβ1H1isTssH+β31H3isTssHexp∫stβ12H2iξH+β13H3iξHdξ+μTt−sds>0,T2it=T2it0exp−∫t0tβ21H1isHds−∫t0tβ23H3isHds−μTt−t0+exp−∫t0tβ21H1isHds−∫t0tβ23H3isHds−μTt−t0∫t0tβ2H2isTssH+β32H3isTssHexp∫stβ21H1iξH+β23H3iξHdξ+μTt−sds>0,H3it=H3it0exp−μHt−t0+exp−μHt−t0∫t0tα3T3isHssH+α12T2isH1isH+α13T3isH1isH+α21T1isH2isH+α23T3isH2isHexpμHt−sds>0.
This implies thatχt∉∂X0 when t>t0. Likewise, we have χt∉∂X0 for t>t0 using condition (2), which contradicts the hypothesis χt∈∂X0. Therefore, it follows that M∂⊆M¯∂, and then M∂=M¯∂.
We know that there exist three equilibriaE0, E1, and E2 in which all are isolated invariant sets in Γ. In addition, we can show that WSE0∩X0=∅, WSE1∩X0=∅, and WSE2∩X0=∅, in which WSEi represents the stable manifold of Ei for i=0,1,2, namely, there is κ>0 such that, for every solution Φtχ0 of system (3) with χ0∈X0, we obtain(48)lim supt⟶+∞dΦtχ0,E0≥κ,lim supt⟶+∞dΦtχ0,E1≥κ,lim supt⟶+∞dΦtχ0,E2≥κ,in which d is a distance function for X0.
Without loss of generality, we only prove thatlimsupt⟶+∞dΦtχ0,E1≥κ. Suppose it does not hold; then, limsupt⟶+∞dΦtχ0,E1<ε for any ε>0, that is, there is t4=tε>0 such that Ts∗−ε≤Ts≤Ts∗+ε, T1i∗−ε≤T1i≤T1i∗+ε, Hs∗−ε≤Hs≤Hs∗+ε, H1i∗−ε≤H1i≤H1i∗+ε, 0≤T2i≤ε, 0≤T3i≤ε, 0≤H2i≤ε, and 0≤H3i≤ε for any t>t4. In system (3) for t>t4, we have(49)dT2itdt≥β2H2iTs∗−εH∗+β32H3iTs∗−εH∗−β21T2iH1i∗+εH∗−β23ε2H∗−μTT2i,dT3itdt≥β3H3iTs∗−εH∗+β12H2i+β13H3iT1i∗−εH∗+β21T2iH1i∗−εH∗−μTT3i,dH2itdt≥α2T2iHs∗−εH∗+α32T3iHs∗−εH∗−α21H2iT1i∗+εH∗−α23ε2H∗−μHH2i,dH3itdt≥α3T3iHs∗−εH∗+α12T2i+α13T3iH1i∗−εH∗+α21H2iT1i∗−εH∗−μHH3i.
Then, we consider the following auxiliary system:(50)dudt=M2εu,in which the vector u=u1,u2,u3,u4T and(51)M2ε=−β21H1i∗+εH∗−μT0β2Ts∗−εH∗β32Ts∗−εH∗β21H1i∗−εH∗−μTβ12T1i∗−εH∗β3Ts∗−ε+β13T1i∗−εH∗α2Hs∗−εH∗α32Hs∗−εH∗−μH−α21T1i∗+εH∗0α12H1i∗−εH∗α3Hs∗−ε+α13H1i∗−εH∗α21T1i∗−εH∗−μH.
It is clear thatM20=F2−V2, where F2 and V2 are given by (16). By referring to Theorem 2 in [14], we get that if ℛ21>1, then sM20>0. As M2ε is continuous for small ε, there exists sufficiently small ε such that sM2ε>0. Consequently, we obtain a positive eigenvalue of M2ε with a positive eigenvector. Let ut=u1t,u2t,u3t,u4tT denote a positive solution for auxiliary system (50), which is strictly increasing, namely, uit⟶+∞ as t⟶+∞ for i=1,2,3,4. Using the comparison principle, we have(52)limt⟶+∞T2it=+∞,limt⟶+∞T3it=+∞,limt⟶+∞H2it=+∞,limt⟶+∞H3it=+∞,and this conflicts with the above assumption.
Thus,E1 is not only an isolated invariant set for Γ, but it also satisfies WSE1∩X0=∅. Analogously, we can also verify that E2 is not only an isolated invariant set for Γ, but it also satisfies WSE2∩X0=∅ when ℛ12>1. In addition, when ℛ1>1 and ℛ2>1, E0 is not only an isolated invariant set for Γ, but it also satisfies WSE0∩X0=∅. According to Thieme [20], system (3) is uniformly persistent if ℛ1>1, ℛ2>1, ℛ21>1, and ℛ12>1.


## 5. Numerical Illustrations and Discussion

In this section, we utilize specific numerical simulations to illustrate the impact of coinfection on the dynamics of two tick-borne pathogens and ensure the effectiveness of the previous theories. Besides, we also propose and verify the uniform persistence of the two diseases and the presence of a coexistence equilibrium in the case in whichℛ1>1 and ℛ2<1<ℛ21 or ℛ2>1 and ℛ1<1<ℛ12 in the following.By referring to some early works [2, 11, 21–23], we fix the parameter values as follows: ΛT=5000day−1, μT=1day−1, ΛH=200day−1, μH=0.1day−1, β1=0.4, α1=0.3, β2=0.25, and α2=0.15. First, we consider subsystems (33) and (37) and set the initial values of these subsystems to Ts0,T1i0,Hs0,H1i0=Ts0,T2i0,Hs0,H2i0=100,250,75,200. According to Figure 2(a), it is clearly shown that the number of infectious individuals with disease 1 consisting of infected ticks and hosts will eventually reach the equilibrium E11=4105,952.5,853.5,1146, which means that disease 1 would remain endemic since ℛ1=1.025>1. On the contrary, it follows from Figure 2(b) that disease 2 will not persist in the relative tick and host populations as expected since ℛ2=0.943<1. Thus, we can conclude that disease 2 cannot survive alone in tick and host populations, but disease 1 can.Figure 2
Solutions of submodels (33) and (37). The numbers of susceptible and infectious individuals are represented by blue and red curves, respectively, where solid lines denote the number of susceptible or infectious ticks, and dotted lines denote the number of susceptible or infectious hosts. (a) Subsystem (33) stabilizes at equilibrium E11 as ℛ1=1.025>1 for this parameter choice, which indicates that disease 1 can be endemic in tick and host populations; (b) susceptible ticks and hosts eventually approach and stabilize at 4990 and 1985, respectively, but infected ticks and hosts die out as ℛ2=0.943<1, which implies that disease 2 cannot spread in the tick-host population.
(a)(b)Next, we introduce disease 2, which cannot spread if taken singularly, into tick and host populations where disease 1 already exists in order to study the influence of coinfection. We choose the following as parameter values:β32=β31=β21=β23=α32=α31=α23=α21=0.1, β3=α3=0.25, and β12=β13=α12=α13=0.12. The initial values of T2i,T3i,H2i,H3i are set to T2i0,T3i0,H2i0,H3i0=100,200,150,100, and we calculate the invasion reproduction number ℛ21=1.01>1. From Figure 3, we see that infected ticks and hosts with only disease 2 approach and stabilize at 76.4 and 65.9, and ones with both pathogens stabilize at 153.7 and 302.9, respectively, which imply that disease 2 can invade tick and host populations in which disease 1 is already endemic if ℛ21>1. In order to visualize theoretical results, we choose β3 and α3 as parameters and plot Figure 4 to show their effects on the invasion reproduction number ℛ21.Figure 3
Solutions of model (3) reflecting the influence of the invasion of disease 2. The numbers of infectious ticks and hosts are represented by blue and red curves, respectively, where dotted lines denote the number of infected individuals with only disease 2, and solid lines denote the number of infected individuals with both diseases. It is shown that infectious individuals related to both diseases tend to persist, and then disease 2 can invade the tick-host population as ℛ21=1.01>1.Figure 4
Effects ofβ3 and α3 on the invasion reproduction number ℛ21.In this paper, we have developed a coinfection model involving two tick-borne diseases to explore their interaction and spread. Based primarily on the works by Van den Driessche [14], we have calculated the basic reproduction numbers concerning the disease-free equilibrium E0 and two boundary equilibria E1 and E2. Considering the invasion of a disease, we have then proposed a new threshold, namely, the invasion reproduction number, to describe the number of secondary cases infected by a disease later introduced into a population where another disease is already endemic, and we have derived the expressions for these two reproduction numbers.According to the comparison principle and the theory of asymptotically semiflow, we have proved the global stability at the disease-free equilibriumE0 if ℛ0<1 and the local stability at the two boundary equilibria E1 and E2 if ℛ1>1ℛ2>1 and ℛ21<1ℛ12<1, as well as the coexistence conditions for these two diseases, which are ℛ1>1, ℛ2>1, ℛ21>1, and ℛ12>1. We have also showed that the two diseases are uniformly persistent, and the transmission system of two diseases admits a coexistence equilibrium when ℛ1>1 (or ℛ2>1) and ℛ2<1<ℛ21 (or ℛ1<1<ℛ12), which implies that the existence of one disease can facilitate the invasion and transmission of another disease. These results have been verified by numerical simulations.The highlight of this paper is the introduction of the invasion reproduction number, which can directly reflect whether one disease can invade a population where another disease is already persistent, in order to investigate the effect of coinfection. We conclude that coinfection has a positive effect on the transmission of two tick-borne diseases in model (3), that is, one disease may invade a tick-host population where another disease already exists, although it cannot remain endemic alone in the population.However, theoretically, pathogens or diseases interacting in a population may not only promote each other but also compete or be independent one from the other [2]. Thus, in future research, we will focus on the negative effects of coinfection on diseases transmission. Furthermore, other factors such as temporal and spatial heterogeneity will also be worth discussing.



---

*Source: 5615173-2020-05-28.xml*