Citation: Zongwei Ma, Hongying Shu. Viral infection dynamics in a spatial heterogeneous environment with cell-free and cell-to-cell transmissions[J]. Mathematical Biosciences and Engineering, 2020, 17(3): 2569-2591. doi: 10.3934/mbe.2020141
[1] | Ting Guo, Zhipeng Qiu . The effects of CTL immune response on HIV infection model with potent therapy, latently infected cells and cell-to-cell viral transmission. Mathematical Biosciences and Engineering, 2019, 16(6): 6822-6841. doi: 10.3934/mbe.2019341 |
[2] | Jiawei Deng, Ping Jiang, Hongying Shu . Viral infection dynamics with mitosis, intracellular delays and immune response. Mathematical Biosciences and Engineering, 2023, 20(2): 2937-2963. doi: 10.3934/mbe.2023139 |
[3] | Jianquan Li, Xiaoyu Huo, Yuming Chen . Threshold dynamics of a viral infection model with defectively infected cells. Mathematical Biosciences and Engineering, 2022, 19(7): 6489-6503. doi: 10.3934/mbe.2022305 |
[4] | Shaoli Wang, Jianhong Wu, Libin Rong . A note on the global properties of an age-structured viral dynamic model with multiple target cell populations. Mathematical Biosciences and Engineering, 2017, 14(3): 805-820. doi: 10.3934/mbe.2017044 |
[5] | Danfeng Pang, Yanni Xiao . The SIS model with diffusion of virus in the environment. Mathematical Biosciences and Engineering, 2019, 16(4): 2852-2874. doi: 10.3934/mbe.2019141 |
[6] | Jinliang Wang, Xiu Dong . Analysis of an HIV infection model incorporating latency age and infection age. Mathematical Biosciences and Engineering, 2018, 15(3): 569-594. doi: 10.3934/mbe.2018026 |
[7] | Chunyang Qin, Yuming Chen, Xia Wang . Global dynamics of a delayed diffusive virus infection model with cell-mediated immunity and cell-to-cell transmission. Mathematical Biosciences and Engineering, 2020, 17(5): 4678-4705. doi: 10.3934/mbe.2020257 |
[8] | Shengqiang Liu, Lin Wang . Global stability of an HIV-1 model with distributed intracellular delays and a combination therapy. Mathematical Biosciences and Engineering, 2010, 7(3): 675-685. doi: 10.3934/mbe.2010.7.675 |
[9] | Cameron J. Browne, Chang-Yuan Cheng . Age-structured viral dynamics in a host with multiple compartments. Mathematical Biosciences and Engineering, 2020, 17(1): 538-574. doi: 10.3934/mbe.2020029 |
[10] | Pengfei Liu, Yantao Luo, Zhidong Teng . Role of media coverage in a SVEIR-I epidemic model with nonlinear incidence and spatial heterogeneous environment. Mathematical Biosciences and Engineering, 2023, 20(9): 15641-15671. doi: 10.3934/mbe.2023698 |
The population dynamics of in-host viral infection models has been studied intensively in the literature [1,2,3,4,5,6,7,8]. Through rigorous mathematical analysis, numerical explonation, and data fitting, the greatly enhanced understanding of viral dynamics can provide us with guidance and support for proposing feasible and effective control strategies to clear viral infections [4,5,9,10]. Much of the existing mathematical modelling has been focused on the cell-free infection modes only [4,5,10]. In cell-free infection, only newly released free virions could infect susceptible target cells. On the other hand, most of the existing works are grounded on ordinary or functional differential equations with constant parameters, and do not consider the spatial heterogeneity, which may induce deficient understanding of the spatial spread of viral infection. So far as we know, only very few works; see for example, [11,12], have taken into account spatial heterogeneity in viral infection modelling.
Assume that cells and virus particles live in a spatially heterogeneous but continuous environment. Let Ω be the spatial habitat with smooth boundary ∂Ω. Denote by u1(x,t), u2(x,t) and u3(x,t) the populations of susceptible target cells, infected target cells and virus particles at location x and time t, respectively. Wu et al. [12] considered a diffusive viral infection model with heterogeneous parameters and distinct dispersal rates for the susceptible and infected target cells:
∂u1∂t=d1Δu1+a(x)−β1(x)u1u3−μ1(x)u1, x∈Ω,t>0,∂u2∂t=d2Δu2+β1(x)u1u3−μ2(x)u2, x∈Ω,t>0,∂u3∂t=k(x)u2−μ3(x)u3, x∈Ω,t>0, | (1.1) |
with nonnegative initial conditions and the homogeneous Neumann boundary condition. Here, d1,d2>0 are the diffusion coefficients of susceptible target cells and infected target cells, respectively; Δ is the Laplacian operator; cell-free infection mode is modelled by the mass action mechanism with β1(x) being the cell-free transmission rate; a(x) is the recruitment rate of susceptible target cells; μ1(x),μ2(x) and μ3(x) are the death rates of susceptible target cells, infected target cells and virus particles, respectively; k(x) is the rate of virus production due to the lysis of infected cells. All these parameters are positive and continuous functions on ˉΩ. In [12], the authors showed that model (1.1) possesses a global attractor, and identified the basic reproduction number R0 and proved its threshold role.
Note that in [11,12], only the cell-free infection mode was considered for the viral infection. It has been recognized that there is another major viral infection mode, namely, the cell-to-cell infection mode [13,14], which allows viral particles to be transferred directly from an infected source cell to a susceptible target cell through the formation of virological synapses [15]. It has been revealed that more than half of viral infections are due to cell-to-cell transmission [15], and even during an antiretroviral therapy, viral particles can be transferred from infected target cells to uninfected ones through virological synapses, and the direct cell-to-cell infection affects the mechanism of HIV-1 transmission in vivo.
Motivated by the previous works, we consider the following general viral infection model incorporating spatial heterogeneity and two infection modes:
∂u1∂t=∇⋅(d1(x)∇u1)+a(x)−f(u1,u2)−g(u1,u3)−μ1(x)u1,∂u2∂t=∇⋅(d2(x)∇u2)+f(u1,u2)+g(u1,u3)−μ2(x)u2,∂u3∂t=k(x)u2−μ3(x)u3, | (1.2) |
for x∈Ω,t>0, with nonnegative initial conditions
ui(x,0)=ϕi(x)≥0 for x∈Ω, i=1,2,3, |
where ∇⋅(di(x)∇ui) describes the divergence of di(x)∇ui and di(x) is the diffusion rate; f(u1,u2) is the cell-to-cell transmission function; and g(u1,u3) is the cell-free transmission function. Here, we consider an isolated habitat Ω, revealed by the Neumann boundary condition
∇ui⋅ν=0,i=1,2,x∈∂Ω,t>0. | (1.3) |
Throughout this paper, we assume that the diffusion rates di(x) with i=1,2, the recruitment rate a(x), the cell-free transmission rate β1(x), the cell-to-cell transmission rate β2(x), the virus production rate k(x), and the death rates μi(x) with i=1,2,3 are positive and continuous functions on ˉΩ. We also make the following biologically motivated assumption.
(H1) f,g∈C1(R+×R+) are strictly increasing with respect to both variables, and f(v,w)=0 (resp. g(v,w)=0) if and only if vw=0. Moreover, ∂2f(v,w)/∂w2≤0 and ∂2g(v,w)/∂w2≤0.
In this paper, we will define the basic reproduction number R0 with a clear biological meaning, and further prove that R0 is a threshold parameter for the global dynamics of model (1.2). As we shall see later, the main challenge is caused by the different dispersal rates of the susceptible and infected target cells and partial degeneration of the model system.
The rest of this paper is organized as follows. In Section 2, we show that our model system admits a unique solution, which exists globally and is ultimately uniformly bounded. In Section 3, we identify the biologically meaningful basic reproduction number R0 for the model using the standard procedure of next generation operator, and further explore the properties of R0 when the dispersal rate for infected target cells varies from zero to infinity. Section 4 is devoted to the global dynamics of the model for the cases of R0≤1 and R0>1, respectively. In Section 5, we consider a special case when all coefficients are spatial homogeneous, and give the global asymptotic stability of the unique chronic infection steady state when R0>1.
Denote by X:=C(ˉΩ,R3) the Banach space of continuous functions on ˉΩ with the supremum norm. The nonnegative cone of X is denoted by X+=C(ˉΩ,R3+), then (X,X+) is a strongly ordered space [16]. For any nonnegative initial condition
u(x,0)=ϕ(x):=(ϕ1,ϕ2,ϕ3)∈X+, |
we define T3(t)ϕ3=e−μ3(⋅)tϕ3. For each i=1,2, let Ti(t) be the C0 semigroups generated by the second-order linear differential operator ∇⋅(di∇)−μi with Neumann boundary condition. It then follows from [16,Corollary 7.2.3] that Ti(t) is compact and strongly positive for all t>0 and i=1,2. Moreover, T(t):=(T1(t),T2(t),T3(t)) is a C0 semigroup on X with an infinitesimal generator A0 [17]. Then the system (1.2) can be written as an abstract differential equation
u′(t)=A0u(t)+F(u(t)) |
with nonnegative initial condition u(0)=ϕ∈X+, where the nonlinear operator F=(F1,F2,F3):X+→X is defined by
F1(φ)(x)=a(x)−f(φ1(x),φ2(x))−g(φ1(x),φ3(x)),F2(φ)(x)=f(φ1(x),φ2(x))+g(φ1(x),φ3(x)),F3(φ)(x)=k(x)φ2(x), |
for any φ=(φ1,φ2,φ3)∈X+. On account of (H1), there exists c>0 such that f(φ1(x),φ2(x))+g(φ1(x),φ3(x))≤cφ1(x) for all x∈ˉΩ. It is easily seen that
φ(x)+ϵF(φ)(x)≥(φ1(x)(1−cϵ),φ2(x),φ3(x))T for x∈Ω. |
By choosing ϵ>0 sufficiently small, we have 1>ϵc and φ+ϵF(φ)∈X+. Particularly,
limϵ→0+1ϵdist(φ+ϵF(φ),X+)=0. |
Thus, by using [16,Theorem 7.3.1] or [18,Corollary 4], we establish the existence of the solution to the system (1.2). Note that the nonlinear operator F is mixed quasimontone, then all solutions are nonnegative due to the comparison principle. To summarize, we obtain the following lemma on the existence and nonnegativity of the solution to (1.2).
Lemma 2.1. For every initial condition ϕ∈X+, system (1.2) with Neumann boundary condition (1.3) has a unique solution u(x,t) on a maximal interval of existence [0,tmax). If tmax<∞, then lim supt→tmax‖u(⋅,t)‖X=∞. Moreover, u(x,t)≥0 for all (x,t)∈Ω×[0,tmax).
To prove that tmax=∞, we need to show that the boundedness of solutions for system (1.2). Before stating this result, we need the following lemma.
Lemma 2.2. For any positive and continuous functions d(x), l(x) and μ(x) on ˉΩ, the scalar reaction-diffusion equation
∂w(x,t)∂t=∇⋅(d(x)∇w(x,t))+l(x)−μ(x)w(x,t), x∈Ω, t>0,∇w(x,t)⋅ν=0, x∈∂Ω, t>0 | (2.1) |
admits a unique and strictly positive steady state w∗(x), which is globally asymptotically stable in C(ˉΩ,R+). Moreover, if d(x)≡d,l(x)≡l and μ(x)≡μ for all x∈Ω, then w∗(x)≡l/μ for all x∈Ω.
Proof. In view of the standard theory of parabolic equations [19], we obtain the existence of a compact semiflow Ψt for (2.1) in C(ˉΩ,R+). Denote
ˉl=maxx∈ˉΩl(x),l_=minx∈ˉΩl(x),ˉμ=maxx∈ˉΩμ(x) and μ_=minx∈ˉΩμ(x). |
It then follows from the comparison theorem and maximum principle [19] that Ψt has a global compact attractor K⊂(l_/ˉμ,ˉl/μ_). This implies that K contains a positive steady state w∗(x) due to Theorem 3.1 in [20]. By using strong maximal principle [21] and the monotonicity of l(x)−μ(x)w(x,t) w.r.t w, we can easily obtain that the positive steady state of (2.1) is unique. According to [20,Theorem 3.2], w∗(x) attracts all solutions of (2.1) with nontrivial initial condition ϕ∈C(ˉΩ,R+). This ends the proof.
Theorem 2.3. For every initial condition ϕ∈X+, system (1.2) has a unique global solution u(x,t)≥0 for t≥0. Moreover, there exists a constant M>0 independent of ϕ such that lim supt→∞ui(x,t)≤M for all x∈Ω and i=1,2,3.
Proof. To establish the solutions of (1.2) exist globally on [0,∞), it suffices to show that the boundedness of the solutions. For any initial condition ϕ∈X+, it follows from comparison principle and Lemma 2.2 that u1(x,t)≤w(x,t) for all t∈[0,tmax), where w(x,t) is the solution of (2.1) with l(x)≡a(x), μ(x)≡μ1(x) and initial condition w(x,0)=ϕ1(x). Note that w(x,t)→w∗(x) as t→∞, which implies that
lim supt→∞u1(x,t)≤w∗(x) uniformly for x∈ˉΩ. | (2.2) |
Thus, there exists K1>maxx∈ˉΩw∗(x), depending on ϕ, such that ‖u1(⋅,t)‖≤K1 for all t≥0.
From the last two equations of (1.2) and the definition of Ti(t) with i=2,3, we have
u2(⋅,t)=T2(t)ϕ2(⋅)+∫t0T2(t−s)(f(u1(⋅,s),u2(⋅,s))+g(u1(⋅,s),u3(⋅,s)))ds,u3(⋅,t)=T3(t)ϕ3(⋅)+∫t0T3(t−s)ku2(⋅,s)ds. |
Let −λ2<0 denote the principal eigenvalue of ∇⋅(d2∇)−μ2 with Neumann boundary condition, and λ3=min{minx∈ˉΩμ3(x),λ2/2}>0. We have ‖T2(t)‖≤e−λ2t and ‖T3(t)‖≤e−λ3t. By (H1) and the boundedness of u1(x,t), there exists m1>0 such that
f(u1(⋅,s),u2(⋅,s))+g(u1(⋅,s),u3(⋅,s))≤m1(‖u2(⋅,s)‖+‖u3(⋅,s)‖) |
for all s∈[0,tmax). It then follows that
‖u2(⋅,t)‖≤e−λ2t‖ϕ2‖+m1∫t0e−λ2(t−s)(‖u2(⋅,s)‖+‖u3(⋅,s)‖)ds,‖u3(⋅,t)‖≤e−λ3t‖ϕ3‖+ˉk∫t0e−λ3(t−s)‖u2(⋅,s)‖ds, | (2.3) |
where ˉk=maxx∈ˉΩk(x). Substituting the second inequality into the first one gives
‖u2(⋅,t)‖≤e−λ2t‖ϕ2‖+m1∫t0e−λ2(t−s)‖u2(⋅,s)‖ds+m1∫t0e−λ2(t−s)(e−λ3s‖ϕ3‖+ˉk∫s0e−λ3(s−r)‖u2(⋅,r)‖dr)ds≤‖ϕ2‖+m1∫t0‖u2(⋅,s)‖ds+m1‖ϕ3‖∫t0e−λ3sds+m1ˉke−λ2t∫t0eλ3r‖u2(⋅,r)‖∫tre(λ2−λ3)sdsdr≤C1+C2∫t0‖u2(⋅,s)‖ds, |
where C1=‖ϕ2‖+m1‖ϕ3‖/λ3>0 and C2=m1+m1ˉk/(λ2−λ3)>0. Thus, Gronwall's inequality implies that
‖u2(⋅,t)‖≤C1eC2t for t∈[0,tmax). |
This together with the second inequality in (2.3) yields
‖u3(⋅,t)‖≤‖ϕ3‖+ˉkC1C2eC2t for t∈[0,tmax). |
On account of Lemma 2.1, tmax=∞ and the solution u(x,t) exists for all t≥0.
Next, we will prove that the solution is ultimately bounded with the bound independent of initial conditions. It follows from (2.2) that there exist a constant M11>0, independent of ϕ, and t1>0 such that u1(x,t)≤M11 for t≥t1. This together with (H1) implies that there exists m2>0 such that
f(u1(x,t),u2(x,t))+g(u1(x,t),u3(x,t))≤m2(u2(x,t)+u3(x,t)), x∈Ω,t≥t1. | (2.4) |
Denote ˉa=maxx∈ˉΩa(x), ˜μ_=minx∈ˉΩ{μi(x):i=1,2,3} and |Ω| is the volume of Ω. By integrating (1.2) for u1 and u2 and adding up, we obtain
∂∂t∫Ω(u1+u2)dx≤|Ω|ˉa−˜μ_∫Ω(u1+u2)dx. |
It then follows from comparison principle that lim supt→∞‖u2(⋅,t)‖1≤|Ω|ˉa/˜μ_. Particularly, there exist t2>t1 and M12>0, such that ‖u2(⋅,t)‖1≤M12 for t≥t2. Similarly, we can easily obtain
∂∂t∫Ωu3dx≤ˉkM12−μ3_∫Ωu3dx for t≥t2, |
where μ3_=minx∈ˉΩμ3(x). Thus, there exist t3>t2 and M13>0, such that ‖u3(⋅,t)‖1≤M13 for t≥t3. Consequently, lim supt→∞(‖u2(⋅,t)‖1+‖u3(⋅,t)‖1)≤M1, where M1=M12+M13 is independent of initial conditions.
Assume that t>t3, we now estimate the upper bound of ‖u2(⋅,t)‖2+‖u3(⋅,t)‖2. By multiplying the equation for u2 (resp. u3) of (1.2) by u2 (resp. u3), and integrating on Ω, it then follows from (2.4) that
12∂∂t∫Ωu22dx≤−d2_∫Ω|∇u2|2dx+m2∫Ω(u22+u2u3)dx−μ_∫Ωu22dx,12∂∂t∫Ωu23dx≤ˉk∫Ωu2u3dx−μ_∫Ωu23dx, |
where d2_=minx∈Ωd2(x). Adding the above two inequalities, together with Young's inequality
u2u3≤μ_4(m2+ˉk)u23+m2+ˉkμ_u22, |
we have
12∂∂t∫Ω(u22+u23)dx≤−d2_∫Ω|∇u2|2dx+C22∫Ωu22dx−μ_∫Ωu22dx−34μ_∫Ωu23dx, |
where C22=m2+(m2+ˉk)2μ_. Making use of the Gagliardo-Nirenberg interpolation inequality: there exists c>0 such that ‖w‖22≤ε‖∇w‖22+cε−n/2‖w‖21 for any w∈W1,2(Ω) and small ε>0, we obtain
12∂∂t∫Ω(u22+u23)dx≤B2M21−δ2∫Ω(u22+u23)dx, |
where B2=C2cε−n/2, δ2=3μ_/4 and ε∈(0,d2_/C22). Therefore,
lim supt→∞(‖u2(⋅,t)‖22+‖u3(⋅,t)‖22)≤B2M21/δ2. |
Especially, there exist t4>t3 and M2>0 such that ‖u2(⋅,t)‖22+‖u3(⋅,t)‖22≤M2 for t≥t4.
Denote Lp=lim supt→∞(‖u2(⋅,t)‖pp+‖u3(⋅,t)‖pp). We multiple the equation for u2 (resp. u3) of (1.2) by 2ku2k−12 (resp. 2ku2k−13) and integrate on Ω, using a similar argument as in the estimation of ‖u2(⋅,t)‖22+‖u3(⋅,t)‖22 to obtain that
12k∂∂t(‖u2(⋅,t)‖2k2k+‖u3(⋅,t)‖2k2k)≤2n2(k−1)B‖u2(⋅,t)‖2k+1−22k−1−δ(‖u2(⋅,t)‖2k2k+‖u3(⋅,t)‖2k2k), |
where B and δ are constants independent of of k and ϕ. Since L2k−1≤M2k−1, there exist t2k−1>0 such that ‖u2(⋅,t)‖2k−12k−1+‖u3(⋅,t)‖2k−12k−1≤L2k−1+1 for all t≥t2k−1. By comparison principle, we obtain L2k≤2n2(k−1)C(L2k−1+1)2, where C is a constant independent of k and ϕ.
Finally, according to the method of induction, we prove that L2k≤∞ for all k=0,1,2,⋯. Define an infinite sequence ak+1=(C+1)2−k−12kn2−k−2ak with a0=L1+1 for nonnegative integer k. It is easily seen that L2k≤a2kk and limk→∞lnak=lnC(L1+1)+nln2/2. Therefore, we have
lim supk→∞2k√L2k≤limk→∞ak=C(L1+1)2n/2. |
Thus we obtain lim supt→∞ui(x,t)≤M:=C(M2+1)2n/2+M1 for all x∈Ω and i=1,2,3. That is, the solution semiflow associated with (1.2) Θ(t) for t≥0 is point dissipative. This completes the proof.
We are now in the position to address the persistence of u1(x,t).
Proposition 2.4. Let u(x,t) be the solution of (1.2) with initial condition ϕ∈X+.
(i) u1(x,t)>0 for all t>0 and x∈Ω. Furthermore, there exists a positive constant m0 independent of ϕ such that
lim inft→∞u1(x,t)≥m0 uniformly for x∈ˉΩ. |
(ii) If there exist some x0∈Ω and t0≥0 such that either u2(x0,t0)>0 or u3(x0,t0)>0, then ui(x,t)>0 for all i=2,3, t≥t0 and x∈Ω.
Proof. (ⅰ) By using the strong maximum principle [21], it is easily seen the positivity of u1(x,t) for t>0 and x∈Ω. We then prove the persistence of u1(x,t). From Theorem 2.3, there exist t0>0 and M>0 such that ui(x,t)<M for all t>t0, i=1,2,3 and x∈Ω. Then the first equation of (1.2) and (H1) imply that
∂u1(x,t)∂t≥∇⋅(d1(x)∇u1(x,t))+a(x)−μ1(x)u1(x,t)−c0u1(x,t) |
for all t≥t0 and some positive constant c0. Thus, Lemma 2.2 and comparison principle yield that u1(x,t) is ultimately bounded below by a unique and strictly positive steady state ˉw∗(x) of (2.1) with d(x)=d1(x), l(x)≡a(x) and μ(x)≡μ1(x)+c0. Denote m0=minx∈ˉΩˉw∗(x), which is a positive constant. Then lim supt→∞u1(x,t)≥m0 for all x∈Ω.
(ⅱ) Assume that either u2(x0,t0)>0 or u3(x0,t0)>0 for some x0∈Ω and t0≥0. Then from the third equation of (1.2), we have
u3(x,t)=e−μ3(x)(t−t0)u3(x,t0)+∫tt0e−μ3(x)(t−s)k(x)u2(x,s)ds>0 |
for all x∈Ω and t>t0. We then apply strong maximum principle to the second equation of (1.2) and obtain u2(x,t)>0 for all t>t0 and x∈Ω.
Note that Lemma 2.2 implies that (2.1) with d(x)=d1(x), l(x)≡a(x) and μ(x)≡μ1(x) has a unique and strictly positive steady state w∗(x). Thus, system (1.2) has a unique infection-free steady state (w∗(x),0,0). For simplicity, we denote
βd(x)=∂f(w∗(x),0)∂u2, βi(x)=∂g(w∗(x),0)∂u3. | (3.1) |
Linearizing system (1.2) for (u2(x,t),u3(x,t)) at (w∗(x),0,0) gives the following cooperative system for the infected cells and free virus,
∂u2∂t=∇⋅(d2(x)∇u2)+βd(x)u2+βi(x)u3−μ2(x)u2, x∈Ω, t>0,∂u3∂t=k(x)u2−μ3(x)u3, x∈Ω, t>0,∇u2⋅ν=0, x∈∂Ω, t>0. | (3.2) |
The suitable functional space for the above system is Y:=C(ˉΩ,R2). The associated linear operator of system (3.2) can be decomposed as A=F+B, where
F=(βd(⋅)βi(⋅)00),B=(∇⋅(d2∇)−μ2(⋅)0k(⋅)−μ3(⋅)). |
It then follows from [22] that the basic reproduction number R0 is defined as the spectral radius of −FB−1, that is, R0=r(−FB−1). We can easily check that B is resolvent-positive with s(B)<0, F is positive and A is also resolvent-positive. Then it follows from [22,Theorem 3.5] that R0−1 has the same sign as s(A), where s(A)=sup{Reλ, λ∈σ(A)} is the spectral bound of A.
Let eBt be the semigroup generated by B. Then the next generation operator is −FB−1=∫∞0FeBtdt. Wang and Zhao [23] proved local asymptotic stability of infection-free steady state when R0<1. Here, we shall prove global asymptotic stability of infection-free steady state when R0≤1. To derive an equivalent formula for R0 such that the direct and indirect transmission mechanisms are clearly separated in the expression, we need to make use of the following result.
Lemma 3.1. Let F=(F11F1200) be a positive operator, B=(∇⋅(d2∇)−V110−V21−V22) be a resolvent-positive operator with s(B)<0. Then we have
r(−FB−1)=r(F11(V11−∇⋅(d2∇))−1−F12V−122V21(V11−∇⋅(d2∇))−1). | (3.3) |
Proof. Note that B is lower triangular and s(B)<0. This implies that both V11−∇⋅(d2∇) and V22 are invertible. Moreover, we can calculate that
−B−1=((V11−∇⋅(d2∇))−10−V−122V21(V11−∇⋅(d2∇))−1V−122). |
Consequently, we obtain
−FB−1=(F11(V11−∇⋅(d2∇))−1−F12V−122V21(V11−∇⋅(d2∇))−1F12V−12200), |
which implies that (3.3) holds. This ends the proof.
By using Lemma 3.1 and a standard variational method, we have an equivalent formula for the basic reproduction number
R0=r(Ad+Ai)=supψ∈H1(Ω),ψ≠0∫Ω(βd(x)+βi(x)μ−13(x)k(x))ψ2(x)dx∫Ω(d2(x)|∇ψ(x)|2+μ2(x)ψ2(x))dx, | (3.4) |
where Ad=βd(μ2−∇⋅(d2∇))−1 is the next generation operator for cell-to-cell transmission, and Ai=βiμ−13k(μ2−∇⋅(d2∇))−1 is the next generation operator for cell-free transmission. In the absence of cell-free transmission, the basic reproduction number for cell-to-cell transmission is
Rd0=r(Ad)=supψ∈H1(Ω),ψ≠0∫Ωβd(x)ψ2(x)dx∫Ω(d2(x)|∇ψ(x)|2+μ2(x)ψ2(x))dx. |
On the other hand, if only cell-free transmission is taken into consideration, the basic reproduction number for cell-free transmission is given by
Ri0=r(Ai)=supψ∈H1(Ω),ψ≠0∫Ωβi(x)μ−13(x)k(x)ψ2(x)dx∫Ω(d2(x)|∇ψ(x)|2+μ2(x)ψ2(x))dx. |
Clearly, R0≤Rd0+Ri0. We then study the dependence of R0 on the diffusion coefficient d2.
Theorem 3.2. (i) R0 is a principal eigenvalue of Ad+Ai associated with a positive eigenfunction.
(ii) Assume that d2 is a constant on ˉΩ, then R0 is a monotone decreasing function of d2, Moreover, we have
R0→¯R0:=maxx∈ˉΩ{βd(x)μ2(x)+βi(x)k(x)μ2(x)μ3(x)} as d2→0,R0→R_0:=∫Ω(βd(x)+βi(x)μ−13(x)k(x))dx∫Ωμ2(x)dx as d2→∞. |
Proof. (ⅰ) Since Ad and Ai are compact and positive, it then follows from Krein-Rutman theorem that R0 is a principal eigenvalue of Ad+Ai with a positive eigenfunction, denoted by ϕ∗(x); namely,
βd(μ2−∇⋅(d2∇))−1ϕ∗+βi(x)k(x)μ3(x)(μ2−∇⋅(d2∇))−1ϕ∗=R0ϕ∗, x∈Ω,∇ϕ∗(x)⋅ν=0, x∈∂Ω. |
Denote ψ∗=ϕ∗/(βd+βiμ−13k), then the above eigenvalue problem can be rewritten as
∇⋅(d2(x)∇ψ(x))−μ2(x)ψ(x)+βd(x)+βi(x)μ−13(x)k(x)R0ψ(x)=0, x∈Ω,∇ψ⋅ν=0, x∈∂Ω, | (3.5) |
(ⅱ) Assume that d2 is a constant on ˉΩ. It is easily seen that R0 is a decreasing function of d2, and the eigenvalue problem (3.5) can be reduced as
d2Δψ∗−μ2ψ∗+βd+βiμ−13kR0ψ∗=0, x∈Ω,∇ψ∗(x)⋅ν=0, x∈∂Ω. | (3.6) |
We first claim that R0≤¯R0, otherwise we have −μ2+(βd+βiμ−13k)/R0<0 and the principal eigenvalue of d2Δ−μ2+(βd+βiμ−13k)/R0 is negative. This contradicts to (3.6). Thus limd2→0R0 exists. We next prove that limd2→0R0=¯R0. Assume to the contrary, then there exists ϵ0>0 such that R0<¯R0−ϵ0 for all positive d2. It follows from the continuity of coefficient functions that there exists a x0∈Ω and a δ>0 such that
βd(x)μ2(x)+βi(x)k(x)μ2(x)μ3(x)>¯R0−ϵ02>R0+ϵ02 for all x∈Bδ(x0), |
which implies that the positivity of βd(x)/μ2(x)+βi(x)k(x)/[μ2(x)μ3(x)] on Bδ(x0). Due to compactness of continuous functions on a bounded domain, there exists ϵ>0 such that
−μ2(x)+βd(x)+βi(x)μ−13(x)k(x)R0>ϵ for all x∈Bδ(x0). |
The above inequality together with (3.6) yields −Δψ∗>ϵψ∗/d2. Denote ψ+(x)=ψ∗(x)/minx∈Bδ(x0)ψ∗(x). Then we have −Δψ+(x)>ϵψ+(x)/d2 and ψ+(x)≥1 on Bδ(x0). Let η>0 be the principal eigenvalue of −Δ on Bδ(x0) under Neumann boundary condition and ψ−(x) the corresponding eigenfunction, we can further normalize ψ−(x) such that ψ−(x)≤1 on Bδ(x0). Then we have −Δψ−(x)=ηψ−(x)<ϵψ−(x)/d2. Thus, ψ+(x) and ψ−(x) are the super- and sub-solutions of −Δφ=ϵφ/d2 with Neumann boundary condition. Thus, ϵ/d2 is an eigenvalue of −Δ on Bδ(x0) with Neumann boundary condition, which contradicts the facts ϵ/d2>η and η is the principal eigenvalue of −Δ. Therefore, R0→¯R0 as d2→0.
It is easily seen from (3.4) that R0≥R_0 for all d2>0. Thus, R0 is uniformly bounded for d2>0 and limd2→∞R0 exists. Then we divide both sides of (3.6) by d2 to obtain
Δψ∗+βd+βiμ−13k−R0μ2R0d2ψ∗=0, x∈Ω. |
It then follows from elliptic regularity [24] that, there exists a positive constant ˉψ such that ψ∗→ˉψ in C(Ω) as d2→∞. Integrating (3.6) by parts over Ω yields
∫Ωμ2ψ∗dx=∫Ωβd+βiμ−13kR0ψ∗dx. |
Letting d2→∞, we obtain R0→R_0. This completes the proof.
From the above theorem, we have a direct application on basic reproduction number.
Proposition 3.3. (i) If βd(x)/μ2(x)+βi(x)k(x)/(μ2(x)μ3(x))≤1 for all x∈Ω, then R0<1 for all d2>0 and Ω is an infection-free environment.
(ii) If ∫Ω(βd(x)+βi(x)μ−13(x)k(x)dx≥∫Ωμ2(x)dx, , then R0>1 for all d2>0 and Ω is a favorable environment for the viral infection.
(iii) If ∫Ω(βd(x)+βi(x)μ−13(x)k(x)dx<∫Ωμ2(x)dx and βd(x)/μ2(x)+βi(x)k(x)/(μ2(x)μ3(x))>1 for some x∈Ω, then there exists a d∗2>0 such that R0≤1 if d2≥d∗2 and R0>1 if d2<d∗2.
Define the continuous semiflow {Θ(t)}t≥0:X+→X+ for the system (1.2) by
Θ(t)ϕ(⋅):=u(⋅,t,ϕ), t≥0. |
It follows from Theorem 2.3 that the semiflow Θ(t) of system (1.2) is point dissipative and the orbit γ+(U)=⋃ϕ∈Uγ+(ϕ) is bounded for any bounded set U⊂X+. To apply the theory in [25], we have to show that Θ(t) is asymptotically smooth. Since Θ(t) is not compact, we introduce the weak compactness condition called κ-contraction, and Kuratowski measure of the noncompactness defined by [25]
κ(U):=inf{r≥0:U has a finite cover of diameter less than r} | (4.1) |
for any bounded set U⊂X+. Clearly, κ(U)=0 if and only if U is precompact. We need to show that Θ(t) is a κ-contraction, that is, there exists a continuous function q(t)∈[0,1):R+→R+ such that κ(Θ(t)U)≤q(t)κ(U) for any bounded set U⊂X+ and t>0. To achieve this, we need the following lemma. The proof is similar to that in [12,Lemma 2.5] with some minor modifications.
Lemma 4.1. For any bounded set U⊂X+ and t>0, {ui(⋅,t,ϕ):ϕ∈U} and {∫t0e−μ3(⋅)(t−s)k(⋅)u2(⋅,s,ϕ)ds:ϕ∈U} with i=1,2 are precompact in C(ˉΩ).
Theorem 4.2. The semiflow Θ(t) is a κ-contraction and asymptotically smooth. Moreover, system (1.2) admits a connected global attractor in X+.
Proof. For any initial condition ϕ=(ϕ1,ϕ2,ϕ3)∈X+, we have Θ(t)ϕ=Θ1(t)ϕ+Θ2(t)ϕ for all t≥0, where
Θ1(t)ϕ=(u1(⋅,t,ϕ),u2(⋅,t,ϕ),∫t0e−μ3(⋅)(t−s)k(⋅)u2(⋅,s,ϕ)ds),Θ2(t)ϕ=(0,0,e−μ3(⋅)tϕ3) |
For any bounded set U⊂X+, it follows from (4.1) that
κ(Θ2(t)U)≤‖e−μ3(⋅)t‖κ(U)≤e−μ_3tκ(U) for all t≥0, |
where μ_3=minx∈ˉΩμ3(x). Note that Lemma 4.1 implies that Θ1(t)U is precompact in C(ˉΩ) for any t>0, that is, κ(Θ1(t)U)=0. Hence, for any t>0, we have
κ(Θ(t)U)≤κ(Θ1(t)U)+κ(Θ2(t)U)≤e−μ_3tκ(U). |
Therefore, Θ(t) is a κ-contraction. It then follows from [25,Lemma 2.3.4] that Θ(t) is asymptotically smooth. Therefore, by Theorem 2.4.6 in [25], system (1.2) admits a connected global attractor in X+.
It follows from Theorem 3.1 in [23] that the infection-free steady state (w∗(x),0,0) is locally asymptotically stable when R0<1. To establish global asymptotic stability of infection-free steady state when R0≤1, we shall first develop the following approach to show local asymptotic stability of infection-free steady state not only when R0<1, but also for the critical case R0=1.
Denote A as the linear operator of (3.2) and eAt the semigroup generated by A. The exponential growth bound of eAt is defined as
ω(eAt):=limt→∞ln‖eAt‖t. |
Lemma 4.3. Assume that R0≤1, then s(A)≤0, ω(eAt)≤0, and there exists a constant Ma>0 such that ‖eAt‖≤Ma.
Proof. By Theorem 3.5 in [22], s(A)≤0 if R0≤1. It follows from [26] that
ω(eAt)=max{s(A),ωess(eAt)}, |
where ωess(eAt) is the essential growth bound of eAt defined by
ωess(eAt):=limt→∞lnσ(eAt)t. |
Here, σ(eAt) denotes the distance of eAt from the set of compact linear operators in Y=C(ˉΩ,R2). To prove that ω(eAt)≤0, it is sufficient to show that ωess(eAt)≤0. For any ˆϕ:=(ϕ2,ϕ3)∈Y, the solution of the linear system (3.2) is eAtˆϕ=Ψ2(t)ˆϕ+Ψ3(t)ˆϕ, where Ψ2(t)ˆϕ=(u2(⋅,t,ˆϕ),∫t0e−μ3(⋅)(t−s)k(⋅)u2(⋅,s,ˆϕ)ds) and Ψ3(t)ˆϕ=(0,e−μ3(⋅)tϕ3). Note that Lemma 4.1 implies that Ψ2(t) is a compact linear operator, that is, σ(Ψ2(t))=0. Thus, we have σ(eAt)=σ(Ψ2(t)+Ψ3(t))=σ(Ψ3(t))≤‖Ψ3(t)‖≤e−μ_3t. Therefore, we compute
ωess(eAt)≤−μ_3<0. |
This implies that the there exists a constant Ma>0 such that ‖eAt‖≤Ma.
Theorem 4.4. Assume that R0≤1, then the infection-free steady state (w∗(x),0,0) of(1.2) is locally asymptotically stable.
Proof. Given any small δ>0, let u(x,t) be any solution of (1.2) with initial condition satisfies ‖u1(x,0)−w∗(x)‖+‖u2(x,0)‖+‖u3(x,0)‖<δ. Denote w1(x,t)=u1(x,t)−w∗(x) and μ_1=minx∈Ωμ1(x)>0 which satisfies
∂w1∂t=∇⋅(d1∇w1)−μ1w1−f(w1+w∗,u2)−g(w1+w∗,u3)≤∇⋅(d1∇w1)−μ_1w1. |
Let −˜λ1<0 be the principle eigenvalue of ˜T1(t), where ˜T1(t) is the C0 semigroup generated by ∇⋅(d1∇)−μ_1 with Neumann boundary condition. It then follows from comparison principle that
w1(x,t)≤‖˜T1(t)w1(x,0)‖≤e−˜λ1t‖u1(x,0)−w∗(x)‖≤δe−˜λ1t |
for all x∈Ω and t≥0. Thus, we have u1(x,t)≤˜u1(x,t):=w∗(x)+δe−˜λ1t on Ω×[0,∞). By (H1), we obtain
f(u1,u2)≤f(˜u1,u2)≤∂f(˜u1,0)∂u2u2 and g(u1,u3)≤g(˜u1,u3)≤∂g(˜u1,0)∂u3u3. |
We obtain from the definitions of βd and βi in (3.1) and the second equation of (1.2) that
∂u2∂t≤∇⋅(d2(x)∇u2)+βd(x)u2+βi(x)u3−μ2(x)u2+p(x,t) |
for x∈Ω and t>0, where
p(x,t)=(∂f(˜u1,0)∂u2−βd)u2+(∂g(˜u1,0)∂u3−βd)u3. |
It follows from system (1.2) and comparison principle that
(u2(⋅,t)u3(⋅,t))≤eAt(u2(⋅,0)u3(⋅,0))+∫t0eA(t−s)(p(⋅,s)0)ds, | (4.2) |
Recall K1>maxx∈ˉΩw∗(x) and ‖u1(x,t)‖≤K1 for all t≥0 in the proof of Theorem 2.3. Denote
ˉf=maxu1∈[0,K1]|∂2f(u1,0)∂u1∂u2|, ˉg=maxu1∈[0,K1]|∂2g(u1,0)∂u1∂u3|. |
We then have p(x,t)≤δe−˜λ1t(ˉfu2+ˉgu3). Set E(t)=max{maxx∈Ωu2(x,t),maxx∈Ωu3(x,t)}. By Lemma 4.3 and inequality (4.2), we obtain
E(t)≤δMa+δMa(ˉf+ˉg)∫t0e−˜λ1sE(s)ds. |
Then Gronwall's inequality yields
E(t)≤δMae∫t0δMa(ˉf+ˉg)e−˜λ1sds≤δMaeδMa(ˉf+ˉg)˜λ1 for all t≥0. |
Thus ‖u2(⋅,t)‖+‖u3(⋅,t)‖=O(δ) as δ→0. We next show that ‖u1(⋅,t)−w∗(x)‖=O(δ) as δ→0. Note that (H1) implies that
f(u1,u2)≤f(K1,u2)≤∂f(K1,0)∂u2u2≤∂f(K1,0)∂u2δMaeδMa(ˉf+ˉg)˜λ1,g(u1,u3)≤g(K1,u3)≤∂g(K1,0)∂u3u3≤∂g(K1,0)∂u3δMaeδMa(ˉf+ˉg)˜λ1. |
It then follows from the above inequalities and the first equation of (1.2) that
∂u1∂t≥∇⋅(d1(x)∇u1)+a(x)−qδ−μ1(x)u1, | (4.3) |
where q=(∂f(K1,0)/∂u2+∂g(K1,0)/∂u3)MaeδMa(ˉf+ˉg)/˜λ1 is positive and finite. By Lemma 2.1, for any small δ>0, the following reaction-diffusion equation
∂w∂t=∇⋅(d1(x)∇w)+a(x)−qδ−μ1(x)w, x∈Ω, t>0,∇w(x,t)⋅ν=0, x∈∂Ω, t>0 |
admits a unique and strictly positive steady state wδ(x), which is globally asymptotically stable in C(ˉΩ,R+). Moreover, ‖w∗(x)−wδ(x)‖=O(δ) as δ→0. Thus, it follows from (4.3) and comparison principle that
u1(x,t)≥wδ(x)≥w∗(x)+O(δ) as δ→0. |
Recall that u1(x,t)≤w∗(x)+δe−˜λ1t on Ω×[0,∞). Therefore, ‖u1(⋅,t)−w∗(⋅)‖+‖u2(⋅,t)‖+‖u2(⋅,t)‖=O(δ) as δ→0, thus proving local stability of (w∗(x),0,0) if R0≤1.
We are now in the position to establish global attractivity of infection-free steady state by constructing a suitable Lyapunov functional and LaSalle invariance principle.
Theorem 4.5. If R0≤1, then the infection-free steady state (w∗(x),0,0) of (1.2) is globally asymptotically stable.
Proof. We establish the global asymptotic stability of (w∗(x),0,0) by proving the following two claims. Define a region D={ϕ∈X+:ϕ(x)≤w∗(x)}.
Claim 1. For any initial data ϕ∈X+, the omega limit set of ϕ is contained in D.
Clearly, for any x∈Ω, if u1(x,t0)≤w∗(x) for some t0≥0, then u1(x,t)≤w∗(x) for all t≥t0. Then we divide the domain Ω into two sub-domains Ω1:={x∈Ω:u1(x,t)>w∗(x)for allt≥0} and Ω2:={x∈Ω:u1(x,t)≤w∗(x)for somet≥0}. Here, Ω2 is closed in Ω, and there exists t0≥0 that u1(x,t)≤w∗(x) for all x∈Ω2. Without loss of generality, we assume t0=0.
For any x∈Ω1, Lemma 2.1 and the first equation of (1.2) imply that ∂u1(x,t)/∂t≤0, that is, u1(x,t) is a decreasing function in t. It then follows from u1(x,t)≥w∗(x) for x∈Ω1 that limt→∞u1(x,t) exists, and limt→∞u1(x,t)≥w∗(x). Moreover, if limt→∞u1(x,t)>w∗(x), then we obtain from the first equation of (1.2) that 0=limt→∞∂u1(x,t)/∂t<0. This is a contradiction. Therefore, limt→∞u1(x,t)=w∗(x), which implies that the omega limit set of ϕ is contained in D.
Claim 2. The infection-free steady state (w∗(x),0,0) attracts all initial profiles in D.
We consider the solution semiflow restricted on the invariant set D and construct a Lyapunov functional V1:D→R given by
V1(u1,u2,u3)=∫D(u22(x,t)+βi(x)k(x)u23(x,t))dx. |
Taking the derivative of V1 along the solution, we obtain
dV1dt=∫D(u2[∇⋅(d2∇u2)+f(u1,u2)+g(u1,u3)−μ2u2]+βiku3(ku2−μ3u3))dx. |
Note that u1(x,t)≤w∗(x) in D, it is readily seen from (H1) that f(u1,u2)≤f(w∗(x),u2)≤βdu2 and g(u1,u3)≤g(w∗(x),u3)≤βiu3. These inequalities and Neumann boundary condition yield that
dV1dt≤∫D(−d2|∇u2|2−(μ2−βd)u22+2βiu2u3−βiμ3ku23)dx≤∫D(−d2|∇u2|2−(μ2−βd)u22+βikμ3u22−βi(√μ3ku3−√kμ3u2)2)dx≤∫D(−d2|∇u2|2−(μ2−βd)u22+βikμ3u22)dx. |
We next prove
∫Ωβikμ3ψ2dx≤∫Ω(d2|∇ψ|2+(μ2−βd)ψ2)dx. | (4.4) |
holds for any ψ∈H1(Ω) if R0≤1. We make another decomposition of the linear operator A=F1+B1 associated with the linear system (3.2), where
F1=(0βi(⋅)00),B1=(∇⋅(d2∇)−(μ2(⋅)−βd(⋅))0k(⋅)−μ3(⋅)). | (4.5) |
Note that Theorem 3.2 implies that μ2>βd when R0≤1. Thus the operator B1 is resolvent-positive with s(B1)<0. Then it follows from [22,Theorem 3.5] and R0≤1 that s(A)≤0 and
r(−F1B−11)=r(βikμ3(μ2−βd−∇⋅(d2∇))−1)=supψ∈H1(Ω),ψ≠0∫Ωβi(x)μ−13(x)k(x)ψ2(x)dx∫Ω(d2(x)|∇ψ(x)|2+(μ2(x)−βd(x))ψ2(x))dx≤1. |
Hence, we obtain (4.4) for any ψ∈H1(Ω) if R0≤1. This implies that dV1/dt≤0 if R0≤1. Moreover, K={(ˉw∗(x),0,0)}, where K is an invariant set on which dV1/dt=0. Note that (ˉw∗(x),0,0) is the unique point in the largest invariant set on which dV1/dt=0. By the LaSalle invariance principal, (ˉw∗(x),0,0) is globally attractive in D.
Finally, it follows from Lemma 1.2.1 in [27] that the omega limit set of any initial data ϕ∈X+ is internally chain transitive. The above two claims and [27,Theorem 1.2.1] yield (w∗(x),0,0) is globally attractive in X+. This, together with the local stability result in Theorem 4.4, implies the global asymptotic stability of (w∗(x),0,0) in X+ when R0≤1. This ends the proof.
By using the same idea in [12,Lemma 3.7], we show that s(A) is actually the principal eigenvalue of A when R0≥1.
Lemma 4.6. If R0≥1, then s(A) is the principal eigenvalue of A with a strongly positive eigenfunction.
Proof. The eigenvalue problem of A is given by
λφ2(x)=∇⋅(d2(x)∇φ2(x))+βd(x)φ2(x)+βi(x)φ3(x)−μ2(x)φ2(x), x∈Ω,λφ3(x)=k(x)φ2(x)−μ3(x)φ3(x), x∈Ω,∇φ2(x)⋅ν=0, x∈∂Ω. | (4.6) |
Then we define a one-parameter family of linear operators with Neumann boundary condition on C(ˉΩ):
Lλ=∇⋅(d2∇)+βd+βikλ+μ3−μ2. |
Let Tλ(t) be the semigroup generated by Lλ. Since βd+βikλ+μ3−μ2 is cooperative and irreducible for all x∈Ω, it then follows from Theorem 7.5.1 in [16] that Tλ(t) is a compact and strongly positive operator for all t>0. By Krein-Rutman theorem, s(Lλ) is a principal eigenvalue of Lλ with positive eigenfunction φ∗2(x). Clearly, s(Lλ) is decreasing and continuously with respect to λ, and s(Lλ) is finite when λ is large.
According to Lemma 2.3(d) in [28], we obtain that R0−1 and s(A) have the same sign as λ0, where λ0 is the principal eigenvalue of L0. This yields that s(L0)=λ0≥0. Thus, there exists a unique λ∗≥0 such that s(Lλ∗)=λ∗. Note that the problem (4.6) can be written as Lλφ2(x)=λφ2(x). Therefore, if R0≥1, then s(Lλ∗)=λ∗>0 is a principal eigenvalue of A with a strongly positive eigenfunction (φ∗2(x),φ∗3(x)), where φ∗3(x)=k(x)s(Lλ)+μ3(x)φ∗2(x). Finally, we can further obtain λ∗=s(A) by using Theorem 2.3 in [23].
To establish the existence of the chronic infection steady state when R0>1, we now apply the permanence theorem in [29,Theorem 3] and use an argument similar to that in the proof of [11,Theorem 2.2] to obtain the following persistence result.
Theorem 4.7. If R0>1, then system (1.2) is uniformly persistent in X+, that is, there exists a η>0 such that for any ϕ∈X0, we have
lim inft→∞ui(x,t,ϕ)≥η, (i=1,2,3) uniformly for all x∈ˉΩ. |
Moreover, system (1.2) admits at least one chronic infection steady state (u∗1(x),u∗2(x),u∗3(x)).
Proof. Denote X0:={(ϕ1,ϕ2,ϕ3)∈X+: ϕ2(⋅)≢ and
\partial X_0: = X^+\backslash X_0 = \{(\phi_1,\phi_2,\phi_3)\in X^+: \ \phi_{2}(\cdot)\equiv0 \ \text{or} \ \phi_{3}(\cdot)\equiv0 \}. |
Obviously, X_0\cap\partial X_0 = \emptyset , X_0\cup\partial X_0 = X^+ , X_0 is open and dense in X^+ , and \Theta(t)\partial X_0\subseteq \partial X_0 . Note that Proposition 2.4(ii) implies that \Theta(t)X_0\subseteq X_0 for all t\ge0 . Denote M_\partial as the largest positively invariant set in \partial X_0 . It follows from Proposition 2.4(ii) that
M_\partial = \{(\phi_1,\phi_2,\phi_3)\in X^+: \ \phi_{2}\equiv0 \ \text{and} \ \phi_{3}\equiv0 \}. |
For any initial data \phi\in M_\partial , we can easily obtain that u_i(x, t, \phi)\equiv0 for all i = 2, 3 , x\in\Omega and t\ge0 . Then in view of Lemma 2.1, the limiting system when u_i\equiv0 for i = 2, 3 has a unique globally asymptotically stable steady state u_1(x, t) = w^*(x) . We then obtain from [30,Theorem 4.1] that (w^*(x), 0, 0) is globally attractive in M_\partial . We now define a generalized distance function \rho: X^+\rightarrow [0, \infty) by
\rho(\phi) = \min\limits_{x\in\bar\Omega}\{\phi_2(x), \phi_3(x)\} \ \ \text{for any} \ \ \phi\in X^+. |
From strong maximum principle, we have \rho(\Theta(t)\phi) > 0 for all \phi\in X_0 . Since \rho^{-1}(0, \infty)\subset X_0 , the condition (P) in [29,Section 3] is satisfied.
Denote W^{s}((w^*(x), 0, 0)) as the stable manifold of (w^*(x), 0, 0) . We next verify that W^{s}((w^*(x), 0, 0))\cap \rho^{-1}(0, \infty) = \emptyset . It suffices to show that there exists a \eta_0 > 0 such that
\limsup\limits_{t\to\infty}\|\Theta_t\phi-(w^*(x),0,0)\|\geq\eta_0 \ \ \text{for any} \ \ \phi\in\rho^{-1}(0,\infty). |
Suppose, to the contrary, for any \eta_0 > 0 there exists \tilde\phi\in\rho^{-1}(0, \infty) such that
\begin{equation} \limsup\limits_{t\to\infty}\|\Theta_t\tilde\phi-(w^*(x),0,0)\| \lt \eta_0. \end{equation} | (4.7) |
In view of Lemma 4.6, \lambda_0 = s(A) > 0 is the principal eigenvalue of A = F+B with a strongly positive eigenfunction if R_0 > 1 . For any sufficiently small \varepsilon > 0 , we consider a small perturbation of F :
F_\varepsilon = \begin{pmatrix} \beta_d-\varepsilon & \beta_i-\varepsilon\\ 0&0 \end{pmatrix}. |
Similar as in the proof of Lemma 4.6, one can show that the eigenvalue problem
\begin{aligned} &\lambda\varphi_2 = \nabla\cdot(d_2\nabla\varphi_2)+(\beta_d-\varepsilon)\varphi_2+(\beta_i-\varepsilon)\varphi_3-\mu_2\varphi_2, \ x\in\Omega,\\ &\lambda\varphi_3 = k\varphi_2-\mu_3\varphi_3, \ x\in\Omega, \\ &\nabla\varphi_2\cdot\nu = 0, \ x\in\partial\Omega. \end{aligned} |
has a principle eigenvalue \lambda_{\varepsilon} with strongly positive eigenfunction (\varphi_2^{\varepsilon}, \varphi_3^{\varepsilon}) . By continuity of the operator, we have \lambda_{\varepsilon}\to\lambda_0 > 0 as \varepsilon\to 0^+ . We then choose a small \varepsilon > 0 such that \lambda_{\varepsilon} > 0 . It follows from (4.7) and (\textbf{H}_1) that there exists a \tilde t > 0 such that
f(u_1,u_2)\ge (\beta_d-\varepsilon) u_2 \ \ \text{and} \ \ g(u_1, u_3)\ge (\beta_i-\varepsilon) u_3 \ \ \text{for all} \ \ t\geq \tilde t. |
Thus, for all t\geq \tilde t , (u_2(x, t, \tilde\phi), u_3(x, t, \tilde\phi)) satisfies
\begin{equation} \begin{array}{lll} &\frac{\partial u_2}{\partial t}\ge\nabla\cdot(d_2(x)\nabla u_2)+(\beta_d-\varepsilon) u_2+(\beta_i-\varepsilon) u_3-\mu_2u_2, \ \ & x\in\Omega, \ t \gt \tilde t,\\ &\frac{\partial u_3}{\partial t} = k(x)u_2-\mu_3(x)u_3, \ \ & x\in\Omega, \ t \gt \tilde t,\\ &\nabla u_2\cdot\nu = 0, \ \ & x\in\partial\Omega, \ t \gt \tilde t. \end{array} \end{equation} | (4.8) |
Since u_i(x, t, \tilde\phi) > 0 for all x\in\bar\Omega , t > 0 and i = 2, 3 , there exists \delta > 0 such that u_2(x, \tilde t, \tilde\phi)\ge\delta\varphi_2^{\varepsilon} and u_3(x, \tilde t, \tilde\phi)\ge\delta\varphi_3^{\varepsilon} . It then follows from (4.8) and comparison principle that
(u_2(x,t,\tilde\phi),u_3(x,t,\tilde\phi))\geq (\delta e^{\lambda_{\varepsilon}(t-\tilde t)}\varphi_2^{\varepsilon}, \ \delta e^{\lambda_{\varepsilon}(t-\tilde t)}\varphi_3^{\varepsilon}) \ \ \text{for} \ \ x\in \bar\Omega, \ t\geq\tilde t. |
Therefore, u_i(x, t, \tilde\phi)\to \infty as t\to \infty for i = 2, 3 , which contradicts to Theorem 2.3. Thus, we prove W^{s}((w^*(x), 0, 0))\cap \rho^{-1}(0, \infty) = \emptyset . Then by applying [29,Theorem 3], there exists \eta_0 > 0 such that \liminf\limits_{t\rightarrow\infty}\rho(\Theta(t)\phi)\geq \eta_0 for any \phi\in X^+ . This, together with Proposition 2.4 implies that \liminf\limits_{t\rightarrow\infty}u_i(x, t)\geq \eta for all i = 1, 2, 3 and x\in\bar\Omega , where \eta = \min\{\eta_0, m_0\} .
Furthermore, in view of [31,Theorem 4.7] and Theorem 4.2, system (1.2) admits at least one positive steady state. This ends the proof.
In this section, we consider the special case where all the coefficients in (1.2) are independent of the variable x , that is, d_1(x) = d_1, \; a(x) = a, \; \mu_j(x) = \mu_j \; (j = 1, 2, 3), \; d_2(x) = d_2, \; k(x) = k for all x\in\bar\Omega . We further assume that
(\textbf{C}) f(u_1, u_2) = h(u_1)f_1(u_2) and g(u_1, u_3) = h(u_1)g_1(u_3) , where h, f_1, g_1\in C^1(\mathbb{R}_+\times\mathbb{R}_+) are increasing functions, and h(v) = 0 (resp. f_1(v) = 0 , g_1(v) = 0 ) if and only if v = 0 . Moreover, d^2 f_1(v)/dv^2\le 0 and d^2 g_1(v)/dv^2\le 0 .
System (1.2) becomes homogeneous, that is,
\begin{equation} \begin{aligned} \frac{\partial u_1}{\partial t} & = d_1\Delta u_1+a-h(u_1)f_1(u_2)-h(u_1)g_1(u_3)-\mu_1u_1,\\ \frac{\partial u_2}{\partial t} & = d_2\Delta u_2+h(u_1)f_1(u_2)+h(u_1)g_1(u_3)-\mu_2u_2,\\ \frac{\partial u_3}{\partial t} & = ku_2-\mu_3u_3, \end{aligned} \end{equation} | (5.1) |
for x\in\Omega, t > 0 with the homogeneous Neumann boundary condition and nonnegative initial conditions. It then follows that w^*(x) = a/\mu_1 . By applying Krein-Rutman theorem, A_d+A_i is a compact and positive operator with a positive eigenfunction 1 corresponding to a positive principle eigenvalue
\begin{equation*} \label{R0-h} R_0 = {\beta_d\over\mu_2}+{\beta_i k\over\mu_2\mu_3}, \end{equation*} |
where \beta_d = h(a/\mu_1)f'_1(0) and \beta_i = h(a/\mu_1)g'_1(0) are constants. This implies that the basic reproduction numbers for system (5.1) and the corresponding diffusive-free ( d_1 = d_2 = 0 ) system are same. Denote (u_1^*, u_2^*, u_3^*) as the positive constant steady state, which satisfy the following equilibrium equations
\begin{equation} a-\mu_1 u_1^* = h(u_1^*)\left(f_1(u_2^*)+g_1(u_3^*)\right) = \mu_2 u_2^* = {\mu_2\mu_3\over k}u_3^*. \end{equation} | (5.2) |
Since the existence of constant steady state for system (5.1) same as for the corresponding ODE system. This, together with Theorem 3.1 in [8], yields the following lemma.
Lemma 5.1. If R_0 > 1 , then system (5.1) has a unique positive constant steady state (u_1^*, u_2^*, u_3^*) .
We next establish that R_0 is a threshold role for the global dynamics of system (5.1), and further give the global stability of the positive constant steady state.
Theorem 5.2. (i) If R_0\le1 , then the infection-free steady state (a/\mu_1, 0, 0) for system (5.1) is globally asymptotically stable in X^+ .
(ii) If R_0 > 1 , then system (5.1) admits a unique chronic infection steady state (u_1^*, u_2^*, u_3^*) , which is also homogeneous and globally asymptotically stable in X_0 .
Proof. Theorem 4.5 implies that (i) holds. We next prove the local asymptotic stability of the positive constant steady state (u_1^*, u_2^*, u_3^*) when R_0 > 1 . Linearizing system (5.1) at (u_1^*, u_2^*, u_3^*) , we obtain
{dU(t)\over dt} = d\Delta U(t)+L(U(t)), |
where U(t) = (u_1(x, t), u_2(x, t), u_3(x, t))^T , b = h'(u_1^*)(f_1(u_2^*)+g_1(u_3^*)) > 0 ,
d\Delta = \begin{pmatrix} d_1\Delta & 0 & 0\\ 0 & d_2\Delta & 0\\ 0 & 0 & 0 \end{pmatrix}, \ L(\phi) = \begin{pmatrix} -b-\mu_1 & -h(u_1^*)f'_1(u_2^*) & -h(u_1^*)g'_1(u_3^*)\\ b & h(u_1^*)f'_1(u_2^*)-\mu_2 & h(u_1^*)g'_1(u_3^*)\\ 0 & k & -\mu_3 \end{pmatrix}, |
and \text{dom}(d\Delta) = \{(u_1, u_2)^T: u_i\in W^{2, 2}(\Omega), {\partial u_i\over \partial\nu} = 0 \ \text{for} \ i = 1, 2.\} . Then the characteristic equation for the above linear system is
\lambda y-d\Delta y-L(y) = 0 \ \ \text{for} \ y\in\text{dom}(d\Delta), \ y\neq0. |
It is well known that the eigenvalue problem
\begin{aligned} -\Delta \psi& = \zeta \psi, \ \ &x\in\Omega,\\ {\partial\psi\over\partial\nu}& = 0, \ \ &x\in\partial\Omega, \end{aligned} |
has eigenvalues 0 = \zeta_0 < \zeta_1\le\zeta_2\le\cdots\le\zeta_n\le\zeta_{n+1}\le\cdots , with the corresponding eigenfunctions \hat\psi_n(x) . Substituting y = \sum\limits_{n = 0}^{\infty}y_n\hat\psi_n(x) into the characteristic equation gives
(\lambda+b+\mu_1+d_1\zeta_n)(\lambda+\mu_2+d_2\zeta_n)(\lambda+\mu_3) = (\lambda+\mu_1+d_1\zeta_n)\Phi_1(\lambda) |
for n = 0, 1, 2, \cdots , where \Phi_1(\lambda) = (\lambda+\mu_3)h(u_1^*)f'_1(u_2^*)+kh(u_1^*)g'_1(u_3^*) . The above characteristic equation is equivalent to
\begin{equation} (\lambda+b+\mu_1+d_1\zeta_n)({\lambda\over\mu_2}+1+{d_2\zeta_n\over\mu_2})(\lambda+\mu_3) = (\lambda+\mu_1+d_1\zeta_n)\Phi_2(\lambda) \end{equation} | (5.3) |
where
\Phi_2(\lambda) = \left({1\over 1+kg'_1(u_3^*)/ (\mu_3f'_1(u_2^*))}\lambda+\mu_3\right) \left({h(u_1^*)f'_1(u_2^*)\over \mu_2}+{kh(u_1^*)g'_1(u_3^*)\over\mu_2\mu_3}\right). |
We claim that all eigenvalues of (5.3) have negative real parts. Otherwise, suppose that \lambda = \sigma+\omega i is an eigenvalue satisfying \sigma\ge0 . Then for any nonnegative integer n , we have
|\lambda+b+\mu_1+d_1\zeta_n| \gt |\lambda+\mu_1+d_1\zeta_n|, \ |{\lambda\over\mu_2}+1+{d_2\zeta_n\over\mu_2}|\ge1. |
It follows from (\textbf{C}) and (5.2) that
{h(u_1^*)f'_1(u_2^*)\over \mu_2}+{kh(u_1^*)g'_1(u_3^*)\over\mu_2\mu_3}\le{h(u_1^*)f_1(u_2^*)\over \mu_2u_2^*} +{kh(u_1^*)g_1(u_3^*)\over\mu_2\mu_3u_3^*} = 1, |
which implies that |\Phi_2(\lambda)|\le|\lambda+\mu_3| . Therefore, we obtain
|(\lambda+b+\mu_1+d_1\zeta_n)({\lambda\over\mu_2}+1+{d_2\zeta_n\over\mu_2})(\lambda+\mu_3)| \gt |(\lambda+\mu_1+d_1\zeta_n)\Phi_2(\lambda)| |
for all integer n\ge0 . This is a contradiction. Hence we proved the claim, and (u_1^*, u_2^*, u_3^*) is locally asymptotically stable when R_0 > 1 .
Denote q(z) = z-1-\ln z . Clearly, q(z)\ge0 for z > 0 , and q(z) = 0 if and only if z = 1 . We next prove global attractiveness of (u_1^*, u_2^*, u_3^*) in X_0 by constructing a Lyapunov functional V_2: X_0\rightarrow\mathbb{R} as follows.
V_2(u_1,u_2,u_3) = \int_{\Omega}W(u_1,u_2,u_3)dx, |
where
W(u_1,u_2,u_3) = u_1-\int_{u_1^*}^{u_1}{h(u_1^*)\over h(s)}ds+u_2^*\; q\left({u_2\over u_2^*}\right) +{h(u_1^*)g_1(u_3^*)u_3^*\over ku_2^*}q\left({u_3\over u_3^*}\right). |
It follows from Theorems 2.3 and 4.7 that the solutions of system (5.1) are bounded and uniform persistent, which implies that V_2 and W are well-defined. Making use of the equilibrium equations (5.2), the time derivative of W along a positive solution of system (5.1) after a tedious calculation, is given by
\begin{align*} {dW\over dt} = &d_1(1-{h(u_1^*)\over h(u_1)})\Delta u_1+d_2(1-{u_2^*\over u_2})\Delta u_2-\mu_1(u_1-u_1^*)\left(1-{h(u_1^*)\over h(u_1)}\right)\\ &-h(u_1^*)g_1(u_3^*)\left[q\left({h(u_1^*)\over h(u_1)}\right)+q\left({u_2u_3^*\over u_2^*u_3}\right) -q\left({u_3g_1(u_3^*)\over u_3^*g_1(u_3)}\right)-q\left({u_2^*h(u_1)g_1(u_3)\over u_2h(u_1^*)g_1(u_3^*)}\right)\right]\\ &-h(u_1^*)f_1(u_2^*)\left[q\left({h(u_1^*)\over h(u_1)}\right)+q\left({u_2f_1(u_2^*)\over u_2^*f_1(u_2)}\right) +q\left({u_2^*h(u_1)f_1(u_2)\over u_2 h(u_1^*)f_1(u_2^*)}\right)\right]\\ &+h(u_1^*)g_1(u_3^*){u_3\over u_3^*}\left({g_1(u_3)\over g_1(u_3^*)}-1\right)\left({u_3^*\over u_3}-{g_1(u_3^*)\over g_1(u_3)}\right)\\ &+h(u_1^*)f_1(u_2^*)\left({u_2\over u_2^*}-{f_1(u_2)\over f_1(u_2^*)}\right)\left({f_1(u_2^*)\over f_1(u_2)}-1\right). \end{align*} |
Note from the Green's identity and Neumann boundary condition that
\begin{align*} \int_{\Omega}d_1\left(1-{h(u_1^*)\over h(u_1)}\right)\Delta u_1dx& = -d_1\int_{\Omega}{h(u_1^*)h'(u_1)\over h^2(u_1)}|\nabla u_1|^2 dx\le0,\\ \int_{\Omega}d_1\left(1-{u_2^*\over u_2}\right)\Delta u_2dx& = -d_2\int_{\Omega}{u_2^*\over u_2^2}|\nabla u_2|^2 dx\le0. \end{align*} |
Since h, f_1 and g_1 are increasing functions, f_1 and g_1 are concave down, then we have (u_1-u_1^*)(1-h(u_1^*)/ h(u_1))\ge0 and
\left({g_1(u_3)\over g_1(u_3^*)}-1\right)\left({u_3^*\over u_3}-{g_1(u_3^*)\over g_1(u_3)}\right)\le0, \ \ \left({u_2\over u_2^*}-{f_1(u_2)\over f_1(u_2^*)}\right)\left({f_1(u_2^*)\over f_1(u_2)}-1\right)\le0. |
Thus, dV_2/dt = \int_\Omega (dW/dt) dx\le0 . The largest invariant subset of dV_2/dt = 0 is the singleton (u_1^*, u_2^*, u_3^*) . By LaSalle-Lyapunov invariance principle, the positive constant steady state (u_1^*, u_2^*, u_3^*) is globally attractive in X_0 . The uniqueness of chronic infection steady state follows immediately from the global attractivity. This, together with the local asymptotic stability, yields that the global asymptotic stability of the positive constant steady state (u_1^*, u_2^*, u_3^*) in X_0 if R_0 > 1 .
H. Shu was partially supported by the National Natural Science Foundation of China (No.11971285, No.11601392), and the Fundamental Research Funds for the Central Universities.
The authors declare there is no conflicts of interest.
[1] | S. Bonhoeffer, R. M. May, G. M. Shaw, M. A. Nowak, Virus dynamics and drug therapy, Proc. Natl. Acad. Sci. USA, 94 (1997), 6971-6976. |
[2] | M. Y. Li, H. Shu, Impact of intracellular delays and target-cell dynamics on in vivo viral infections, SIAM J. Appl. Math., 70 (2010), 2434-2448. |
[3] | A. S. Perelson, D. E. Kirschner, R. de Boer, Dynamics of HIV infection of CD4 T cells, Math. Biosci., 114 (1993), 81-125. |
[4] | A. S. Perelson, P. W. Nelson, Mathematical analysis of HIV-I dynamics in vivo, SIAM Rev., 41 (1999), 3-44. |
[5] | A. S. Perelson, A. U. Neumann, M. Markowitz, M. J. Leonard, D. D. Ho, HIV-1 dynamics in vivo: virion clearance rate, infected cell life-span, and viral generation time, Science, 271 (1996), 1582-1586. |
[6] | H. Shu, L. Wang, J. Watmough, Global stability of a nonlinear viral infection model with infinitely distributed intracellular delays and CTL immune responses, SIAM J. Appl. Math., 73 (2013), 1280-1302. |
[7] | H. Shu, L. Wang, J. Watmough, Sustained and transient oscillations and chaos induced by delayed antiviral immune response in an immunosuppressive infection model, J. Math. Biol., 68 (2014), 477-503. |
[8] | H. Shu, Y. Chen, L. Wang, Impacts of the cell-free and cell-to-cell infection modes on viral dynamics, J. Dyn. Diff. Equat., 30 (2018), 1817-1836. |
[9] | N. M. Dixit, M. Markowitz, D. D. Ho, A. S. Perelson, Estimates of intracellular delay and average drug efficacy from viral load data of HIV-infected individuals under antiretroviral therapy, Antivir. Ther., 9 (2004), 237-246. |
[10] | M. A. Nowak, S. Bonhoeffer, A. M. Hill, R. Boehme, H. C. Thomas, Viral dynamics in hepatitis B virusinfection, Proc. Natl. Acad. Sci. USA, 93 (1996), 4398-4402. |
[11] | F. Wang, Y. Huang, X. Zou, Global dynamics of a PDE in-host viral model, Appl. Anal., 93 (2014), 2312-2329. |
[12] | Y. Wu, X. Zou, Dynamics and profiles of a diffusive host-pathogen system with distinct dispersal rates, J. Differential Equations, 264 (2018), 4989-5024. |
[13] | N. Martin, Q. Sattentau, Cell-to-cell HIV-1 spread and its implications for immune evasion, Curr. Opin. HIV AIDS, 4 (2009), 143-149. |
[14] | Q. Sattentau, Avoiding the void: cell-to-cell spread of human viruses, Nat. Rev. Microbiol., 6 (2008), 28-41. |
[15] | W. Hübner, G. P. McNerney, P. Chen, B. M. Dale, R. E. Gordan, F. Y. S. Chuang, et al., Quantitative 3D video microscopy of HIV transfer across T cell virological synapses, Science, 323 (2009), 1743-1747. |
[16] | H. L. Smith, Monotone Dynamical Systems: an introduction to the theory of competitive and cooperative systems, American Mathematical Society, Providence, 1995. |
[17] | A. Pazy, Semigroups of linear operators and application to partial differential equations, Springer, Berlin, 1983. |
[18] | R. Jr. Martin, H. L. Smith, Abstract functional differential equations and reaction-diffusion systems, Trans. AMS, 321 (1990), 1-44. |
[19] | C. V. Pao, Nonlinear parabolic and elliptic equations, Plenum, New York, 1992. |
[20] | M. W. Hirsch, The dynamical systems approach to differential equations, Bull. Am. Math. Soc., 11 (1984), 1-64. |
[21] | M. H. Protter, H. F. Weinberger, Maximum Principles in Differential Equations, Springer-Verlag, New York, 1984. |
[22] | H. R. Thieme, Spectral bound and reproduction number for infinite-dimensional population structure and time heterogeneity, SIAM J. Appl. Math., 70 (2009), 188-211. |
[23] | W. Wang, X-Q. Zhao, Basic reproduction numbers for reaction-diffusion epidemic models, SIAM J. Appl. Dyn. Syst., 11 (2012), 1652-1673. |
[24] | I. D. Gilbarg, N. S. Trudinger, Elliptic Partial Differential Equations of Second Order, 2nd ed., Springer, Berlin, 1983. |
[25] | J. K. Hale, Asymptotic Behavior of Dissipative Systems, American Mathematical Society, Providence, 1988. |
[26] | K.-J. Engel, R. Nagel, One-parameter Semigroups for Linear Evolution Equations, Graduate Texts in Mathematics, 194, Springer-Verlag, New York, 2000. |
[27] | X.-Q. Zhao, Dynamical systems in population biology, Second edition, CMS Books in Mathematics, Springer, Cham, 2017. |
[28] | L. J. S. Allen, B. M. Bolker, Y. Lou, A. L. Nevai, Asymptotic profiles of the steady states for an SIS epidemic reaction-diffusion model, Discrete Contin. Dyn. Syst., 21 (2008), 1-20. |
[29] | H. L. Smith, X-Q. Zhao, Robust persistence for semidynamical systems, Nonlinear Anal., 47 (2001), 6169-6179. |
[30] | H. R. Thieme, Convergence results and Poincaré-Bendixson trichotomy for asymptotically autonomous differential equations, J. Math. Biol., 30 (1992), 755-763. |
[31] | P. Magal, X-Q. Zhao, Global attractors and steady states for uniformly persistent dynamical systems, SIAM. J. Math. Anal., 37 (2005), 251-275. |
1. | Juping Ji, Genghong Lin, Lin Wang, Ali Mai, Spatiotemporal dynamics induced by intraguild predator diffusion in an intraguild predation model, 2022, 85, 0303-6812, 10.1007/s00285-022-01772-w | |
2. | Shu-Min Liu, Zhenguo Bai, Gui-Quan Sun, Global dynamics of a reaction-diffusion brucellosis model with spatiotemporal heterogeneity and nonlocal delay, 2023, 36, 0951-7715, 5699, 10.1088/1361-6544/acf6a5 | |
3. | Khellaf Ould Melha, Medjahed Djilali, Vaijanath L. Chinchane, Asha B. Nale, Sabri T. M. Thabet, Imed Kedim, Uniqueness and stability results of the abstract fractional spatial heterogeneous viral infection model, 2025, 2025, 1687-2770, 10.1186/s13661-025-02077-9 |