Processing math: 48%
Research article Special Issues

Identification and validation of biomarkers for epithelial-mesenchymal transition-related cells to estimate the prognosis and immune microenvironment in primary gastric cancer by the integrated analysis of single-cell and bulk RNA sequencing data


  • Received: 08 April 2023 Revised: 08 May 2023 Accepted: 23 May 2023 Published: 16 June 2023
  • Background: The epithelial-mesenchymal transition (EMT) is associated with gastric cancer (GC) progression and immune microenvironment. To better understand the heterogeneity underlying EMT, we integrated single-cell RNA-sequencing (scRNA-seq) data and bulk sequencing data from GC patients to evaluate the prognostic utility of biomarkers for EMT-related cells (ERCs), namely, cancer-associated fibroblasts (CAFs) and epithelial cells (ECs). Methods: scRNA-seq data from primary GC tumor samples were obtained from the Gene Expression Omnibus (GEO) database to identify ERC marker genes. Bulk GC datasets from the Cancer Genome Atlas (TCGA) and GEO were used as training and validation sets, respectively. Differentially expressed markers were identified from the TCGA database. Univariate Cox, least-absolute shrinkage, and selection operator regression analyses were performed to identify EMT-related cell-prognostic genes (ERCPGs). Kaplan-Meier, Cox regression, and receiver-operating characteristic (ROC) curve analyses were adopted to evaluate the prognostic utility of the ERCPG signature. An ERCPG-based nomogram was constructed by integrating independent prognostic factors. Finally, we evaluated the correlations between the ERCPG signature and immune-cell infiltration and verified the expression of ERCPG prognostic signature genes by in vitro cellular assays. Results: The ERCPG signature was comprised of seven genes (COL4A1, F2R, MMP11, CAV1, VCAN, FKBP10, and APOD). Patients were divided into high- and low-risk groups based on the ERCPG risk scores. Patients in the high-risk group showed a poor prognosis. ROC and calibration curves suggested that the ERCPG signature and nomogram had a good prognostic utility. An immune cell-infiltration analysis suggested that the abnormal expression of ERCPGs induced the formation of an unfavorable tumor immune microenvironment. In vitro cellular assays showed that ERCPGs were more abundantly expressed in GC cell lines compared to normal gastric tissue cell lines. Conclusions: We constructed and validated an ERCPG signature using scRNA-seq and bulk sequencing data from ERCs of GC patients. Our findings support the estimation of patient prognosis and tumor treatment in future clinical practice.

    Citation: Kaiyu Shen, Shuaiyi Ke, Binyu Chen, Tiantian Zhang, Hongtai Wang, Jianhui Lv, Wencang Gao. Identification and validation of biomarkers for epithelial-mesenchymal transition-related cells to estimate the prognosis and immune microenvironment in primary gastric cancer by the integrated analysis of single-cell and bulk RNA sequencing data[J]. Mathematical Biosciences and Engineering, 2023, 20(8): 13798-13823. doi: 10.3934/mbe.2023614

    Related Papers:

    [1] Frédéric Mazenc, Michael Malisoff, Patrick D. Leenheer . On the stability of periodic solutions in the perturbed chemostat. Mathematical Biosciences and Engineering, 2007, 4(2): 319-338. doi: 10.3934/mbe.2007.4.319
    [2] Gonzalo Robledo . Feedback stabilization for a chemostat with delayed output. Mathematical Biosciences and Engineering, 2009, 6(3): 629-647. doi: 10.3934/mbe.2009.6.629
    [3] Xue Liu, Xin You Meng . Dynamics of Bacterial white spot disease spreads in Litopenaeus Vannamei with time-varying delay. Mathematical Biosciences and Engineering, 2023, 20(12): 20748-20769. doi: 10.3934/mbe.2023918
    [4] Paolo Fergola, Marianna Cerasuolo, Edoardo Beretta . An allelopathic competition model with quorum sensing and delayed toxicant production. Mathematical Biosciences and Engineering, 2006, 3(1): 37-50. doi: 10.3934/mbe.2006.3.37
    [5] Xiaomeng Ma, Zhanbing Bai, Sujing Sun . Stability and bifurcation control for a fractional-order chemostat model with time delays and incommensurate orders. Mathematical Biosciences and Engineering, 2023, 20(1): 437-455. doi: 10.3934/mbe.2023020
    [6] Zhichao Jiang, Xiaohua Bi, Tongqian Zhang, B.G. Sampath Aruna Pradeep . Global Hopf bifurcation of a delayed phytoplankton-zooplankton system considering toxin producing effect and delay dependent coefficient. Mathematical Biosciences and Engineering, 2019, 16(5): 3807-3829. doi: 10.3934/mbe.2019188
    [7] Ying He, Yuting Wei, Junlong Tao, Bo Bi . Stationary distribution and probability density function analysis of a stochastic Microcystins degradation model with distributed delay. Mathematical Biosciences and Engineering, 2024, 21(1): 602-626. doi: 10.3934/mbe.2024026
    [8] Yan Wang, Minmin Lu, Daqing Jiang . Viral dynamics of a latent HIV infection model with Beddington-DeAngelis incidence function, B-cell immune response and multiple delays. Mathematical Biosciences and Engineering, 2021, 18(1): 274-299. doi: 10.3934/mbe.2021014
    [9] Rundong Zhao, Qiming Liu, Huazong Zhang . Dynamical behaviors of a vector-borne diseases model with two time delays on bipartite networks. Mathematical Biosciences and Engineering, 2021, 18(4): 3073-3091. doi: 10.3934/mbe.2021154
    [10] Hong Yang . Global dynamics of a diffusive phytoplankton-zooplankton model with toxic substances effect and delay. Mathematical Biosciences and Engineering, 2022, 19(7): 6712-6730. doi: 10.3934/mbe.2022316
  • Background: The epithelial-mesenchymal transition (EMT) is associated with gastric cancer (GC) progression and immune microenvironment. To better understand the heterogeneity underlying EMT, we integrated single-cell RNA-sequencing (scRNA-seq) data and bulk sequencing data from GC patients to evaluate the prognostic utility of biomarkers for EMT-related cells (ERCs), namely, cancer-associated fibroblasts (CAFs) and epithelial cells (ECs). Methods: scRNA-seq data from primary GC tumor samples were obtained from the Gene Expression Omnibus (GEO) database to identify ERC marker genes. Bulk GC datasets from the Cancer Genome Atlas (TCGA) and GEO were used as training and validation sets, respectively. Differentially expressed markers were identified from the TCGA database. Univariate Cox, least-absolute shrinkage, and selection operator regression analyses were performed to identify EMT-related cell-prognostic genes (ERCPGs). Kaplan-Meier, Cox regression, and receiver-operating characteristic (ROC) curve analyses were adopted to evaluate the prognostic utility of the ERCPG signature. An ERCPG-based nomogram was constructed by integrating independent prognostic factors. Finally, we evaluated the correlations between the ERCPG signature and immune-cell infiltration and verified the expression of ERCPG prognostic signature genes by in vitro cellular assays. Results: The ERCPG signature was comprised of seven genes (COL4A1, F2R, MMP11, CAV1, VCAN, FKBP10, and APOD). Patients were divided into high- and low-risk groups based on the ERCPG risk scores. Patients in the high-risk group showed a poor prognosis. ROC and calibration curves suggested that the ERCPG signature and nomogram had a good prognostic utility. An immune cell-infiltration analysis suggested that the abnormal expression of ERCPGs induced the formation of an unfavorable tumor immune microenvironment. In vitro cellular assays showed that ERCPGs were more abundantly expressed in GC cell lines compared to normal gastric tissue cell lines. Conclusions: We constructed and validated an ERCPG signature using scRNA-seq and bulk sequencing data from ERCs of GC patients. Our findings support the estimation of patient prognosis and tumor treatment in future clinical practice.



    Chemostat models have been extensively studied in the literature. However, many issues still remain elusive and understudied. In particular, in this paper we will focus on the global stability of periodic solutions for chemostats described by the following family of delay differential equations with periodic inputs

    {˙s(t)=Ds0(t)Ds(t)μ(s(t))x(t)˙x(t)=x(t)[μ(s(tτ))D]whent0, (1.1)
    s(θ)=φ(θ)andx(0)=x0ifθ[τ,0] (1.2)

    whose detailed description will be given below.

    The existence and uniqueness of periodic solutions for some cases of (1.1) and (1.2) has been recently addressed in [1] and generalized in [2,3]. Nevertheless, to the best of our knowledge, there are no results about the global stability for the periodic solutions mentioned above. The main contribution of this paper is to provide a set of sufficient conditions for global stability.

    The chemostat is a continuous bioreactor where one or several microbial species are cultivated in an homogeneous liquid medium, which contains all the nutrients ensuring its growth with the exception of a specific one named the limiting substrate, or just substrate. There exists a vast literature devoted to both the modeling of the chemostat and its applications and the study of the arising dynamical issues. We would like to refer the reader to the seminal works [4,5,6] and the monographs [7,8,9] where more information on the subject can be found.

    More specifically, the system of differential delay equations (1.1) describes the interaction between a limiting substrate whose concentration at time t is denoted by s(t) and one species of microbial biomass whose concentration at time t is denoted by x(t). Moreover, the initial conditions of (1.1) are described by the continuous function φ:[τ,0][0,+) and x00 in (1.2).

    In the first equation of (1.1), the term Ds0(t) means that the limiting substrate is pumped inside the chemostat at a rate D>0, with concentration described by a positive continuous and ω–periodic function ts0(t). On the other hand, the terms Ds(t) and Dx(t) in the first and second equation respectively, indicate that the mixture of microbial biomass and substrate is pumped outside at similar rate D>0, which ensures a constant volume.

    The consumption of the substrate is described by the term μ(s(t))x(t) in the first equation of (1.1), while the per–capita growth of the microbial biomass is described by the term ˙x(t)x(t)=μ(s(tτ)) in the second equation. The term μ() denotes the uptake function, which makes a link between the process of consumption of substrate and its consequences on the microbial growth. In this context, the case τ=0 means that the consumption of substrate has an immediate effect on the microbial growth, while the case τ>0 refers to the existence of an interval of time [0,τ], which is necessary in order to metabolize the substrate.

    We will assume that the uptake function μ:[0,+)[0,+) is the Monod function

    μ(s)=μmsks+swith ks>0 and μm>0. (1.3)

    Note that μ(0)=0 means that there is no microbial growth in the absence of substrate, while μ(s)>0 states that the per–capita growth rate is an increasing function with respect to the substrate, which is upper bounded by the constant μm>0, namely, the maximal growth rate, and verifies

    μ(s)<μmfor any s0 and lims+μ(s)=μm. (1.4)

    On the other hand, the constant ks is named as the half–saturation constant because μ(ks)=μm/2. We refer the reader to [9] for further details.

    Finally, observe for later use that

    μ(a)μ(b)=μmks(ab)(ks+a)(ks+b), (1.5)

    for all a0 and b0.

    System (1.1) encompasses two relevant particular cases: On one hand, when considering a constant input of nutrient s0>0, system (1.1) becomes the autonomous DDE system

    {˙s(t)=Ds0Ds(t)μ(s(t))x(t)˙x(t)=x(t)[μ(s(tτ))D]whent0 (1.6)

    which was introduced by Wangersky & Cunningham in [10], and initially studied by Caperon [11] and Thingstad [12]. The local asymptotic stability results for the equilibria are obtained in [13,14,15] by studying the roots of the characteristic quasipolynomial equation associated with the linear approximation of (1.6) around its positive equilibrium. The existence of periodic equations of (1.6) and related models has been addressed in [16,17], and we also refer to [9, Ch.10]. Sufficient conditions for global asymptotic stability have been obtained in [13,18,19] by constructing Lyapunov–Krasovskii functions.

    When considering τ=0, system (1.1) becomes the nonautonomous ODE system

    {˙s(t)=Ds0(t)Ds(t)μ(s(t))x(t)˙x(t)=x(t)[μ(s(t))D], (1.7)

    which has been employed, among other things, to emulate oscillations in marine ecosystems and to study optimization issues in periodically operating bioprocesses. We refer the reader to Section 2 of [1] and references therein for additional details.

    The mathematical study of (1.7) has been carried out by a numerical approach [20], and also by the construction of Poincaré maps [21,22,23]. In order to study the existence of ω–periodic solutions of (1.1) and (1.7), a pivotal tool is the study of the properties of the scalar equation

    ˙ξ=Ds0(t)Dξ (1.8)

    which corresponds to (1.1) and (1.7) in the absence of microbial biomass. It is well known that (1.8) has a globally attractive, positive, and ω–periodic solution tΣ(t), which is named as the washout solution:

    Σ(t)=tDeD(tr)s0(r)dr. (1.9)

    In [22,23], Wolkowicz and Zhao proved that the following average inequality involving the uptake function and the washout solution

    1ωω0μ(Σ(r))dr>D (1.10)

    implies the existence, uniqueness, and attractiveness of a nontrivial ω–periodic solution of the undelayed system (1.7), which will be denoted by t(sp(t),xp(t)) throughout this article. Moreover, it is useful to recall that a direct consequence of the average inequality (1.10) is

    D<μm. (1.11)

    The ω–periodic DDE system (1.1) can be seen as the coupling of the autonomous DDE system (1.6) with the ω–periodic ODE system (1.7). Nevertheless, in comparison with the extensive knowledge we have on systems (1.6) and (1.7), there exist fewer results about the qualitative properties of (1.1). A first step to fill this gap is given in [1] and [2,3], where the results of [22,23] are partially generalized:

    Proposition 1. [1, Th.2 and Th.3] System (1.1) has a non-trivial and positive ω–periodic solution if and only if the average condition (1.10) is satisfied.

    Furthermore:

    a) If t(s(t),x(t)) is a non-trivial ω–periodic solution, then

    0<s(t)<Σ(t)and0<x(t)for any t0.

    b) The non-trivial ω–periodic solution is unique when the delay τ is sufficiently small and will be denoted by t(s,τ(t),x,τ(t)).

    As pointed out in the discussion of article [1], an open question arises from this last result, which is to obtain a set of sufficient conditions ensuring the global attractiveness of the unique non–trivial ω–periodic solution t(s,τ(t),x,τ(t)) for small enough delays. The main contribution of this article is to provide a partial answer by considering a particular uptake function μ(), namely Monod's function defined by (1.3) and a family of functions s0(). The approach to this problem has already been proposed in [1, Section 6], which is the construction of a Lyapunov like function.

    Section 2 describes the three main results of this article: Theorem 1 is a result of comparison between the solutions of (1.1) with the ω–periodic solution t(sp(t),xp(t)) of the undelayed system (1.7), and Theorem 2 gives sufficient conditions for the local asymptotic stability of the ω–periodic solution t(s,τ(t),x,τ(t)), while Theorem 3 furnishes conditions that ensure global asymptotic stability. Both conditions are obtained by constructing Lyapunov–like functions fashioned along ideas and techniques of [24]. Sections 3–5 are devoted to the proof of Theorems 1–3, respectively. Section 6 provides an example based on the culture of the microalgae Dunaliella tertiolecta, and several numerical simulations are presented to illustrate our results. Section 7 contains a brief discussion about the scope and limitations of the results of this work.

    Throughout this article, we will work under the assumptions that the average condition (1.10) is satisfied, and the delay τ is small enough such that statement b) of the Proposition 1 ensures the existence and uniqueness of a non-trivial and positive ω–periodic solution t(s,τ(t),x,τ(t)) of system (1.1).

    The main results are described in terms of an upper bound for the limiting substrate and a system equivalent to (1.1). In consequence, a necessary first step in order to state the main results is to obtain an explicit upper bound for the substrate and to introduce a couple of transformations to system (1.1).

    In this subsection, we will establish that the solutions of the system (1.1) are bounded for sufficiently small delays. First, by (1.9) it is easy to see that, for any tR,

    0<sm=s0min:=mint[0,ω]{s0(t)}Σ(t)maxt[0,ω]{s0(t)}:=s0M. (2.1)

    It can be easily verified that the positive orthant of system (1.1) is positively invariant. In addition, by (2.1) we deduce that s0(t)s0M for any t0.

    Since the average condition (1.10) is satisfied, Proposition 1 then states that system (1.1) admits a positive ω–periodic solution. Moreover, we define a fixed number s such that s>s0M. Then, we notice for later use that, by using the fact that μ is increasing and Σ(t)s0M, we obtain the inequalities

    Dμ(s0M)μ(s). (2.2)

    The following result states the uniform boundedness for the limiting substrate.

    Lemma 1. If t(s(t),x(t)) is a solution of (1.1), then there exists t0 dependent on the initial conditions such that, for all tt, the following inequality is satisfied:

    s(t)s. (2.3)

    Proof. By using the positiveness of the solution, we can see that the substrate satisfies the differential inequality

    ˙sDs0(t)Ds(t).

    Now, the result is a consequence of s0(t)s0M<s for any t0, combined with comparison results for scalar differential equations such as [25, Ch.4] and [26, Lemma 3.4].

    In order to state our main results, it is useful to introduce the variable a(t)=x(t+τ), which transforms (1.1) into

    {˙s(t)=Ds0(t)Ds(t)μ(s(t))a(tτ)˙a(t)=[μ(s(t))D]a(t). (2.4)

    Furthermore, by using the identity

    a(tτ)=a(t)ettτ(Dμ(s(r)))dr (2.5)

    we obtain the alternative formulation

    {˙s(t)=Ds0(t)Ds(t)μ(s(t))ettτ(Dμ(s(r)))dra(t)˙a(t)=[μ(s(t))D]a(t). (2.6)

    The study of (2.6) has technical advantages with respect to (1.1). It is important to emphasize that the existence and uniqueness of the non-trivial and positive ω–periodic solution t(s,τ(t),x,τ(t)) of system (1.1) is equivalent to the existence and uniqueness of the non-trivial and positive ω–periodic solution t(s,τ(t),a,τ(t)) of system (2.6).

    The first result considers the solutions of the transformed system (2.6) when the delay is sufficiently small, and compares them with the non–trivial ω–periodic solution t(sp(t),xp(t)) of the undelayed system (1.7):

    Theorem 1. Consider the DDE system (1.1) and (1.2) with an uptake function μ() described by (1.3), and input nutrient ts0(t) such that the average condition (1.10) is satisfied. Let τ>0 be small enough such that

    τμ(s)(s)eτ(s)D4 (2.7)

    with

    (s)=max{μ(s)D,D} (2.8)

    is satisfied. Then, there exist positive constants (p1,p2) and τμ(0,τ) such that, when τ[0,τμ], there is t0 such that, for all tt, the solutions of (2.6) satisfy

    |s(t)sp(t)|5p13μ(s)s(eτμ(s)1),|a(t)xp(t)|5p23μ(s)s(eτμ(s)1). (2.9)

    The second result is concerned with the local asymptotic stability of the above mentioned non-trivial and positive ω–periodic solution t(s,τ(t),a,τ(t)), and it also assumes that the delay must be upper bounded as follows:

    ττ:=1Dμmln(1D8μm) (2.10)

    which is well defined due to D<μm, as it was stated in (1.11).

    Theorem 2. Consider the DDE system (1.1) and (1.2) with an uptake function μ() described by (1.3) and input nutrient ts0(t) such that the average condition (1.10) is satisfied. Then, there exists a small enough delay τ(0,τ] such that, when τ[0,τ], the unique nontrivial ω–periodic solution t(s,τ(t),a,τ(t)) of system (2.4) is a locally exponentially stable solution.

    The last result faces an unsolved problem from [1], namely, to obtain sufficient conditions ensuring the global asymptotic stability of the unique non–trivial ω–periodic solution t(s,τ(t),x,τ(t)) when considering small delays. In this context, the main contribution of this article is to provide a partial answer, by considering a nutrient input described by the ω–periodic continuous function

    s0(t)=Sc+εγ(t) (2.11)

    where |γ(t)|¯γ for all t0 and the positive constants Sc and ε are such that s0() is positive.

    Let us observe that, in this case, the washout solution tΣ(t), defined in (1.9), becomes

    Σ(t)=Sc+εϕ(t,γ)withϕ(t,γ):=tDeD(tr)γ(r)dr.

    Theorem 3. Consider the DDE system (1.1) and (1.2) with an uptake function μ() and input nutrient ts0(t), respectively described by (1.3) and (2.11). If μ, Sc, ε, and γ are such that

    1ωω0μ(Sc+εϕ(t,γ))dt>D

    then the ω–periodic solution t(s,τ(t),x,τ(t)) is globally attractive when the delay is sufficiently small.

    This theorem shows that the unique positive ω–periodic solution of the DDE system (1.1) and the solutions of the ODE system (1.7) are close when the time is large and the delay τ is small enough. To prove this property, it will be useful to recall that system (1.1) can be transformed into (2.6). In addition, notice that (2.6) can be rewritten as follows:

    {˙s(t)=Ds0(t)Ds(t)μ(s(t))a(t)+ϖ(t)˙a(t)=[μ(s(t))D]a(t) (3.1)

    where the term ϖ(t) is given by

    ϖ(t)=μ(s(t))[1ettτ(Dμ(s(r)))dr]a(t).

    Notice that (3.1) can be seen as a particular case of the ODE system (1.7) with an additive perturbation tδ(t) described by

    {˙s(t)=Ds0(t)Ds(t)μ(s(t))x(t)+δ(t)˙x(t)=x(t)[μ(s(t))D]. (3.2)

    It will be useful to quantify the effect caused by an input δ(), which is assumed to be a bounded function verifying additional properties which will be described later.

    Under the above perspective, we can see that systems (3.1) and (3.2) have a similar structure. This sheds light on a way to prove Theorem 1: to study the asymptotic behavior of system (3.2) by considering perturbations tδ(t) having asymptotic behavior emulating the properties of the map tϖ(t) when τ is small enough.

    The following result concludes the uniform boundedness for the microbial biomass of system (1.1).

    Lemma 2. If τ[0,τ] with τ satisfying inequality (2.7), then there is tt+τ dependent on the initial conditions, such that, for all tt+τ, the biomass component of any positive solution of (1.1) satisfies

    x(t)53s. (3.3)

    Proof. Let t(s(t),a(t)) be a positive solution of (2.6). Let us note that the change of variable b(t)=s(t)+a(t) leads to the equation with distributed delay

    ˙b(t)=Ds0(t)Db(t)+[1ettτ(Dμ(s(r)))dr]μ(s(t))a(t).

    By Lemma 1, combined with the inequality |1eu||u|e|u|, we deduce that, when tt,

    ˙b(t)Ds0(t)Db(t)+|1ettτ(Dμ(s(r)))dr|μ(s)a(t)Ds0(t)Db(t)+|ttτη(r)dr|e|ttτη(r)dr|μ(s)a(t)

    where η(r)=Dμ(s(r)).

    Lemma 1 ensures that s(r)<s when tt+τ. In addition, by using the fact that μ() is increasing, from inequality (2.2), it follows that μ(s(r))μ(s) and Dμ(s). Then, we obtain that

    |ttτη(r)dr|ttτmax{μ(s(r))D,Dμ(s(r))}dr

    which implies that

    |ttτη(r)dr|(s)τ (3.4)

    with defined in (2.8). By using (2.2) combined with a(t)b(t), we deduce that, when tt+τ,

    ˙b(t)DsDb(t)+τ(s)μ(s)eτ(s)a(t)Ds+[τμ(s)(s)eτ(s)D]b(t).

    Then, it follows that, when τ[0,τ], inequality (2.7) implies that

    ˙b(t)Ds34Db(t).

    By using the aforementioned comparison results for scalar differential inequalities, we deduce the existence of tt+τ, dependent on the initial conditions, such that the inequality

    b(t)53s

    is satisfied for all tt. This inequality together with b(t)>a(t)=x(t+τ) allows us to conclude the proof.

    The concept of Input to State Stability (ISS) was proposed for the first time by E. Sontag in [27], and we refer the reader to [28] for additional details. It characterizes the behavior of the solutions of a system in terms of the external input. More specifically, a system ˙u=f(t,u,δ) is said to be Input to State Stable with restriction if there are ¯δ>0 and functions γK and βKL such that, for any bounded and piecewise continuous input δ() such that |δ(t)|¯δ for all t0, the corresponding solution tϕ(t,t0,u0,δ) passing through u0 at t=t0 satisfies

    |ϕ(t,t0,u0,δ)|β(|u0|,tt0)+γ(maxs[t0,t]|δ(s)|),for alltt0. (3.5)

    Let us recall that γK if it is continuous, increasing, and such that γ(0)=0 and βKL if sβ(s,t)K for any t0 and, for any r>0, the function tβ(r,t) strictly decreases to zero, see [26] for details.

    Now, we show that system (3.2) is ISS with restriction. First, as t(sp(t),xp(t)) is the non–trivial ω–periodic solution of the ODE system (1.7) and the change of variables b(t)=s(t)+x(t) leads to the system

    {˙b(t)=Ds0(t)Db(t)˙x(t)=x(t)[μ(b(t)x(t))D] (3.6)

    and allows us to define the ω–periodic functions

    bp(t)=sp(t)+xp(t),χp(t)=ln(xp(t))andχminp=mint[0,ω]{χp(t)}.

    It is clear that t(bp(t),xp(t)) is the unique ω–periodic solution of the system above.

    Now, we are ready to state the following result:

    Lemma 3. Let us consider perturbations tδ(t) of system (3.2) with the following property: there exist a positive constant ¯δ and a threshold T0:=T0(δ)>0 such that

    suptT0|δ(t)|¯δ. (3.7)

    Furthermore, if ¯δ is small enough and there exists T1>T0 dependent on the initial conditions such that any solution of (3.2) verifies s(t)s for any t>T1, then:

    a) There exists taT0, dependent on the initial conditions, such that

    |b(t)bp(t)|2¯δDfor all tta. (3.8)

    b) For any ε small enough there exists tb>ta, dependent on the initial conditions such that

    e(1+ε)<x(t)xp(t)eln(2e1)+εfor any ttb. (3.9)

    c) There exists tc>tb, dependent on the initial conditions, and a positive constant p such that

    |x(t)xp(t)|xp(t)(ep¯δ1)for any ttc. (3.10)

    d) There exist two positive constants, p1 and p2 such that any solution t(s(t),x(t)) of (3.2) satisfies

    |s(t)sp(t)|p1¯δand|x(t)xp(t)|p2¯δfor all ttc. (3.11)

    Proof. See Appendix.

    A careful reading of the proof shows that:

    i) The constants p, p1, and p2 mentioned in statements c) and d) depend on the non–trivial ω–periodic solution t(sp(t),xp(t)) of the undelayed system (1.7) and the set of parameters ks, D, and s. In fact, it will be proved that

    p=3(ks+s)2k2sDe(1+ε)eχminpwhereχminp:=mint[0,ω]ln(xp(t))

    and ε is a constant from (3.9).

    ii) The statements a), b), c), and d) implicitly impose different smallness conditions for ¯δ. In fact, inequality (3.8) is valid for any positive ¯δ, while inequalities (3.9) and (3.10) are obtained for ¯δ having an explicit upper bound. Finally, estimation (3.11) is deduced for ¯δ small enough such that ln(1+p¯δ) can be approximated by its first–order MacLaurin expansion.

    Remark 1. From estimations (3.8) and (3.10), we deduce that if t>tb, then

    2¯δD+bp(t)b(t)andxp(t)(ep¯δ1)+xp(t)x(t).

    Thus, when ¯δ is sufficiently small,

    bb(t),andxx(t)

    with b=12mint[0,ω]bp(t) and x=12mint[0,ω]xp(t) for all t>tb.

    Now, we go back to study system (3.1) which, as we know, is equivalent to (2.6). First, notice that the input tϖ(t) verifies property (3.7). In fact, as t(s(t),a(t)) is solution of (2.6), Lemmas 1 and 2 ensure that if t>t>t+τ>t with τ[0,τ] satisfying inequality (2.7), then s(t)<s and a(t)53s for any t>t. These facts, combined with the inequalities |1ex|e|x|1 and D<μ(s), imply that for any t>t,

    |ϖ(t)|=|μ(s(t))||1ettτ(Dμ(s(r)))dr||a(t)|μ(s)|1ettτ(Dμ(s(r)))dr||x(t+τ)|53sμ(s)|1ettτ(Dμ(s(r)))dr|53sμ(s)(e|ttτ[μ(s(r))D]dr|1)53μ(s)s(eτμ(s)1)=¯ϖ(τ):=¯ϖ.

    We observe that, by an appropriate choice of τ, the constant ¯ϖ can be rendered as small as desired.

    Second, as we know that systems (3.1) and (3.2) have the same structure, statement c) of Lemma 3, with ϖ(t) instead of δ(t), implies the existence of a constant p and a time tγ>t, dependent on the initial conditions, such that

    |x(t)xp(t)|xp(t)(ep¯ϖ1)for any ttγ.

    Finally, when considering delays small enough such that p¯ϖ is a good approximation of ep¯ϖ1, statement d) of Lemma 3 implies the existence of a couple of positive constants (p1,p2), such that

    |s(t)sp(t)|p1¯ϖand|x(t)xp(t)|p2¯ϖfor any ttγ.

    This concludes the proof.

    Let us recall that, under the assumptions of Proposition 1, we know that system (1.1) has a unique positive ω–periodic solution (s,τ(t),x,τ(t)). Moreover, by performing a change of variables on the original system (1.1), we obtained system (2.4) introduced in Section 2:

    {˙s(t)=Ds0(t)Ds(t)μ(s(t))a(tτ)˙a(t)=[μ(s(t))D]a(t).

    We denoted its unique positive ω–periodic solution by (s,τ(t),a,τ(t)). Nevertheless, to simplify the notation, we simply write (s(t),a(t)).

    The proof of Theorem 2 is decomposed into several intermediate results. Moreover, we point out that we obtain an explicit value for τ.

    Lemma 4. Let inequality (2.10) be satisfied. Then, the unique ω–periodic and positive solution of system (2.4) t(s(t),a(t)) satisfies the inequalities

    a(t)43s0Mands(t)s0Mfor any t0. (4.1)

    Proof. Let us observe that

    ˙s(t)+˙a(t)=Ds0(t)D[s(t)+a(t)]+μ(s(t))[a(t)a(tτ)]Ds0MD[s(t)+a(t)]+μ(s(t))[a(t)a(tτ)].

    By (2.5) we know that a(t)=ettτ[μ(s(r))D]dra(tτ), which gives

    ˙s(t)+˙a(t)Ds0MD[s(t)+a(t)]+μ(s(t))[1ettτ[Dμ(s(r))]dr]a(t).

    Now, by (1.4), we obtain

    ttτ[Dμ(s(r))]drτ(Dμm)

    for all t0, and we deduce that

    1ettτ[Dμ(s(r))]dr1eτ(Dμm)

    which allows us to obtain the estimations

    ˙s(t)+˙a(t)Ds0MD[s(t)+a(t)]+μm(1eτ(Dμm))a(t)Ds0M+[μm(1eτ(Dμm))D][s(t)+a(t)].

    Bearing in mind (1.11), we note that (2.10) is equivalent to the inequality

    μm(1eτ(Dμm))D8

    which ensures that

    ˙s(t)+˙a(t)Ds0M78D[s(t)+a(t)].

    Then, by the positiveness of the periodic solutions together with the comparison result for scalar differential inequalities, we deduce the existence of t0 such that

    a(t)87s0M<43s0Mfor alltt.

    On the other hand, s(t)s0M is a consequence of Proposition 1 and (2.1).

    In order to simplify the local stability analysis of system (2.4), let us perform the change of variable ξ=ln(a). This gives

    {˙s(t)=Ds0(t)Ds(t)μ(s(t))eξ(tτ)˙ξ(t)=μ(s(t))D. (4.2)

    Let us observe that (s(t),ξ(t)), where ξ(t)=ln(a(t)) is a ω–periodic solution of (4.2). Moreover, it is interesting to note that the washout solution tΣ(t) defined in (1.9) can be written as follows:

    Σ(t)=D1eDωttωeD(rt)s0(r)dr. (4.3)

    Lemma 5. The unique ω–periodic and positive solution of system (2.4), namely, t(s(t),a(t)), verifies the following inequalities for any t[0,ω]:

    Σ(t)μmˇaDs(t)<Σ(t)with ˇa=maxt[0,ω]a(t) (4.4)

    and

    a(t)a_:=DeDωμmmint[0,ω]{Σ(t)s(t)}>0. (4.5)

    Proof. The right inequality of (4.4) follows directly from statement a) of Proposition 1. In order to prove the left inequality, we use (1.4) and the constant ˇa defined in (4.4) to deduce that

    ˙s(t)Ds0(t)Ds(t)μmˇafor all tR.

    The integration of the above inequality combined with the ω–periodicity of s() implies that

    eDts(t)eD(tω)s(t)+ttωeDr[Ds0(r)μmˇa]dr.

    We deduce that

    (1eDω)s(t)DttωeD(tr)s0(r)dr(1eDω)μmˇaD.

    By using (4.3), we obtain that

    s(t)Σ(t)μmˇaD.

    Thus, the left inequality of (4.4) holds.

    A direct consequence of inequalities (4.4) together with the ω–periodicity of Σ() and s() is that

    ˇaDμm{Σ(t)s(t)}Dμmmint[0,ω]{Σ(t)s(t)}>0. (4.6)

    Now, for all t1R and t2R, we have

    a(t1)=a(t2)et1t2[μ(s(r))D]dr.

    Thus, if we choose t2[0,ω] such that a(t2)=ˇa, then we obtain

    a(t)=ˇaett2[μ(s(r))D]drˇaeDω

    for any t[0,ω]. From (4.6), we deduce that

    a(t)DeDωμmmint[0,ω]{Σ(t)s(t)}:=a_

    for all tR. This concludes the proof.

    The next result provides the first order approximations of the solutions of (4.2) around t(s(t),ξ(t)):

    Lemma 6. The linear approximation of system (4.2) is described by the system

    {˙¯s(t)=D¯s(t)μ(s(t))eξ(t)+α2(t)¯s(t)μ(s(t))eξ(t)+α2(t)¯ξ(t)+α1(t,¯st)˙¯ξ(t)=μ(s(t))¯s(t) (4.7)

    where

    α1(t,¯st)=μ(s(t))eξ(tτ)ttτμ(s(r))¯s(r)dr (4.8)

    and

    α2(t)=ttτ[Dμ(s(r))]dr. (4.9)

    Proof. One can check easily that the linear approximation around the solution (s(t),ξ(t)) of system (4.2) is

    {˙¯s(t)=D¯s(t)μ(s(t))eξ(tτ)¯s(t)μ(s(t))eξ(tτ)¯ξ(tτ)˙¯ξ(t)=μ(s(t))¯s(t)

    which also can be written as the linear ω–periodic system

    ˙ζ(t)=A(t)ζ(t)+B(t)ζ(tτ) (4.10)

    where

    ζ(t)=(¯s(t)¯ξ(t))T,A(t)=[{D+μ(s(t))eξ(tτ)}0μ(s(t))0]andB(t)=[0μ(s(t))eξ(tτ)}00].

    Now, observing that

    ¯ξ(t)¯ξ(tτ)=ttτμ(s(r))¯s(r)dr

    the above system becomes

    {˙¯s(t)=D¯s(t)μ(s(t))eξ(tτ)¯s(t)μ(s(t))eξ(tτ)¯ξ(t)+α1(t,¯st)˙¯ξ(t)=μ(s(t))¯s(t),

    with α1(t,) defined in (4.8). Now, by noticing that

    ξ(tτ)=ξ(t)ttτ[μ(s(r))D]dr

    the above system becomes (4.7). This allows us to conclude the proof.

    Now, we start to study the stability properties of (4.7). To ease this analysis, we adopt the simplified notation:

    {˙¯s=D¯sμ(s)eξ+α2(t)¯sμ(s)eξ+α2(t)¯ξ+α1(t,¯st)˙¯ξ=μ(s)¯s. (4.11)

    Let us introduce a time-varying change of coordinates:

    ¯z(t)=¯s(t)+eξ(t)¯ξ(t).

    Then, by using the simplified notation, we obtain

    {˙¯z=D(¯zeξ¯ξ)μ(s)eξ+α2(t)(¯zeξ¯ξ)μ(s)eξ+α2(t)¯ξ+α1(t,¯st)+eξ˙ξ(t)¯ξ+eξμ(s)(¯zeξ¯ξ)˙¯ξ=μ(s)(¯zeξ¯ξ). (4.12)

    By recalling that ˙ξ(t)=μ(s(t))D, we obtain

    ˙¯z=D¯zμ(s)eξ+α2(t)(¯zeξ¯ξ)μ(s)eξ+α2(t)¯ξ+eξμ(s)¯ξ+eξμ(s)(¯zeξ¯ξ)+α1(t,¯st)

    and, as an immediate consequence, we obtain the equivalent version

    ˙¯z=D¯z+μ(s)eξ(1eα2(t))(¯zeξ¯ξ)+eξμ(s)(1eα2(t))¯ξ+α1(t,¯st).

    Then, by grouping the terms again, we obtain

    {˙¯z=D¯z+α3(t)¯z+α4(t)¯ξ+α1(t,¯st)˙¯ξ=μ(s)¯zμ(s)eξ¯ξ (4.13)

    with

    α3(t)=μ(s(t))eξ(t)(1eα2(t))

    and

    α4(t)=(μ(s(t))μ(s(t))eξ(t))eξ(t)(1eα2(t)).

    In order to study the stability of system (4.13), we adopt a Lyapunov approach. Let us introduce the positive definite function

    V1(t):=V1(¯ξ(t))=12¯ξ(t)2.

    Its derivative along the trajectories of (4.13) satisfies

    ˙V1(t)=μ(s)¯ξ¯zμ(s)eξ¯ξ2.

    Let us recall that Lemma 5 ensures that a(t)a_, which leads to

    ξ(t)ξ_

    since ξ_=ln(a_). Now, by (1.3) and the right estimation of (4.1) we know that μ(s)μ(0) and μ(s)μ(s0M). In addition, by using Young's inequality

    μ(0)|¯z||¯ξ|δ2¯ξ2+12δμ(0)2¯z2

    with δ=μ(S0M)eξ_, we deduce that

    ˙V1(t)μ(0)|¯ξ||¯z|μ(s0M)eξ_¯ξ212μ(s0M)eξ_¯ξ2+12μ(s0M)eξ_μ(0)2¯z2μ(s0M)eξ_¯ξ2=κ¯z212μ(s0M)eξ_¯ξ2,with κ=μ(0)22μ(s0M)eξ_. (4.14)

    Now, let us introduce the positive definite function

    V2(t):=V2(¯z(t))=1D¯z2(t).

    Its derivative along the trajectories of (4.13) satisfies

    ˙V2(t)2(1+max{0,α3(t)}D)¯z2+2|α4(t)|D|¯z¯ξ|+2D¯zα1(t,¯st).

    By applying Young's inequality to the products |α4(t)¯ξ|D|¯z| and |¯z||α1(t,¯st)| with δ=D2, we obtain

    ˙V2(t)(1+2max{0,α3(t)}D)¯z2+2(α4(t)D¯ξ)2+2α1(t,¯st)2D2. (4.15)

    Now, we consider the candidate Lyapunov function

    V3(t):=V3(¯z(t),¯ξ(t))=12κV1(¯ξ(t))+V2(¯z(t))

    which is well–defined because κ>0. From (4.14) and (4.15), we deduce that its derivative along the trajectories of system (4.13) satisfies

    ˙V3(t)μ(s0M)4κeξ_¯ξ2+(12+2max{0,α3(t)}D)¯z2+2α4(t)2D2¯ξ2+2α1(t,¯st)2D2. (4.16)

    We can estimate the function α2 defined by (4.9) as done in (3.4):

    |α2(t)|μ(s)τμmτfor all t0

    where the last inequality follows from (1.4). Notice that, as ts(t) is ω–periodic, the above estimation is valid for any t0, while in (3.4), is verified after a finite time.

    A direct consequence of the above estimation is

    |1eα2(t)|eτμm1for all t0 (4.17)

    which follows from the inequality |1ex|e|x|1 for any xR.

    The inequality (4.17) combined with (1.4) and μ(s)μ(0) for any s0 allows to estimate α3 and α4 as follows:

    |α3(t)|μ(0)eξ(t)(eτμm1)and|α4(t)|(μm+μ(0)eξ(t))eξ(t)(eτμm1).

    By using (4.1) we deduce that

    a(t)=eξ(t)43s0Mfor all t0

    which makes it possible to obtain new estimations

    |α3(t)|43μ(0)s0M(eτμm1)

    and

    |α4(t)|43(μm+43μ(0)s0M)s0M(eτμm1)

    for all t0. Moreover, using again (4.1), we obtain

    α1(t,¯st)2(43μms0Mμ(0)ttτ|¯s(r)|dr)2βτttτ¯s(r)2dr (4.18)

    with β=169(μms0Mμ(0))2, where the last inequality is a consequence of the Cauchy-Schwarz inequality.

    By gathering the above estimations for α1(t,¯st), α3(t), and α4(t), we obtain a sharper estimation for ˙V3:

    ˙V3(t)μ(s0M)4κeξ_¯ξ2+[12+2D43μ(0)s0M(eτμm1)]¯z2+2D2[43(μm+43μ(0)s0M)s0M(eτμm1)]2¯ξ2+2βτD2ttτ¯s(r)2dr=[μ(s0M)4κeξ_+λ1(eτμm1)2]¯ξ2+[12+λ2(eτμm1)]¯z2+2βτD2ttτ¯s(r)2dr

    with \lambda_1 = \frac{32}{9 D^2} \left[\left(\mu_{m} + \frac{4}{3} \mu'(0) s_{\mathcal{M}}^{0}\right) s_{\mathcal{M}}^{0}\right]^{2} and \lambda_2 = \frac{8}{3 D} \mu'(0) s_{\mathcal{M}}^{0} .

    Bearing in mind that \overline{z}(t) = \overline{s}(t) + e^{\xi_{\star}(t)} \overline{\xi}(t) , let us introduce a Lyapunov-Krasovskii functional

    \begin{equation} V_4(\overline{z}_t, \overline{\xi}_t) = V_3(\overline{z}(t), \overline{\xi}(t)) + \frac{2 \beta \tau}{D^2} \int_{t - \tau}^{t} \int_{\ell}^{t} \overline{s}(r)^2 \, dr\, d\ell. \end{equation} (4.19)

    A simple calculation gives

    \frac{d}{dt}\left(\int_{t-\tau}^{t}\int_{\ell}^{t}\overline{s}^{2}(r)\, dr \, d\ell\right) = \overline{s}^{2}(t)\tau-\int_{t-\tau}^{t}\overline{s}^{2}(r)\, dr

    which leads to

    \begin{array}{rcl} \dot{V}_4(t) & \leq & \left[- \frac{\mu'\left(s_{\mathcal{M}}^{0}\right) }{4\kappa} e^{\underline{\xi}} + \lambda_1 \left(e^{\tau \mu_{m}} - 1\right)^2\right] \overline{\xi}^2 + \left[- \frac{1}{2} + \lambda_2\left(e^{\tau \mu_{m}} - 1\right)\right] \overline{z}^2 \\\\ && + \frac{2 \beta \tau^2}{D^2} \overline{s}^2. \end{array}

    By using the definition of \overline{s} = \overline{z}-e^{\xi_{\star}(t)}\overline{\xi} and noticing that \overline{s}^2\leq 2(\overline{z}^2+e^{2\xi_{\star}(t)}\overline{\xi}^2) , we obtain

    \begin{array}{rcl} \dot{V}_4(t) & \leq & \left[- \frac{\mu'\left(s_{\mathcal{M}}^{0}\right) }{4\kappa} e^{\underline{\xi}} + \lambda_1 \left(e^{\tau \mu_{m}} - 1\right)^2\right] \overline{\xi}^2 + \left[- \frac{1}{2} + \lambda_2\left(e^{\tau \mu_{m}} - 1\right)\right] \overline{z}^2 \\\\ && + \frac{4 \beta \tau^2}{D^2} \left(\overline{z}^2 + e^{2\xi_{\star}(t)} \overline{\xi}^2\right) \\\\ & \leq & \left[- \frac{\mu'\left(s_{\mathcal{M}}^{0}\right) }{4\kappa} e^{\underline{\xi}} + \lambda_1 \left(e^{\tau \mu_{m}} - 1\right)^2 + \frac{64 \beta \tau^2}{9 D^2} (s_{\mathcal{M}}^{0})^2\right] \overline{\xi}^2 + \left[- \frac{1}{2} + \lambda_2\left(e^{\tau \mu_{m}} - 1\right) + \frac{4 \beta \tau^2}{D^2} \right] \overline{z}^2 \end{array}

    where the last inequality uses (4.1) with a_{\star}(t) = e^{\xi_{\star}}(t)\leq \frac{4}{3}s_{\mathcal{M}}^{0} .

    Employing the inequality \lambda_2\left(e^{\tau \mu_{m}} - 1\right) \leq \frac{1}{4} + \lambda_2^2\left(e^{\tau \mu_{m}} - 1\right)^2 , we obtain the estimation

    \begin{array}{rcl} \dot{V}_4(t) & \leq & \left[- \frac{\mu'\left(s_{\mathcal{M}}^{0}\right) }{4\kappa} e^{\underline{\xi}} + \lambda_1 \left(e^{\tau \mu_{m}} - 1\right)^2 + \frac{64 \beta \tau^2}{9 D^2} (s_{\mathcal{M}}^0)^2\right] \overline{\xi}^2 + \left[- \frac{1}{4} + \lambda_2^2\left(e^{\tau \mu_{m}} - 1\right)^2 + \frac{4 \beta}{D^2} \tau^2\right] \overline{z}^2. \end{array}

    Given that \tau \mu_{m} \leq e^{\mu_{m}\tau} - 1 holds for all \tau \geq 0 , and the quadratic function is monotonically increasing on [0, + \infty) , it follows that

    \begin{array}{rcl} \dot{V}_4(t) & \leq & \left[- \frac{\mu'\left(s_{\mathcal{M}}^{0}\right) }{4\kappa} e^{\underline{\xi}} + \lambda_1 \left(e^{\tau \mu_{m}} - 1\right)^2 + \frac{64 \beta}{9 (D\mu_{m})^2} \left(e^{\tau \mu_{m}} - 1\right)^2 (s_{\mathcal{M}}^0)^2\right] \overline{\xi}^2 \\\\ && + \left[- \frac{1}{4} + \lambda_2^2\left(e^{\tau \mu_{m}} - 1\right)^2 + \frac{4 \beta}{(D \mu_{m})^2} \left(e^{\tau \mu_{m}} - 1\right)^2\right] \overline{z}^2 \\\\ & = & \underbrace{\left[- \frac{\mu'\left(s_{\mathcal{M}}^0\right)}{4\kappa} e^{\underline{\xi}} + \left(\lambda_1 + \frac{64 \beta}{9 (D\mu_{m})^2} (s_{\mathcal{M}}^{0})^2 \right) \left(e^{\tau\mu_{m}} - 1\right)^2 \right]}_{ = \zeta_{1}} \overline{\xi}^2 \\\\ &&+ \underbrace{\left[- \frac{1}{4} + \left(\lambda_2^2 + \frac{4 \beta}{(D\mu_{m})^2}\right) \left(e^{\tau \mu_{m}} - 1\right)^2\right]}_{ = \zeta_{2}} \overline{z}^2. \end{array}

    The above inequality implies that, if \tau \in (0, \tau_{0}^{*}) with \tau_{0}^{*} given by

    \tau_{0}^{*} = \frac{1}{\mu_{m}}\min\left\{ \ln\left(1 + \sqrt{\frac{\mu'\left(s_{\mathcal{M}}^{0}\right)}{4\kappa \left(\lambda_1 + \frac{64 \beta}{9 (D\mu_{m})^2} (s_{\mathcal{M}}^{0})^2 \right)} e^{\underline{\xi}} } \right), \ln\left(1 + \frac{1}{2 \sqrt{\lambda_2^2 + \frac{4 \beta}{(D\mu_{m})^2}}}\right) \right\}

    then \zeta_{1} < 0 and \zeta_{2} < 0 . By integrating the above inequality between 0 and t \geq 0 , we obtain

    |\zeta_{1}| \int_{0}^{t}\overline{\xi}^{2}(r)\, dr + |\zeta_{2}|\int_{0}^{t}\overline{z}^{2}(r)\, dr \leq V_{4}(0) - V_{4}(t) \leq V_{4}(0)

    which implies that \overline{\xi} and \overline{z} belong to L^{2}(\mathbb{R}_{0}^{+}) .

    In addition, as \overline{\xi} and \overline{z} , as well as their derivatives, are bounded on [0, + \infty) , it is easy to deduce that \overline{\xi}^{2} and \overline{z}^{2} are uniformly continuous. Then, Barbalat's Lemma [26] ensures that

    \lim\limits_{t \to +\infty}\overline{\xi}(t) = \lim\limits_{t\to +\infty}\overline{z}(t) = 0 .

    Thus, the equilibrium (\overline{z}, \overline{\xi}) = (0, 0) is a globally asymptotically stable solution of system (4.12).

    Finally, from our analysis, we can deduce that if \tau \in (0, \tau_{\ddagger}) with

    \begin{equation} \tau_{\ddagger} = \min\left\{\frac{1}{D - \mu_{m}}\ln\left(1 - \frac{D}{8 \overline{\mu}}\right), \tau_{0}^{*}\right\} \end{equation} (4.20)

    then the equilibrium (\overline{s}, \overline{\xi}) = (0, 0) is a globally exponentially stable solution of system (4.7). This concludes the proof.

    Let us consider the family of systems (2.6) with input substrate s^{0}(t) = S_{c}+\varepsilon\gamma(t) , namely

    \begin{equation} \left\{ \begin{array}{rcl} \dot{s}(t)& = & D[S_c+\varepsilon \gamma(t)]-Ds(t)-\mu(s(t))e^{\int_{t-\tau}^t[D-\mu(s(r))]dr}a(t) \\ \dot{a}(t)& = & a(t)[\mu(s(t))-D] \end{array} \right. \end{equation} (5.1)

    where \gamma(\cdot) is a function of class C^{1} and \omega –periodic such that \max\limits_{t\in [0, \omega]}|\gamma(t)|\leq \overline{\gamma} . In addition, the positive constants S_{c} and \varepsilon are such that \varepsilon \in (0, S_{c}\overline{\gamma}^{-1}) , which implies the positiveness of the washout solution. Finally, by Proposition 1, it follows that system (5.1) admits an \omega –periodic solution (s_{\varepsilon}, a_{\varepsilon}) .

    The goal of this subsection is, roughly speaking, to prove that \dot{s}_{\varepsilon}(t) and \dot{a}_{\varepsilon}(t) are small when the parameters \varepsilon and \tau are small.

    Lemma 7. There is a delay \tau_{\sharp} > 0 such that if \tau \in [0, \tau_{\sharp}] and \varepsilon is small enough, then the derivatives of s_{\varepsilon}(\cdot) and a_{\varepsilon}(\cdot) satisfy

    \begin{equation} \max\limits_{t\in [0, \omega]}|\dot{s}_{\varepsilon}(t)|\leq L_{1}\varepsilon \quad and \quad \max\limits_{t\in [0, \omega]}|\dot{a}_{\varepsilon}(t)|\leq L_{2}\varepsilon, \end{equation} (5.2)

    for some positive constants L_{1} and L_{2} .

    Proof. First, notice that since s_{\varepsilon}(\cdot) and a_{\varepsilon}(\cdot) are \omega –periodic, it follows that their derivatives are also \omega –periodic. Now, let us introduce the auxiliary function

    \begin{equation} \phi_{\varepsilon}(t) = \varepsilon D\gamma(t)+\mu(s_{\varepsilon}(t))[1-e^{\int_{t-\tau}^t[D-\mu(s_{\varepsilon}(r))]dr}]a_{\varepsilon}(t). \end{equation} (5.3)

    By using Lemma 2 combined with the inequality |1-e^{x}|\leq e^{|x|}-1 and inequalities (2.3) and (3.4), we can see that

    \begin{equation} \begin{array}{rcl} |\phi_{\varepsilon}(t)|&\leq & \varepsilon D|\gamma(t)|+|\mu(s_{\varepsilon}(t))|\, \left|1-e^{\int_{t-\tau}^t[D-\mu(s_{\varepsilon}(r))]dr}\right | \, |a_{\varepsilon}(t)| \\ & \leq &\varepsilon D\overline{\gamma}+\frac{5}{3}\mu(s_{\triangle})s_{\triangle}[e^{\tau\mu(s_{\triangle})}-1] \end{array} \end{equation} (5.4)

    for t\geq t_{\sharp}+\tau . Then, when \tau is sufficiently small, it follows that

    \begin{equation} |\phi_{\varepsilon}(t)|\leq 2\varepsilon D\overline{\gamma}. \end{equation} (5.5)

    By introducing the change of variables z_{\varepsilon} = s_{\varepsilon} + a_{\varepsilon} and considering (5.1), we have

    \begin{equation} \begin{array}{rcl} \dot{z}_{\varepsilon} & = &D[S_c+\varepsilon \gamma(t)]-Ds_{\varepsilon}(t)-\mu(s_{\varepsilon}(t))e^{\int_{t-\tau}^t[D-\mu(s_{\varepsilon}(r))]dr}a_{\varepsilon}(t) \\ && + a_{\varepsilon}(t)[\mu(s_{\varepsilon}(t))-D]\\ & = &DS_c+\varepsilon D\gamma(t)-Dz_{\varepsilon}(t)+\mu(s_{\varepsilon}(t))\left[1-e^{\int_{t-\tau}^t[D-\mu(s_{\varepsilon}(r))]dr}\right]a_{\varepsilon}(t) \\ & = &DS_c+\phi_{\varepsilon}(t)-Dz_{\varepsilon}(t) \end{array} \end{equation} (5.6)

    where \phi_{\varepsilon} is defined in (5.3). By using the fact that z_{\varepsilon} is \omega –periodic, we can prove that

    \begin{equation} z_{\varepsilon}(t) = \dfrac{1}{1-e^{-D\omega}}\int_{t-\omega}^{t}e^{D(r-t)}[DS_c+\phi_{\varepsilon}(r)]\, dr \end{equation} (5.7)

    and (5.6) can be written as

    \begin{array}{rcl} \dot{z}_{\varepsilon} & = & DS_c+\phi_{\varepsilon}(t)-\dfrac{D}{1-e^{-D\omega}}\int_{t-\omega}^{t}e^{D(r-t)}[DS_c+\phi_{\varepsilon}(r)]\, dr \\ & = & \phi_{\varepsilon}(t)-\dfrac{D}{1-e^{-D\omega}}\int_{t-\omega}^{t}e^{D(r-t)}\phi_{\varepsilon}(r)\, dr \end{array}

    which, combined with (5.5), allows us to deduce that

    \begin{equation} |\dot{z}_{\varepsilon}(t)|\leq 2|\phi_{\varepsilon}|_{\infty} \leq 4\varepsilon D \overline{\gamma} \quad \text{for any}\ t\in \mathbb{R} . \end{equation} (5.8)

    Now, the change of variables \xi_{\varepsilon} = \ln(a_{\varepsilon}) gives

    \dot{\xi}_{\varepsilon} = \mu(s_{\varepsilon}(t))-D \quad \text{and} \quad \ddot{\xi}_{\varepsilon} = \mu'(s_{\varepsilon}(t))\dot{s}_{\varepsilon}(t) .

    Then, the Eq (5.1) combined with the identity

    \begin{equation} \dot{a}_{\varepsilon}(t) = a_{\varepsilon}(t)\dot{\xi}_{\varepsilon}(t) \end{equation} (5.9)

    and the definition of z_{\varepsilon} allows us to conclude that

    \begin{array}{rcl} \ddot{\xi}_{\varepsilon}& = &\mu'(s_{\varepsilon}(t))[\dot{z}_{\varepsilon}(t)-\dot{a}_{\varepsilon}(t)] = -\mu'(s_{\varepsilon}(t))a_{\varepsilon}(t)\dot{\xi}_{\varepsilon}(t)+\mu'(s_{\varepsilon}(t))\dot{z}_{\varepsilon}(t). \end{array}

    It is interesting to observe that \dot{\xi}_{\varepsilon} is an \omega –periodic solution of the \omega –periodic time-varying equation

    \dot{u} = -\mu'(s_{\varepsilon}(t))a_{\varepsilon}(t)u+\mu'(s_{\epsilon}(t))\dot{z}_{\epsilon}(t).

    Moreover, since \int_{0}^{\omega}\mu'(s_{\epsilon}(\ell))a_{\epsilon}(\ell)d\ell > 0 , we deduce similarly as in (5.7), the identity

    \dot{\xi}_{\varepsilon}(t) = \dfrac{1}{1-e^{-\int_{0}^{\omega}\mu'(s_{\varepsilon}(\ell))a_{\epsilon}(\ell)d\ell}} \int_{t-\omega}^{t}e^{-\int_r^{t}\mu'(s_{\varepsilon}(\ell))a_{\epsilon}(\ell)d\ell}\mu'(s_{\varepsilon}(r))\dot{z}_{\varepsilon}(r)\, dr.

    Now, by using \mu''(s) < 0 for any s\geq 0 , Lemma 1 implies the inequality

    |\dot{\xi}_{\varepsilon}(t)|\leq \dfrac{\mu'(0)}{1-e^{-\mu'(s_{\triangle})\int_{0}^{\omega}a_{\varepsilon}(\ell)d\ell}} \int_{t-\omega}^{t}|\dot{z}_{\varepsilon}(r)|\, dr

    when t is sufficiently large. Now we recall that, as a consequence of Theorem 1, there are constants \overline{s} > \underline{s} > 0 , \overline{\xi} > \underline{\xi} , \tau_s > 0 , and \varepsilon_s > 0 such that for \tau\in [0, \tau_s] and \varepsilon\in [0, \varepsilon_s] , when t is sufficiently large, it follows that

    \begin{equation} \underline{s}\leq s_{\varepsilon}(t)\leq \overline{s}\quad \text{and} \quad \underline{\xi}\leq \xi_{\varepsilon}(t)\leq \overline{\xi} \end{equation} (5.10)

    which leads to

    |\dot{\xi}_{\varepsilon}(t)|\leq \dfrac{\mu'(0)}{1-e^{-\mu'(s_{\triangle})\, \omega\, e^{\underline{\xi}}}} \int_{t-\omega}^{t}|\dot{z}_{\varepsilon}(r)|\, dr.

    As an immediate consequence of (5.8), we have that

    |\dot{\xi}_{\varepsilon}(t)|\leq \dfrac{2\mu'(0)\omega}{1-e^{-\mu'(s_{\triangle})\, \omega\, e^{\underline{\xi}}}} |\phi_{\varepsilon}|_{\infty}.

    By (5.9) and the above estimation, we can conclude

    |\dot{a}_{\varepsilon}(t)|\leq |a_{\varepsilon}(t)||\dot{\xi}_{\varepsilon}(t)| \quad \text{with} \quad |\dot{\xi}_{\varepsilon}(t)|\leq \dfrac{4\mu'(0)\omega}{1-e^{-\mu'(s_{\triangle})\omega e^{\underline{\xi}}}} D\overline{\gamma}\varepsilon.

    Then, inequalities (5.2) can be deduced by using (5.10) combined with the definitions of z_{\varepsilon} and the above estimations.

    By using the transformation \xi = \ln(a) , system (1.1) becomes

    \begin{equation} \left\{ \begin{array}{rcl} \dot{s}(t)& = & Ds^0(t)-Ds(t)-\mu(s(t))e^{\xi(t)} \\ \dot{\xi}(t)& = & \mu(s(t-\tau))-D. \end{array} \right. \end{equation} (5.11)

    Now, let us introduce the variables

    \begin{equation} \tilde{s}(t) = s(t)-s_{\star}(t), \quad \tilde{\xi}(t) = \xi(t)-\xi_{\star}(t) \end{equation} (5.12)

    where (s_{\star}, \xi_{\star}) denotes a periodic solution of (5.11). Then, by using (1.5) and (5.12), we can deduce the new representation

    \begin{equation} \left\{ \begin{array}{rcl} \dot{\tilde{s}}(t)& = & -D\tilde{s}(t)-\frac{k_s\mu_m\tilde{s}(t)e^{\xi(t)}}{(k_s+s_{\star}(t))(k_s+s(t))}+\mu(s_{\star}(t))e^{\xi_{\star}(t)}[1-e^{\tilde{\xi}(t)}] \\ \dot{\tilde{\xi}}(t)& = & k_s\mu_m\frac{\tilde{s}(t-\tau)}{(k_s+s_{\star}(t-\tau))(k_s+s(t-\tau))}. \end{array} \right. \end{equation} (5.13)

    It is useful to point out that \tilde{s} is bounded on [0, +\infty) since it is a difference of two bounded functions. This implies that \dot{\tilde{s}} is bounded on [0, +\infty) . Consequently, \tilde{s} is uniformly continuous on the same interval.

    Let us define

    V_1(\tilde{\xi}(t)) = e^{\tilde{\xi}(t)}-\tilde{\xi}(t)-1

    and note that V_{1}(\cdot) is nonnegative with V_{1}(0) = 0 . Moreover, we have

    \dot{V}_{1}(\tilde{\xi}(t+\tau)): = \frac{d}{dt}V_1(\tilde{\xi}(t+\tau)) = k_s\mu_m\frac{\tilde{s}(t)}{(k_s+s_{\star}(t))(k_s+s(t))}[e^{\tilde{\xi}(t+\tau)}-1].

    On the other hand, by (5.13) we obtain the identity

    \tilde{\xi}(t+\tau) = \tilde{\xi}(t)+ k_s\mu_m\int_{t-\tau}^{t}\frac{\tilde{s}(r)}{(k_s+s_{\star}(r))(k_s+s(r))}dr .

    Hence, it follows that

    \dot{V}_1(\tilde{\xi}(t+\tau)) = \frac{k_s\mu_m\tilde{s}(t)}{(k_s+s_{\star}(t))(k_s+s(t))}[e^{\tilde{\xi}(t)+ \int_{t-\tau}^{t}\frac{k_s\mu_m\tilde{s}(r)}{(k_s+s_{\star}(r))(k_s+s(r))}dr}-1].

    Now, we define

    V_2(t, \tilde{s}(t)) = \frac{k_s\mu_m}{\mu(s_{\star}(t))e^{\xi_{\star}(t)}} \int_0^{\tilde{s}(t)} \frac{r}{(k_s+s_{\star}(t))(k_s+s_{\star}(t)+r)}dr.

    Notice that, by (1.3) and \xi_{\star}(t) = \ln(x_{\star}(t)) , we can deduce

    V_2(t, \tilde{s}(t)) = \frac{k_s}{x_{\star}(t)s_{\star}(t)} \int_0^{\tilde{s}(t)} \frac{r}{(k_s+s_{\star}(t)+r)}dr.

    It is easy to show that V_{2}(t, \tilde{s}(t)) is nonnegative. Moreover, its derivative along the trajectories of system (5.13) satisfies

    \dot{V}_2(t, \tilde{s}(t)): = \frac{d}{dt}V_2(t, \tilde{s}(t)) = \frac{k_s}{s_{\star}(t)e^{\xi_{\star}(t)}}\frac{\tilde{s}(t)}{(k_s+s_{\star}(t)+\tilde{s}(t))}\dot{\tilde{s}}(t)+\kappa(t)

    with \kappa(t) defined by

    \begin{array}{rcl} \kappa(t) & = & -\frac{k_s \dot{s}_{\star}(t)}{x_{\star}(t)s_{\star}(t)} \int_0^{\tilde{s}(t)} \frac{r}{(k_s+s_{\star}(t)+r)^2}\, dr \\\\ && -\frac{k_s\left(\dot{x}_{\star}(t)s_{\star}(t)+x_{\star}(t)\dot{s}_{\star}(t)\right)}{(x_{\star}(t)s_{\star}(t))^2} \int_0^{\tilde{s}(t)} \frac{r}{(k_s+s_{\star}(t)+r)}dr . \end{array}

    From (5.13), it follows that

    \begin{array}{rcl} \dot{V}_2(t, \tilde{s}(t))& = & -\frac{k_s}{s_{\star}(t)e^{\xi_{\star}(t)}}\frac{\tilde{s}(t)}{(k_s+s_{\star}(t)+\tilde{s}(t))}\left[D\tilde{s}(t)+\frac{k_s\mu_m\tilde{s}(t)e^{\xi(t)}}{(k_s+s_{\star}(t))(k_s+s(t))}-\mu(s_{\star}(t))e^{\xi_{\star}(t)}[1-e^{\tilde{\xi}(t)}] \right]\\\\ && +\kappa(t)\\ & = & -\frac{k_s}{s_{\star}(t)e^{\xi_{\star}(t)}}\frac{1}{(k_s+s_{\star}(t)+\tilde{s}(t))}\left[D+k_s\mu_m\frac{e^{\xi(t)}}{(k_s+s_{\star}(t))(k_s+s(t))}\right]\tilde{s}(t)^2 \\ & & +\frac{k_s \mu_m \tilde{s}(t)}{(k_s+s_{\star}(t))(k_s+s_{\star}(t)+\tilde{s}(t))}[1-e^{\tilde{\xi}(t)}] +\kappa(t). \end{array}

    Now, we introduce

    V_{3}(t): = V_3(t, \tilde{s}(t), \tilde{\xi}(t+\tau)) = V_1(\tilde{\xi}(t+\tau))+V_2(t, \tilde{s}(t)) .

    This function is nonnegative and

    \begin{array}{rcl} \dot{V}_3(t)& = & -\frac{k_s}{s_{\star}(t)e^{\xi_{\star}(t)}}\frac{1}{(k_s+s_{\star}(t)+\tilde{s}(t))}\left[D+k_s\mu_m\frac{e^{\xi(t)}}{(k_s+s_{\star}(t))(k_s+s(t))}\right]\tilde{s}(t)^2\\ & & +\underbrace{\frac{k_s \mu_m \tilde{s}(t)}{(k_s+s_{\star}(t))(k_s+s(t))} \left[e^{k_s \mu_m \int_{t-\tau}^{t}\frac{\tilde{s}(r)}{(k_s+s_{\star}(r))(k_s+s(r))}dr}-1\right]e^{\tilde{\xi}(t)}}_{ = \dot{V}_{31}(t)}+\kappa(t). \end{array}

    Now, we have

    |\kappa(t)|\leq \frac{1}{2}\left[ \frac{|\dot{s}_{\star}(t)|}{k_s \underline{x_{\star}s_{\star}}}+\frac{(|\dot{x}_{\star}(t)|s_{\star}(t)+x_{\star}(t)|\dot{s}_{\star}(t)|)}{(\underline{x_{\star}s_{\star}})^2} \right]\tilde{s}(t)^2.

    In addition, by using the mean value theorem combined with Young and Jensen's inequalities, we have that

    \begin{array}{rcl} \dot{V}_{31}(t) & = & \frac{k_s \mu_m \tilde{s}(t)e^{\tilde{\xi}(t)}} {(k_s+s_{\star}(t))(k_s+s(t))} \left[k_{s}\mu_{m}\int_{t-\tau}^{t}\frac{\tilde{s}(r)}{(k_s+s_{\star}(r))(k_s+s(r))}\, dr \, e^{\eta}\right] \\\\ &\leq& \frac{\mu_m \tilde{s}(t)e^{\tilde{\xi}(t)}} {k_s} \left[\frac{\mu_{m}}{k_{s}} \int_{t-\tau}^{t}\tilde{s}(r)\, dr \, e^{\eta}\right] \\\\ &\leq& \left(\frac{\mu_m \tilde{s}(t)e^{\tilde{\xi}(t)}} {k_s}\right)^{2}\frac{\tau}{2} +\frac{1}{2\tau}\left(\frac{\mu_{m}e^{\eta}}{k_{s}}\right)^{2} \left[ \int_{t-\tau}^{t}\tilde{s}(r)\, dr\right]^{2} \\\\ & = & \left(\frac{\mu_m \tilde{s}(t)e^{\tilde{\xi}(t)}} {k_s}\right)^{2}\frac{\tau}{2} +\frac{\tau}{2}\left(\frac{\mu_{m}e^{\eta}}{k_{s}}\right)^{2} \left[\frac{1}{\tau} \int_{t-\tau}^{t}\tilde{s}(r)\, dr\right]^{2} \\\\ &\leq &\left(\frac{\mu_m \tilde{s}(t)e^{\tilde{\xi}(t)}} {k_s}\right)^{2}\frac{\tau}{2} +\frac{\tau}{2}\left(\frac{\mu_{m}e^{\eta}}{k_{s}}\right)^{2} \int_{t-\tau}^{t}\tilde{s}^{2}(r)\, dr \end{array}

    where \eta is a number between 0 and k_{s}\mu_{m}\int_{t-\tau}^{t}\frac{\widetilde{s}(r)}{(k_{s}+s_{\star}(r))(k_{s}+s(r))}\, dr .

    Now, we deduce that, when |\dot{s}_{\star}(t)| and |\dot{x}_{\star}(t)| are sufficiently small, we can find constants \nu_1 > 0 , \nu_2\geq 0 , and \nu_3\geq 0 such that

    \dot{V}_3(t)\leq -\nu_1 \tilde{s}(t)^2+\tau \nu_2\int_{t-\tau}^t\tilde{s}(\ell)^2 d\ell +\nu_3\tau\tilde{s}(t)^2 .

    Thus, when \tau is sufficiently small, by Lemma 7 we obtain the inequality

    \begin{equation} \dot{V}_3(t)\leq -\frac{3\nu_1}{4} \tilde{s}(t)^2+\tau \nu_2\int_{t-\tau}^t\tilde{s}(\ell)^2 d\ell. \end{equation} (5.14)

    Let us introduce the functional

    V_{4}(t): = V_4(t, \tilde{s}_{t}, \tilde{\xi}(t)) = V_3(t, \tilde{s}(t), \tilde{\xi}(t))+\tau \nu_2 \int_{t-\tau}^{t}\int_{r}^{t}\tilde{s}(\ell) d\ell\, dr.

    Since t\mapsto\tilde{s}(t) is bounded, there exists L > 0 such that \tilde{s}(t) > -L for any t\geq 0 . It follows that

    \int_{t-\tau}^{t}\int_{r}^{t}\tilde{s}(\ell) d\ell\, dr\geq -\frac{L}{2}\tau^{2}

    and the nonnegativeness of V_{3}(\cdot) imply that V_{4}(\cdot) is lower bounded.

    By using (5.14) combined with the integral term of V_{4}(\cdot) , and considering \tau small enough, we can see that \dot{V}_4(t) = \frac{d}{dt}V_4(t, \tilde{s}_{t}, \tilde{\xi}(t)) satisfies

    \dot{V}_4(t) \leq -\frac{3\nu_1}{4} \tilde{s}(t)^2+\tau^2 \nu_2\tilde{s}(t)^2.

    Thus, when \tau is sufficiently small,

    \dot{V}_4(t)\leq -\frac{\nu_1}{2} \tilde{s}(t)^2

    which implies that t\mapsto V_{4}(t) is not increasing. By integrating the above inequality, we obtain

    \frac{\nu_1}{2} \int_0^t \tilde{s}(m)^2 dm \leq V_4(0)-V_{4}(t)\leq V_{4}(0) \quad \text{for any}\ t\geq 0 ,

    which implies that \tilde{s}\in L^{2}(\mathbb{R}_{0}^{+}) .

    Since \tilde{s} and \dot{\tilde{s}} are bounded on [0, + \infty) , it follows that \tilde{s}^{2} is uniformly continuous on [0, +\infty) . By Barbalat's Lemma [26], it follows that \lim\limits_{t\to+\infty}\tilde{s}(t) = 0 .

    In addition, by considering the uniform continuity of \tilde{s} , we deduce again from Barbalat's lemma that \lim\limits_{t\to+\infty}\dot{\tilde{s}}(t) = 0 . Finally, by studying the \tilde{s} -dynamics from (5.13) and recalling that \lim\limits_{t\to+\infty}\tilde{s}(t) = 0 , we can deduce that \lim\limits_{t\to+\infty}\tilde{\xi}(t) = 0 and the result follows.

    In order to illustrate our results and their applications, we will examine the culture of the microalgae Dunaliella tertiolecta by considering nitrate as the limiting nutrient. We will revisit the numerical simulations carried out in [1] for the solutions of (1.1) where the supply of nitrate into the chemostat is described by the \omega –periodic function

    s^0(t) = S_c+\varepsilon\cos\left(\frac{2\pi t}{\omega}\right).

    As pointed out in [29,30], considering that the parameters s^{0} and/or D vary periodically over time allows us to mimic environmental oscillations in aquatic ecosystems and study the effects of these variations in the microalgae physiology.

    As mentioned in [1], the numerical simulations rely on a previous work of Vatcheva, Bernard, and Mars [31] which analyzes the merits and drawbacks of several alternative mathematical models based on a set of available experimental data from [32] describing the growth of the microalgae under a limitation of nitrate. In particular, when considering a growth described by Monod's function (1.3), the parameters involved in the model (1.1) are summarized in Table 1.

    Table 1.  Dunaliella tertiolecta culture conditions growing with nitrate as a limiting substrate.
    Parameter Values Units Meaning Source
    \mu_{\max} [1.2\, , \, 1.6] d^{-1} maximum growth rate [31, p.491]
    k_s [0.01\, , \, 0.2] \mu Mol L^{-1} half-saturation constant [31, p.491]
    D arbitrary d^{-1} dilution rate [31, p.491]
    S_c [80\, , \, 120] \mu Mol L^{-1} average input nutrient concentration [31, p.491]

     | Show Table
    DownLoad: CSV

    We carried out our numerical simulations* in the R software version 4.0.2 by using the libraries PBSDDESOLVE, GGPLOT2, RESHAPE2, LATEX2EXP, and GRIDEXTRA to build our graphics.

    *The code is available in https://github.com/Oehninger/Time-varying-chemostat-stability

    Figures 1 and 2 illustrate numerical approximations for some solutions of system (1.1) by considering a set of initial conditions with colors black, red, and green, respectively, and different delays. Furthermore, in each case, the left graph illustrates the global dynamics (with transient phase) while the asymptotic dynamics (without transient phase) are illustrated by the right graph. The parameters considered are

    S_{c} = 90, \quad \varepsilon = 10, \quad \mu_{\max} = 1.4, \quad k_{s} = 0.1, \quad \text{and} \quad D = 1.

    One can prove that

    \Sigma^{*}(t) = S_{c}+\frac{25 D\varepsilon}{4\pi^{2}+25D^{2}}\left\{\frac{2\pi}{5}\sin\left(\frac{2\pi t}{5}\right)+D \cos\left(\frac{2\pi t}{5}\right)\right\} .

    Moreover, a numerical integration by Riemann sums shows that condition (1.10) is satisfied since

    \frac{1}{5}\int_{0}^{5}\mu(\Sigma^{*}(t))\, dt\approx 1.305277.
    Figure 1.  Dynamics of the total biomass and substrate. The initial conditions are: (s, x) = (20, 3) (black curve), (s, x) = (35, 1.5) (red curve) and (s, x) = (50, 0.8) (green curve).
    Figure 2.  Dynamics of the total biomass and substrate. The initial conditions are: (s, x) = (20, 3) (black curve), (s, x) = (35, 1.5) (red curve) and (s, x) = (50, 0.8) (green curve).

    Figures 1 and 2 give the numerical solution for three initial conditions (black, green, and red) by considering several delays from \tau = 1.0\, d to \tau = 2.0\, d . The left side of each sub-figure presents the graphic of the solutions (nitrate and phytoplankton) for the first 50 days while the right side only considers from the 25th to the 50th days in order to have a better description of the asymptotic behavior.

    Figure 1ac illustrates our theoretical results since the convergence towards an \omega –periodic solution is observed, even for delays not so small.

    Figure 2ac considers bigger delays, leading to a longer transient phase. It seems that the solutions are not convergent to an \omega –periodic solution and this issue will deserve more attention in a future work.

    In this subsection, we will consider the parameters

    S_{c} = 90, \quad \varepsilon = 10, \quad \mu_{\max} = 1.4, \quad k_{s} = 0.1, \quad \text{and} \quad \tau = 1,

    combined with several values for the dilution rate D .

    The numerical simulations show that, if we increase D but the average condition (1.10) is still verified, the global attractivity is preserved but the convergence towards the periodic solution is slower.

    Figures 3 and 4 illustrate numerical approximations for some solutions of system (1.1) by considering a set of initial conditions with colors black, blue, and orange, respectively, and different dilution rates D taking values between 1.25 d^{-1} and 1.45 d^{-1} . As before, in each case, the left graph illustrates the global dynamics (with transient phase) while the asymptotic dynamics (without transient phase) are illustrated by the right graph.

    Figure 3.  Dynamics of the total biomass and substrate. The initial conditions are: (s, x) = (20, 3) (black curve), (s, x) = (35, 1.5) (blue curve) and (s, x) = (50, 0.8) (orange curve).
    Figure 4.  Dynamics of the total biomass and substrate. The initial conditions are: (s, x) = (20, 3) (black curve), (s, x) = (35, 1.5) (blue curve) and (s, x) = (50, 0.8) (orange curve).

    Figure 3ac illustrates the global attractivity of the unique \omega –periodic positive solution, and how the rate of convergence is decreasing as D increases. In order to contextualize these simulations, we can refer to [1] where the chemostat can be seen as a device to produce microbial biomass and the stabilization of biomass is strongly dependent of the dilution rate. In fact, note that an increasing of the dilution rate from 1.25 d^{-1} to 1.35 d^{-1} led to stabilization times from 30 to 80 days.

    Figure 4a, b considers bigger dilution rates such that the average assumption (1.10) does not hold, and therefore the washout solution is globally stable and the microbial biomass is driven to extinction.

    We have studied the \omega –periodic system of delay differential equations (1.1), which describes the dynamics of a limiting substrate and microbial biomass in a well–stirred chemostat. This system was considered in [1], where a result of existence and uniqueness of a nontrivial \omega –periodic solution is proved for small enough delays. In this article, we showed that this \omega –periodic solution is, in fact, globally asymptotically stable for Monod's uptake functions (1.3) and \omega –periodic inputs s^{0}(\cdot) of type (2.11).

    Although Monod's uptake functions are ubiquitous in chemostat modeling, an interesting open question is whether our global asymptotic stability results can be extended to general uptake functions having the same qualitative properties as Monod ones, namely, smoothness, fixed point at the origin, and the fact that they are increasing. This problem is certainly elusive, but the mathematical study of chemostat models shows a clear trend: a myriad of results initially obtained for Monod's uptake functions have subsequently been extended to wider families of uptake functions. We expect to cope with this problem in future research.

    Numerical simulations suggest the possible existence of a threshold for the delay ensuring asymptotic stability. This shows some conservativeness of our results, which is due the existence of several technical conditions involved in the estimations deduced in the construction of the Lyapunov–like function. This is also true for the proof of Theorem 2 and it will be extremely interesting to study the stability of the linear periodic delayed system (4.10) by alternative ways such as Floquet's theory and Bohl exponents.

    Last but not least, another important remark is that our delayed chemostat model is inspired by the autonomous system (1.6). Nevertheless, current research has given priority in consideration to another delayed model: the one introduced by Ellermeyer and Freedman in [33,34,35,36]. We also have been working on this approach and, in [37], we recently proposed a version of the Ellermeyer & Freedman delayed model adapted to the periodic case described by the differential delay equation

    \begin{equation} \left\{\begin{array}{rclcl} \dot{s}(t) & = & D(t)s^{0}(t)-D(t)s(t)-\gamma^{-1}\mu(s(t))x(t) & \text{if} & t\geq 0\\\\ \dot{x}(t) & = &e^{-\int_{t-\tau}^{t}D(s)}\mu(s(t-\tau))x(t-\tau)-D(t)x(t) & \text{if} & t\geq 0 \end{array}\right. \end{equation} (7.1)

    which recovers the model presented in [33,34,35,36] when s^{0}(\cdot) and D(\cdot) are positive constants. It is important to emphasize that the existence result was obtained by following a completely different method to the one carried out in [1], while the uniqueness for small delays followed along similar lines. Additional results for (7.1) have been obtained in [38], and a stochastic version has been considered in [39]. We expect to emulate the ideas presented in this work and to start a stability study for (7.1).

    The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

    This article is partially supported by MATHAMSUD regional cooperation program (Grant MATH2020006). The second author is supported by FONDECYT Regular 1210733.

    The authors declare there is no conflict of interest.

    The lemma is composed by four statements whose proof will be divided into several steps: Step 1 is devoted to preliminary estimations leading to the proof of statement a). The proof of statement b) will be a consequence of Steps 2–4. Finally, statements c) and d) are respectively proved in Steps 5 and 6.

    Step 1: Preliminaries. Given any positive solution of the undelayed perturbed system (3.2) denoted by t \mapsto (s(t), x(t)) , its error with respect to the unique positive \omega –periodic solution t\mapsto (s_{p}(t), x_{p}(t)) of the undelayed system (1.7) is defined by

    \begin{equation} \tilde{s}(t) = s(t)-s_p(t), \quad \text{and} \quad \tilde{x}(t) = x(t)-x_p(t). \end{equation} (A.2)

    In addition, the change of variables s(t) = b(t)-x(t) transforms (3.2) into

    \begin{equation} \left\{ \begin{array}{rcl} \dot{b}(t) & = & Ds^{0}(t)-Db(t)+\delta(t)\\ \dot{x}(t) & = & x(t)[\mu(s(t))-D]. \end{array} \right. \end{equation} (A.3)

    Let b_{p}(t) : = s_{p}(t) + x_{p}(t) . The error \tilde{b}(t) = b(t)-b_p(t) satisfies the scalar equation \dot{\tilde{b}}(t) = -D\tilde{b}(t)+\delta(t) and, consequently, the differential inequality

    \dot{\tilde{b}}(t) \leq - D\tilde{b}(t) + \overline{\delta}

    since \delta(t)\leq \bar{\delta} for any t\geq T_{0} . Note that \overline{\delta}/D is a globally attractive solution of the equation \dot{\nu}(t) = -D\nu(t)+\overline{\delta} . It follows, by the previously mentioned comparison results for scalar differential inequalities, the existence of t_{a}\geq T_0 , dependent on the initial conditions, such that (3.8) is verified, that is

    |\widetilde{b}(t)| = |b(t)-b_{p}(t)|\leq 2\frac{\overline{\delta}}{D} \quad \text{for all}\ t\geq t_{a}

    and statement a) of the lemma has been proved.

    Step 2: Logarithmic estimation for the biomass. We introduce the transformation \chi = \ln(x) , which yields

    \begin{equation} \dot{\chi}(t) = \frac{\mu_m s(t)}{k_s+s(t)}-D. \end{equation} (A.4)

    We also define \widetilde{\chi}(t) = \chi(t)-\chi_p(t) , where \chi_p = \ln(x_p) . Let us note that

    \dot{\widetilde{\chi}}(t) = \frac{\mu_m s(t)}{k_s+s(t)}-\frac{\mu_m s_p(t)}{k_s+s_p(t)}.

    From (1.5) combined with (A.2) and

    \widetilde{s}(t) = s(t)-s_{p}(t) = \widetilde{b}(t)-e^{\chi(t)}+e^{\chi_{p}(t)}

    we deduce that

    \begin{equation} \begin{array}{rcl} \dot{\tilde{\chi}}(t)& = & \frac{\mu_m k_s(s(t)-s_p(t))}{(k_s+s(t))(k_s+s_p(t))}\\ & = & \frac{\mu_m k_s(\tilde{b}(t)-e^{\chi(t)}+e^{\chi_p(t)})}{(k_s+s(t))(k_s+s_p(t))}\\ & = & \frac{\mu_m k_s(-e^{\chi(t)}+e^{\chi_p(t)})}{(k_s+s(t))(k_s+s_p(t))}+\frac{\mu_m k_s\tilde{b}(t)}{(k_s+s(t))(k_s+s_p(t))}\\ & = & \frac{\mu_m k_s(-e^{\chi(t)-\chi_p(t)}+1)e^{\chi_p(t)}}{(k_s+s(t))(k_s+s_p(t))}+\frac{\mu_m k_s\tilde{b}(t)}{(k_s+s(t))(k_s+s_p(t))}\\ & = & \frac{\mu_m k_s(-e^{\tilde{\chi}(t)}+1)e^{\chi_p(t)}}{(k_s+s(t))(k_s+s_p(t))}+\frac{\mu_m k_s\tilde{b}(t)}{(k_s+s(t))(k_s+s_p(t))}. \end{array} \end{equation} (A.5)

    Let us introduce the Lyapunov function

    U_1(\widetilde{\chi}(t)) = \frac{1}{2}\widetilde{\chi}^2(t) .

    Its derivative with respect to t is

    \begin{equation} \begin{array}{rcl} \dot{U}_1(t) & = & \left\{ \frac{\mu_m k_s(-e^{\widetilde{\chi}(t)}+1)e^{\chi_p(t)}}{(k_s+s(t))(k_s+s_p(t))} +\frac{\mu_m k_s\tilde{b}(t)}{(k_s+s(t))(k_s+s_p(t))}\right\}\widetilde{\chi}(t). \end{array} \end{equation} (A.6)

    Furthermore, inequalities (3.8) and s(t)\leq s_{\triangle} when t > T_{1} , and the positiveness and boundedness of s(t) and s_{p}(t) combined with

    z (-e^{z}+1) \leq 0 \quad \text{and} \quad -|z||e^{z}-1| = z(-e^{z}+1) \quad \text{for any}\ z\in \mathbb{R}

    imply that for any t\geq \max\{t_{a}, T_{1}\}

    \begin{array}{rcl} \dot{U}_1(t) & = & \left[\frac{\mu_m k_s e^{\chi_p(t)}}{(k_s+s(t))(k_s+s_p(t))}(-e^{\widetilde{\chi}(t)}+1)+\frac{\mu_m k_s \tilde{b}(t)}{(k_s+s(t))(k_s+s_p(t))}\right]\tilde{\chi}(t)\\ &\leq & -\frac{\mu_m k_s e^{\chi_p(t)}}{(k_s+s(t))(k_s+s_p(t))}|\tilde{\chi}(t)||e^{\widetilde{\chi}(t)}-1|+\frac{\mu_m \widetilde{\chi}(t) \tilde{b}(t)}{k_s}\\ &\leq & -\frac{\mu_m k_s e^{\chi_p(t)}}{(k_s+s_{\triangle})^2}|\tilde{\chi}(t)||e^{\widetilde{\chi}(t)}-1|+2\frac{\mu_m}{k_s}\frac{\overline{\delta}}{D}|\widetilde{\chi}(t)|. \end{array}

    The map t \mapsto \chi_{p}(t) = \ln(x_{p}(t)) is \omega –periodic and \chi_{p}^{\min}: = \min\limits\limits_{t\in [0, \omega]}\chi_p(t) is well defined. Moreover, we can deduce that

    \begin{equation} \begin{array}{rcl} \dot{U}_{1}(t) & \leq & \left[-c_{0}|e^{\widetilde{\chi}(t)}-1|+d_{0}\overline{\delta}\right] |\tilde{\chi}(t)| \end{array} \end{equation} (A.7)

    with

    c_{0} : = \frac{\mu_{m} k_s e^{\chi_{p}^{\min}}}{(k_s+s_{\triangle})^2} \quad \text{and} \quad d_{0}: = 2\frac{\mu_m}{k_sD}.

    Now, the sufficiently small constant \overline{\delta} from Lemma 3.2 will be rewritten as

    \begin{equation} \overline{\delta}: = \frac{\Delta_{0}\, k^2_s e^{\chi_{p}^{\min}}D}{2(k_s+s_{\triangle})^2} = \frac{c_0}{d_0}\Delta_0 \quad \text{where}\ \Delta_{0} < \ln(2-e^{-1})e^{-3/2} . \end{equation} (A.8)

    Then, by using \ln(2-e^{-1})e^{-3/2} < 1 - e^{-1} and d_{0}\, \overline{\delta} = c_{0}\, \Delta_{0} , it follows that

    \begin{array}{rcl} \dot{U}_1(t) & < & c_{0}\left[-|e^{\tilde{\chi}(t)}-1|+\Delta_0 \right]|\tilde{\chi}(t)| \\ & < & c_{0}\left[-|e^{\tilde{\chi}(t)}-1|+(1-e^{-1}) \right]|\tilde{\chi}(t)| \\ & = & c_{0}\, \mathcal{F}(\tilde{\chi}(t))|\, \tilde{\chi}(t)| \end{array}

    where \mathcal{F}(u) = (1-e^{-1})-|e^{u}-1| .

    It is easy to see that, with u_{\min} = - 1 and u_{\max} = \ln(2 - e^{-1}) , \mathcal{F}(u)\geq 0 for any u\in [u_{\min}, u_{\max}] and \mathcal{F}(u) < 0 otherwise, and \mathcal{F}(u_{\min}) = \mathcal{F}(u_{\max}) = 0 . Now, let \varepsilon \in (0, \frac{1}{2}) , and define u_{\min}(\varepsilon) = u_{\min} - \varepsilon and u_{\max}(\varepsilon) = u_{\max} + \varepsilon .

    Step 3: Lower asymptotic bound for \tilde{\chi} . We will prove the existence of T_{*}: = T_{*}(\tilde{\chi}(0)) > \max\{t_{a}, T_{1}\} such that \tilde{\chi}(t) > u_{\min}(\varepsilon) for any t > T_{*} .

    This will be proved by contradiction: If we assume that \tilde{\chi}(t)\leq u_{\min}(\varepsilon) = -1 for any t\geq 0 , it follows that \mathcal{F}(\tilde{\chi}(t)) < 0 , which in turns implies that

    \dot{U}_{1}(t) = \tilde{\chi}(t)\dot{\widetilde{\chi}}(t) < c_{0}\, \mathcal{F}(\tilde{\chi}(t)) < 0 \quad \text{for any} \quad t\geq \max\{t_{a}, T_{1}\}.

    The above inequality combined with \tilde{\chi}(t)\leq -1 for any t\in J: = [\max\{t_{a}, T_{1}\}, +\infty) implies that \dot{\tilde{\chi}}(t) > 0 for any t\in J . Then, we conclude that t\mapsto \tilde{\chi}(t) is strictly increasing on J and upper bounded. Furthermore, there exists \chi^{*}\leq u_{\min}(\varepsilon) such that

    \lim\limits_{t\to \infty}\tilde{\chi}(t) = \chi^{*}\leq u_{\min}(\varepsilon) \quad \text{and} \quad \lim\limits_{t\to +\infty}\dot{\widetilde{\chi}}(t) = 0 .

    Now, letting t\to \infty and noticing that \mathcal{F}(\chi^{*}) < 0 , it follows that

    0 = \lim\limits_{t\to +\infty}\tilde{\chi}(t)\dot{\tilde{\chi}}(t) \leq \lim\limits_{t\to +\infty}c_{0}\, \mathcal{F}(\tilde{\chi}(t))|\tilde{\chi}(t)| = c_{0}\, \mathcal{F}(\chi^{*})|\chi^{*}| < 0 .

    This leads to a contradiction, which allows us to deduce the existence of T_{*} > 0 such that two mutually exclusive behaviors could take place: either \tilde{\chi}(t) > u_{\min}(\varepsilon) for any t > T_{*} , or there exists T_{**} > T_{*} such that

    \tilde{\chi}(t) > u_{\min}(\varepsilon) \, \, \, \text{for any} \, \, \, t\in (T_{*}, T_{**}) \, \, \, \text{with} \, \, \, \tilde{\chi}(T_{**}) = u_{\min}(\varepsilon)\, \, \, \text{and} \, \, \, \dot{\tilde{\chi}}(T_{**}) < 0.

    Now, it is important to emphasize that the last behavior is not possible. Indeed, otherwise we would have

    \dot{U}_{1}(\widetilde{\chi}(T_{**})) = u_{\min}(\varepsilon)\dot{\widetilde{\chi}}(T_{**}) > 0 \, \, \text{and} \, \, \dot{U}_{1}(\tilde{\chi}(T_{**}))\leq c_{0}\, \mathcal{F}(u_{\min}(\varepsilon))|\, u_{\min}(\varepsilon)| < 0

    which gives a contradiction. Then, it follows that \tilde{\chi}(t) > u_{\min}(\varepsilon) for any t > T_{*} .

    Step 4: Upper asymptotic bound for \tilde{\chi} . It can be proved similarly as in the previous step that if \tilde{\chi}(0) > u_{\max}(\varepsilon) , there exists T^{*}: = T^{*}(\tilde{\chi}(0)) > \max\{t_{a}, T_{1}\} such that \tilde{\chi}(t) < u_{\max}(\varepsilon) for any t > T^{*} .

    Next, we fix \varepsilon\in (0, \frac{1}{2}) , and observe that a direct consequence of steps 3 and 4 is the existence of t_{b}\geq \max\{T_{*}, T^{*}\} such that if t\geq t_{b} , then

    -(1+\varepsilon) = u_{\min}(\varepsilon)\leq \tilde{\chi}(t) = \ln(x(t))-\ln(x_{p}(t)) < u_{\max} = (\varepsilon) = \ln(2-e^{-1})+\varepsilon

    and statement b) of the lemma has been proved.

    Step 5: A result of asymptotic invariance. By using the mean value theorem, we can see that e^{\tilde{\chi}}-1 = e^{\theta}\tilde{\chi} with \theta between \tilde{\chi} and zero. Note that if \tilde{\chi}\in[u_{\min}(\varepsilon), u_{\max}(\varepsilon)] , then |e^{\tilde{\chi}}-1| = |e^{\theta}||\tilde{\chi}|\geq e^{u_{\min}(\varepsilon)} |\tilde{\chi}| . As a consequence, when t\geq t_b , it follows from (A.7) that

    \begin{equation} \dot{U}_1(t)\leq \left[-c_{0}e^{u_{\min}(\varepsilon)}|\widetilde{\chi}(t)|+d_{0}\overline{\delta}\right] |\tilde{\chi}(t)| \end{equation} (A.9)

    and, in order to obtain an asymptotic invariance interval for \tilde{\chi}(t) , we will consider three possible cases:

    \bullet Case a) This case corresponds to \tilde{\chi}(t) > 0 after a finite time T > t_{b} . By recalling that \dot{U}_{1}(\tilde{\chi}(t)) = \tilde{\chi}(t)\dot{\tilde{\chi}}(t) , for any t > T , we deduce from (A.9) that

    \begin{array}{rcl} \dot{\tilde{\chi}}(t)&\leq& \left[-c_0e^{u_{\min}(\varepsilon)}\tilde{\chi}(t)+d_0\overline{\delta} \right]. \end{array}

    Since \tilde{\chi}(t)\in (0, u_{\max}(\varepsilon)] , by using comparison results combined with (A.8), we deduce that

    \limsup\limits_{t\to +\infty}\tilde{\chi}(t)\leq x_+: = \frac{d_0\, \overline{\delta}}{c_0 e^{u_{\min}(\varepsilon)}} .

    Recalling that u_{\max} = \ln(2-e^{-1}) and u_{\min}(\varepsilon) = - 1 - \varepsilon , and by using again (A.8) together with \varepsilon\in (0, 1/2) , it follows that

    x_+ < \frac{\Delta_{0}}{e^{u_{\min}(\varepsilon)}} = \Delta_{0}e^{1+\varepsilon} < \ln(2-e^{-1})e^{\varepsilon-1/2} < u_{\max}.

    \bullet Case b) This case corresponds to \tilde{\chi}(t) < 0 after a finite time T > t_{b} . Then, we deduce from (A.9) that, for any t > T ,

    \begin{array}{rcl} \dot{\tilde{\chi}}(t)&\geq& -\left[c_0 e^{u_{\min}(\varepsilon)}\tilde{\chi}(t)+d_0\overline{\delta} \right]. \end{array}

    As \tilde{\chi}(t)\in [u_{\min}(\varepsilon), 0) , using argument analogous to those of the previous case, we have that

    \liminf\limits_{t\to +\infty}\tilde{\chi}(t)\geq x_-: = -\frac{d_0 \overline{\delta}}{c_0 e^{u_{\min}(\varepsilon)}}

    and, as before, from (A.8), d_{0}\, \overline{\delta} = c_{0}\, \Delta_{0} , and \varepsilon \in (0, 1/2) , we obtain

    x_{-} > - \frac{\Delta_{0}}{e^{u_{\min}(\varepsilon)}} = -\Delta_{0}e^{1+\varepsilon} > -\ln(2-e^{-1})e^{\varepsilon-1/2} > -\ln(2-e^{-1}) > u_{\min}.

    \bullet Case c) We assume that \tilde{\chi}(t) \in [u_{\min}(\varepsilon), u_{\max}(\varepsilon)] , but with infinite changes of sign. By using the fluctuations lemma [40, Lemma A.1], there exist divergent sequences \{s_{n}\}_{n} and \{t_{n}\}_{n} of minimum and maximum values of \tilde{\chi}(t) satisfying \dot{\tilde{\chi}}(s_{n}) = \dot{\tilde{\chi}}(t_{n}) = 0 such that

    \lim\limits_{n\to +\infty}\tilde{\chi}(s_{n}) = \liminf\limits_{t\to +\infty}\tilde{\chi}(t) \quad \text{and} \quad \lim\limits_{n\to +\infty}\tilde{\chi}(t_{n}) = \limsup\limits_{t\to +\infty}\tilde{\chi}(t).

    By using the above differential inequalities, it can be deduced that

    \tilde{\chi}(t_{n})\leq x_{+} \quad \text{and} \quad x_{-}\leq \tilde{\chi}(s_{n}) \quad \text{for}\ n \ \text{large enough}.

    Letting n\to +\infty , we obtain

    \limsup\limits_{t\to +\infty}\tilde{\chi}(t)\leq x_+ \quad \text{and} \quad \liminf\limits_{t\to +\infty}\tilde{\chi}(t)\geq x_- .

    Summarizing the above three cases leads to the following property: for any \eta > 0 , there exists t_{c}(\eta)\geq 0 such that

    x_{-}-\eta\leq \tilde{\chi}(t)\leq x_{+}+\eta \quad \text{for}\ t\geq t_{c}(\eta) .

    By noticing that |x_{-}| = x_{+} , the previous inequalities can be written as

    \left|\widetilde{\chi}(t)\right|\leq x_{+}+\eta \quad \text{for}\ t\geq t_{c}(\eta)

    and, after choosing \eta = x_{+}/2 and recalling the definition of x_{+} , we can deduce the existence of t_c\geq t_b such that, for all t\geq t_{c} , we have the estimation

    \begin{equation} |\tilde{\chi}(t)| = |\chi(t)-\chi_{p}(t)| = \left|\ln\left(\frac{x(t)}{x_{p}(t)}\right)\right|\leq p\overline{\delta} \quad \text{with} \quad p = \frac{3d_{0}}{2c_{0}}e^{-u_{\min}(\varepsilon)} \end{equation} (A.10)

    which implies that

    \begin{array}{ccccc} x_p(t)e^{-p\overline{\delta}}-x_p(t)& \leq& x(t)-x_p(t) &\leq& x_p(t)e^{p\overline{\delta}}-x_p(t) \end{array}

    and also leads to

    |x(t) - x_{p}(t)| \leq x_p(t)\left(e^{p\overline{\delta}} - 1\right) .

    Then, inequality (3.10) is obtained and statement c) of the lemma has been proved.

    Step 6: End of proof.

    When \bar{\delta} is sufficiently small, e^{p\overline{\delta}}-1 \leq 2 p\bar{\delta} . Then, from (3.10), we deduce that the second inequality of (3.11) is satisfied with p_{2} = 2 p \max\{t\in [0, \omega]\colon |x_{p}(t)|\} . This fact combined with (3.8) allows us to deduce that the first inequality of (3.11) is satisfied with p_1 = \left(\frac{2}{D} + p_2\right) \overline{\delta} . This concludes the proof.



    [1] H. Sung, J. Ferlay, R. L. Siegel, M. Laversanne, I. Soerjomataram, A. Jemal, et al., Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries, CA: Cancer J. Clin., 71 (2021), 209–249. https://doi.org/10.3322/caac.21660 doi: 10.3322/caac.21660
    [2] M. Das, Neoadjuvant chemotherapy: survival benefit in gastric cancer, Lancet Oncol., 18 (2017), 307. https://doi.org/10.1016/s1470-2045(17)30321-2 doi: 10.1016/s1470-2045(17)30321-2
    [3] National Health Commission of the People's Republic of China, Chinese guidelines for diagnosis and treatment of gastric cancer 2018 (English version), Chin. J. Cancer Res., 31 (2019), 707–737. https://doi.org/10.21147/j.issn.1000-9604.2019.05.01
    [4] M. Orditura, G. Galizia, V. Sforza, V. Gambardella, A. Fabozzi, M. M. Laterza, et al., Treatment of gastric cancer, World J. Gastroenterol., 20 (2014), 1635–1649. https://doi.org/10.3748/wjg.v20.i7.1635 doi: 10.3748/wjg.v20.i7.1635
    [5] M. P. Lutz, J. R. Zalcberg, M. Ducreux, J. A. Ajani, W. Allum, D. Aust, et al., Highlights of the EORTC St. Gallen International Expert Consensus on the primary therapy of gastric, gastroesophageal and oesophageal cancer-Differential treatment strategies for subtypes of early gastroesophageal cancer, Eur. J. Cancer, 48 (2012), 2941–2953. https://doi.org/10.1016/j.ejca.2012.07.029 doi: 10.1016/j.ejca.2012.07.029
    [6] I. Thomassen, Y. R. van Gestel, B. van Ramshorst, M. D. Luyer, K. Bosscha, S. W. Nienhuijs, et al., Peritoneal carcinomatosis of gastric origin: A population-based study on incidence, survival and risk factors, Int. J. Cancer, 134 (2014), 622–628. https://doi.org/10.1002/ijc.28373 doi: 10.1002/ijc.28373
    [7] H. Nakagawa, M. Fujita, Whole genome sequencing analysis for cancer genomics and precision medicine, Cancer Sci., 109 (2018), 513–522. https://doi.org/10.1111/cas.13505 doi: 10.1111/cas.13505
    [8] A. Dongre, R. A. Weinberg, New insights into the mechanisms of epithelial-mesenchymal transition and implications for cancer, Nat. Rev. Mol. Cell. Biol., 20 (2019), 69–84. https://doi.org/10.1038/s41580-018-0080-4 doi: 10.1038/s41580-018-0080-4
    [9] M. Izumiya, A. Kabashima, H. Higuchi, T. Igarashi, G. Sakai, H. Iizuka, et al., Chemoresistance is associated with cancer stem cell-like properties and epithelial-to-mesenchymal transition in pancreatic cancer cells, Anticancer Res., 32 (2012), 3847–3853.
    [10] R. B. Hazan, R. Qiao, R. Keren, I. Badano, K. Suyama, Cadherin switch in tumor progression, Ann. N. Y. Acad. Sci., 1014 (2004), 155–163. https://doi.org/10.1196/annals.1294.016 doi: 10.1196/annals.1294.016
    [11] C. Cai, J. Yu, J. Wu, R. Lu, X. Ni, S. Wang, et al., CD133 promotes the invasion and metastasis of gastric cancer via epithelial-mesenchymal transition, Chin. J. Gastrointest. Surg., 16 (2013), 662–667.
    [12] C. Zeltz, I. Primac, P. Erusappan, J. Alam, A. Noel, D. Gullberg, Cancer-associated fibroblasts in desmoplastic tumors: emerging role of integrins, Semin. Cancer Biol., 62 (2020), 166–181. https://doi.org/10.1016/j.semcancer.2019.08.004 doi: 10.1016/j.semcancer.2019.08.004
    [13] D. F. Quail, J. A. Joyce, Microenvironmental regulation of tumor progression and metastasis, Nat. Med., 19 (2013), 1423–1437. https://doi.org/10.1038/nm.3394 doi: 10.1038/nm.3394
    [14] N. Kemi, M. Eskuri, A. Herva, J. Leppanen, H. Huhta, O. Helminen, et al., Tumour-stroma ratio and prognosis in gastric adenocarcinoma, Br. J. Cancer, 119 (2018), 435–439. https://doi.org/10.1038/s41416-018-0202-y doi: 10.1038/s41416-018-0202-y
    [15] L. Huang, A. M. Xu, S. Liu, W. Liu, T. J. Li, Cancer-associated fibroblasts in digestive tumors, World. J. Gastroenterol., 20 (2014), 17804–17818. https://doi.org/10.3748/wjg.v20.i47.17804 doi: 10.3748/wjg.v20.i47.17804
    [16] A. C. Johansson, A. Ansell, F. Jerhammar, M. B. Lindh, R. Grenman, E. Munck-Wikland, et al., Cancer-associated fibroblasts induce matrix metalloproteinase–mediated cetuximab resistance in head and neck squamous cell carcinoma cells, Mol. Cancer Res., 10 (2012), 1158–1168. https://doi.org/10.1158/1541-7786.Mcr-12-0030 doi: 10.1158/1541-7786.Mcr-12-0030
    [17] R. A. Saito, P. Micke, J. Paulsson, M. Augsten, C. Pena, P. Jonsson, et al., Forkhead box F1 regulates tumor-promoting properties of cancer-associated fibroblasts in lung cancer, Cancer Res., 70 (2010), 2644–2654. https://doi.org/10.1158/0008-5472.Can-09-3644 doi: 10.1158/0008-5472.Can-09-3644
    [18] K. Pietras, K. Rubin, T. Sjoblom, E. Buchdunger, M. Sjoquist, C. Heldin, et al., Inhibition of PDGF receptor signaling in tumor stroma enhances antitumor effect of chemotherapy, Cancer Res., 62 (2002), 5476–5484.
    [19] X. Liu, K. M. Chu, E-Cadherin and gastric cancer: Cause, consequence, and applications, BioMed Res. Int., 2014 (2014), 637308. https://doi.org/10.1155/2014/637308 doi: 10.1155/2014/637308
    [20] S. Herbertz, J. S. Sawyer, A. J. Stauber, I. Gueorguieva, K. E. Driscoll, S. T. Estrem, et al., Clinical development of galunisertib (LY2157299 monohydrate), a small molecule inhibitor of transforming growth factor-beta signaling pathway, Drug Des., Dev. Ther., 9 (2015), 4479–4499. https://doi.org/10.2147/dddt.S86621 doi: 10.2147/dddt.S86621
    [21] M. Singh, N. Yelle, C. Venugopal, S. K. Singh, EMT: Mechanisms and therapeutic implications, Pharmacol. Ther., 182 (2018), 80–94. https://doi.org/10.1016/j.pharmthera.2017.08.009 doi: 10.1016/j.pharmthera.2017.08.009
    [22] X. Zhang, S. L. Marjani, Z. Hu, S. M. Weissman, X. Pan, S. Wu, Single-cell sequencing for precise cancer research: Progress and prospects, Cancer Res., 76 (2016), 1305–1312. https://doi.org/10.1158/0008-5472.Can-15-1907 doi: 10.1158/0008-5472.Can-15-1907
    [23] G. Sun, Z. Li, D. Rong, H. Zhang, X. Shi, W. Yang, et al., Single-cell RNA sequencing in cancer: Applications, advances, and emerging challenges, Mol. Ther Oncolytics., 21 (2021), 183–206. https://doi.org/10.1016/j.omto.2021.04.001 doi: 10.1016/j.omto.2021.04.001
    [24] A. Sathe, S. M. Grimes, B. T. Lau, J. Chen, C. Suarez, R. J. Huang, et al., Single-cell genomic characterization reveals the cellular reprogramming of the gastric tumor microenvironment, Clin. Cancer Res., 26 (2020), 2640–2653. https://doi.org/10.1158/1078-0432.Ccr-19-3231 doi: 10.1158/1078-0432.Ccr-19-3231
    [25] M. Zhang, S. Hu, M. Min, Y. Ni, Z. Lu, X. Sun, et al., Dissecting transcriptional heterogeneity in primary gastric adenocarcinoma by single cell RNA sequencing, Gut, 70 (2021), 464–475. https://doi.org/10.1136/gutjnl-2019-320368 doi: 10.1136/gutjnl-2019-320368
    [26] Y. Li, X. Hu, R. Lin, G. Zhou, L. Zhao, D. Zhao, et al., Single-cell landscape reveals active cell subtypes and their interaction in the tumor microenvironment of gastric cancer, Theranostics, 12 (2022), 3818–3833. https://doi.org/10.7150/thno.71833 doi: 10.7150/thno.71833
    [27] B. Wang, Y. Zhang, T. Qing, K. Xing, J. Li, T. Zhen, et al., Comprehensive analysis of metastatic gastric cancer tumour cells using single-cell RNA-seq, Sci. Rep., 11 (2021), 10. https://doi.org/10.1038/s41598-020-80881-2 doi: 10.1038/s41598-020-80881-2
    [28] A. C. Obenauf, J. Massague, Surviving at a distance: Organ-specific metastasis, Trends Cancer, 1 (2015), 76–91. https://doi.org/10.1016/j.trecan.2015.07.009 doi: 10.1016/j.trecan.2015.07.009
    [29] D. X. Nguyen, P. D. Bos, J. Massague, Metastasis: from dissemination to organ-specific colonization, Nat. Rev. Cancer, 9 (2009), 274–284. https://doi.org/10.1038/nrc2622 doi: 10.1038/nrc2622
    [30] X. Zhang, Y. Lan, J. Xu, F. Quan, E. Zhao, C. Deng, et al., CellMarker: a manually curated resource of cell markers in human and mouse, Nucleic. Acids. Res., 47 (2019), 721–728. https://doi.org/10.1093/nar/gky900 doi: 10.1093/nar/gky900
    [31] T. Yan, W. Qiu, J. Song, Y. Fan, Z. Yang, ARHGAP36 regulates proliferation and migration in papillary thyroid carcinoma cells, J. Mol. Endocrinol., 66 (2021), 1–10. https://doi.org/10.1530/jme-20-0230 doi: 10.1530/jme-20-0230
    [32] C. Han, T. Liu, R. Yin, Biomarkers for cancer-associated fibroblasts, Biomark. Res., 8 (2020). https://doi.org/10.1186/s40364-020-00245-w doi: 10.1186/s40364-020-00245-w
    [33] S. Togo, U. M. Polanska, Y. Horimoto, A. Orimo, Carcinoma-associated fibroblasts are a promising therapeutic target, Cancers, 5 (2013), 149–169. https://doi.org/10.3390/cancers5010149 doi: 10.3390/cancers5010149
    [34] G. Corso, J. Figueiredo, S. P. De Angelis, F. Corso, A. Girardi, J. Pereira, et al., E-cadherin deregulation in breast cancer, J. Cell. Mol. Med., 24 (2020), 5930–5936. https://doi.org/10.1111/jcmm.15140 doi: 10.1111/jcmm.15140
    [35] Y. A. Lyons, S. Y. Wu, W. W. Overwijk, K. A. Baggerly, A. K. Sood, Immune cell profiling in cancer: molecular approaches to cell-specific identification, npj Precision Onc., 1 (2017). https://doi.org/10.1038/s41698-017-0031-0 doi: 10.1038/s41698-017-0031-0
    [36] Z. Chen, Z. Han, H. Nan, J. Fan, J. Zhan, Y. Zhang, et al., A novel pyroptosis-related gene signature for predicting the prognosis and the associated immune infiltration in colon adenocarcinoma, Front. Oncol., 12 (2022). https://doi.org/10.3389/fonc.2022.904464 doi: 10.3389/fonc.2022.904464
    [37] S. Han, K. Huang, Z. Gu, J. Wu, Tumor immune microenvironment modulation-based drug delivery strategies for cancer immunotherapy, Nanoscale, 12 (2020), 413–436. https://doi.org/10.1039/c9nr08086d doi: 10.1039/c9nr08086d
    [38] X. Geng, H. Chen, L. Zhao, J. Hu, W. Yang, G. Li, et al., Cancer-Associated Fibroblast (CAF) heterogeneity and targeting therapy of CAFs in pancreatic cancer, Front. Cell Dev. Biol., 9 (2021). https://doi.org/10.3389/fcell.2021.655152 doi: 10.3389/fcell.2021.655152
    [39] L. A. Aparicio, M. Blanco, R. Castosa, A. Concha, M. Valladares, L. Calvo, et al., Clinical implications of epithelial cell plasticity in cancer progression, Cancer Lett., 366 (2015), 1–10. https://doi.org/10.1016/j.canlet.2015.06.007 doi: 10.1016/j.canlet.2015.06.007
    [40] T. Baslan, J. Hicks, Unravelling biology and shifting paradigms in cancer with single-cell sequencing, Nat. Rev. Cancer, 17 (2017), 557–569. https://doi.org/10.1038/nrc.2017.58 doi: 10.1038/nrc.2017.58
    [41] Z. Zhang, S. Zheng, Y. Lin, J. Sun, N. Ding, J. Chen, et al., Genomics and prognosis analysis of epithelial-mesenchymal transition in colorectal cancer patients, BMC Cancer, 20 (2020). https://doi.org/10.1186/s12885-020-07615-5 doi: 10.1186/s12885-020-07615-5
    [42] C. Xiong, G. Wang, D. Bai, A novel prognostic models for identifying the risk of hepatocellular carcinoma based on epithelial-mesenchymal transition-associated genes, Bioengineered, 11 (2020), 1034–1046. https://doi.org/10.1080/21655979.2020.1822715 doi: 10.1080/21655979.2020.1822715
    [43] W. Dai, Y. Xiao, W. Tang, J. Li, L. Hong, J. Zhang, et al., Identification of an EMT-related gene signature for predicting overall survival in gastric cancer, Front. Genet., 12 (2021). https://doi.org/10.3389/fgene.2021.661306 doi: 10.3389/fgene.2021.661306
    [44] Y. Hu, Z. Hu, T. Liao, Y. Li, Y. Pan, LncRNA SND1-IT1 facilitates TGF-beta 1-induced epithelialto-mesenchymal transition via miR-124/COL4A1 axis in gastric cancer, Cell Death Discov., 8 (2022). https://doi.org/10.1038/s41420-021-00793-6 doi: 10.1038/s41420-021-00793-6
    [45] T. Wang, H. Jin, J. Hu, X. Li, H. Ruan, H. Xu, et al., COL4A1 promotes the growth and metastasis of hepatocellular carcinoma cells by activating FAK-Src signaling, J. Exp. Clin. Cancer Res., 39 (2020). https://doi.org/10.1186/s13046-020-01650-7 doi: 10.1186/s13046-020-01650-7
    [46] X. Xie, H. He, N. Zhang, X. Wang, W. Rui, D. Xu, et al., Overexpression of DDR1 promotes migration, invasion, though EMT-Related molecule expression and COL4A1/DDR1/MMP-2 signaling axis, Technol. Cancer Res. Treat., 19 (2020). https://doi.org/10.1177/1533033820973277 doi: 10.1177/1533033820973277
    [47] M. Miyake, Y. Morizawa, S. Hori, Y. Tatsumi, S. Onishi, T. Owari, et al., Diagnostic and prognostic role of urinary collagens in primary human bladder cancer, Cancer Sci., 108 (2017), 2221–2228. https://doi.org/10.1111/cas.13384 doi: 10.1111/cas.13384
    [48] Y. Zhang, X. Qu, C. Li, Y. Fan, X. Che, X. Wang, et al., miR-103/107 modulates multidrug resistance in human gastric carcinoma by downregulating Cav-1, Tumor Biol., 36 (2015), 2277–2285. https://doi.org/10.1007/s13277-014-2835-7 doi: 10.1007/s13277-014-2835-7
    [49] D. S. Sun, S. A. Hong, H. S. Won, S. H. Yoo, H. H. Lee, O. Kim, et al., Prognostic value of metastatic tumoral caveolin-1 expression in patients with resected gastric cancer, Gastroenterol. Res. Pract., 2017 (2017). https://doi.org/10.1155/2017/5905173 doi: 10.1155/2017/5905173
    [50] W. Liu, N. Yin, H. Liu, K. Nan, Cav-1 promote lung cancer cell proliferation and invasion through lncRNA HOTAIR, Gene, 641 (2018), 335–340. https://doi.org/10.1016/j.gene.2017.10.070 doi: 10.1016/j.gene.2017.10.070
    [51] D. Fujimoto, Y. Hirono, T. Goi, K. Katayama, S. Matsukawa, A. Yamaguchi, The activation of Proteinase-Activated Receptor-1 (PAR1) mediates gastric cancer cell proliferation and invasion, BMC Cancer, 10 (2010). https://doi.org/10.1186/1471-2407-10-443 doi: 10.1186/1471-2407-10-443
    [52] A. K. S. Arakaki, W. Pan, H. Wedegaertner, I. Roca-Mercado, L. Chinn, T. S. Gujral, et al., α-Arrestin ARRDC3 tumor suppressor function is linked to GPCR-induced TAZ activation and breast cancer metastasis, J. Cell Sci., 134 (2021). https://doi.org/10.1242/jcs.254888 doi: 10.1242/jcs.254888
    [53] N. Smoktunowicz, M. Plate, A. O. Stern, V. D. Antongiovanni, E. Robinson, V. Chudasama, et al., TGF beta upregulates PAR-1 expression and signalling responses in A549 lung adenocarcinoma cells, Oncotarget, 7 (2016), 65471–65484. https://doi.org/10.18632/oncotarget.11472 doi: 10.18632/oncotarget.11472
    [54] H. Deng, R. Guo, W. Li, M. Zhao, Y. Lu, Matrix metalloproteinase 11 depletion inhibits cell proliferation in gastric cancer cells, Biochem. Biophys. Res. Commun., 326 (2005), 274–281. https://doi.org/10.1016/j.bbrc.2004.11.027 doi: 10.1016/j.bbrc.2004.11.027
    [55] G. Xu, B. Zhang, J. Ye, S. Cao, J. Shi, Y. Zhao, et al., Exosomal miRNA-139 in cancer-associated fibroblasts inhibits gastric cancer progression by repressing MMP11 expression, Int. J. Biol. Sci., 15 (2019), 2320–2329. https://doi.org/10.7150/ijbs.33750 doi: 10.7150/ijbs.33750
    [56] H. B. Han, J. Gu, H. Zuo, Z. Chen, W. Zhao, M. Li, et al., Let-7c functions as a metastasis suppressor by targeting MMP11 and PBX3 in colorectal cancer, J. Pathol., 226 (2012), 544–555. https://doi.org/10.1002/path.3014 doi: 10.1002/path.3014
    [57] Y. Zhuang, X. Li, P. Zhan, G. Pi, G. Wen, MMP11 promotes the proliferation and progression of breast cancer through stabilizing Smad2 protein, Oncol. Rep., 45 (2021). https://doi.org/10.3892/or.2021.7967 doi: 10.3892/or.2021.7967
    [58] L. Zhai, W. Chen, B. Cui, B. Yu, Y. Wang, H. Liu, Overexpressed versican promoted cell multiplication, migration and invasion in gastric cancer, Tissue Cell, 73 (2021). https://doi.org/10.1016/j.tice.2021.101611 doi: 10.1016/j.tice.2021.101611
    [59] S. P. Evanko, P. Y. Johnson, K. R. Braun, C. B. Underhill, J. Dudhia, T. N. Wight, Platelet-derived growth factor stimulates the formation of versican-hyaluronan aggregates and pericellular matrix expansion in arterial smooth muscle cells, Arch. Biochem. Biophys., 394 (2001), 29–38. https://doi.org/10.1006/abbi.2001.2507 doi: 10.1006/abbi.2001.2507
    [60] Y. Zhang, X. Zou, W. Qian, X. Weng, L. Zhang, L. Zhang, et al., Enhanced PAPSS2/VCAN sulfation axis is essential for Snail-mediated breast cancer cell migration and metastasis, Cell Death Differ., 26 (2019), 565–579. https://doi.org/10.1038/s41418-018-0147-y doi: 10.1038/s41418-018-0147-y
    [61] R. Wang, D. Zhang, C. Zhao, Q. Wang, H. Qu, Q. He, FKBP10 functioned as a cancer-promoting factor mediates cell proliferation, invasion, and migration via regulating PI3K signaling pathway in stomach adenocarcinoma, Kaohsiung J. Med. Sci., 36 (2020), 311–317. https://doi.org/10.1002/kjm2.12174 doi: 10.1002/kjm2.12174
    [62] G. Ramadori, R. M. Ioris, Z. Villanyi, R. Firnkes, O. O. Panasenko, G. Allen, et al., FKBP10 regulates protein translation to sustain lung cancer growth, Cell Rep., 30 (2020), 3851–3863. https://doi.org/10.1016/j.celrep.2020.02.082 doi: 10.1016/j.celrep.2020.02.082
    [63] H. Cai, M. Zhang, Z. Cheng, J. Yu, Q. Yuan, J. Zhang, et al., FKBP10 promotes proliferation of glioma cells via activating AKT-CREB-PCNA axis, J. Biomed. Sci., 28 (2021). https://doi.org/10.1186/s12929-020-00705-3 doi: 10.1186/s12929-020-00705-3
    [64] R. Murad, A. Avanes, X. Ma, S. Geng, A. Mortazavi, J. Momand, Transcriptome and chromatin landscape changes associated with trastuzumab resistance in HER2+breast cancer cells, Gene, 799 (2021), 145808. https://doi.org/10.1016/j.gene.2021.145808 doi: 10.1016/j.gene.2021.145808
    [65] S. Ashida, H. Nakagawa, T. Katagiri, M. Furihata, M. Iiizumi, Y. Anazawa, et al., Molecular features of the transition from prostatic intraepithelial neoplasia (PIN) to prostate cancer: genome-wide gene-expression profiles of prostate cancers and PINs, Cancer Res., 64 (2004), 5963–5972. https://doi.org/10.1158/0008-5472.Can-04-0020 doi: 10.1158/0008-5472.Can-04-0020
    [66] J. Vazquez, L. Gonzalez, A. Merino, F. Vizoso, Expression and clinical significance of apolipoprotein D in epithelial ovarian carcinomas, Gynecol. Oncol., 76 (2000), 340–347. https://doi.org/10.1006/gyno.1999.5678 doi: 10.1006/gyno.1999.5678
    [67] J. Huo, L. Wu, Y. Zang, Construction and validation of a universal applicable prognostic signature for gastric cancer based on seven immune-related gene correlated with tumor associated macrophages, Front. Oncol., 11 (2021). https://doi.org/10.3389/fonc.2021.635324 doi: 10.3389/fonc.2021.635324
    [68] X. Guo, X. Liang, Y. Wang, A. Cheng, H. Zhang, C. Qin, et al., Significance of tumor mutation burden combined with immune infiltrates in the progression and prognosis of advanced gastric cancer, Front. Genet., 12 (2021). https://doi.org/10.3389/fgene.2021.642608 doi: 10.3389/fgene.2021.642608
    [69] Z. Peng, C. Wang, E. Fang, G. Wang, Q. Tong, Role of epithelial-mesenchymal transition in gastric cancer initiation and progression, World J. Gastroenterol., 20 (2014), 5403–5410. https://doi.org/10.3748/wjg.v20.i18.5403 doi: 10.3748/wjg.v20.i18.5403
    [70] W. Li, X. Zhang, F. Wu, Y. Zhou, Z. Bao, H. Li, et al., Gastric cancer-derived mesenchymal stromal cells trigger M2 macrophage polarization that promotes metastasis and EMT in gastric cancer, Cell Death Dis., 10 (2019). https://doi.org/10.1038/s41419-019-2131-y doi: 10.1038/s41419-019-2131-y
    [71] S. Su, Q. Liu, J. Chen, J. Chen, F. Chen, C. He, et al., A positive feedback loop between mesenchymal-like cancer cells and macrophages is essential to breast cancer metastasis, Cancer Cell, 25 (2014), 605–620. https://doi.org/10.1016/j.ccr.2014.03.021 doi: 10.1016/j.ccr.2014.03.021
    [72] P. Xiao, X. Long, L. Zhang, Y. Ye, J. Guo, P. Liu, et al., Neurotensin/IL-8 pathway orchestrates local inflammatory response and tumor invasion by inducing M2 polarization of Tumor-Associated macrophages and epithelial-mesenchymal transition of hepatocellular carcinoma cells, OncoImmunology, 7 (2018). https://doi.org/10.1080/2162402x.2018.1440166 doi: 10.1080/2162402x.2018.1440166
    [73] M. C. A. Wouters, B. H. Nelson, Prognostic significance of tumor-infiltrating B cells and plasma cells in human cancer, Clin. Cancer Res., 24 (2018), 6125–6135. https://doi.org/10.1158/1078-0432.Ccr-18-1481 doi: 10.1158/1078-0432.Ccr-18-1481
    [74] C. Gu-Trantien, S. Loi, S. Garaud, C. Equeter, M. Libin, A. de Wind, et al., CD4(+) follicular helper T cell infiltration predicts breast cancer survival, J. Clin. lnvestigation, 123 (2013), 2873–2892. https://doi.org/10.1172/jci67428 doi: 10.1172/jci67428
    [75] G. R. Gunassekaran, S. M. P. Vadevoo, M. Baek, B. Lee, M1 macrophage exosomes engineered to foster M1 polarization and target the IL-4 receptor inhibit tumor growth by reprogramming tumor-associated macrophages into M1-like macrophages, Biomaterials, 278 (2021), 121137. https://doi.org/10.1016/j.biomaterials.2021.121137 doi: 10.1016/j.biomaterials.2021.121137
  • mbe-20-08-614 supplementary.zip
  • This article has been cited by:

    1. Pablo Amster, Gonzalo Robledo, Daniel Sepúlveda, Dynamics of a delayed discrete size-structured chemostat with periodic nutrient supply, 2024, 132, 10075704, 107904, 10.1016/j.cnsns.2024.107904
    2. Mauro Rodriguez Cartabia, Daniel Sepulveda Oehninger, Uniform persistence criteria for a variable inputs chemostat model with delayed response in growth and complete analysis of the periodic case, 2024, 10075704, 108505, 10.1016/j.cnsns.2024.108505
    3. Pablo Amster, 2024, Chapter 3, 978-3-031-73273-7, 15, 10.1007/978-3-031-73274-4_3
  • Reader Comments
  • © 2023 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(3905) PDF downloads(223) Cited by(1)

Figures and Tables

Figures(12)

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog