Research article Special Issues

Linking state-and-transition simulation and timber supply models for forest biomass production scenarios

  • We linked state-and-transition simulation models (STSMs) with an economics-based timber supply model to examine landscape dynamics in North Carolina through 2050 for three scenarios of forest biomass production. Forest biomass could be an important source of renewable energy in the future, but there is currently much uncertainty about how biomass production would impact landscapes. In the southeastern US, if forests become important sources of biomass for bioenergy, we expect increased land-use change and forest management. STSMs are ideal for simulating these landscape changes, but the amounts of change will depend on drivers such as timber prices and demand for forest land, which are best captured with forest economic models. We first developed state-and-transition model pathways in the ST-Sim software platform for 49 vegetation and land-use types that incorporated each expected type of landscape change. Next, for the three biomass production scenarios, the SubRegional Timber Supply Model (SRTS) was used to determine the annual areas of thinning and harvest in five broad forest types, as well as annual areas converted among those forest types, agricultural, and urban lands. The SRTS output was used to define area targets for STSMs in ST-Sim under two scenarios of biomass production and one baseline, business-as-usual scenario. We show that ST-Sim output matched SRTS targets in most cases. Landscape dynamics results indicate that, compared with the baseline scenario, forest biomass production leads to more forest and, specifically, more intensively managed forest on the landscape by 2050. Thus, the STSMs, informed by forest economics models, provide important information about potential landscape effects of bioenergy production.

    Citation: Jennifer K. Costanza, Robert C. Abt, Alexa J. McKerrow, Jaime A. Collazo. Linking state-and-transition simulation and timber supply models for forest biomass production scenarios[J]. AIMS Environmental Science, 2015, 2(2): 180-202. doi: 10.3934/environsci.2015.2.180

    Related Papers:

    [1] Zhenguo Bai . Threshold dynamics of a periodic SIR model with delay in an infected compartment. Mathematical Biosciences and Engineering, 2015, 12(3): 555-564. doi: 10.3934/mbe.2015.12.555
    [2] Curtis L. Wesley, Linda J. S. Allen, Michel Langlais . Models for the spread and persistence of hantavirus infection in rodents with direct and indirect transmission. Mathematical Biosciences and Engineering, 2010, 7(1): 195-211. doi: 10.3934/mbe.2010.7.195
    [3] Jie He, Zhenguo Bai . Global Hopf bifurcation of a cholera model with media coverage. Mathematical Biosciences and Engineering, 2023, 20(10): 18468-18490. doi: 10.3934/mbe.2023820
    [4] Cuicui Jiang, Kaifa Wang, Lijuan Song . Global dynamics of a delay virus model with recruitment and saturation effects of immune responses. Mathematical Biosciences and Engineering, 2017, 14(5&6): 1233-1246. doi: 10.3934/mbe.2017063
    [5] Lin Zhao, Zhi-Cheng Wang, Liang Zhang . Threshold dynamics of a time periodic and two–group epidemic model with distributed delay. Mathematical Biosciences and Engineering, 2017, 14(5&6): 1535-1563. doi: 10.3934/mbe.2017080
    [6] Yilei Tang, Dongmei Xiao, Weinian Zhang, Di Zhu . Dynamics of epidemic models with asymptomatic infection and seasonal succession. Mathematical Biosciences and Engineering, 2017, 14(5&6): 1407-1424. doi: 10.3934/mbe.2017073
    [7] Yan Zhang, Sanyang Liu, Zhenguo Bai . Global dynamics of a di usive single species model with periodic delay. Mathematical Biosciences and Engineering, 2019, 16(4): 2293-2304. doi: 10.3934/mbe.2019114
    [8] Shengqiang Liu, Lin Wang . Global stability of an HIV-1 model with distributed intracellular delays and a combination therapy. Mathematical Biosciences and Engineering, 2010, 7(3): 675-685. doi: 10.3934/mbe.2010.7.675
    [9] Qian Ding, Jian Liu, Zhiming Guo . Dynamics of a malaria infection model with time delay. Mathematical Biosciences and Engineering, 2019, 16(5): 4885-4907. doi: 10.3934/mbe.2019246
    [10] Rong Ming, Xiao Yu . Global dynamics of an impulsive vector-borne disease model with time delays. Mathematical Biosciences and Engineering, 2023, 20(12): 20939-20958. doi: 10.3934/mbe.2023926
  • We linked state-and-transition simulation models (STSMs) with an economics-based timber supply model to examine landscape dynamics in North Carolina through 2050 for three scenarios of forest biomass production. Forest biomass could be an important source of renewable energy in the future, but there is currently much uncertainty about how biomass production would impact landscapes. In the southeastern US, if forests become important sources of biomass for bioenergy, we expect increased land-use change and forest management. STSMs are ideal for simulating these landscape changes, but the amounts of change will depend on drivers such as timber prices and demand for forest land, which are best captured with forest economic models. We first developed state-and-transition model pathways in the ST-Sim software platform for 49 vegetation and land-use types that incorporated each expected type of landscape change. Next, for the three biomass production scenarios, the SubRegional Timber Supply Model (SRTS) was used to determine the annual areas of thinning and harvest in five broad forest types, as well as annual areas converted among those forest types, agricultural, and urban lands. The SRTS output was used to define area targets for STSMs in ST-Sim under two scenarios of biomass production and one baseline, business-as-usual scenario. We show that ST-Sim output matched SRTS targets in most cases. Landscape dynamics results indicate that, compared with the baseline scenario, forest biomass production leads to more forest and, specifically, more intensively managed forest on the landscape by 2050. Thus, the STSMs, informed by forest economics models, provide important information about potential landscape effects of bioenergy production.


    The term hantavirus represents several groups of RNA-containing viruses (that are members of the virus family of Bunyaviridae) that are carried by rodents, particularly deer mice. The virus is found in their urine and feces, but it does not make the animal sick. The hantavirus is a rare but potentially very serious disease that affects a small number of people every year. People become infected through contact with hantavirus-infected rodents or their urine and droppings. Infection with hantavirus can progress to Hantavirus Pulmonary Syndrome (HPS), which can be fatal. The Sin Nombre hantavirus, first recognized in 1993, is one of several New World hantaviruses circulating in the US. Old World hantaviruses, including Seoul virus, are found across the world and can cause Hemorrhagic Fever with Renal Syndrome (HFRS). Hantavirus does not seem to spread from human to human. The humans get the virus from the mice but have no feedback effects on the mice in the infection process. Furthermore, the transmission of the infection is horizontal, i.e., no mice are born infected [1], infection may only be contracted in adulthood from other mice mainly through agammaessive encounters, such as fights among them, or through inhalation of the aerosolized virus [2,3,4]. The infection and persistence of hantavirus in its rodent host has little or no effect on survival [5]. In order to have a better understanding of the spread of the disease in humans, it is necessary to understand the transmission dynamics of hantavirus in rodent populations.

    Several epidemic models that have been applied to hantavirus infection in rodents are available. For instance, Abramson and Kenkre [6] and Abramson and Kenkre et al. [7] formulated a two-equation reaction diffusion model for susceptible and infected deer mice. Logistic growth is assumed with carrying capacity $ K $. $ K $ is chosen as a control parameter of the dynamics. They analyzed the traveling waves of a model of the hantavirus infection. Allen et al. [8] developed an SI epidemic model for a host with two viral infections. The model was applied to a hantavirus and an arenavirus that infect cotton rats. The first virus is transmitted horizontally whereas the second is transmitted vertically. Considering the movement characteristics of the mice that carry the infection, Kenkre et al. [9] considered two types of mice, stationary adult mice and itinerant juvenile mice. Gedeon et al. [10] applied their model to hantavirus infection in deer mice. Their goal was to compare the relative importance of direct (from infected to healthy individuals) and indirect (by the contaminated environment) transmission in sylvan and peridomestic environments. Sauvage et al. [11] formulated a model which was applied to Puumala virus infection in bank voles. Host population was divided into two age classes, juveniles and adults. It was assumed that individuals can live in favourable or unfavourable patches. Their study showed that indirect transmission significantly increased the probability for the virus to persist in the host population. These two transmission modes also have been discussed by Wolf [12] and Wolf et al. [13]. Abramson and Kenkre [1] gave a stage-dependent model with maturation delay, in their model, a virus-free young mice variable is introduced, the adult population was subdivided into susceptibles and infectives. The spatial version of [1] was presented in [14]. Based on the model in [6], Buceta et al. [15] studied the impact of seasonality on hantavirus, they have shown that the alternation of seasons may cause outbreaks of the disease even if neither season by itself satisfies the environmental requirements for propagation of the disease. Allen et al. [16] constructed two gender-structured SEIR model, the first model is a system of ordinary differential equations, while the second model extends the first model to a system of stochastic differential equations. These models are studied mainly from the numerical simulation point of view. A spatio-temporal SEIR compartmental model was proposed in Burger et al. [17].

    It has been noted that environmental conditions are directly connected to outbreaks of Hanta. For example, the Four Corners Region, where an important number of cases of HPS have occurred, has a desert climate. The largest climate variations within this region come from periods of rain and of drought [15]. The influence of the environmental conditions plays a fundamental role in the evolution of the population. Resource availability would fluctuate as seed and fruit production vary over a 3-year period [18].

    As in most infectious diseases, there is a lag between exposure and infectivity, which is usually called the incubation period. Because the life expectancy of rodents is relatively short, then the incubation period cannot be neglected [16]. Since infected mice that survive the incubation period will remain infectious for the rest of their lives. Thus, the incubation period directly influences the number of infectious mice. However, explicit delay effects related to finite incubation periods have received little attention.

    Motivated by the works of [10,11,15,16,17], in this paper we formulate a periodic time-delayed model by taking into account the seasonality. The model contains a time delay accounting for the incubation period of the virus.

    The rest of this paper is organized as follows. In the next section, based on the work of Gedeon et al. [10] we present the model, and study its well-posedness, also we introduce the basic reproduction number $ R_0 $. In section 3, we establish the threshold dynamics in terms of $ R_0 $. In section 4, we study the autonomous case of the periodic model, and prove the global stability of the virus-free equilibrium and the global attractivity of the endemic equilibrium. In section 5, we perform numerical simulations to illustrate our analytic results. A brief discussion section completes the paper.

    Our model was built on the framework of [10]. In [10], $ S(t) $, $ I(t) $ and $ G(t) $ denote the susceptible, infectious adult mice at time $ t $, and the number of contaminated sites in the environment. In our model, we add the exposed class, we let $ E(t) $ denote the exposed adult mice at time $ t $. We make the same assumptions as those in [10], that is, we assume that there is a discrete number of sites that are visited by mice. Each site is small enough to be infected by a single mouse. Also we assume that the total number of sites $ M $ is large and we represent it by a continuous variable.

    We let $ N(t) $ be the total population of mice, then $ N(t) = S(t)+E(t)+I(t) $. The direct contact rate $ c(t) $ of mice is the average number of direct contacts between mice per mouse per unit time at time $ t $. This rate depends on a number of factors, and in particular, climatic ones, but for simplicity in this paper we assume $ c(t) $ to be periodic. These contacts may involve biting and scratching and are thought to occur predominantly between sexually active males [2,19]. Suppose the transmission probability that given an active contact between susceptible and infected mice is denoted by $ \beta $, which is called a direct transmission rate.

    For the indirect transmission route, we let $ \bar{c}(t) $ be the number of contacts between a mouse and all the potentially contaminated sites per susceptible mouse per unit time at time $ t $, which is also assumed to be periodic. The probability that given a contact between susceptible mouse and contaminated site is denoted by $ \alpha $, this is called an indirect transmission rate.

    Indirect transmission is due to the fact that infected individuals can excrete viruses in their feces, vomit, urine or others. To model the process of site contamination by infected mice, let $ d(t) $ be the number of contacts between the mice and the site that can lead to transmission of the infection per uncontaminated site per unit time at time $ t $. The probability of site contamination is called the site contamination rate, which is denoted by $ \gamma $, given a contact between an uncontaminated site and an infected mouse.

    Assume that $ B(t) $ is the recruitment rate for mice and $ \mu(t) $ is the death rate of the mice. Since hantavirus is not lethal to mice, then we assume the death rates are the same for the infected and the susceptible classes [3]. The mice are infected for the rest of their lives [3,20], so there is no recovery. We assume that the environment eliminates viruses with time at a rate $ \delta(t) $, which is called the disinfection rate. Suppose that $ \tau $ is the average incubation period. Then we obtain the following system

    $ {dS(t)dt=B(t)μ(t)S(t)β1(t)S(t)I(t)N(t)β2(t)S(t)G(t)M,dE(t)dt=β1(t)S(t)I(t)N(t)+β2(t)S(t)G(t)Mμ(t)E(t)[β1(tτ)S(tτ)I(tτ)N(tτ)+β2(tτ)S(tτ)G(tτ)M]ettτμ(r)dr,dI(t)dt=[β1(tτ)S(tτ)I(tτ)N(tτ)+β2(tτ)S(tτ)G(tτ)M]ettτμ(r)drμ(t)I(t),dG(t)dt=a(t)I(t)N(t)(MG(t))δ(t)G(t),
    $
    (2.1)

    where $ \beta_1(t) = \beta c(t), \beta_2(t) = \alpha\bar{c}(t), a(t) = \gamma d(t). $ All parameters are positive, continuous, and $ \omega $-periodic functions for some $ \omega > 0 $. It is easy to see that the equation for $ E(t) $ can be rewritten as one integral equation

    $ E(t)=ttτ[β1(ξ)S(ξ)I(ξ)N(ξ)+β2(ξ)S(ξ)G(ξ)M]etξμ(r)drdξ.
    $
    (2.2)

    The dynamics of the mouse population is governed by the following equation

    $ dN(t)dt=B(t)μ(t)N(t).
    $
    (2.3)

    It is easy to see that (2.3) has a unique positive $ \omega $-periodic solution

    $ N(t)=[t0B(r)er0μ(s)dsdr+ω0B(r)er0μ(s)dsdreω0μ(s)ds1]et0μ(s)ds,
    $

    which is globally asymptotically stable.

    Let $ C = C([-\tau, 0], \mathbb{R}^4) $, $ C^+ = C([-\tau, 0], \mathbb{R}^4_{+}) $. Then $ (C, C^{+}) $ is an ordered Banach space equipped with the maximum norm and the partial order induced by the positive cone $ C^{+} $. For any given continuous function $ x: [-\tau, \sigma)\rightarrow \mathbb{R}^4 $ with $ \sigma > 0 $, we can define $ x_t\in C $ as $ x_t(\theta) = x(t+\theta) $, $ \forall \theta\in[-\tau, 0] $ for any $ t\in[0, \sigma) $.

    For a given continuous $ \omega $-periodic function $ g(t) $, let

    $ \widehat{g} = \max\limits_{t\in[0, \omega]}g(t), \quad \overline{g} = \min\limits_{t\in[0, \omega]}g(t). $

    Let

    $ W: = C([-\tau, 0], \mathbb{R}^3_{+})\times C([-\tau, 0], [0, M]). $

    In view of (2.2), we choose the initial data for system (2.1) in $ \mathcal{X}_{\delta_0} $, which is defined as

    $ Xδ0={ϕW:3i=1ϕi(s)δ0,s[τ,0],ϕ2(0)=0τ[β1(ξ)ϕ1(ξ)ϕ3(ξ)3i=1ϕi(ξ)+β2(ξ)ϕ1(ξ)ϕ4(ξ)M]e0ξμ(r)drdξ}
    $

    for small $ \delta_0\in\left(0, \frac{\overline{B}}{\widehat{\mu}}\right) $.

    Lemma 2.1. For any $ \phi\in \mathcal{X}_{\delta_0} $, system (2.1) has a unique nonnegative solution $ u(t, \phi) $ with $ u_0 = \phi $ for all $ t\geq0 $, and solutions are ultimately bounded and uniformly bounded.

    Proof. Given $ \phi\in \mathcal{X}_{\delta_0} $, define $ f(t, \phi) = (f_1(t, \phi), f_2(t, \phi), f_3(t, \phi), f_4(t, \phi)) $ with

    $ f1(t,ϕ)=B(t)μ(t)ϕ1(0)β1(t)ϕ1(0)ϕ3(0)3i=1ϕi(0)β2(t)ϕ1(0)ϕ4(0)M,f2(t,ϕ)=β1(t)ϕ1(0)ϕ3(0)3i=1ϕi(0)+β2(t)ϕ1(0)ϕ4(0)Mμ(t)ϕ2(0)[β1(tτ)ϕ1(τ)ϕ3(τ)3i=1ϕi(τ)+β2(tτ)ϕ1(τ)ϕ4(τ)M]ettτμ(r)dr,f3(t,ϕ)=[β1(tτ)ϕ1(τ)ϕ3(τ)3i=1ϕi(τ)+β2(tτ)ϕ1(τ)ϕ4(τ)M]ettτμ(r)drμ(t)ϕ3(0),f4(t,ϕ)=a(t)ϕ3(0)3i=1ϕi(0)(Mϕ4(0))δ(t)ϕ4(0).
    $

    Since $ f(t, \phi) $ is continuous in $ (t, \phi)\in \mathbb{R}_{+}\times \mathcal{X}_{\delta_0} $, and $ f(t, \phi) $ is Lipschitz in $ \phi $ on each compact subset of $ \mathcal{X}_{\delta_0} $, it then follows that system (1) has a unique solution $ u(t, \phi) $ on its maximal interval $ [0, \sigma_\phi) $ of existence with $ u_0 = \phi $ (see [21,Theorems 2.2.1 and 2.2.3]).

    In view of the second equation of system (2.1), we have

    $ (dE(t)dt+μ(t)E(t))et0μ(s)ds=[β1(t)S(t)I(t)N(t)+β2(t)S(t)G(t)M(β1(tτ)S(tτ)I(tτ)N(tτ)+β2(tτ)S(tτ)G(tτ)M)ettτμ(r)dr]et0μ(s)ds.
    $

    By integrating on both sides, we obtain

    $ E(t)et0μ(s)dsE(0)=t0[β1(s)S(s)I(s)N(s)+β2(s)S(s)G(s)M]es0μ(ρ)dρdst0[β1(sτ)S(sτ)I(sτ)N(sτ)+β2(sτ)S(sτ)G(sτ)M]esτ0μ(ρ)dρds=t0[β1(s)S(s)I(s)N(s)+β2(s)S(s)G(s)M]es0μ(ρ)dρdstττ[β1(s)S(s)I(s)N(s)+β2(s)S(s)G(s)M]es0μ(ρ)dρds=ttτ[β1(s)S(s)I(s)N(s)+β2(s)S(s)G(s)M]es0μ(ρ)dρds0τ[β1(s)S(s)I(s)N(s)+β2(s)S(s)G(s)M]es0μ(ρ)dρds.
    $

    Hence, if $ E(0) = \int^0_{-\tau} [\frac{\beta_1(s)S(s)I(s)}{N(s)} +\frac{\beta_2(s)S(s)G(s)}{M}]e^{\int^{s}_0\mu(\rho)d\rho}ds $ is satisfied, we then have

    $ E(t)=ttτ[β1(s)S(s)I(s)N(s)+β2(s)S(s)G(s)M]etsμ(ρ)dρds.
    $
    (2.4)

    We see that $ u_2(t)\geq0 $, $ \forall t\in[0, m] $ whenever $ u_i(t)\geq0 $ for all $ i\neq2 $ and $ t\in[0, m]\subseteq[0, \sigma_\phi) $. Let $ \phi = (\phi_1, \phi_2, \phi_3, \phi_4)\in \mathcal{X}_{\delta_0} $ be given. If $ \phi_i(0) = 0 $ for some $ i\in\{1, 3, 4\} $, then $ f_i(t, \phi)\geq0 $. If $ \phi_4(0) = M $, then $ f_4(t, \phi)\leq0 $. By [22,Theorem 5.2.1] and its proof, it follows that for any $ \phi\in \mathcal{X}_{\delta_0} $, $ u_i(t, \phi)\geq0 $ for $ i = 1, 3, 4 $ for all $ t\in [0, \sigma_\phi) $. By equation (2.4), we have $ E(t)\geq0 $ for all $ t\in [0, \sigma_\phi) $. Therefore, it follows that for any $ \phi\in \mathcal{X}_{\delta_0} $, system (2.1) has a unique nonnegative solution $ u(t, \phi) $ with $ u_0 = \phi $ satisfies $ u(t, \phi)\in W $ for all $ t\in[0, \sigma_\phi). $

    Note that $ \frac{dN(t)}{dt}\geq\overline{B}-\widehat{\mu}N(t) $. For the system

    $ \frac{dy}{dt} = \overline{B}-\widehat{\mu}y(t), $

    the equilibrium $ \frac{\overline{B}}{\widehat{\mu}} $ is globally asymptotically stable. For any $ 0 < \delta_0 < \frac{\overline{B}}{\widehat{\mu}} $, $ \frac{dy}{dt}|_{y = \delta_0} = \overline{B}-\widehat{\mu}\delta_0 > 0 $. So if $ y(0)\geq\delta_0 $, then $ y(t)\geq \delta_0 $ for any $ t\geq0 $. By the comparison principle,

    $ N(t)\geq \delta_0\, \, \mbox{if}\, \, N(0) = \sum\limits^3_{i = 1}\phi_i(0)\geq \delta_0. $

    This implies that $ u(t, \phi)\in \mathcal{X}_{\delta_0} $ for all $ t\in[0, \sigma_\phi). $ Also we have

    $ \dfrac{dG(t)}{dt}\leq a(t)(M-G(t))-\delta(t)G(t) = a(t)M-(a(t)+\delta(t))G(t) $

    for all $ t\in[0, \sigma_\phi). $ Thus, both $ N(t) $ and $ G(t) $ are bounded on $ [0, \sigma_\phi) $, it follows that $ \sigma_\phi = \infty $ (see [21,Theorem 2.3.1]), and that all solutions are ultimately bounded. Moreover, when $ N(t) > \max\{\frac{\widehat{B}}{\overline{\mu}}, \frac{\widehat{a}M}{\overline{a}+\overline{\delta}}\} $ and $ G(t) > \max\{\frac{\widehat{B}}{\overline{\mu}}, \frac{\widehat{a}M}{\overline{a}+\overline{\delta}}\} $, we have

    $ \frac{dN(t)}{dt} \lt 0 \, \, \mbox{and}\, \, \frac{dG(t)}{dt} \lt 0. $

    This implies that all solutions are uniformly bounded.

    It is easy to see that system (2.1) has a unique virus-free periodic solution $ E_0(t) = (N^*(t), 0, 0, 0) $, where $ N^*(t) $ is the positive periodic solution of (2.3). Linearizing system (2.1) at its virus-free periodic solution $ E_0(t) = (N^*(t), 0, 0, 0) $, we then obtain the following system of periodic linear equations for the infective variables $ E, I, $ and $ G $

    $ {dE(t)dt=β1(t)I(t)+β2(t)N(t)MG(t)μ(t)E(t)[β1(tτ)I(tτ)+β2(tτ)N(tτ)MG(tτ)]ettτμ(r)dr,dI(t)dt=[β1(tτ)I(tτ)+β2(tτ)N(tτ)MG(tτ)]ettτμ(r)drμ(t)I(t),dG(t)dt=a(t)MN(t)I(t)δ(t)G(t).
    $
    (2.5)

    Since the first equation of system (2.5) is decoupled from the second and third equations of system (2.5), it suffices to use the following system to define the basic reproduction number

    $ {dI(t)dt=[β1(tτ)I(tτ)+β2(tτ)N(tτ)MG(tτ)]ettτμ(r)drμ(t)I(t),dG(t)dt=a(t)MN(t)I(t)δ(t)G(t).
    $
    (2.6)

    Let $ F: \mathbb{R}\rightarrow \mathcal{L}(C, \mathbb{R}^2) $ be a map and $ V(t) $ be a continuous $ 2\times2 $ matrix function on $ \mathbb{R} $ defined as follows

    $ F(t)ϕ=([β1(tτ)ϕ1(τ)+β2(tτ)N(tτ)Mϕ2(τ)]ettτμ(r)dra(t)MN(t)ϕ1(0)),
    $

    and

    $ V(t)=(μ(t)00δ(t)).
    $

    Let $ \Phi(t, s) $ be the evolution matrices of the linear ordinary differential system

    $ \frac{dy}{dt} = -V(t)y, $

    that is

    $ \frac{\partial \Phi(t, s)}{\partial t} = -V(t)\Phi(t, s), \forall t\geqslant s \, \, \mbox{and} \, \, \Phi(s, s) = I, \forall s\in \mathbb{R}, $

    where $ I $ is the $ 2\times2 $ identity matrix. It then easily follows that

    $ Φ(t,s)=(etsμ(r)dr00etsδ(r)dr).
    $

    Let $ C_\omega $ be the ordered Banach space of all continuous and $ \omega $-periodic functions from $ \mathbb{R} $ to $ \mathbb{R}^2 $, which is equipped with the maximum norm and the positive cone $ C^+_\omega = \{v\in C_\omega: v(t)\geqslant0, \forall t\in \mathbb{ R}\} $. Then we can define one linear operator on $ C_\omega $ by

    $ [Lv](t) = \int_0^\infty\Phi(t, t-s)F(t-s)v(t-s+\cdot)ds, \forall t\in \mathbb{R}, v\in C_\omega. $

    According to [23], the basic reproduction number is defined as $ R_0 = r(L) $, the spectral radius of $ L $.

    Let

    $ \mathcal{Y} = C([-\tau, 0], \mathbb{R}^2), \quad \mathcal{Y}_{+} = C([-\tau, 0], \mathbb{R}_{+}^2). $

    Then $ (\mathcal{Y}, \mathcal{Y}_{+}) $ is an ordered Banach space.

    Let $ P(t) $ be the solution maps of (2.6), that is, $ P(t)\phi = u_t(\phi) $, where $ u(t, \phi) $ is the unique solution of (2.6) with $ u_0 = \phi\in \mathcal{Y} $. Then $ P: = P(\omega) $ is the Poincar$ \acute{e} $ map associated with (2.6). Let $ r(P) $ be the spectral radius of $ P $. By [23,Theorem 2.1], we have the following result:

    Lemma 2.2. $ R_0-1 $ has the same sign as $ r(P)-1 $.

    The following lemma indicates that the periodic semiflow $ P(t) $ is eventually strongly monotone.

    Lemma 2.3. For any $ \varphi $ and $ \psi $ in $ \mathcal{Y}_{+} $ with $ \varphi > \psi $ (that is, $ \varphi\geq\psi $ but $ \varphi\neq\psi $), the solutions $ v(t, \varphi) $ and $ v(t, \psi) $ of system (2.6) with $ v_0 = \varphi $ and $ v_0 = \psi $, respectively, satisfy $ v_i(t, \varphi) > v_i(t, \psi) $ for all $ t > \tau $, $ i = 1, 2 $, and hence, $ P(t)\varphi \gg P(t)\psi $ in $ \mathcal{Y}_{+} $ for all $ t > 2\tau $.

    Proof. For any $ \varphi, \psi\in \mathcal{Y}_{+} $ with $ \varphi\geq \psi $, let $ v(t, \varphi) $ and $ v(t, \psi) $ be the unique solutions of system (3.1) satisfying $ v_0 = \varphi $ and $ v_0 = \psi $, respectively. By [22,Theorem 5.1.1], we have $ v(t, \varphi)\geq v(t, \psi) $ for all $ t\geq0 $; that is, $ P(t):\mathcal{Y}_{+}\rightarrow\mathcal{Y}_{+} $ is monotone. Next we prove that $ P(t):\mathcal{Y}_{+}\rightarrow\mathcal{Y}_{+} $ is eventually strongly monotone. Let $ \varphi, \psi\in \mathcal{Y_{+}} $ satisfy $ \varphi > \psi $. Denote $ v(t, \varphi) = (\overline{y}_1(t), \overline{y}_2(t)) $ and $ v(t, \psi) = (y_1(t), y_2(t)) $.

    $ Claim\, 1. \, \, \, There \, \, exists\, \, t_0\in[0, \tau]\, \, such\, \, that\, \, \overline{y}_1(t) > y_1(t)\, \, for\, \, all\, \, t\geq t_0. $

    We first prove that $ \overline{y}_1(t_0) > y_1(t_0) $ for some $ t_0\in[0, \tau] $. Otherwise, we have $ \overline{y}_1(t) = y_1(t) $, $ \forall t\in [0, \tau] $, and hence $ \frac{d\overline{y}_1(t)}{dt} = \frac{dy_1(t)}{dt} $ for all $ t\in (0, \tau) $. Thus, we have

    $ \beta_1(t-\tau)[\overline{y}_1(t-\tau)-y_1(t-\tau)]+\frac{\beta_2(t-\tau)N^*(t-\tau)}{M}[\overline{y}_2(t-\tau)-y_2(t-\tau)] = 0, \forall t\in[0, \tau]. $

    Since $ P(t) $ is monotone, it follows that $ \overline{y}_1(t-\tau) = y_1(t-\tau) $ and $ \overline{y}_2(t-\tau) = y_2(t-\tau) $ for all $ t\in [0, \tau] $, that is, $ \varphi_1(\theta) = \psi_1(\theta) $ and $ \varphi_2(\theta) = \psi_2(\theta) $ for all $ \theta\in [-\tau, 0] $, a contradiction to the assumption that $ \varphi > \psi $.

    Let $ g_1(t, y) = \left[\beta_1(t-\tau)y_1(t-\tau)+\frac{\beta_2(t-\tau)N^*(t-\tau)}{M}y_2(t-\tau)\right]e^{-\int^t_{t-\tau}\mu(r)dr}-\mu(t)y $. Since

    $ d¯y1(t)dt=[β1(tτ)¯y1(tτ)+β2(tτ)N(tτ)M¯y2(tτ)]ettτμ(r)drμ(t)¯y1(t)[β1(tτ)y1(tτ)+β2(tτ)N(tτ)My2(tτ)]ettτμ(r)drμ(t)¯y1(t)=g1(t,¯y1(t)),
    $

    we have $ \frac{d\overline{y}_1(t)}{dt}-g_1(t, \overline{y}_1(t))\geq0 = \frac{dy_1(t)}{dt}-g_1(t, y_1(t)), \forall t\geq t_0 $. Since $ \overline{y}_1(t_0) > y_1(t_0) $, the comparison theorem for ordinary differential equations [24,Theorem 4] implies that $ \overline{y}_1(t) > y_1(t) $ for all $ t\geq t_0 $.

    $ Claim\, 2. \, \, \, \overline{y}_2(t) > y_2(t)\, \, for\, \, all\, \, t > t_0. $

    Let $ g_2(t, y) = \frac{a(t)M}{N^*(t)}y_1(t)-\delta(t)y $. Then we have

    $ d¯y2(t)dt=a(t)MN(t)¯y1(t)δ(t)¯y2(t)>a(t)MN(t)y1(t)δ(t)¯y2(t)=g2(t,¯y2(t))t>t0,
    $

    and hence, $ \frac{d\overline{y}_2(t)}{dt}-g_2(t, \overline{y}_2(t)) > 0 = \frac{dy_2(t)}{dt}-g_2(t, y_2(t)) $ $ \forall t > t_0 $. Since $ \overline{y}_2(t_0)\geq y_2(t_0) $, it follows from [24,Theorem 4] that $ \overline{y}_2(t) > y_2(t) $ for all $ t > t_0 $.

    In view of the above two claims, we obtain

    $ (\overline{y}_1(t), \overline{y}_2(t))\gg (y_1(t), y_2(t)), \forall t \gt t_0. $

    Since $ t_0\in [0, \tau] $, it follows that $ (\overline{y}_{1t}, \overline{y}_{2t})\gg (y_{1t}, y_{2t}) $, $ \forall t > 2\tau $, that is, $ v_t(\varphi)\gg v_t(\psi) $ for all $ t > 2\tau $. This shows $ P(t): \mathcal{Y}_{+}\rightarrow \mathcal{Y}_{+} $ is strongly monotone for any $ t > 2\tau $.

    First we show that the virus will be endemic if $ R_0 > 1 $. Let

    $ X0={ϕ=(ϕ1,ϕ2,ϕ3,ϕ4)Xδ0:ϕ3(0)>0,ϕ4(0)>0},X0=Xδ0X0={ϕXδ0:ϕ3(0)=0orϕ4(0)=0}.
    $

    Theorem 3.1. Assume that $ R_0 > 1 $, then there exists a positive constant $ \eta > 0 $ such that any solution $ (S(t, \phi), E(t, \phi), I(t, \phi), G(t, \phi)) $ of system (2.1) with $ \phi \in X_0 $ satisfies

    $ \liminf\limits_{t\rightarrow\infty}(I(t, \phi), G(t, \phi))\geq(\eta, \eta). $

    Proof. Let $ Q(t) $ be the solution maps of (2.1) on $ \mathcal{X}_{\delta_0} $, that is, $ Q(t)\phi = u_t(\phi) $, $ t\geq0 $, where $ u(t, \phi) $ is the unique solution of (2.1) satisfying $ u_0 = \phi\in \mathcal{X}_{\delta_0} $. Then $ Q: = Q(\omega) $ is the Poincar$ \acute{e} $ map associated with (2.1).

    From the third, fourth equations of system (2.1) and equation (2.2), it is easy to see that $ Q(t)X_0\subseteq X_0 $ for all $ t\geq0 $. Lemma 2.1 implies that $ \{Q^n: \mathcal{X}_{\delta_0}\rightarrow \mathcal{X}_{\delta_0}\}_{n\geq0} $ is point dissipative and $ Q^n $ is compact for sufficiently large $ n $. It then follows from [25,Theorem 2.9] that $ Q $ admits a global attractor in $ \mathcal{X}_{\delta_0} $. Now we prove that $ Q $ is uniformly persistent with respect to $ (X_0, \partial X_0) $.

    Let $ M_1 = (N^*_0, 0, 0, 0) $, where $ N^*_0(\theta) = N^*(\theta) $ for all $ \theta\in[-\tau, 0]. $ Since $ \lim_{\phi\rightarrow M_1}\|Q(t)\phi-Q(t)M_1\| = 0 $ uniformly for $ t\in[0, \omega] $, for any given $ \varepsilon > 0 $, there exists a $ \eta_1 > 0 $ such that for any $ \phi\in X_0 $ with $ \|\phi-M_1\| < \eta_1 $, we have $ \|Q(t)\phi-Q(t)M_1\| < \varepsilon $ for all $ t\in[0, \omega] $.

    $ Claim \, 3. \quad\limsup_{n\to\infty}\|Q(n\omega)\phi-M_1\|\geq\eta_1\, \, for\, \, all\, \, \phi\in X_0. $

    Suppose, by contradiction, then $ \limsup_{n\to\infty}\|Q(n\omega)\psi-M_1\| < \eta_1 $ for some $ \psi\in X_0 $. Then there is an integer $ N_1\geq1 $ such that $ \|Q(n\omega)\psi-M_1\| < \eta_1 $ for all $ n\geq N_1 $. For $ t\geq N_1\omega $, we have $ t = n\omega+t_1 $, with $ n\geq N_1 $, $ t_1\in[0, \omega) $ and $ \|Q(t)\psi-Q(t)M_1\| = \|Q(t_1)Q(n\omega)\psi-Q(t_1)Q(n\omega)M_1\| = \|Q(t_1)Q(n\omega)\psi-Q(t_1)M_1\| < \varepsilon $ for all $ t\geq N_1\omega $. Therefore, $ N^*(t)-\varepsilon < S(t) < N^*(t)+\varepsilon $, $ 0\leq E(t) < \varepsilon $, $ 0 < I(t) < \varepsilon $, $ 0 < G(t) < \varepsilon $ for all $ t\geq N_1\omega $.

    let $ P_\varepsilon(t) $ be the solution maps of the following system on $ \mathcal{Y}_{+} $:

    $ {d˘I(t)dt=[β1(tτ)(N(tτ)ε)N(tτ)+3ε˘I(tτ)+β2(tτ)(N(tτ)ε)M˘G(tτ)]ettτμ(r)drμ(t)˘I(t),d˘G(t)dt=a(t)(Mε)N(t)+3ε˘I(t)δ(t)˘G(t).
    $
    (3.1)

    and $ P_\varepsilon: = P_\varepsilon(\omega) $. Since $ R_0 > 1 $, then $ \lim_{\varepsilon\rightarrow0^+}r(P_\varepsilon) = r(P) > 1 $. Then we can fix a sufficiently small $ \varepsilon > 0 $, such that $ r(P_\varepsilon) > 1 $, and $ N^*(t)-\varepsilon > 0 $, $ M-\varepsilon > 0 $ for all $ t\geq 0 $. By the arguments similar to system (2.6), it is easy to verify that $ P_\varepsilon(t) $ is strongly monotone on $ \mathcal{Y}_{+} $ for $ t > 2\tau $. It follows from [21,Theorem 3.6.1] that the linear operator $ P_\varepsilon(t) $ is compact on $ \mathcal{Y}_{+} $. Choose an integer $ n_0 > 0 $ such that $ n_0\omega > 2\tau $. Since $ P^{n_0}_\varepsilon = P_\varepsilon(n_0\omega) $, [26,Lemma 3.1] implies that $ r(P_\varepsilon) $ is a simple eigenvalue of $ P_\varepsilon $ having a strongly positive eigenvector. It then follows from [27,Lemma 1] that there is a positive $ \omega $-periodic function $ v^*(t) = (v^*_1(t), v^*_2(t)) $ such that $ u^*_\varepsilon(t) = e^{\frac{\ln r(P_\varepsilon)}{\omega}t}v^*(t) $ is a positive solution of system (3.1).

    For all $ t\geq N_1\omega+\tau $, by system (2.1), we have

    $ {dI(t)dt[β1(tτ)(N(tτ)ε)N(tτ)+3εI(tτ)+β2(tτ)(N(tτ)ε)MG(tτ)]ettτμ(r)drμ(t)I(t),dG(t)dta(t)(Mε)N(t)+3εI(t)δ(t)G(t).
    $

    Since $ I(t, \psi) > 0 $, $ G(t, \psi) > 0 $ for all $ t\geq0 $, we can choose a sufficiently small $ k > 0 $ such that

    $ I(t, \psi), G(t, \psi)\geq ku^*_\varepsilon(t), \forall t\in [N_1\omega+\tau, N_1\omega+2\tau]. $

    By [22,Theorem 5.1.1] it follows that $ I(t, \psi), G(t, \psi)\geq ku^*_\varepsilon(t), \forall t\geq N_1\omega+2\tau. $ Thus we have $ I(t, \psi)\rightarrow \infty $, $ G(t, \psi)\rightarrow \infty $, a contradiction.

    Define

    $ M_\partial: = \{\phi\in \partial X_0: Q^n(\phi)\in \partial X_0, \forall n\geq0\}. $

    For any given $ \psi\in \partial X_0 $, we have $ \psi_3(0) = 0 $ or $ \psi_4(0) = 0 $.

    If $ \psi_3(0) = 0 $, we have the following two cases.

    Case 1. $ I(t, \psi) = 0 $ for all $ t\geq0 $.

    From the third equation of (2.1), we have $ G(t-\tau, \psi) = 0 $ for all $ t\geq \tau $. Then from equation (2.2), we get $ E(t, \psi) = 0 $ for all $ t\geq \tau $. In this case, $ Q^n(\psi)\rightarrow M_1 $ as $ n\rightarrow\infty $.

    Case 2. $ I(t_2, \psi) > 0 $ for some $ t_2 > 0 $.

    From the third and fourth equations of (2.1), we have $ I(t, \psi) > 0 $ and $ G(t, \psi) > 0 $ for all $ t\geq t_2 $. Thus $ Q^n(\psi)\in X_0 $ for $ n\omega > t_2 $.

    For the case $ \psi_4(0) = 0 $, we can do similar analysis. Finally, we see that for any $ \psi\in M_\partial $, $ Q^n(\psi)\rightarrow M_1 $ as $ n\rightarrow\infty $. Thus $ \omega(\psi) = M_1 $ for any $ \psi\in M_\partial $, and $ M_1 $ can not form a cycle in $ \partial X_0 $.

    By Claim 3, we see that $ M_1 $ is an isolated invariant set for $ Q $ in $ X $, and $ W^s(M_1)\cap X_0 = \emptyset $, where $ W^s(M_1) $ is the stable set of $ M_1 $ for $ Q $. By the acyclicity theorem on uniform persistence for maps (see [28,Theorem 1.3.1 and Remark 1.3.1], it follows that $ Q : X\rightarrow X $ is uniformly persistent with respect to $ X_0 $. Thus, [28,Theorem 3.1.1] implies that the periodic semiflow $ Q(t) : X\rightarrow X $ is also uniformly persistent with respect to $ X_0 $.

    It remains to prove the practical uniform persistence. By [25,Theorem 4.5], with $ \rho(x) = d(x, \partial X_0) $, it then follows that $ Q: X_0\rightarrow X_0 $ has a compact global attractor $ A_0 $. For any $ \phi\in A_0 $, we have $ \phi_i(0) > 0 $ for all $ i = 3, 4 $. Let $ B_0 = \cup_{t\in[0, \omega]}Q(t)A_0 $. Then $ \psi_i(0) > 0 $, $ i = 3, 4 $, for all $ \psi\in B_0 $. It is easy to see that $ \lim_{t\rightarrow\infty}d(Q(t)\phi, B_0) = 0 $ for all $ \phi\in X_0 $. Let $ X_{+} = C([-\tau, 0], \mathbb{R}^{4}_{+}) $, and define a continuous function $ p: X_{+}\rightarrow \mathbb{R}_{+} $ by

    $ p(\phi) = \min\{\phi_3(0), \phi_4(0)\}, \forall \phi \in X_{+}. $

    Clearly, $ p(\phi) > 0 $ for all $ \phi\in B_0 $. Since $ B_0 $ is a compact subset of $ X_+ $, we have $ \inf_{\phi\in B_0}p(\phi) = \min_{\phi\in B_0}p(\phi) > 0 $. By the attractiveness of $ B_0 $, it then follows that there exists $ \eta > 0 $ such that

    $ \liminf\limits_{t\rightarrow\infty}\min(I(t, \phi), G(t, \phi)) = \liminf\limits_{t\rightarrow\infty}p(Q(t)\phi)\geq \eta, \forall\phi\in X_0. $

    This completes the proof.

    The following theorem shows that the virus will be cleared from the population if $ R_0 < 1 $.

    Theorem 3.2. If $ R_0 < 1 $, then the virus-free periodic solution $ E_0(t) = (N^*(t), 0, 0, 0) $ is globally attractive for system (2.1) in $ \mathcal{X}_{\delta_0} $.

    Proof. By the aforementioned analysis for (2.3), we know that (2.3) has a positive $ \omega $-periodic solution $ N^*(t) $ which is globally asymptotically stable. It then follows that for any $ \epsilon > 0 $, we can choose a sufficiently large integer $ n_1 > 0 $ such that $ n_1\omega > \tau $ and $ N^*(t)-\epsilon < N(t) = S(t)+E(t)+I(t) < N^*(t)+\epsilon $ for $ t > n_1\omega-\tau $. Thus, for $ t > n_1\omega $, we have

    $ {dI(t)dt[β1(tτ)I(tτ)+β2(tτ)(N(tτ)+ϵ)MG(tτ)]ettτμ(r)drμ(t)I(t),dG(t)dta(t)MN(t)ϵI(t)δ(t)G(t).
    $

    It then suffices to show that positive solutions of the auxiliary system

    $ {d˜I(t)dt=[β1(tτ)˜I(tτ)+β2(tτ)(N(tτ)+ϵ)M˜G(tτ)]ettτμ(r)drμ(t)˜I(t),d˜G(t)dt=a(t)MN(t)ϵ˜I(t)δ(t)˜G(t)
    $
    (3.2)

    tend to zero when $ t $ tends to infinity. Let $ P_\epsilon(t) $ be the solution maps of system (3.2) on $ \mathcal{Y} $, and $ P_\epsilon: = P_\epsilon(\omega) $. Since $ R_0 < 1 $, then by Lemma 2.2, we have $ \lim_{\epsilon\rightarrow 0^+} r(P_\epsilon) = r(P) < 1 $. Thus we can fix $ \epsilon > 0 $ small enough such that $ r(P_\epsilon) < 1 $. According to [27,Lemma 1], there is a positive $ \omega $-periodic function $ \varsigma(t) = (\varsigma_1(t), \varsigma_2(t)) $, such that $ u_\epsilon(t) = e^{\frac{\ln r(P_\epsilon)}{\omega}t}\varsigma(t) $ is a positive solution of system (3.2).

    Let $ l > 0 $ large enough, such that $ (I(t), G(t))\leq lu_\epsilon(t) $ for $ t\in[n_1\omega, n_1\omega+\tau] $. Then [22,Theorem 5.1.1] implies that $ (I(t), G(t))\leq lu_\epsilon(t) $ for all $ t\geq n_1\omega+\tau $. Hence, $ I(t)\rightarrow0 $, $ G(t)\rightarrow0 $ as $ t\rightarrow\infty $. It then follows from the theory of asymptotically periodic semiflow (see [28,Theorem 3.2.1]) that

    $ \lim\limits_{t\rightarrow\infty}(S(t)-N^*(t)) = 0, \, \lim\limits_{t\rightarrow\infty}E(t) = 0. $

    This completes the proof.

    In this section, we study the corresponding autonomous system of system (2.1), that is, we suppose that all coefficients in system (2.1) are time independent and positive. In this case, system (2.1) becomes

    $ {dS(t)dt=BμS(t)β1S(t)I(t)N(t)β2S(t)G(t)M,dE(t)dt=β1S(t)I(t)N(t)+β2S(t)G(t)MμE(t)[β1S(tτ)I(tτ)N(tτ)+β2S(tτ)G(tτ)M]eμτ,dI(t)dt=[β1S(tτ)I(tτ)N(tτ)+β2S(tτ)G(tτ)M]eμτμI(t),dG(t)dt=aI(t)N(t)(MG(t))δG(t).
    $
    (4.1)

    For system (4.1), there is always the virus-free equilibrium $ P_0 = (N^*, 0, 0, 0) $ with $ N^* = \frac{B}{\mu} $.

    Let

    $ F1=(β1eμτβ2NeμτM00),F2=(00aMN0)andV=(μ00δ).
    $

    It then follows that $ F(\phi) = F_1\phi(-\tau)+F_2\phi(0) $. By [23,Corollary 2.1], we can define the basic reproduction number for system (4.1) as $ \mathbf{R_0} = r(\hat{F}V^{-1}) $, where $ \hat{F} = F_1+F_2 $. It then easily follows that

    $ ˆFV1=(β1eμτμβ2NeμτδMaMμN0)
    $

    and

    $ \mathbf{R_0} = \frac{1}{2}\left(C_1+\sqrt{C_1^2+4C_0}\right), $

    where $ C_0 = \frac{a\beta_2e^{-\mu\tau}}{\mu\delta} $, $ C_1 = \frac{\beta_1e^{-\mu\tau}}{\mu} $. We define $ \Theta = \frac{\beta_1\delta+\beta_2a}{\mu\delta}e^{-\mu\tau}, $ it is easy to verify that $ \mathbf{R_0}-1 $ has the same sign as $ \Theta-1 $.

    Now consider endemic equilibria with $ I = I^* > 0 $. Let the right-hand sides be zero, then system (4.1) admits another equilibrium: $ P^* = (S^*, E^*, I^*, G^*) $, where

    $ S=NeμτI,E=(eμτ1)I,G=aMIaI+δN,
    $

    and $ I^* $ is a positive real root of equation $ g_1(I) = g_2(I) $, where

    $ g1(I)=β1N+β2aaI+δN,g2(I)=μeμτNeμτI.
    $

    Note that $ g_1(0) = \frac{\beta_1}{N^*}+\frac{\beta_2a }{\delta N^*} $, $ g_1(I) $ is decreasing with respect to $ I > 0 $ and $ \lim_{I\rightarrow+\infty}g_1(I) = \frac{\beta_1}{N^*} $. Also $ g_2(0) = \frac{\mu e^{\mu\tau}}{N^*} $, $ g'_2(I) > 0 $ for $ I\in(0, N^* e^{-\mu\tau}) $, and $ \lim_{x\rightarrow (N^* e^{-\mu\tau})^-} g_2(x) = +\infty $. Then $ g_1(I) = g_2(I) $ has a unique positive root in $ (0, N^* e^{-\mu\tau}) $ if and only if $ g_1(0) > g_2(0) $, which is equivalent to $ \mathbf{R_0} > 1 $. Therefore, if $ \mathbf{R_0} > 1 $, system (4.1) has a unique endemic equilibrium $ P^* $.

    Define

    $ Ωδ0={ϕW:3i=1ϕi(s)δ0,s[τ,0],ϕ2(0)=0τ[β1ϕ1(ξ)ϕ3(ξ)3i=1ϕi(ξ)+β2ϕ1(ξ)ϕ4(ξ)M]eμξdξ}.
    $

    The following theorem provides the global stability of the virus-free equilibrium.

    Theorem 4.1. If $ \mathbf{R_0} < 1 $, then the virus-free equilibrium $ P_0 $ is globally asymptotically stable for system (4.1) in $ \Omega_{\delta_0} $, and is unstable if $ \mathbf{R_0} > 1 $.

    Proof. Linearizing system (4.1) at the virus-free equilibrium $ P_0 $, we obtain the characteristic equation

    $ (λ+μ)2g(λ,τ)=0,
    $

    with

    $ g(\lambda, \tau) = \lambda^2+(\mu+\delta)\lambda-e^{-\mu\tau}e^{-\lambda\tau}(\beta_1\lambda+\beta_1\delta+\beta_2a)+\mu\delta. $

    It is clear that the local stability of $ P_0 $ is determined by the roots of $ g(\lambda, \tau) = 0. $ If $ \tau = 0, $ then $ g(\lambda, 0) = \lambda^2+(\mu+\delta-\beta_1)\lambda+\mu\delta-(\beta_1\delta+\beta_2a). $ If $ \mathbf{R_0} < 1 $, then $ \Theta = \frac{\beta_1\delta+\beta_2a}{\mu\delta} < 1, $ which implies that $ \mu+\delta-\beta_1 > 0 $ and $ \mu\delta-(\beta_1\delta+\beta_2a) > 0 $. Then by the Routh-Hurwitz criterion, all the roots of $ g(\lambda, 0) = 0 $ have negative real parts, while if $ \mathbf{R_0} > 1 $, then $ \mu\delta-(\beta_1\delta+\beta_2a) < 0 $, and $ g(\lambda, 0) = 0 $ has one positive real root. That is, for $ \tau = 0 $, $ P_0 $ is locally asymptotically stable for $ \mathbf{R_0} < 1 $, and is unstable if $ \mathbf{R_0} > 1 $.

    Next we will show that $ g(\lambda, \tau) = 0 $ has no pure imaginary roots if $ \mathbf{R_0} < 1 $. Obviously, $ i\rho $ ($ \rho\in \mathbb{R} $) is a root of $ g(\lambda, \tau) = 0 $ if and only if $ \rho $ satisfies

    $ -\rho^{2}+(\mu+\delta)\rho i-e^{-\mu\tau}(\cos\rho\tau-i\sin\rho\tau)(\beta_1\rho i+\beta_1\delta+\beta_2a)+\mu\delta = 0. $

    Separating the real and the imaginary parts, we have

    $ {μδρ2=β1eμτρsinρτ+(β1δ+β2a)eμτcosρτ,(μ+δ)ρ=β1eμτρcosρτ(β1δ+β2a)eμτsinρτ.
    $

    Squaring the two equations and adding them gives

    $ z2+(μ2+δ2β21e2μτ)z+μ2δ2(1Θ2)=0,
    $
    (4.2)

    where $ z = \rho^2 $. If $ \mathbf{R_0} < 1 $, we have $ \beta_1e^{-\mu\tau} < \mu $, hence, (4.2) has no nonnegative root, which implies that $ g(\lambda, \tau) = 0 $ has no zero root and pure imaginary roots. Therefore, $ P_0 $ is locally asymptotically stable for all $ \tau\geq0 $ if $ \mathbf{R_0} < 1 $. Together with Theorem 3.2, we can show that $ P_0 $ is globally asymptotically stable in $ \Omega_{\delta_0} $ when $ \mathbf{R_0} < 1 $.

    We let

    $ S^\infty = \limsup\limits_{t\rightarrow\infty}S(t), \quad S_\infty = \liminf\limits_{t\rightarrow\infty} S(t). $

    We can define $ I^\infty, I_\infty, G^\infty, G_\infty $ in a similar way. Now we will use the method of fluctuations (see, e.g., [29]) to show the global attractivity of the endemic equilibrium $ P^* $. To do this, we need the following assumption:

    (H) $ (e^{-\mu\tau}a-\delta)\beta_2a > \beta_1(a+\delta)^2. $

    The condition (H) is technically needed in the arguments supporting the global attractivity of the endemic equilibrium $ P^* $ for system (4.1).

    Theorem 4.2. Let $ \mathrm{(H)} $ hold. If $ \mathbf{R_0} > 1 $, then (4.1) has a unique endemic equilibrium $ P^* = (S^*, E^*, I^*, G^*) $ such that $ \lim_{t\rightarrow\infty}u(t, \phi) = P^* $ for any $ \phi\in \Omega_{\delta_0} $ with $ \phi_3(0) > 0 $ and $ \phi_4(0) > 0 $.

    Proof. The whole mouse population satisfies the following equation:

    $ \frac{dN(t)}{dt} = B-\mu N(t). $

    Then $ N^* $ is globally asymptotically stable. Hence, we have the following limiting system:

    $ {dS(t)dt=BμS(t)˜β1S(t)I(t)~β2S(t)G(t),dE(t)dt=˜β1S(t)I(t)+~β2S(t)G(t)μE(t)[˜β1S(tτ)I(tτ)+˜β2S(tτ)G(tτ)]eμτ,dI(t)dt=[˜β1S(tτ)I(tτ)+˜β2S(tτ)G(tτ)]eμτμI(t),dG(t)dt=˜aI(t)(MG(t))δG(t),
    $
    (4.3)

    where $ \widetilde{\beta}_1 = \frac{\beta_1}{N^*} $, $ \widetilde{\beta}_2 = \frac{\beta_2}{M} $, $ \widetilde{a} = \frac{a}{N^*} $. Since the second equation of system (4.3) is decoupled from the other three equations of system (4.3), we then consider the following system:

    $ {dS(t)dt=BμS(t)˜β1S(t)I(t)~β2S(t)G(t),dI(t)dt=[˜β1S(tτ)I(tτ)+˜β2S(tτ)G(tτ)]eμτμI(t),dG(t)dt=˜aI(t)(MG(t))δG(t).
    $
    (4.4)

    From the proof of Lemma 2.1, we can show that

    $ \Gamma = C\left([-\tau, 0], \left[0, \frac{B}{\mu}\right]^2\times[0, M]\right) $

    is positively invariant for system (4.4). By the arguments similar to those in Theorem 3.1, it is easy to verify that system (4.4) is uniformly persistent in the sense that there exists a $ \eta_2 > 0 $ such that for any given $ \psi = (\psi_1, \psi_2, \psi_3)\in \Gamma $ with $ \psi_2(0) > 0 $, $ \psi_3(0) > 0 $, the solution $ (S(t, \psi), I(t, \psi), G(t, \psi)) $ of (4.4) satisfies $ \liminf_{t\rightarrow\infty}(I(t, \psi), G(t, \psi))\geq (\eta_2, \eta_2). $

    For any given $ \psi\in \Gamma $ with $ \psi_2(0) > 0 $, $ \psi_3(0) > 0 $, let $ (S(t), I(t), G(t)) = (S(t, \psi), I(t, \psi), G(t, \psi)) $. It is clear that $ S^\infty\geq S_\infty $, $ \frac{B}{\mu}\geq I^\infty\geq I_\infty\geq \eta_2 > 0 $, $ G^\infty\geq G_\infty\geq\eta_2 > 0 $. Moreover, there exist sequences $ t^i_n $ and $ \sigma^i_n $, $ i = 1, 2, 3 $, such that

    $ \lim\limits_{n\rightarrow \infty} S(t^1_n) = S^\infty, \, S'(t^1_n) = 0, \, \lim\limits_{n\rightarrow \infty}S(\sigma^1_n) = S_\infty, \, S'(\sigma^1_n) = 0, \forall n\geq1; $
    $ \lim\limits_{n\rightarrow \infty} I(t^2_n) = I^\infty, \, I'(t^2_n) = 0, \, \lim\limits_{n\rightarrow \infty}I(\sigma^2_n) = I_\infty, \, I'(\sigma^2_n) = 0, \forall n\geq1; $
    $ \lim\limits_{n\rightarrow \infty} G(t^3_n) = G^\infty, \, G'(t^3_n) = 0, \, \lim\limits_{n\rightarrow \infty}G(\sigma^3_n) = G_\infty, \, G'(\sigma^3_n) = 0, \forall n\geq1. $

    It then follows from the first equation of (4.4), we have

    $ B-\mu S^\infty-\widetilde{\beta}_1S^\infty I_\infty-\widetilde{\beta_2}S^\infty G_\infty\geq0 \geq B-\mu S^\infty-\widetilde{\beta}_1S^\infty I^\infty-\widetilde{\beta_2}S^\infty G^\infty, $

    and

    $ B-\mu S_\infty-\widetilde{\beta}_1S_\infty I_\infty-\widetilde{\beta_2}S_\infty G_\infty\geq0 \geq B-\mu S_\infty-\widetilde{\beta}_1S_\infty I^\infty-\widetilde{\beta_2}S_\infty G^\infty, $

    we see that

    $ Bμ+˜β1I+˜β2GSSBμ+˜β1I+˜β2G.
    $
    (4.5)

    By the second and third equations of (4.4), by a similar argument, it then follows that

    $ (˜β1SI+~β2SG)eμτμII(˜β1SI+~β2SG)eμτμ,
    $
    (4.6)

    and

    $ ˜aMI˜aI+δGG˜aMI˜aI+δ.
    $
    (4.7)

    Using (4.5) and (4.6), we get

    $ (˜β1I+~β2G)Beμτμ(μ+˜β1I+˜β2G)II(˜β1I+~β2G)Beμτμ(μ+˜β1I+˜β2G).
    $
    (4.8)

    Then (4.7) and (4.8) imply that

    $ \frac{Be^{-\mu\tau}}{\mu}(\widetilde{\beta}_1+\frac{\widetilde{\beta}_2\widetilde{a}M}{\widetilde{a}I^\infty+\delta}) \geq\mu+\widetilde{\beta}_1 I_\infty+\frac{\widetilde{\beta}_2\widetilde{a}MI_\infty}{\widetilde{a}I_\infty+\delta}, $
    $ \frac{Be^{-\mu\tau}}{\mu}(\widetilde{\beta}_1+\frac{\widetilde{\beta}_2\widetilde{a}M}{\widetilde{a}I_\infty+\delta}) \leq\mu+\widetilde{\beta}_1 I^\infty+\frac{\widetilde{\beta}_2\widetilde{a}MI^\infty}{\widetilde{a}I^\infty+\delta}. $

    Simplifying the above two inequalities, we get

    $ \widetilde{\beta}_1(\widetilde{a}I^\infty+\delta)(\widetilde{a}I_\infty+\delta)(I^\infty-I_\infty)\geq \widetilde{\beta}_2\widetilde{a}M\left(\frac{Be^{-\mu\tau}}{\mu}\widetilde{a}-\delta\right)(I^\infty-I_\infty). $

    Since $ \widetilde{\beta}_1(\widetilde{a}I^\infty+\delta)(\widetilde{a}I_\infty+\delta)\leq \widetilde{\beta}_1(\widetilde{a}\frac{B}{\mu}+\delta)^2 $, then if condition (H) holds, we have $ I^\infty = I_\infty $. By (4.5) and (4.7), we get $ S^\infty = S_\infty $ and $ G^\infty = G_\infty $. Then we have $ \lim_{t\rightarrow\infty}(S(t), I(t), G(t)) = (S^*, I^*, G^*) $ for any $ \psi\in \Gamma $ with $ \psi_2(0) > 0 $ and $ \psi_3(0) > 0 $. By the theory of chain transitive sets [28], we can lift the global attractivity for system (4.4) to system (4.1). It follows that $ \lim_{t\rightarrow\infty} u(t, \phi) = P^* $, for any $ \phi\in \Omega_{\delta_0} $ with $ \phi_3(0) > 0 $ and $ \phi_4(0) > 0 $.

    In this section, we give numerical simulations that support the theory presented in the previous sections. To compute the basic reproduction number $ R_0 $ numerically, we first write the operator $ L $ into the integral form in [30]. For $ v\in C_\omega $, we get

    $ [Lv](t)=0Φ(t,ts)F(ts)v(ts+)ds=0Φ(t,ts)×([β1(tsτ)v1(tsτ)+β2(tsτ)N(tsτ)Mv2(tsτ)]etstsτμ(r)dra(ts)MN(ts)v1(ts))ds=0([β1(tsτ)v1(tsτ)+β2(tsτ)N(tsτ)Mv2(tsτ)]ettsτμ(r)dra(ts)MN(ts)ettsδ(r)drv1(ts))ds=(τ[β1(ts)v1(ts)+β2(ts)N(ts)Mv2(ts)]ettsμ(r)drds0a(ts)MN(ts)ettsδ(r)drv1(ts)ds)=0K(t,s)v(ts)ds,
    $

    where

    $ K(t,s)=(β1(ts)ettsμ(r)drβ2(ts)N(ts)Mettsμ(r)dra(ts)MN(ts)ettsδ(r)dr0)ifsτ,
    $

    and

    $ K(t,s)=(00a(ts)MN(ts)ettsδ(r)dr0)ifs<τ.
    $

    By the $ \omega $-periodicity of $ v $, we obtain

    $ [Lv](t)=ω0G(t,s)v(ts)ds,
    $

    where $ G(t, s) = \sum^\infty_{n = 0} K(t, s+n\omega) $. Consequently, we can use the numerical method in [30] to compute $ R_0 $.

    The time unit is taken as month. Baseline parameters are $ \mu = 0.246 $, $ \beta = 0.0024 $, $ \alpha = 0.0016 $, $ \gamma = 0.5 $, $ d = 60 $, $ \delta = 9.21 $, $ M = 50 $. These parameters values are taken from [10]. To reflect the seasonality, we suppose that the birth rate of mice is $ B(t) = 5.85(1+0.8\sin(\pi t/6)) $, the direct contact rate and the indirect contact rate are equal, that is, $ c(t) = \bar{c}(t) = 60(1+0.8\sin(\pi t/6)) $. Since the incubation period is about one week, then we take $ \tau = \frac{1}{3} $.

    With this set of parameters, we have $ R_0 = 1.3586 > 1 $, the virus will be endemic and the infection is persistent in the mouse population (see Figure 1(a)). If we can increase the incubation period to $ \tau = 2 $, then $ R_0 = 0.9876 < 1 $. In this case, the long-term behavior of the infectious mice is shown in Figure 1(b), which implies that the virus will be eradicated from the mouse population. These results are coincident with Theorems 3.1 and 3.2.

    Figure 1.  Long-term behavior of the infectious mouse population (a) $ R_0 > 1 $ and (b) $ R_0 < 1 $.

    First, we discuss the influence of the incubation period. Let $ \tau $ varies and keep other parameters as above. Our numerical computations demonstrate that $ R_0 $ is a decreasing function (see Figure 2). Hence, we may try to prolong the incubation period via medical drugs or control measures to control hantavirus infection. For example, to eradicate the virus, we should keep $ \tau > 1.923 $.

    Figure 2.  The graph of the basic reproduction numbers $ R_0 $ and $ \mathbf{R}_0 $ when $ \tau $ varies.

    The dependence of $ \mathbf{R_0} $ on $ \tau $ is also illustrated in Figure 2. One can observe that for each fixed $ \tau $, $ \mathbf{R_0} $ is greater than or equal to $ R_0 $. Thus, the autonomous model may overestimate the value of $ R_0 $.

    Since contact rate plays a very important role in the spread of hantavirus. To explore the influence of direct contact rate and indirect contact rate on hantavirus transmission, we replace $ c(t) $ with $ \hat{c}(t) = (1-k)c(t) $ and $ \bar{c}(t) $ with $ \hat{\bar{c}}(t) = (1-l)\bar{c}(t) $ in our model. Then, Figure 3 shows $ R_0 $ is a decreasing function with respect to $ k $ and $ l $. Moreover, Figure 3(a) shows that the virus will be endemic even if $ k = 1 $ (that is, $ \hat{c}(t) = 0 $) for $ l = 0 $, Figure 3(b) shows that if we keep $ l > 0.63 $, then the virus can be eradicated for $ k = 0 $. For the parameters we have chosen, $ R_0 $ is more sensitive to the indirect transmission than the direct transmission. Thus, indirect transmission can not be neglected.

    Figure 3.  $ R_0 $ vs $ k $ and $ l $. (a) Relationship between $ R_0 $ and $ k $. (b) Relationship between $ R_0 $ and $ l $.

    In this paper, we formulate and analyze a compartmental model for hantavirus infection that incorporating the seasonality, incubation period of the mice, direct and indirect transmission. Using the theory developed in [23], we first introduce the basic reproduction number $ R_0 $ for the model. Then we show that $ R_0 $ is a threshold parameter for the persistence and extinction of the virus. More precisely, the virus will be endemic when $ R_0 > 1 $, and the virus will be cleared if $ R_0 < 1 $. For the corresponding autonomous system, we obtain an explicit expression of $ \mathbf{R_0} $ and establish a threshold result on the global stability in terms of $ \mathbf{R_0} $.

    Numerical simulations are performed to illustrate our analytic results. Figure 2 shows that an increase of the incubation period could reduce $ R_0 $, which implies that the virus infection can be relieved by prolonging the length of the incubation period. Moreover, our numerical results show that completely abolishing direct transmission would not eradicate virus (see Figure 3(a)), indirect transmission may play an important role in the virus transmission. On the other hand, Figure 3(b) shows that lowering the indirect transmission to $ 37\% $ of the original value would lead to eradication of the virus. Hence, outdoor disinfection is very important to control the transmission of hantavirus.

    In this paper, we ignore the age-structure of the mouse population, since juvenile mice must leave to find and establish their own home ranges, while the adults do not move, so it should be more reasonable to consider two types of mice: the adult and juvenile mice as noted in [9]. Also the diffusion of the juvenile mice plays an important role in hantavirus infection. We leave these problems for further investigation.

    This work was supported by the National Natural Science Foundation of China (No. 11801431), the Natural Science Basic Research Plan in Shaanxi Province of China (No. 2018JM1011), and China Scholarship Council.

    The author declares there is no conflict of interest.

    [1] Sedjo RA, Sohngen B (2013) Wood as a Major Feedstock for Biofuel Production in the United States: Impacts on Forests and International Trade. J Sustain For 32: 195-211. doi: 10.1080/10549811.2011.652049
    [2] Goh CS, Junginger M, Cocchi M, et al. (2013) Wood pellet market and trade: a global perspective. Biofuels Bioprod Bioref 7: 24-42. doi: 10.1002/bbb.1366
    [3] Immerzeel DJ, Verweij PA, van der Hilst F, et al. (2014) Biodiversity impacts of bioenergy crop production: A state-of-the-art review. GCB Bioenergy 6: 183-209. doi: 10.1111/gcbb.12067
    [4] McDonald RI, Fargione J, Kiesecker J, et al. (2009) Energy sprawl or energy efficiency: Climate policy impacts on natural habitat for the United States of America. PLoS One 4: e6802. doi: 10.1371/journal.pone.0006802
    [5] Stoms DM, Davis FW, Jenner MW, et al. (2012) Modeling wildlife and other trade-offs with biofuel crop production. GCB Bioenergy 4: 330-341. doi: 10.1111/j.1757-1707.2011.01130.x
    [6] Dale VH, Kline KL, Wiens J, et al. (2010) Biofuels: Implications for Land Use and Biodiversity. The Ecological Society of America Biofuels and Sustainability Reports. Available from: http://www.esa.org/biofuelsreports/files/ESA%20Biofuels%20Report_VH%20Dale%20et%20al.pdf.
    [7] Wiens J, Fargione J, Hill J (2011) Biofuels and biodiversity. Ecol Appl 21: 1085-1095. doi: 10.1890/09-0673.1
    [8] Daystar J (2014) Environmental Impacts of Cellulosic Biofuels Made in the South East: Implications of Impact Assessment Methods and Study Assumptions. North Carolina State University: 264 pages.
    [9] Wear D, Abt R, Alavalapati J, et al. (2010) The South's Outlook for Sustainable Forest Bioenergy and Biofuels Production. The Pinchot Institute Report. Available from: http://www.pinchot.org/uploads/download?fileId=512.
    [10] Fletcher RJ, Robertson BA, Evans J, et al. (2011) Biodiversity conservation in the era of biofuels: risks and opportunities. Front Ecol Environ 9: 161-168. doi: 10.1890/090091
    [11] Riffell S, Verschuyl J, Miller D, et al. (2011) Biofuel harvests, coarse woody debris, and biodiversity – A meta-analysis. For Ecol Manage 261: 878-887. doi: 10.1016/j.foreco.2010.12.021
    [12] Dale VH, Lowrance R, Mulholland P, et al. (2010) Bioenergy Sustainability at the Regional Scale. Ecol Soc 15: 23.
    [13] Wear DN, Huggett R, Li R, et al. (2013) Forecasts of Forest Conditions in U.S. Regions under Future Scenarios: A Technical Document Supporting the Forest Service 2010 RPA Assessment. Gen Tech Rep SRS-170.
    [14] Lubowski RN, Plantinga AJ, Stavins RN (2008) What Drives Land-Use Change in the United States? A National Analysis of Landowner Decisions. Land Econ 84: 529-550.
    [15] Daniel CJ, Frid L (2012) Predicting Landscape Vegetation Dynamics Using State-and-Transition Simulation Models. Proc First Landsc State-and-Transition Simul Model Conf June 14-16 2011: 5-22.
    [16] Bestelmeyer BT, Herrick JE, Brown JR, et al. (2004) Land management in the American southwest: a state-and-transition approach to ecosystem complexity. Environ Manage 34: 38-51.
    [17] Costanza JK, Hulcr J, Koch FH, et al. (2012) Simulating the effects of the southern pine beetle on regional dynamics 60 years into the future. Ecol Modell 244: 93-103. doi: 10.1016/j.ecolmodel.2012.06.037
    [18] Wilson T, Costanza J, Smith J, et al. (2014) Second State-and-Transition Simulation Modeling Conference. Bull Ecol Soc Am 96: 174-175.
    [19] Halofsky J, Halofsky J, Burscu T, et al. (2014) Dry forest resilience varies under simulated climate-management scenarios in a central Oregon, USA landscape. Ecol Appl 24: 1908-1925. doi: 10.1890/13-1653.1
    [20] Provencher L, Forbis TA, Frid L, et al. (2007) Comparing alternative management strategies of fire, grazing, and weed control using spatial modeling. Ecol Modell 209: 249-263. doi: 10.1016/j.ecolmodel.2007.06.030
    [21] Abt R, Cubbage F, Abt K (2009) Projecting southern timber supply for multiple products by subregion. For Prod J 59: 7-16.
    [22] Abt KL, Abt RC, Galik CS, et al. (2014) Effect of Policies on Pellet Production and Forests in the U.S. South: A Technical Document Supporting the Forest Service Update of the 2010 RPA Assessment. Gen Tech Rep GTR-SRS-202.
    [23] U.S. Geological Survey National Gap Analysis Program (2013) Protected Areas Database-US (PAD-US), Version 1.3. Available from: http://gapanalysis.usgs.gov/padus/.
    [24] Terando A, Costanza JK, Belyea C, et al. (2014) The southern megalopolis: using the past to predict the future of urban sprawl in the Southeast U.S. PLoS One 9: e102261. doi: 10.1371/journal.pone.0102261
    [25] Noss RF, Platt WJ, Sorrie BA, et al. (2015) How global biodiversity hotspots may go unrecognized: lessons from the North American Coastal Plain. Richardson D, ed. Divers Distrib 21: 236-244. doi: 10.1111/ddi.12278
    [26] Southeast Gap Analysis Project (SEGAP) (2008) Southeast GAP regional land cover [digital data]. Available from: http://www.basic.ncsu.edu/segap/.
    [27] Burke S, Hall BR, Shahbazi G, et al. (2007) North Carolina's Strategic Plan for Biofuels Leadership. Available from: http://www.ces.ncsu.edu/fletcher/mcilab/publications/NC_Strategic_Plan_for_Biofuels_Leadership.pdf.
    [28] Forisk Consulting LLC (2014) Wood bioenergy US database 2013. Available by subscription.
    [29] Lal P, Alavalapati JRR, Marinescu M, et al. (2011) Developing Sustainability Indicators for Woody Biomass Harvesting in the United States. J Sustain For 30: 736-755. doi: 10.1080/10549811.2011.571581
    [30] Evans A, Perschel R, Kittler B, et al. (2010) Revised assessment of biomass harvesting and retention guidelines. For Guild, St Fe, NM, USA: 33.
    [31] Janowiak MK, Webster CR (2010) Promoting Ecological Sustainability in Woody Biomass Harvesting. J For 108: 16-23.
    [32] Apex Resource Management Solutions (2014) ST-Sim state-and-transition simulation model software. Available from: http//www.apexrms.com/stsm.
    [33] Rollins MG (2009) LANDFIRE: a nationally consistent vegetation, wildland fire, and fuel assessment. Int J Wildl Fire 18: 235-249. doi: 10.1071/WF08088
    [34] Comer P, Faber-Langendoen D, Evans R, et al. (2003) Ecological Systems of the United States: A Working Classification of U.S. Terrestrial Systems. Arlington, VA, USA. NatureServe, 82 pages.
    [35] Costanza JK, Terando AJ, McKerrow AJ, et al. (2015) Modeling climate change, urbanization, and fire effects on Pinus palustris ecosystems of the southeastern U.S. J Environ Manage 151: 186-199. doi: 10.1016/j.jenvman.2014.12.032
    [36] LANDFIRE (2014) LANDFIRE 2008 (version 1.1.0) Succession Class (S-Class) Layer. U.S. Department of Interior, Geological Survey. Available from: Http://landfire.cr.usgs.gov/viewer.
    [37] Multi-Resolution Land Characteristics Consortium (MRLC) (2011) National Land Cover Database, USFS Tree Canopy Cartographic, 2014. Available from: http://www.mrlc.gov/nlcd11_data.php.
    [38] Mackie R, Mason J, Curcio G (2007) LANDFIRE biophysical setting model for Southern Piedmont Dry Oak(-Pine) Forest. Available from: http://www.landfire.gov/national_veg_models_op2.php.
    [39] USDA Forest Service (2012) Forest Inventory and Analysis Data. Available from: http://apps.fs.fed.us/fiadb-downloads/datamart.html.
    [40] Young T, Wang Y, Guess F, et al. (2015) Understanding the Characteristics of Non-industrial Private Forest Landowners Who Harvest Trees. Small-scale For 1-13.
    [41] Hardie I, Parks P, Gottleib P, et al. (2000) Responsiveness of Rural and Urban Land Uses to Land Rent Determinants in the U.S. South. Land Econ 76: 659. doi: 10.2307/3146958
    [42] USDA Natural Resources Conservation Service (2000) 1997 National Resources Inventory Data, Revised December 2000.
    [43] Dale VH, Kline KL, Wright LL, et al. (2011) Interactions among bioenergy feedstock choices, landscape dynamics, and land use. Ecol Appl 21: 1039-1054. doi: 10.1890/09-0501.1
    [44] Evans JM, Fletcher RJ, Alavalapati JRR, et al. (2013) Forestry Bioenergy in the Southeast United States: Implications for Wildlife Habitat and Biodiversity. Availbale from: http://www.nwf.org/News-and-Magazines/Media-Center/Reports/Archive/2013/12-05-13-Forestry-Bioenergy-in-the-Southeast.aspx.
    [45] Owens AK, Moseley KR, McCay TS, et al. (2008) Amphibian and reptile community response to coarse woody debris manipulations in upland loblolly pine (Pinus taeda) forests. For Ecol Manage 256: 2078-2083. doi: 10.1016/j.foreco.2008.07.030
    [46] Otto CR V, Kroll AJ, McKenny HC (2013) Amphibian response to downed wood retention in managed forests: A prospectus for future biomass harvest in North America. For Ecol Manage 304: 275-285. doi: 10.1016/j.foreco.2013.04.023
    [47] Davis JC, Castleberry SB, Kilgo JC (2010) Influence of coarse woody debris on herpetofaunal communities in upland pine stands of the southeastern Coastal Plain. For Ecol Manage 259: 1111-1117. doi: 10.1016/j.foreco.2009.12.024
    [48] Wood P, Sheehan J, Keyser P, et al. (2013) Cerulean Warbler: Management Guidelines for Enhancing Breeding Habitat in Appalachian Hardwood Forests. American Bird Conservancy. The Plains, VA, USA. 28 Pages.
    [49] Perry RW, Thill RE (2013) Long-term responses of disturbance-associated birds after different timber harvests. For Ecol Manage 307: 274-283. doi: 10.1016/j.foreco.2013.07.026
    [50] Wilson MD, Watts BD (2000) Breeding bird communities in pine plantations on the coastal plain of North Carolina. Chat 64: 1-14.
    [51] Peet RK, Allard DJ (1993) Longleaf Pine Vegetation of the Southern Atlantic an Eastern Gulf Coast Regions: A Preliminary Classification. In: Hermann SM, ed. Proceedings of the Tall Timbers Fire Ecology Conference, No. 18, The Longleaf Pine Ecosystem: Ecology, Restoration and Management. Tallahassee, FL, USA: Tall Timbers Research Station.; 1993: 45-81.
  • This article has been cited by:

    1. Mahmoud A. Ibrahim, Attila Dénes, A mathematical model for Lassa fever transmission dynamics in a seasonal environment with a view to the 2017–20 epidemic in Nigeria, 2021, 60, 14681218, 103310, 10.1016/j.nonrwa.2021.103310
    2. Robert R. Parmenter, Gregory E. Glass, Hantavirus outbreaks in the American Southwest: Propagation and retraction of rodent and virus diffusion waves from sky-island refugia, 2022, 36, 0217-9792, 10.1142/S021797922140052X
    3. Juan Pablo Gutiérrez Jaraa, María Teresa Quezada, Modeling of hantavirus cardiopulmonary syndrome, 2022, 22, 07176384, e002526, 10.5867/medwave.2022.03.002526
    4. 惠振阳 Hui Zhenyang, 李卓宣 Li Zhuoxuan, 程朋根 Cheng Penggen, 蔡诏晨 Cai Zhaochen, 郭先春 Guo Xianchun, 基于多约束图形分割的点云对象基元获取方法, 2024, 61, 1006-4125, 1037001, 10.3788/LOP231575
    5. Juan Pablo Gutiérrez-Jara, María Teresa Muñoz-Quezada, Fernando Córdova-Lepe, Alex Silva-Guzmán, Mathematical Model of the Spread of Hantavirus Infection, 2023, 12, 2076-0817, 1147, 10.3390/pathogens12091147
    6. Yuhang Li, Yanni Xiao, Effects of nonlinear impulsive controls and seasonality on hantavirus infection, 2025, 00255564, 109378, 10.1016/j.mbs.2025.109378
  • Reader Comments
  • © 2015 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(7405) PDF downloads(1455) Cited by(15)

Figures and Tables

Figures(7)  /  Tables(2)

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog