# Global Dynamics and Applications of an Epidemiological Model for Hepatitis C Virus Transmission in China

Authors: Mingwang Shen; Yanni Xiao; Weike Zhou; Zhen Li
Journal: Discrete Dynamics in Nature and Society
 (2015)

---

## Abstract

An epidemiological model is proposed and studied to understand the
transmission dynamics and prevalence of HCV infection in China.
Theoretical analysis indicates that the basic reproduction numberR0 provides a threshold value determining whether the disease
dies out or not. Two Lyapunov functions are constructed to prove the
global asymptotic stability of the disease-free and the endemic
equilibria, respectively. Based on data reported by the National Health and
Family Planning Commission of China, the basic reproduction number
is estimated as approximately R0=1.9897, which is much less than
that for the model when a treatment strategy is not considered. An
ever-increasing HCV infection is predicted in the near future.
Numerical simulations, performed to investigate the potential effect
of antiviral treatment, show that increasing the treatment cure rate
and enlarging the treatment rate for patients at the chronic
stage remain effective in reducing the number of new infections and the
equilibrium prevalence. The finding suggests that treatment measures
are significantly beneficial for disease control in terms of reducing
new infections and, in particular, more attention should be paid to
treatment for patients at the chronic stage.


---

## Body

## 1. Introduction

Hepatitis C virus (HCV) has been considered as a leading cause of liver cirrhosis and hepatocellular carcinoma and is becoming a major and growing global health problem [1, 2]. HCV is an enveloped single-stranded RNA virus in the Flaviviridae family and mutates so rapidly that no vaccine is currently available [3]. The spread of HCV mainly results from blood-to-blood contact through blood transfusions, intravenous drug use (IDU), and the use of inadequately sterilized or unsterilized medical equipment. According to the World Health Organization (WHO) estimates [4], nearly 3% of the world’s population (more than 170 million) has been infected with HCV. In a recent cross-sectional study [5], for which 8,762 Chinese subjects from six areas of China were randomly selected, the overall average prevalence of anti-HCV was estimated to be 0.58% in China, which is much lower than the WHO estimates. Surveillance data show that the annual numbers of newly reported HCV cases increased sharply from 21,145 in 2003 to 223,094 in 2013 [6], which indicates that HCV infection is becoming a serious threat to public health in China. Therefore, it is urgently necessary to understand the present epidemic situation and to provide suggestions on how to control HCV infections.Mathematical models have been used to analyze the spread and control of HCV infection in [7–17] and provide some insights into the disease’s transmission. Martcheva and Castillo-Chavez [7] considered an epidemiological model of hepatitis C in a varying population, divided into susceptible and infected individuals with acute and chronic hepatitis C, and studied the role of the chronic infectious stage on the long-term dynamics of HCV. This model was extended by Das et al. [8] by incorporating the immune class and was also extended by Yuan and Yang [9] by incorporating the latent period. The impacts of HCV treatment on prevalence among active injecting drug users (IDUs) have been studied in [10–12] where the treated populations are assumed not to infect the susceptible populations. Imran et al. [14] and Khan et al. [17] formulated epidemic models of HCV containing an isolation class and analyzed the effects of the isolation class on the transmission dynamics of the disease. Few models have been proposed to understand the transmission dynamics of HCV in mainland China. In [18], Zhang and Zhou extended the model in [9] by considering recovery from the acute infected stage and simulated HCV transmission in the near future based on the available HCV epidemic data from 2003 to 2010 in China. Note that their model does not involve the treatment class although a treatment strategy has been available in mainland China for several years [19]. Antiviral treatment, especially at the early stages of HCV recurrence, is an extremely effective strategy. How effective antiviral treatment remains is therefore an issue of great importance for HCV control, and quantifying this impact through a mathematical modeling framework falls within the scope of this study.The purpose of this study is to formulate a mathematical model which involves almost all possible stages during HCV infection with the aim of accessing the potential impact of antiviral treatment on HCV prevalence of China and then providing reliable quantitative information on controlling the HCV epidemic in China. On the basis of the model in [18] we incorporate the treated and recovered classes with partial immunity. Note that we consider that the treated individuals can also infect the susceptible individuals [16], which adds to the complexity of our model and makes it considerably more insightful from an epidemiological perspective than previous models. We analyze the global dynamical behavior of the proposed system by developing suitable Lyapunov functions [20] and suggest some measures to control HCV infection in China.The paper is organized as follows. In Section2, we introduce the model and derive the basic reproduction number. The global stability of the disease-free equilibrium and the endemic equilibrium are studied in Sections 3 and 4, respectively. In Section 5, we apply the model to simulate the HCV data in China and investigate various control strategies with numerical simulations. Finally, we conclude the paper by a discussion in Section 6.


## 2. Model Formulation

We propose a mathematical model to understand the transmission dynamics and prevalence of HCV in mainland China using a system of ordinary differential equations. The population are divided into six classes:S (susceptible), E (exposed), A (infected with acute hepatitis C), C (infected with chronic hepatitis C), T (treated population), and R (recovered population with partial immunity). Let N denote the total population; that is, N=S+E+A+C+T+R. New susceptible individuals enter into the S class at a fixed rate Λ. Let m be the natural death rate of the population. Susceptible individuals are infected by patients in the A,C, and T classes at rates of β1, β2, and β3, respectively. It is generally thought that the acute phase is more infectious than the the chronic one and that the treated individuals have the lowest infectiousness [15, 18], so it is reasonable to assume that β1>β2>β3. Once infected, the individuals move into the exposed class E and then progress to the acute stage at a rate of ϵ. An individual at the acute stage who can spontaneously clear the virus recovers at rate ξ; otherwise he/she will progress to the chronic stage at rate σ. Acutely and chronically infected individuals are treated at rates δ and μ, respectively. After the treatment, some patients succeed in clearing HCV and move to the class R at rate η, while the others who have not responded to the treatment move back to the chronic state C at rate θ. The disease-induced death rate in the C class is denoted as α. Individuals in the class R lose their immunity and eventually return to the susceptible class S at rate γ. A flow diagram for the model is shown in Figure 1 and the variables and parameters are described in Table 1. The complete dynamical model is as follows:
(1)dSdt=Λ-(β1A+β2C+β3T)S-mS+γR,dEdt=β1A+β2C+β3TS-m+ϵE,dAdt=ϵE-m+σ+δ+ξA,dCdt=σA+θT-m+α+μC,dTdt=δA+μC-m+θ+ηT,dRdt=ξA+ηT-(m+γ)R.Table 1
Parameters and initial data chosen for the simulation.
Variable and parameter
Description
Initial or default values
Source
Variables
S(t)
Susceptible population
1.287
×
1
0
9
LS
E(t)
Exposed population
4187
LS
A(t)
Acutely infected population
151360
LS
C(t)
Chronically infected population
1.04
×
1
0
6
LS
T(t)
Treated population
50000
LS
R(t)
Recovered population with partial immunity
2
×
1
0
5
LS
Parameters
Λ
Recruitment rate
1.573
×
1
0
7 year−1
See text
m
Natural death rate
0.007 year−1
[18, 24]
β1
Transmission rate of acutely infected population
3.0769
×
1
0
-
11
LS
β2
Transmission rate of chronically infected population
2.5846
×
1
0
-
11
LS
β3
Transmission rate of treated population
2.0846
×
1
0
-
11
LS
ϵ
Rate of progression to acute stage from the exposed
6 year−1
[18]
ξ
Recovery rate for the acute state
0.5  year−1
[15]
σ
Rate of moving from acute stage to chronic stage
4  year−1
[18]
δ
Treatment rate of acutely infected population
0.1545 year−1
see text
μ
Treatment rate of chronically infected population
0.04 year−1
[15]
η
Treatment cure rate
0.67 year−1
[15]
θ
Treatment failure rate
0.82 year−1
[15]
α
HCV induced death rate at the chronic stage
0.001 year−1
[18]
γ
Rate of waning immunity
0.025 year−1
[15]
LS, least square.Figure 1
A schematic flow diagram illustrating the transmission dynamics of the HCV infection with treatment.SincedN/dt=Λ-mN-αC≤Λ-mN, we can study the dynamical behavior of system (1) in the positively invariant set:
(2)D=S,E,A,C,T,R∈R+6:ΛmS+E+A+C+T+R≤Λm.
By using the next generation matrix approach given in [21], we obtain the basic reproduction number:
(3)R0=Λϵm(m+ϵ)(m+σ+δ+ξ)β1+β2δθ+σm+θ+ηm+α+μm+η+m+αθ-1+β3μσ+δm+α+μ·m+α+μm+η+m+αθ-1,
which can be interpreted as the average number of new infections generated by a single infectious individual in the acute (A), chronic (C), and treated (T) classes.To investigate the effect of therapy, we study variation in the basic reproduction numberR0 with treatment rates δ, μ, treatment failure rate θ, and cure rate η. Calculating the derivative of R0 with respect to μ gives
(4)∂R0∂μ=-Λϵm(m+ϵ)(m+σ+δ+ξ)·δθ+σm+θ+ηβ2-β3m+β2η-β3αm+α+μm+η+m+αθ2.
Note that β1>β2>β3; then ∂R0/∂μ<0 if the condition β2η-β3α>0 holds. This condition is likely to be true because the disease-induced death rate α is always relatively low [7] (i.e., η>α). Thus, R0 decreases with increasing treatment uptake for the patients at the chronic stage. Similarly, by finding the derivative of R0 with respect to θ,
(5)∂R0∂θ=Λϵmm+ϵm+σ+δ+ξ·μσ+δm+α+μβ2-β3m+β2η-β3αm+α+μm+η+m+αθ2,
we obtain ∂R0/∂θ>0 provided β2η-β3α>0, which implies that decreasing the treatment failure rate θ results in R0 declining and thus is beneficial to disease control.Regardless,R0 is inversely related to η:
(6)∂R0∂η=-Λϵmm+ϵm+σ+δ+ξ·μσ+δm+α+μβ2θ+β3m+α+μm+α+μm+η+m+αθ2<0.
Therefore, a higher cure rate leads to lower R0 and a lower intensity of the epidemic. Differentiating R0 with respect to δ yields
(7)∂R0∂δ=-Λϵmm+ϵm+α+μm+η+m+αθm+σ+δ+ξ2·m+σ+δ+ξ2-1×(β1-β3)m(m+α+μ)+(β1-β2)mθ+(β2-β3)mσ+(m+α+μ)(β1η-β3ξ)+β2(ση-θξ)+α(β1θ-β3σ).
If (β3/β1)σ<θ<(η/ξ)σ, then ∂R0/∂δ<0. Thus, we get a sufficient condition to make sure that R0 is inversely related to δ, which means that increasing treatment uptake for the patients at the acute stage is beneficial to control the disease for certain conditions of the treatment failure rate.Model (1) always has a disease-free equilibrium P¯=(Λ/m,0,0,0,0,0). There is an endemic equilibrium P*=(S*,E*,A*,C*,T*,R*) if R0>1, where
(8)S*=ΛmR0,E*=(m+σ+δ+ξ)ϵA*,A*=Λ(R0-1)A~R0,C*=δθ+σ(m+θ+η)(m+α+μ)(m+η)+(m+α)θA*,T*=μσ+δ(m+α+μ)(m+α+μ)(m+η)+(m+α)θA*,R*=1m+γξ+ημσ+δ(m+α+μ)(m+α+μ)(m+η)+(m+α)θA*,
with
(9)A~=m+ϵϵm+σ+δ+ξ-γm+γξ+ημσ+δ(m+α+μ)(m+α+μ)(m+η)+(m+α)θ>0.


## 3. Global Stability of the Disease-Free Equilibrium

Theorem 1.
(i) IfR0≤1, then the unique disease-free equilibrium P¯ of system (1) is globally asymptotically stable in D.
(ii) IfR0>1, then P¯ is unstable and system (1) is uniformly persistent; that is, there exists a constant ϵ0>0 such that, for all initial values (S(0),E(0),A(0),C(0),T(0),R(0))∈R+×Int(R+5), the solutions of system (1) satisfy liminf⁡t→∞S(t)>ϵ0, liminf⁡t→∞E(t)>ϵ0, liminf⁡t→∞A(t)>ϵ0, liminf⁡t→∞C(t)>ϵ0, liminf⁡t→∞T(t)>ϵ0, and liminf⁡t→∞R(t)>ϵ0.Proof.
(i) Construct a continuously differentiable and positive definite Lyapunov function(10)L=mΛE+mm+ϵΛϵA+Hβ3μ+β2m+θ+ηC+H(β2θ+β3(m+α+μ))T,
where H is
(11)H=1R0((m+α+μ)(m+η)+(m+α)θ).
Calculating the derivative of L along (1), we obtain
(12)dLdt=mΛdEdt+mm+ϵΛϵdAdt+Hβ3μ+β2m+θ+ηdCdt+Hβ2θ+β3m+α+μdTdt=mΛβ1A+β2C+β3TS-m+ϵE+mm+ϵΛϵϵE-m+σ+δ+ξA+Hβ3μ+β2m+θ+η·σA+θT-m+α+μC+Hβ2θ+β3m+α+μ·δA+μC-m+θ+ηT=mΛβ1A+β2C+β3TS-1R0β1A+β2C+β3T=1R0β1A+β2C+β3TR0mΛS-1.
In the domainD, it is easy to see that S≤Λ/m, which gives (m/Λ)S≤1. Thus, if R0≤1, then we have dL/dt≤0 for all S,A,C,T>0. Every solution of system (1) converges to M, where M is the largest invariant set in {(S,E,A,C,T,R)∈D:dL/dt=0}. When R0<1, the equality dL/dt=0 holds if and only if A=C=T=0, which implies that E=R=0 and S=Λ/m from (1). When R0=1, the equality dL/dt=0 holds if and only if A=C=T=0 or S=Λ/m. Both of these cases indicate that M=P¯. By the LaSalle largest invariant set theorem, the disease-free equilibrium P¯ is globally asymptotically stable if R0≤1.
(ii) The Jacobian matrix of system (1) at P¯ is (13)JP¯=-m0-β1Λm-β2Λm-β3Λmγ0-m+ϵβ1Λmβ2Λmβ3Λm00ϵ-m+σ+δ+ξ00000σ-m+α+μθ000δμ-m+θ+η000ξ0η-m+γ.For convenience of denotation, let h1=m+ϵ, h2=m+σ+δ+ξ, h3=m+α+μ, and h4=m+θ+η. Then we obtain the characteristic equation at P¯ as follows:
(14)(λ+m)(λ+m+γ)·(λ+h1)(λ+h2)(λ+h3)(λ+h4)Λm-μθ(λ+h1)(λ+h2)-β1Λmϵ(λ+h3)(λ+h4)+β1Λmϵμθ-β2Λmϵσ(λ+h4)-β2Λmϵδθ-β3Λmϵμσ-β3Λmϵδ(λ+h3)=λ+mλ+m+γλ4+k1λ3+k2λ2+k3λ+k4=0,
where
(15)k1=h1+h2+h3+h4,k2=h1h2+h1h3+h1h4+h2h3+h2h4+h3h4-μθ-β1Λmϵ,k3=h1h2h3+h1h2h4+h1h3h4+h2h3h4-μθ(h1+h2)-β1Λmϵ(h3+h4)-β2Λmϵσ-β3Λmϵδ,k4=h1h2h3h4-μθh1h2-β1Λmϵh3h4+β1Λmϵμθ-β2Λmϵσh4-β2Λmϵδθ-β3Λmϵμσ-β3Λmϵδh3.
Denote the eigenvalues of the characteristic equation (14) as λi, i=1,2,…,6. Obviously, λ1=-m<0, λ2=-(m+γ)<0, and λ3, λ4, λ5, λ6 satisfy the equation (λ4+k1λ3+k2λ2+k3λ+k4)=0. When R0>1, we get
(16)k4=h1h2(h3h4-μθ)-β1Λmϵ(h3h4-μθ)-β2Λmϵ(σh4+δθ)-β3Λmϵ(δh3+μσ)=h1h2(h3h4-μθ)·1-Λϵmh1h2β1+β2(δθ+σh4)+β3(μσ+δh3)h3h4-μθ=(m+ϵ)(m+σ+δ+ξ)·[(m+α+μ)(m+η)+(m+α)θ]×1-Λϵmm+ϵm+σ+δ+ξ·β1+β2δθ+σm+θ+ηm+α+μm+η+m+αθ-1+β3μσ+δm+α+μΛϵmm+ϵm+σ+δ+ξ·m+α+μm+η+m+αθ-1=(m+ϵ)(m+σ+δ+ξ)·[(m+α+μ)(m+η)+(m+α)θ](1-R0)<0,
which means that λ3λ4λ5λ6=k4<0. This indicates that at least one of λ3,λ4,λ5,λ6 has positive real part. Thus, P¯ is unstable if R0>1.
Define(17)X=S,E,A,C,T,R∣S≥0,E≥0,A≥0,C≥0,T≥0,R≥0,X0=S,E,A,C,T,R∣S≥0,E>0,A>0,C>0,T>0,R>0,∂X0=X∖X0,M∂=S0,E0,A0,C0,T0,R0∣St,Et,At,Ct,Tt,RTsatisfies(1),St,Et,At,Ct,Tt,RT∈∂X0,R0∀t≥0.
It then suffices to show that system (1) is uniformly persistent with respect to (X0,∂X0). First, it is easy to see that both X and X0 are positively invariant. Obviously, ∂X0 is relatively closed in X. Furthermore, it follows from (2) that system (1) is point dissipative. Next we prove that M∂={(S,0,0,0,0,0)∣S≥0}. Assume, by contradiction, that φ0=(S(0),E(0),A(0),C(0),T(0),R(0))∈M∂ and there exists a t0≥0 such that at least one of E(t),A(t),C(t),T(t),R(t) is greater than zero at t=t0; for example, A(t0)>0 and E(t0)=C(t0)=T(t0)=R(t0)=0. Then we get C′(t0)=σA(t0)>0, T′(t0)=δA(t0)>0, R′(t0)=ξA(t0)>0. This indicates that there exists a t1>0 such that C(t)>0, T(t)>0, R(t)>0, and E(t)>0, A(t)>0, for t0<t<t0+t1, which means that (S(t),E(t),A(t),C(t),T(t),R(t))∉∂X0 for t0<t<t0+t1, contradicting the assumption that φ0∈M∂. Other cases can be proved in the same way. Thus we have M∂={(S,0,0,0,0,0)∣S≥0}. Let φ0 be an initial value. There is only one equilibrium P¯ in M∂, so ∪φ0∈M∂ω(φ0)=P¯. Therefore, P¯ is a compact and isolated invariant set for φ0 in ∂X0.
Next we claim that there exists a positive constantδ0 such that any solution u(t,φ0),φ0∈X0, satisfies
(18)lim⁡t→∞sup⁡u(t,φ0)-P¯≥δ0,i.e.,Ws(P¯)∩X0=∅.
Suppose the claim is not true. Then lim⁡t→∞sup⁡u(t,φ0)-P¯<δ0, for any δ0>0; namely, there exists a positive constant η0 such that Λ/m-δ0<S(t)<Λ/m+δ0, E(t)<δ0, A(t)<δ0, C(t)<δ0, T(t)<δ0, and R(t)<δ0, for any t>η0. While t>η0, we have
(19)dEdt>(β1A+β2C+β3T)Λm-δ0-(m+ϵ)E.
Consider the auxiliary system
(20)dx1dt=(β1x2+β2x3+β3x4)Λm-δ0-(m+ϵ)x1,dx2dt=ϵx1-(m+σ+δ+ξ)x2,dx3dt=σx2+θx4-(m+α+μ)x3,dx4dt=δx2+μx3-(m+θ+η)x4,dx5dt=ξx2+ηx4-(m+γ)x5.
Define (21)M1(δ0)=-(m+ϵ)β1Λm-δ0β2Λm-δ0β3Λm-δ00ϵ-(m+σ+δ+ξ)0000σ-(m+α+μ)θ00δμ-(m+θ+η)00ξ0η-(m+γ).Let s(M1(δ0)) be the maximum real part of the eigenvalues of M1(δ0). Since M1(δ0) is irreducible and has nonnegative off-diagonal elements, s(M1(δ0)) is a simple eigenvalue of M1(δ0) with a positive eigenvector. It follows from Lemma 2.1 in [22] that s(M1(0))>0 when R0>1. Since s(M1(δ0)) is continuous for small δ0, there exists a positive constant δ0 small enough such that s(M1(δ0))>0. Therefore, there is a positive eigenvalue of M1(δ0) with a positive eigenvector. It is easy to see xi→∞ as t→∞, i=1,2,…,5. Then according to the comparison principle we have
(22)lim⁡t→∞Et=∞,lim⁡t→∞At=∞,lim⁡t→∞Ct=∞,lim⁡t→∞T(t)=∞,lim⁡t→∞R(t)=∞,
which contradicts our assumption. This completes the proof of the claim, which implies that P¯ is an isolated invariant set in X and Ws(P¯)∩X0=∅. Using the uniform persistence theory (Theorem 4.2) in [23], we obtain that system (1) is uniformly persistent if R0>1. This completes the proof.


## 4. Global Stability of the Endemic Equilibrium

The global stability of the endemic equilibrium is studied under the simplified assumption that the immune loss rate is zero (i.e.,γ=0). Thus, we can omit the decoupled equation for dR/dt and the following theorem holds.Theorem 2.
WhenR0>1, if γ=0, then the unique endemic equilibrium P* is globally asymptotically stable in D.Proof.
Construct a continuously differentiable and nonnegative Lyapunov function:(23)V=S-S*-S*ln⁡SS*+E-E*-E*ln⁡EE*+a1A-A*-A*ln⁡AA*+a2C-C*-C*ln⁡CC*+a3T-T*-T*ln⁡TT*,
where
(24)a1=m+ϵϵ,a2=(β3μ+β2(m+θ+η))S*(m+α+μ)(m+η)+(m+α)θ,a3=(β2θ+β3(m+α+μ))S*(m+α+μ)(m+η)+(m+α)θ.
It can be verified that the global minimum of V occurs at the endemic equilibrium P*, and the function V takes the value V=0 at the endemic equilibrium P*.
DifferentiatingV along the solutions of system (1) and using the equilibrium relations, we obtain
(25)dVdt=1-S*S(Λ-β1SA-β2SC-β3ST-mS)+1-E*E(β1SA+β2SC+β3ST-(m+ϵ)E)+a11-A*A(ϵE-(m+σ+δ+ξ)A)+a21-C*C(σA+θT-(m+α+μ)C)+a31-T*T(δA+μC-(m+θ+η)T)=M-mS-β1E*SAE-β2E*SCE-β3E*STE-ΛS*S-a1ϵA*EA-a2σC*AC-a2θC*TC-a3δT*AT-a3μT*CT=b12-SS*-S*S+b23-S*S-E*S*A*SAE-A*E*EA+b34-S*S-E*S*C*SCE-A*E*EA-C*A*AC+b45-S*S-E*S*C*SCE-A*E*EA-C*T*TC-T*A*AT+b55-S*S-E*S*T*STE-A*E*EA-C*A*AC-T*C*CT+b64-S*S-E*S*T*STE-A*E*EA-T*A*AT+b72-C*T*TC-T*C*CT,
where M=Λ+mS*+(m+ϵ)E*+a1(m+σ+δ+ξ)A*+a2(m+α+μ)C*+a3(m+θ+η)T* and bi (i=1,…,7) satisfy the following equations:
(26)2b1+3b2+4b3+5b4+5b5+4b6+2b7=M,b1=mS*,b1+b2+b3+b4+b5+b6=Λ,b2+b3+b4+b5+b6=a1ϵE*,b3+b5=a2σA*,b4+b7=a2θT*,b4+b6=a3δA*,b5+b7=a3μC*,b2=β1S*A*,b3+b4=β2S*C*,b5+b6=β3S*T*.
Obviously,b1>0 and b2>0. Next we will show that bi≥0(i=3,…,7). From the equations above, we have
(27)b3=a2σA*-a3μC*+b7=β2S*C*-a2θT*+b7,b4=a2θT*-b7,b5=a3μC*-b7,b6=a3δA*-a2θT*+b7=β3S*T*-a3μC*+b7.
At the equilibrium, the equalities σA*+θT*=(m+α+μ)C*, δA*+μC*=(m+θ+η)T* hold. Then by substituting the values of a2 and a3, it is easy to obtain that
(28)a2σA*-a3μC*=β2S*C*-a2θT*,a3δA*-a2θT*=β3S*T*-a3μC*.
In order to assurebi≥0 (i=3,4,5,6), from (27) we know that b7 must satisfy the inequalities:
(29)max⁡⁡0,a3μC*-a2σA*,a3μC*-β3S*T*≤b7≤min⁡a2θT*,a3μC*.
Using (28) yields
(30)a3μC*>0,a3μC*>a3μC*-a2σA*,a3μC*>a3μC*-β3S*T*,a2θT*>0,a2θT*>a2θT*-β2S*C*=a3μC*-a2σA*,a2θT*>a2θT*-a3δA*=a3μC*-β3S*T*,
and then inequalities (29) must have a nonnegative solution, and so has (27). This indicates that bi≥0(i=3,…,7). Because the geometric mean is always less than or equal to the arithmetic mean, we have dV/dt≤0 and the equality holds if and only if (S,E,A,C,T) take the equilibrium values (S*,E*,A*,C*,T*). Therefore, by LaSalle’s Invariance Principle, it follows that the endemic equilibrium P* is globally asymptotically stable in the feasible region D.


## 5. Numerical Results

The annual reported HCV case numbers have been released by the National Health and Family Planning Commission of China [6], shown in Figure 2(a). The birth rate is fixed as 0.0121 in [24] and the total population of China is about 1.3×109, so we obtain that the recruitment rate is Λ=0.0121×1.3×109=1.573×107 per year. Okosun [16] chose the rate of progression for treatment from acute infected individuals as in the interval (0.12, 0.189); thus, we choose the mean treatment rate for the patients at the acute stage as δ=0.1545 per year. By fitting (1) to the annual reported HCV data (Figure 2(a)) we obtain estimates for the transmission rates and the initial population size, which are listed in Table 1 and the goodness of fit is shown in Figure 2(b). Figure 2(b) shows that the estimated reported case numbers will reach 809,970 in 2021 in mainland China if the current surveillance, testing, and interventions are unchanged.(a) Annual newly reported HCV cases for mainland China; (b) goodness of fit and prediction of HCV trends until 2021. Stars represent the reported number of people with HCV by year. The solid line shows the fit based on the current circumstances (μ=0.04). The dashed and dash dot lines show the prediction at higher treatment rates μ=0.12 and μ=0.2 for patients at the chronic stage from 2015, respectively. All the other parameters are as shown in Table 1.
(a)
(b)Using the estimated parameter values we calculated the basic reproduction numberR0 as 1.9897, which is similar to other estimates in the literature [25, 26] but is smaller than the value of 4.0636 estimated in [18]. This is because the treated individuals included in our model have relatively low infectiousness, which results in a smaller basic reproduction number. Moreover, it is interesting to note that the prevalence of HCV infection in China (Figure 3) is estimated as 0.51% in 2013 and 0.58% in 2014, which is in good agreement with the cross-sectional study [5]. Note that this estimated HCV prevalence is much lower than the WHO estimation. This is partial because strict HCV screening in blood began in the 1990s in China, which reduced the transmission rate of HCV by blood transfusion and other sources of iatrogenic infection.Plots of the prevalence against time at varying parameter values: (a) the treatment rateμ for patients at the acutely infected stage; (b) the treatment rate δ for patients at the chronic stage; (c) the treatment cure rate η; (d) the treatment failure rate θ. All the other parameters are as shown in Table 1.
(a)
(b)
(c)
(d)To examine the impact of treatment on HCV transmission dynamics and prevalence and identify the most effective measures to control the transmission of HCV in mainland China we investigated variation in the basic reproduction number and prevalence (i.e.,(E(t)+A(t)+C(t)+T(t))/N(t)) with parameters associated with treatment. If the treatment rate μ for the infected people at the chronic stage increases threefold, then the predicted reported case number will decrease by 13.21% to 702,960 in 2021, while if μ increases fivefold, then the predicted number will decrease by 24.74% to 609,580 in 2021, as shown in Figure 2(b). It follows from Figure 3 that the HCV prevalence in China will continue to rise and reach 1.63% in year 2021 under the current circumstances. Figure 3(a) shows that increasing the treatment rate μ for the patients at the chronic stage by 200% and 400% from the baseline value can decrease the prevalence in 2021 by 19.75% and 47.15%, respectively. Figure 3(b) shows that increasing the treatment rate δ for the patients at the acute stage by 200% and 400% from the baseline value can decrease the prevalence in 2021 by 32.16% and 61.28%, respectively. Figure 3(c) shows that increasing the cure rate η by 200% from η=0.67 can decrease the prevalence in 2021 by 11.01%. It follows from Figure 3(d) that reducing the treatment failure rate θ by 50% and 25% from the baseline value can decrease the prevalence in 2021 by 8.37% and 14.44%, respectively. These results indicate that, to decrease the prevalence effectively, it will be better to enlarge the treatment rate for the patients at the acute and chronic stage and cure rate and decrease the treatment failure rate in the short term.To access the effectiveness of treatment interventions in the long term, we examine the effects of the corresponding treatment parametersμ,δ,η, and θ on the basic reproduction number and prevalence at the endemic equilibrium (i.e., (E*+A*+C*+T*)/N*). Let k (0≤k≤1) be the rate of change of each parameter; then we could increase the treatment rate for the patients at the chronic (or acute) stage (represented by (1+k)μ or ((1+k)δ)), the cure rate ((1+k)η), and decrease the treatment failure rate ((1-k)θ) by increasing parameter k. In particular, increasing μ by 30% or decreasing θ by 30% from baseline values (while keeping other parameters fixed) can reduce R0 by 16.56% or 11.83% and can reduce the equilibrium prevalence by 27.98% or 19.54%, respectively. However, increasing δ or η by 30% from baseline values can only reduce R0 by 0.36% or 9.42% and reduce the equilibrium prevalence by 0.56% or 15.47%, respectively. It follows from Figure 4 that increasing the treatment rate μ for the patients at the chronic stage and the cure rate η and decreasing the treatment failure rate θ are more effective than increasing the treatment rate δ for the patients at the acute stage in terms of reducing both the basic reproduction number and the equilibrium prevalence in the long run.(a) Plots of the basic reproduction numberR0 by varying k (the rates of change of μ; δ; η; and θ); (b) plots of the prevalence at the steady state by varying k (the rates of change of μ; δ; η; and θ). Each column of markers denote that k increases 10% per time. All the other parameters are as shown in Table 1.
(a)
(b)To examine the sensitivity of our results to parameter variation, we used Latin hypercube sampling (LHS) and partial rank correlation coefficients (PRCCs) [27, 28] to examine the dependence of R0 and the equilibrium prevalence on each parameter. Because of limited information on the distributions of each parameter, we chose a uniform distribution as in [28] for all input parameters with ranges listed in Table 2. To know whether the significance of any parameter varies over an entire time interval during model dynamics, PRCC values were calculated for numerous times and plotted versus time. This enables us to assess whether the significance of one parameter changes over an entire time interval during the progression of the model dynamics. Figure 5(a) shows PRCC values plotted from 2003 to 2021 and it indicates that there are five PRCC values that are significantly different from zero and that the PRCC values of the eight examined parameters vary little with time and stabilize at fixed values in about 2006. It follows from Figure 5(b) that the first four parameters with the most significant impact on R0 are the transmission rate for the chronically infected individuals β2 and the transmission rate for the treated population β3, the cure rate η, and the treatment rate μ for the population with chronic infection. This implies that a greater efficacy of treatment, a larger treatment uptake for the population with chronic infection, and lower transmission rates at the chronic stage definitely result in lower new HCV infections. Moreover, these four parameters η,μ,β2,β3 also have significant impact on equilibrium prevalence, as shown in (Figure 5(c)). It is worth mentioning that the rate of waning immunity γ, although it has no effect on R0, greatly affects the equilibrium prevalence because quick waning of immunity significantly increases the number of the susceptible individuals and hence the prevalence. It should be noted that the transmission rate β1 and the treatment rate δ for the population with acute infection have little effect on both the basic reproduction number R0 and the prevalence at the endemic equilibrium because of the relative short duration of the acute stage.Table 2
PRCC values forR0 and prevalence at endemic equilibrium.
Input parameters
Distributions
R
0
Prevalence
PRCC
P value
PRCC
P value
β
1
U
(
7.6923
×
1
0
-
12
,
7.6923
×
1
0
-
10
)
0.0395
0.0781
0.0337
0.1328
β
2
U
(
7.6923
×
1
0
-
12
,
7.6923
×
1
0
-
10
)
0.7609
0
0.4823
0
β
3
U
(
7.6923
×
1
0
-
12
,
7.6923
×
1
0
-
10
)
0.4123
0
0.2284
0
η
U
(
0.01,1
)
−0.8065
0
−0.7549
0
δ
U
(
0.01,1
)
−0.0311
0.1645
−0.0380
0.0902
μ
U
(
0.01,1
)
−0.7644
0
−0.7151
0
θ
U
(
0.01,1
)
0.0386
0.0845
0.2219
0
γ
U
(
0.01,1
)
—
—
0.6550
0(a) PRCCs of the eight parameters for prevalence from 2003 to 2021. (b) PRCCs forR0. (c) PRCCs for the equilibrium prevalence. Sample size is set to 2000. * denotes PRCCs are significant (P value <0.01). Parameter values and ranges are as shown in Tables 1 and 2.
(a)
(b)
(c)


## 6. Discussion

In order to understand the prevalence of HCV infection in China based on the reported data [6] and examine the role that treatment plays in the transmission dynamics, we proposed a mathematical model which includes realistic features of HCV transmission such as treatment and partial immunity. Theoretically, the global dynamics of our model are determined by the basic reproduction number R0. The disease-free equilibrium is globally asymptotically stable if R0≤1, which means that hepatitis C can be entirely eliminated from the population. When R0>1, hepatitis C will persist in the population and the endemic equilibrium is globally asymptotically stable for a special case. It is worth mentioning that although the constructed Lyapunov function V(t) is deterministic when proving the global stability of the positive steady state, the coefficients bk(k=3,4,…,7) of function dV/dt which appeared in (25) were chosen to be a positive solution of (27) and are nonunique, which is more general than the proof of Theorem 3.3 in [15].When the model was applied to HCV transmission in China, we estimated the basic reproduction numberR0 as 1.9897, which is similar to other estimates in the literature [25, 26] but is smaller than the value 4.0636 estimated in [18]. This is because the treated population has relatively low infectiousness which is included in our model, resulting in a smaller basic reproduction number. Goodness of fit and prediction of HCV trends (Figure 2(b)) show a more accurate result than the overestimated simulation in [18] and the prediction indicates that newly reported cases will continue to rise rapidly in the near future. Moreover, we estimate that the prevalence of HCV infection in China (Figure 3) was 0.51% in 2013 and 0.58% in 2014, in good agreement with the cross-sectional study [5], and will reach 1.63% in 2021. These findings lead us to believe that the exact HCV prevalence is much lower than what the WHO estimation indicates. Also, HCV screening in blood, a practice that may significantly reduce the transmission rate of HCV by blood transfusion and other sources of iatrogenic infection, began in the 1990s in China, and thus it is reasonable to believe that strict blood screening and other procedural measures are preventing the spread of HCV and are leading to a lower prevalence.It follows from sensitivity analysis (Figure5) that the transmission rate β3 of the treated population contributes greatly to the transmission of HCV throughout the period of the disease spread. So the prevalence of HCV may be underestimated in [10–12] because in their models the treated population are assumed not to infect the susceptible populations. Figures 4 and 5 show that the basic reproduction number R0 and equilibrium prevalence are not sensitive to the transmission rate β1 and the treatment rate δ for the population with acute infection, indicating that the acute stage does not substantially affect the transmission of HCV in the long run due to its relatively short duration, but it may affect the prevalence at the beginning of the epidemic (Figure 3(b)).It should be acknowledged that one limitation of our results is that the reported national data may not be completely composed of exposed people who enter into the acute stage. The data may contain some cases diagnosed at the acute or chronic stage; though the number of these cases is low, it may still result in a slight overestimate of the prevalence of HCV in China. However, the slightly overestimated results were not caused by our model, but rather by the deficiency of data which did not distinguish which stage the reported cases came from. More realistic models about HCV infection on complex networks [29–31] will be studied in the future work.



---

*Source: 543029-2015-01-28.xml*