# Transmission Dynamics of Hepatitis C with Control Strategies

Authors: Adnan Khan; Sultan Sial; Mudassar Imran
Journal: Journal of Computational Medicine
 (2014)

---

## Abstract

We present a rigorous mathematical analysis of a deterministic model, for the transmission dynamics of hepatitis C, using a standard incidence function. The infected population is divided into three distinct compartments featuring two distinct infection stages (acute and chronic) along with an isolation compartment. It is shown that for basic reproduction numberR0≤1, the disease-free equilibrium is locally and globally asymptotically stable. The model also has an endemic equilibrium for R0>1. Uncertainty and sensitivity analyses are carried out to identify and study the impact of critical parameters on R0. In addition, we have presented the numerical simulations to investigate the influence of different important parameters on R0. Since we have a locally stable endemic equilibrium, optimal control is applied to the deterministic model to reduce the total infected population. Two different optimal control strategies (vaccination and isolation) are designed to control the disease and reduce the infected population. Pontryagin’s Maximum Principle is used to characterize the optimal controls in terms of an optimality system which is solved numerically. Numerical results for the optimal controls are compared against the constant controls and their effectiveness is discussed.


---

## Body

## 1. Introduction

Hepatitis C (HCV) is an important public health problem, as it is the common cause of liver diseases throughout the world [1]. The disease was first recognized in 1975 and its causative agent was identified in 1989. Hepatitis C is characterized by an acute (often asymptotic) stage, which, in most cases, is followed by a chronic stage that can result in cirrhosis and liver cancer. The hepatitis C virus (causative agent) is an enveloped RNA virus, which is further characterized to be a positive-sense single stranded virus belonging to the family Flaviviridae and is considerably small in size. Replication of the RNA-based virus involves the use of the enzyme RNA-dependent RNA polymerase (RdRP), which has a high error rate while going through this process. World Health Organization’s report suggests that around 3% of the world population has been infected with HCV. The population infected with chronic HCV, who are at risk of developing liver cancer or cirrhosis, is estimated to be around 170 million. Furthermore, nearly 350,000 people die annually throughout the globe as a result of HCV-related liver diseases [2].Hepatitis C can be characterized by two distinct stages: an acute stage and a chronic stage. Initially, infection by HCV causes an acute HCV which is usually asymptotic. Only about15% of the cases show mild symptoms like decreased appetite, nausea, fatigue, joint or muscle pains, and weight loss. In 20% of the cases, the infection may resolve spontaneously. And the remaining 80% of the people exposed to HCV progress to the chronic stage of the infection by developing a chronic infection, which can last for decades. During the starting years of infection, most people experience minimal or no symptoms at all. However, HCV becomes the main cause of liver cancer and cirrhosis after several years of living with it. About 1%–5% of chronic HCV patients die from liver cancer or cirrhosis and nearly 5%–20% develop cirrhosis over 30 years. Patients with cirrhosis are 20 times more likely to develop hepatocellular carcinoma, at the rate of 13% each year. Moreover, 27% of cirrhosis cases and 25% of hepatocellular carcinoma cases worldwide are estimated to be caused by HCV [3–5].Depending on the genotype of the HCV, the standard treatment of infected patients includes a combination of pegylated interferon (peg IFN-α) and the robust antiviral drug Ribavirin, for a period of 24 or 48 weeks [6]. The reaction to treatment also differs by genotype and lies between 70% and 80% for genotypes 2 and 3, while it is almost non existent for genotype 6. Recently, there have been promising treatment advances of genotype 1 using directly acting antiviral agents (DAAs). However, to prevent the infection, there is still a need to create effective vaccine strategy. Vaccines, which are used in preclinical and clinical trials, involve DNA-based proteins, recombinant proteins, synthetic peptide vaccines, and so forth. The future design of vaccines, along with the use and success of previously mentioned vaccines, has been discussed here [7]. However, no long-term immunity is granted in recovering from hepatitis C infection. So, this absence of acquired immunity must be shown in any model for hepatitis C. This is done by allowing recovered patients to become susceptible again.Several studies have been carried out [8–13] and they are pertinent to our work. Most of these papers classify individuals in the population into different states and then formulate a system of ordinary differential equations (ODE) to analyze the time-evolution of each of these population states. Reade et al. [12] present an ODE model of infections with acute and chronic stages. Similarly, by Luo and Xiang, [10] a four state system was analyzed with exposed, acute and chronic states. Suna et al. [13] present a study on a SEIRS model where it was assumed that recovered individuals lose their infection-acquired immunity. Martcheva and Castillo-Chavez in [11] have formulated a model for hepatitis C lacking an exposed class and have discussed the stability of the equilibrium states.We will formulate a five-state deterministic model with individuals in the population being classified as susceptible, acute, chronic, isolated, and recovered. Individuals suffering from acute and chronic stages of the infections are represented by acute and chronic states, respectively. The isolated state represents the chronically infected individuals getting isolated. The isolation of those with disease symptoms is probably the first infection control measure in recorded human history [14]. Over the decades, quarantine and isolation have been used to reduce the transmission of numerous emerging and reemerging human diseases such as leprosy, plague, cholera, typhus, yellow fever, smallpox, diphtheria, tuberculosis, measles, Ebola, pandemic influenza, and, more recently, SARS [15–20].Previous models of HCV, particularly the models calculating the threshold quantity basic reproduction numberR0, have included the treatment and/or vaccination and have discussed the control of the disease by looking at the role of disease transmission parameters in the reduction of R0 and the prevalence of the disease. However, these models did not account for time-dependent control strategies since their discussions are based on prevalence of the disease at equilibria. Optimal control theory has been employed to make decisions involving epidemic and biological models. The desired results and performance of the control functions depend on the different situations. Lenhart’s HIV models [21, 22] used optimal control to design the treatment strategies. Jung et al. [23] provide a very good example of deciding how to divide the efforts between two treatment strategies (case holding and case finding) of the two-strain TB model. Yan et al. [24] used an optimal isolation campaign to fight the SARS epidemic. Study control strategies produce valuable theoretical results which can be used to suggest or design epidemic control programs. Depending on a chosen goal (or goals), various objective criteria may be adopted.Our model extends previous work done on modeling the spread of hepatitis C in several key ways. First, we introduce an isolation class and qualitatively assess the effects of this isolation class on the transmission dynamics. Quarantine of individuals suspected of being exposed to a disease and the isolation of those with disease symptoms constitute what probably is the first infection control measure since the beginning of recorded human history [14]. However, almost no analysis of the effects of an isolation class on diseases with a chronic stage has been done and, therefore, our paper will be one of the first attempts to study the effect of isolation on the spread of a disease with a chronic stage. In addition, we will model the force of infection by a proportionate mixing, with the possibility of secondary infections due to contact with individuals who belong to the acute, chronic, or isolation class. Furthermore, we consider the disease-induced death rates for HCV in our model and will also take into account the possibility of recovery at every stage of the disease. These features add to the complexity of our model and make it considerably more insightful from an epidemiological perspective than previous models [11, 25].In the paper, Section2 presents a rigorous mathematical analysis of the deterministic model. It is shown that for basic reproduction number R0≤1, the disease-free equilibrium is locally and globally asymptotically stable. Further, the model has an endemic equilibrium, which exists if R0>1 and persists in this case. The effect of using isolation on population is discussed using a threshold quantity. Sensitivity and uncertainty analyses are carried out to study the impact of crucial parameters on R0. The existence of a locally stable endemic equilibrium, in case of R0>1, encourages us to use a time-dependent optimal control strategy to prevent and control HCV. In Section 3, we designed two different optimal control strategies (vaccination and isolation) and performed numerical simulations to illustrate the effects of an optimal control strategy. Conclusion is presented in the last section.


## 2. Model Formulation and Steady State Analysis

### 2.1. Model Formulation

We formulate a five-state model with individuals classified as susceptible, acute, chronic, isolated and recovered. Hep C has an extremely slow progression that makes it difficult to characterize the natural history of the disease [3]. The following assumptions will therefore be made.(1)
All infected individuals develop the acute form of Hep C first.(2)
Individuals with either the acute or chronic form of Hep C are capable of transmitting the disease.(3)
Individuals with the acute form of the disease either progress to the chronic form or recover naturally. Since the acute form of the disease is largely asymptomatic, there is little chance of treatment at this stage.(4)
There is no permanent immunity against HCV after recovery; thus the recovered individuals move back to the susceptible.In addition, the model will assume that the susceptible populationS has a constant recruitment rate Π and natural death rate μ. Susceptible individuals who get infected suffer from the acute form of hepatitis C and move to the compartment A with the force of infection given by λ. Individuals in A, in addition to the natural death rate μ, die at a disease-induced death rate δa. They also have a natural recovery rate of κ. Individuals with the acute form of the infection progress to the chronic form of the disease at a rate ξ, in which case the individual is shifted to compartment C. Individuals in C, in addition to the natural death rate μ, also die at a disease-induced death rate δc. Furthermore, these individuals recover at a rate ψ and thus move to the recovered compartment S. Also, the individuals in compartment C are isolated and moved to compartment Q at a rate α. Individuals in Q, in addition to the natural death rate μ, also die at a disease-induced death rate δq. Isolated individuals can either recover at a rate γf or become acutely infected with HCV at a rate γ(1-f). The individuals in R are prone to the natural death rate μ, or they become susceptible again and enter the S compartment at a rate of ω. Recovery from HCV does not result in immunity.Mathematically, the model is as follows:(1)dSdt=Π+ωR-λS-μS,dAdt=λS+γ(1-f)Q-(ξ+κ+μ+δa)A,dCdt=ξA-(α+ψ+μ+δc)C,dQdt=αC-(γ+μ+δq)Q,dRdt=γfQ+κA+ψC-(ω+μ)R,
where
(2)λ=β(ηA+C+ζQ)N.
The descriptions of variables and parameters of model (1) are as follows.Variable
N(t): Total populationS(t): Population of susceptible individualsA(t): Population of individuals with acute Hep CC(t): Population of individuals with chronic Hep CQ(t): Population of isolated individualsR(t): Population of recovered individuals.Parameter
Π: Recruitment rateμ: Natural death rateδa: Disease-induced death rate for individuals with acute Hep Cδc: Disease-induced death rate for individuals with chronic Hep Cδq: Disease-induced death rate for isolated individualsγ: Recovery rate of isolationsf: Fraction of isolated that becomes susceptiblesξ: Progression rate from acute to chronicα: Isolation rate of chronic individualsκ: Natural recovery rate for acute individualsψ: Recovery rate for chronic individualsω: Progression rate of recovered individuals to susceptible individualsβ: Effective contact rateη: Modification parameter for reduction in infectiousness of acute individualsζ: Modification parameter for reduction in infectiousness of quarantined individuals.Since (1) is a model for human populations, all the associated parameters are nonnegative. Furthermore, the following result holds and can be proved easily.Lemma 1.
The variables of model (1) are nonnegative for all time t>0. In other words, solutions of the system (1) with positive initial data will exist and remain positive for all t>0. Moreover, the closed set
(3)D={(S,A,C,Q,R)∈R+5:S+A+C+Q+R≤Πμ}
is positively invariant and global attractor.Since the regionD is positively invariant and global attractor, it is sufficient to consider the dynamics of the flow generated by model (1) (Figure 1) in region D, where the model is considered to be epidemiologically and mathematically well posed.Figure 1
Schematic diagram of model (1).


### 2.2. Disease-Free Equilibrium (DFE)

In this section, we discuss the existence and stability of the disease-free equilibrium (DFE).


#### 2.2.1. Local Stability

Model (1) has a disease-free equilibrium DFE, obtained by setting the right-hand sides of the equations in (1) to zero, given by
(4)ℵ0=(S*,A*,C*,Q*,R*)=(Πμ,0,0,0,0).
The local stability property of ℵ0 will be determined using the next generation operator method described in [26]. The nonnegative matrix F, of the new infection terms, and the M-matrix, V, of the transition terms associated with model (1), are given by
(5)F=(βηββζ000000),V=(ξ+κ+μ+δa0-γ(1-f)-ξα+ψ+μ+δc00-αγ+μ+δq).The eigenvalues of matrixFV-1 are
(6){0,0,β[ηk2k3+ξk3+ζαξ]k1k2k3-αξk4}.It follows that the basic reproduction numberR0=ρ(FV-1) is given by
(7)R0=β[ηk2k3+ξk3+ζαξ]k1k2k3-αξk4>0,
where
(8)k1=(ξ+κ+μ+δa),k2=(α+ψ+μ+δc),k3=(γ+μ+δq),k4=γ(1-f).
The basic reproduction number is interpreted as the average number of new infections that one infectious individual can produce if introduced into a population composed of susceptible individuals. Susceptible individuals acquire infection following contact with either an acute (A), chronic (C), or isolated (Q) individual. The number of infections produced by an acutely infected individual (near the DFE) is βη/k1 given by the product of the infection rate of an acute individual (βη) and the average duration in the acute class (1/k1). Furthermore, the number of infections produced by a chronically infected individual (near the DFE) is βξ/k1k2 given by the product of the infection rate of a chronic individual (β), the average duration in the chronic C class (1/k2), and the probability that an acute individual survives and progresses to the chronic stage ξ/k1. Similarly, the number of infections produced by an isolated individual (near the DFE) is βζξα/k1k2k3 given by the product of the infection rate of an isolated individual (βζ), the average duration in the isolated class (1/k3), and the probability that an acute individual survives and progresses to the isolated stage ξα/k1k2. Finally, we observe that a fraction ξαγ(1-f)/k1k2k3 of newly infected individuals will reenter the acute class. Thus, the average number of new infections generated by a single infectious individual is given by
(9)(βηk1+βξk1k2+βζξαk1k2k3)∑n=0∞‍[ξαγ(1-f)k1k2k3]n=(βηk1+βξk1k2+βζξαk1k2k3)[11-(ξαγ(1-f)/k1k2k3)]=R0.
The local stability of the DFE holds due to Theorem 2 of [26].Lemma 2.
The DFE,ℵ0, of model (1) is locally asymptotically stable if R0<1 and unstable if R0>1.Lemma 2.2 implies that withR0<1, a small influx of infectious individuals will not lead to a large outbreak of the disease. To ensure that disease elimination is independent of the initial sizes of subpopulations, it is necessary to show that the DFE is globally asymptotically stable if R0<1. This is explored below.


#### 2.2.2. Global Stability

Theorem 3.
The DFE of model (1), given by (4), is globally asymptotically stable whenever R0≤1.Proof.
Consider the following Lyapunov function:(10)L=aA+cC+qQ,
where
(11)a=k2k3k1k2k3-αξk4,c=βk3+α[βζ+k4]k1k2k3-αξk4,q=k2[βζ+k4]k1k2k3-αξk4.
Clearly L is positive definite. We have
(12)L˙=aA˙+cC˙+qQ˙=aλS+ak4Q-ak1A+cξA-qk3Q+qαC-ck2C=aβ[ηA+C+ζQ]NS+[ak4-qk3]Q+[cξ-ak1]A+[qα-ck2]C≤[aβη+cξ-ak1]A+[aβ+qα-ck2]C+[aβζ+ak4-qk3]Q=[βηk2k3+ξβk3+ξαβζ+ξαβk4-k1k2k3]k1k2k3-αξk4A+[k2k3β+αk2βζ+αk2k4-βk3k2-αk2βζ-αk2k4]k1k2k3-αξk4C+[k2k3βζ+k2k3k4-k2k3βζ-k2k3k4]k1k2k3-αξk4Q=(R0-1)A.
Thus,
(13)L˙≤(R0-1)A≤0,forR0<1.
It follows that L˙≤0 for R0< 1 with L˙=0 if and only if A=C=Q=0 or R0=1. Hence, L is a Lyapunov function on D.
The largest invariant set in{(S,A,C,Q,R)∈D∣L˙=0} is the singleton {ℵ0}. According to the LaSalle Invariance Principle, ℵ0 is globally asymptotically stable in D if R0<1. This means that, with R0<1, every solution to the system (1) with initial conditions in D approaches ℵ0 as t→∞:
(14)(S,A,C,Q,R)⟶ℵ0=(Πμ,0,0,0,0)ast⟶∞.The epidemiological implication of the above result is that the disease can be eliminated from the population if the basic reproduction numberR0 can be brought down to and maintained at a value less than unity (i.e., the condition R0<1 is sufficient and necessary for disease elimination). Figure 2 depicts numerical results by simulating model (1) using various initial conditions with R0<1. It is evident from the simulation that all initial solutions converged to DFE, ℵ0, in-line with Theorem 3.Figure 2
Disease-free equilibrium:R0=0.7994. Simulation shows the total infected population with different initial infected population sizes. The parameter values are given in the appendix.


### 2.3. Endemic Equilibrium

In this section the existence and stability of endemic equilibrium of model (1) will be discussed. We define endemic equilibrium to be those fixed points of the system (1) in which at least one of the infected compartments of the model is nonzero.Letℵ1=(S**,A**,C**,Q**,R**) denote an arbitrary endemic equilibrium of model (1) so that N**=S**+A**+C**+Q**+R**. Solving the equations of model (1) at steady state gives
(15)A**=k2ξC**,Q**=αk3C**,S**=1λ**(k1k2k3-αξk4ξk3)C**,R**=1k5(fγαk3+κk2ξ+ψ)C**,
where
(16)λ**=β[ηA**+C**+ζQ**]N**,k5=(ω+μ).ConsiderS**. Then, using λ** and (15) from above we have the following endemic states:
(17)A**=k2ξ[R0-1Y]S**,C**=[R0-1Y]S**,Q**=αk3[R0-1Y]S**,R**=1k5[k2k3+k3κξ+γfαk2k2k3][R0-1Y]S**,
where
(18)Y=[k2ξ+1+αk3].
Hence, we have the following result.Lemma 4.
Model (1) has endemic equilibria, given by ℵ1 (15), whenever R0>1.Now, we address the question of uniform persistence of the infected population.Theorem 5.
IfR0>1, then the disease is uniformly persistent: there exists an ϵ>0 such that
(19)lim⁡t→∞inf⁡A(t)>ϵ,lim⁡t→∞inf⁡C(t)>ϵ,lim⁡t→∞inf⁡Q(t)>ϵ,
for all solutions (S,A,C,Q,R) of (1) with A(0)>0, C(0)>0, and Q(0)>0.Proof.
LetX={(S,A,C,Q,R)∈R+5:A=C=Q=0}. Thus, X is the set of all disease-free states of (1) and it can be easily verified that X is positively invariant. Let M=D∩X. Since both D and X are positively invariant, M is also positively invariant. Also note that ℵ0∈M and ℵ0 attracts all the solutions in X. So, Ω(M)={ℵ0}. The equations for the infected components of (1) can be written as
(20)x′(t)=Y(x)x(t),
where x(t)=(A(t),C(t),Q(t))T, Y(x)=[(S/N)F-V]. It is clear that Y(ℵ0)=F-V. Also it is easy to check that Y(ℵ0) is irreducible. We will apply Lemma A.4 in [27] to show that M is a uniform weak repeller. Since ℵ0 is a steady state solution, we can consider it to be a periodic orbit of period T=1. P(t,x), the fundamental matrix of the solutions for (20), is etY. Since the spectral radius of Y(ℵ0)=R0>1, the spectral radius of eY(ℵ0)>1. So condition 2 of Lemma A.4 is satisfied. Taking x=ℵ0, we get P(T,ℵ0)=eY(ℵ0) which is a primitive matrix, because Y(ℵ0) is irreducible, as mentioned in Theorem A.12(i) [28]. This satisfies the condition 1 of Lemma A.4. Thus, M is a uniform weak repeller and disease is weakly persistent. M is trivially closed and bounded relative to D and hence, compact. Therefore, by Theorem 1.3 [29], we have that M is a uniform strong repeller and disease is uniformly persistent.The epidemiological implication of Theorem5 is that the disease will persist in the population whenever R0>1. Numerical simulation results, depicted in Figure 3 using different initial conditions, shows convergence of solutions to the ℵ1 in-line with Theorem 5.Figure 3
Endemic equilibrium:R0=1.5545. Simulation shows the total infected population with different initial infected population sizes. The parameter values are given in the appendix.


### 2.4. Sensitivity Analysis

The asymptotic dynamics of the model are completely determined by the threshold quantityR0, which determines the prevalence of the disease. Since we have a deterministic model, the only uncertainty is generated by the input variation and parameters. Therefore, we present parameter-related global uncertainty and sensitivity analyses on R0. Parameter estimates can be uncertain because of many reasons including natural variation, error in measurements, or a lack of measuring techniques. Uncertainty analysis qualitatively decides which parameters are most influential in the model output (R0 in our case) and quantifies the degree of confidence in the existing data and parameter estimation. On the other hand, the sensitivity analysis identifies critical model parameters and quantifies the impact of each input parameter on the value of an output, in the presence of other input parameters.Ideally, uncertainty and sensitivity analyses should be performed simultaneously. Here we use the Latin-hypercube sampling based method to quantify uncertainty and sensitivity ofR0 as a function of 12 model parameters (μ, γ, ξ, α, κ, ψ, β, η, ζ, δa, δc and δq). For the sensitivity analysis, partial rank correlation coefficient (PRCC) measures the impact of the parameters on the output variable. PRCC provides a measure of monotonicity after the removal of the linear effects of all but one variable. PRCC method uses the rank transformation of the data (i.e., replacing the values with their ranks) to reduce the effects of nonlinearity. The Rank Correlation Coefficient (RCC) indicates the degree of monotonicity between the input and output variables. The resultant data are considered partially in some sense, that is, partial rank correlation coefficients (PRCC) are computed that take into account the correlations among other input variables. The basic reproduction number R0 is the output measure in the sensitivity and uncertainty analyses.The assumed distributions of the model parameters used in the two analyses are mentioned in the appendix. Our estimate ofR0 for Hep C from uncertainty analysis is 1.33 with 95% CI (1.1, 1.95) as shown in Figure 4. The probability that R0>1 is 90%. This suggests that Hep C will get endemic under the preset conditions. However, the time taken to reach that state could be large.Figure 4
Uncertainty analysis: the probability thatR0>1 is 90% with 95% confidence interval (1.1, 1.95). This suggests that hepatitis C will get endemic under the present conditions. However, the time taken to reach that state could be large. 10000 values were generated for each parameter according to their distributions and mean values. Values of parameters given in Appendix were used to calculate R0.The sensitivity analysis suggests that the most significant (PRCC values above 0.5 or below −0.5 in Figure5) sensitivity parameters to R0 are α, κ, β, and ζ. This suggests that these parameters need to be estimated with precision to capture the transmission dynamics of the Hep C. The analyses further suggest that the isolation strategy aimed to reduce the infected population yields the desired result.Figure 5
Sensitivity analysis: the proportion of chronically infected being quarantinedα, proportion of acute infections recovering spontaneously κ1, effective contact rate β, and modification parameter for infectiousness of quarantined ζ are the most significant parameters. This means that even a small error in the estimation of these parameters can greatly affect the value of R0 and hence the analysis of our model. Partial rank correlation coefficients (PRCC) are calculated with respect to R0. Parameters with modulus of PRCC values in excess of 0.5 are declared sensitive to R0.Since we are interested in the influence of critical model parameters on the basic reproductive number and hence the prevalence of chronic Hep C, we conduct numerical simulation to investigate it. In order to qualitatively measure the effect of isolation on the transmission dynamics of Hep C, a threshold analysis of the parameter associated with the isolation of chronically infected individuals is discussed (α). We computed the partial derivative of R0 with respect to this parameter.For the case of the isolation of chronically infected individuals, it is easy to see that(21)∂R0∂α=((βηk1+βζξ)(k1k2k3-αξk4)-β(ηk1k2+ξk3+ζαξ)(k1k3-ξk4))×((k1k2k3-αξk4)2)-1
which simplifies to
(22)∂R0∂α=βξk1k3[ζ(k2-α)+ηk4(k2-α)k3+ξk4k1-k3],
where
(23)k2-α=ψ+μ+δc>0.
Now we have
(24)∂R0∂α<0(>0)iffk3′<k3(>k3),
where
(25)k3′=ζ(k2-α)+ηk4(k2-α)k3+ξk4k1>0.
Thus, the isolation of chronically infected individuals will reduce R0 and, therefore, reduce disease burden (new infections, mortality, etc.) if k3′ does not exceed the threshold k3 (which is the rate of transfer of individuals out of the isolation state). This case is presented in Figure 6(a).Effect of isolation on the infected population.
(a)
Positive effect of isolation measures on the infected population withk3′=0.0414<k3=.2017
(b)
Negative effect of isolation measures on the infected population withk3′=0.0803>k3=.0517On the other hand, ifk3′>k3, then the use of isolation (of chronically infected individuals) will increase R0, and, consequently, increase disease burden (hence, the use of isolation is detrimental to the community in this case). This case is presented in Figure 6(b). This important result is summarized below.Lemma 6.
The use of isolation of the chronically infected individuals will have positive (negative) population-level impact ifk3′<(>)k3.Now we present the simulations of the critical parameters (as identified by sensitivity analysis) andR0. Figure 7 presents the dependence of the basic reproductive number on the parameters α and ξ, where α denotes the isolation rate of chronic and ξ denotes the progression rate to chronic from acute. From the contour plot, we see that if ξ is larger, then R0 is always greater than one, which implies that it is important to control the acute Hep C. Figure 7(b) shows that the basic reproductive number may be less than one if α and ξ can be restricted to a range, leading to the potential extinction of the disease.Plots of the basic reproductive numberR0 in terms of the parameters α and ξ, which show the estimated effects of α and ξ on R0. (a) A contour plot of the surface R0 for the values of R0=1,2,3. (b) Two surfaces, R0 and the constant 1, are plotted to show the curve on which R0=1.
(a)
(b)The dependence of basic reproductive numberR0 on the recovery rate κ, isolation rate α, and effective contact rate is explored in Figure 8. From Figure 8(a), it is clear that high isolation rate with low effective contact will result in smaller value of R0. Furthermore, R0 is very sensitive to β and basic reproductive number increases sharply when β is slightly increased. Therefore, keeping the effective contact rate low will result in disease extinction. In Figure 8(b), larger recovery rate of chronic individuals κ results in smaller values of R0. However, still the R0 increases as β increases but smoothly and not sharply as seen in Figure 8(a).Plots of the basic reproductive numberR0 in terms of the parameters β, α, and κ, which show the effects of β, α, and κ on R0. (a) Effect of isolation rate α and effective contact rate β on the R0. (b) Effect of recovery rate of chronic κ and effective contact rate β on the R0.
(a)
(b)


## 2.1. Model Formulation

We formulate a five-state model with individuals classified as susceptible, acute, chronic, isolated and recovered. Hep C has an extremely slow progression that makes it difficult to characterize the natural history of the disease [3]. The following assumptions will therefore be made.(1)
All infected individuals develop the acute form of Hep C first.(2)
Individuals with either the acute or chronic form of Hep C are capable of transmitting the disease.(3)
Individuals with the acute form of the disease either progress to the chronic form or recover naturally. Since the acute form of the disease is largely asymptomatic, there is little chance of treatment at this stage.(4)
There is no permanent immunity against HCV after recovery; thus the recovered individuals move back to the susceptible.In addition, the model will assume that the susceptible populationS has a constant recruitment rate Π and natural death rate μ. Susceptible individuals who get infected suffer from the acute form of hepatitis C and move to the compartment A with the force of infection given by λ. Individuals in A, in addition to the natural death rate μ, die at a disease-induced death rate δa. They also have a natural recovery rate of κ. Individuals with the acute form of the infection progress to the chronic form of the disease at a rate ξ, in which case the individual is shifted to compartment C. Individuals in C, in addition to the natural death rate μ, also die at a disease-induced death rate δc. Furthermore, these individuals recover at a rate ψ and thus move to the recovered compartment S. Also, the individuals in compartment C are isolated and moved to compartment Q at a rate α. Individuals in Q, in addition to the natural death rate μ, also die at a disease-induced death rate δq. Isolated individuals can either recover at a rate γf or become acutely infected with HCV at a rate γ(1-f). The individuals in R are prone to the natural death rate μ, or they become susceptible again and enter the S compartment at a rate of ω. Recovery from HCV does not result in immunity.Mathematically, the model is as follows:(1)dSdt=Π+ωR-λS-μS,dAdt=λS+γ(1-f)Q-(ξ+κ+μ+δa)A,dCdt=ξA-(α+ψ+μ+δc)C,dQdt=αC-(γ+μ+δq)Q,dRdt=γfQ+κA+ψC-(ω+μ)R,
where
(2)λ=β(ηA+C+ζQ)N.
The descriptions of variables and parameters of model (1) are as follows.Variable
N(t): Total populationS(t): Population of susceptible individualsA(t): Population of individuals with acute Hep CC(t): Population of individuals with chronic Hep CQ(t): Population of isolated individualsR(t): Population of recovered individuals.Parameter
Π: Recruitment rateμ: Natural death rateδa: Disease-induced death rate for individuals with acute Hep Cδc: Disease-induced death rate for individuals with chronic Hep Cδq: Disease-induced death rate for isolated individualsγ: Recovery rate of isolationsf: Fraction of isolated that becomes susceptiblesξ: Progression rate from acute to chronicα: Isolation rate of chronic individualsκ: Natural recovery rate for acute individualsψ: Recovery rate for chronic individualsω: Progression rate of recovered individuals to susceptible individualsβ: Effective contact rateη: Modification parameter for reduction in infectiousness of acute individualsζ: Modification parameter for reduction in infectiousness of quarantined individuals.Since (1) is a model for human populations, all the associated parameters are nonnegative. Furthermore, the following result holds and can be proved easily.Lemma 1.
The variables of model (1) are nonnegative for all time t>0. In other words, solutions of the system (1) with positive initial data will exist and remain positive for all t>0. Moreover, the closed set
(3)D={(S,A,C,Q,R)∈R+5:S+A+C+Q+R≤Πμ}
is positively invariant and global attractor.Since the regionD is positively invariant and global attractor, it is sufficient to consider the dynamics of the flow generated by model (1) (Figure 1) in region D, where the model is considered to be epidemiologically and mathematically well posed.Figure 1
Schematic diagram of model (1).


## 2.2. Disease-Free Equilibrium (DFE)

In this section, we discuss the existence and stability of the disease-free equilibrium (DFE).


### 2.2.1. Local Stability

Model (1) has a disease-free equilibrium DFE, obtained by setting the right-hand sides of the equations in (1) to zero, given by
(4)ℵ0=(S*,A*,C*,Q*,R*)=(Πμ,0,0,0,0).
The local stability property of ℵ0 will be determined using the next generation operator method described in [26]. The nonnegative matrix F, of the new infection terms, and the M-matrix, V, of the transition terms associated with model (1), are given by
(5)F=(βηββζ000000),V=(ξ+κ+μ+δa0-γ(1-f)-ξα+ψ+μ+δc00-αγ+μ+δq).The eigenvalues of matrixFV-1 are
(6){0,0,β[ηk2k3+ξk3+ζαξ]k1k2k3-αξk4}.It follows that the basic reproduction numberR0=ρ(FV-1) is given by
(7)R0=β[ηk2k3+ξk3+ζαξ]k1k2k3-αξk4>0,
where
(8)k1=(ξ+κ+μ+δa),k2=(α+ψ+μ+δc),k3=(γ+μ+δq),k4=γ(1-f).
The basic reproduction number is interpreted as the average number of new infections that one infectious individual can produce if introduced into a population composed of susceptible individuals. Susceptible individuals acquire infection following contact with either an acute (A), chronic (C), or isolated (Q) individual. The number of infections produced by an acutely infected individual (near the DFE) is βη/k1 given by the product of the infection rate of an acute individual (βη) and the average duration in the acute class (1/k1). Furthermore, the number of infections produced by a chronically infected individual (near the DFE) is βξ/k1k2 given by the product of the infection rate of a chronic individual (β), the average duration in the chronic C class (1/k2), and the probability that an acute individual survives and progresses to the chronic stage ξ/k1. Similarly, the number of infections produced by an isolated individual (near the DFE) is βζξα/k1k2k3 given by the product of the infection rate of an isolated individual (βζ), the average duration in the isolated class (1/k3), and the probability that an acute individual survives and progresses to the isolated stage ξα/k1k2. Finally, we observe that a fraction ξαγ(1-f)/k1k2k3 of newly infected individuals will reenter the acute class. Thus, the average number of new infections generated by a single infectious individual is given by
(9)(βηk1+βξk1k2+βζξαk1k2k3)∑n=0∞‍[ξαγ(1-f)k1k2k3]n=(βηk1+βξk1k2+βζξαk1k2k3)[11-(ξαγ(1-f)/k1k2k3)]=R0.
The local stability of the DFE holds due to Theorem 2 of [26].Lemma 2.
The DFE,ℵ0, of model (1) is locally asymptotically stable if R0<1 and unstable if R0>1.Lemma 2.2 implies that withR0<1, a small influx of infectious individuals will not lead to a large outbreak of the disease. To ensure that disease elimination is independent of the initial sizes of subpopulations, it is necessary to show that the DFE is globally asymptotically stable if R0<1. This is explored below.


### 2.2.2. Global Stability

Theorem 3.
The DFE of model (1), given by (4), is globally asymptotically stable whenever R0≤1.Proof.
Consider the following Lyapunov function:(10)L=aA+cC+qQ,
where
(11)a=k2k3k1k2k3-αξk4,c=βk3+α[βζ+k4]k1k2k3-αξk4,q=k2[βζ+k4]k1k2k3-αξk4.
Clearly L is positive definite. We have
(12)L˙=aA˙+cC˙+qQ˙=aλS+ak4Q-ak1A+cξA-qk3Q+qαC-ck2C=aβ[ηA+C+ζQ]NS+[ak4-qk3]Q+[cξ-ak1]A+[qα-ck2]C≤[aβη+cξ-ak1]A+[aβ+qα-ck2]C+[aβζ+ak4-qk3]Q=[βηk2k3+ξβk3+ξαβζ+ξαβk4-k1k2k3]k1k2k3-αξk4A+[k2k3β+αk2βζ+αk2k4-βk3k2-αk2βζ-αk2k4]k1k2k3-αξk4C+[k2k3βζ+k2k3k4-k2k3βζ-k2k3k4]k1k2k3-αξk4Q=(R0-1)A.
Thus,
(13)L˙≤(R0-1)A≤0,forR0<1.
It follows that L˙≤0 for R0< 1 with L˙=0 if and only if A=C=Q=0 or R0=1. Hence, L is a Lyapunov function on D.
The largest invariant set in{(S,A,C,Q,R)∈D∣L˙=0} is the singleton {ℵ0}. According to the LaSalle Invariance Principle, ℵ0 is globally asymptotically stable in D if R0<1. This means that, with R0<1, every solution to the system (1) with initial conditions in D approaches ℵ0 as t→∞:
(14)(S,A,C,Q,R)⟶ℵ0=(Πμ,0,0,0,0)ast⟶∞.The epidemiological implication of the above result is that the disease can be eliminated from the population if the basic reproduction numberR0 can be brought down to and maintained at a value less than unity (i.e., the condition R0<1 is sufficient and necessary for disease elimination). Figure 2 depicts numerical results by simulating model (1) using various initial conditions with R0<1. It is evident from the simulation that all initial solutions converged to DFE, ℵ0, in-line with Theorem 3.Figure 2
Disease-free equilibrium:R0=0.7994. Simulation shows the total infected population with different initial infected population sizes. The parameter values are given in the appendix.


## 2.2.1. Local Stability

Model (1) has a disease-free equilibrium DFE, obtained by setting the right-hand sides of the equations in (1) to zero, given by
(4)ℵ0=(S*,A*,C*,Q*,R*)=(Πμ,0,0,0,0).
The local stability property of ℵ0 will be determined using the next generation operator method described in [26]. The nonnegative matrix F, of the new infection terms, and the M-matrix, V, of the transition terms associated with model (1), are given by
(5)F=(βηββζ000000),V=(ξ+κ+μ+δa0-γ(1-f)-ξα+ψ+μ+δc00-αγ+μ+δq).The eigenvalues of matrixFV-1 are
(6){0,0,β[ηk2k3+ξk3+ζαξ]k1k2k3-αξk4}.It follows that the basic reproduction numberR0=ρ(FV-1) is given by
(7)R0=β[ηk2k3+ξk3+ζαξ]k1k2k3-αξk4>0,
where
(8)k1=(ξ+κ+μ+δa),k2=(α+ψ+μ+δc),k3=(γ+μ+δq),k4=γ(1-f).
The basic reproduction number is interpreted as the average number of new infections that one infectious individual can produce if introduced into a population composed of susceptible individuals. Susceptible individuals acquire infection following contact with either an acute (A), chronic (C), or isolated (Q) individual. The number of infections produced by an acutely infected individual (near the DFE) is βη/k1 given by the product of the infection rate of an acute individual (βη) and the average duration in the acute class (1/k1). Furthermore, the number of infections produced by a chronically infected individual (near the DFE) is βξ/k1k2 given by the product of the infection rate of a chronic individual (β), the average duration in the chronic C class (1/k2), and the probability that an acute individual survives and progresses to the chronic stage ξ/k1. Similarly, the number of infections produced by an isolated individual (near the DFE) is βζξα/k1k2k3 given by the product of the infection rate of an isolated individual (βζ), the average duration in the isolated class (1/k3), and the probability that an acute individual survives and progresses to the isolated stage ξα/k1k2. Finally, we observe that a fraction ξαγ(1-f)/k1k2k3 of newly infected individuals will reenter the acute class. Thus, the average number of new infections generated by a single infectious individual is given by
(9)(βηk1+βξk1k2+βζξαk1k2k3)∑n=0∞‍[ξαγ(1-f)k1k2k3]n=(βηk1+βξk1k2+βζξαk1k2k3)[11-(ξαγ(1-f)/k1k2k3)]=R0.
The local stability of the DFE holds due to Theorem 2 of [26].Lemma 2.
The DFE,ℵ0, of model (1) is locally asymptotically stable if R0<1 and unstable if R0>1.Lemma 2.2 implies that withR0<1, a small influx of infectious individuals will not lead to a large outbreak of the disease. To ensure that disease elimination is independent of the initial sizes of subpopulations, it is necessary to show that the DFE is globally asymptotically stable if R0<1. This is explored below.


## 2.2.2. Global Stability

Theorem 3.
The DFE of model (1), given by (4), is globally asymptotically stable whenever R0≤1.Proof.
Consider the following Lyapunov function:(10)L=aA+cC+qQ,
where
(11)a=k2k3k1k2k3-αξk4,c=βk3+α[βζ+k4]k1k2k3-αξk4,q=k2[βζ+k4]k1k2k3-αξk4.
Clearly L is positive definite. We have
(12)L˙=aA˙+cC˙+qQ˙=aλS+ak4Q-ak1A+cξA-qk3Q+qαC-ck2C=aβ[ηA+C+ζQ]NS+[ak4-qk3]Q+[cξ-ak1]A+[qα-ck2]C≤[aβη+cξ-ak1]A+[aβ+qα-ck2]C+[aβζ+ak4-qk3]Q=[βηk2k3+ξβk3+ξαβζ+ξαβk4-k1k2k3]k1k2k3-αξk4A+[k2k3β+αk2βζ+αk2k4-βk3k2-αk2βζ-αk2k4]k1k2k3-αξk4C+[k2k3βζ+k2k3k4-k2k3βζ-k2k3k4]k1k2k3-αξk4Q=(R0-1)A.
Thus,
(13)L˙≤(R0-1)A≤0,forR0<1.
It follows that L˙≤0 for R0< 1 with L˙=0 if and only if A=C=Q=0 or R0=1. Hence, L is a Lyapunov function on D.
The largest invariant set in{(S,A,C,Q,R)∈D∣L˙=0} is the singleton {ℵ0}. According to the LaSalle Invariance Principle, ℵ0 is globally asymptotically stable in D if R0<1. This means that, with R0<1, every solution to the system (1) with initial conditions in D approaches ℵ0 as t→∞:
(14)(S,A,C,Q,R)⟶ℵ0=(Πμ,0,0,0,0)ast⟶∞.The epidemiological implication of the above result is that the disease can be eliminated from the population if the basic reproduction numberR0 can be brought down to and maintained at a value less than unity (i.e., the condition R0<1 is sufficient and necessary for disease elimination). Figure 2 depicts numerical results by simulating model (1) using various initial conditions with R0<1. It is evident from the simulation that all initial solutions converged to DFE, ℵ0, in-line with Theorem 3.Figure 2
Disease-free equilibrium:R0=0.7994. Simulation shows the total infected population with different initial infected population sizes. The parameter values are given in the appendix.


## 2.3. Endemic Equilibrium

In this section the existence and stability of endemic equilibrium of model (1) will be discussed. We define endemic equilibrium to be those fixed points of the system (1) in which at least one of the infected compartments of the model is nonzero.Letℵ1=(S**,A**,C**,Q**,R**) denote an arbitrary endemic equilibrium of model (1) so that N**=S**+A**+C**+Q**+R**. Solving the equations of model (1) at steady state gives
(15)A**=k2ξC**,Q**=αk3C**,S**=1λ**(k1k2k3-αξk4ξk3)C**,R**=1k5(fγαk3+κk2ξ+ψ)C**,
where
(16)λ**=β[ηA**+C**+ζQ**]N**,k5=(ω+μ).ConsiderS**. Then, using λ** and (15) from above we have the following endemic states:
(17)A**=k2ξ[R0-1Y]S**,C**=[R0-1Y]S**,Q**=αk3[R0-1Y]S**,R**=1k5[k2k3+k3κξ+γfαk2k2k3][R0-1Y]S**,
where
(18)Y=[k2ξ+1+αk3].
Hence, we have the following result.Lemma 4.
Model (1) has endemic equilibria, given by ℵ1 (15), whenever R0>1.Now, we address the question of uniform persistence of the infected population.Theorem 5.
IfR0>1, then the disease is uniformly persistent: there exists an ϵ>0 such that
(19)lim⁡t→∞inf⁡A(t)>ϵ,lim⁡t→∞inf⁡C(t)>ϵ,lim⁡t→∞inf⁡Q(t)>ϵ,
for all solutions (S,A,C,Q,R) of (1) with A(0)>0, C(0)>0, and Q(0)>0.Proof.
LetX={(S,A,C,Q,R)∈R+5:A=C=Q=0}. Thus, X is the set of all disease-free states of (1) and it can be easily verified that X is positively invariant. Let M=D∩X. Since both D and X are positively invariant, M is also positively invariant. Also note that ℵ0∈M and ℵ0 attracts all the solutions in X. So, Ω(M)={ℵ0}. The equations for the infected components of (1) can be written as
(20)x′(t)=Y(x)x(t),
where x(t)=(A(t),C(t),Q(t))T, Y(x)=[(S/N)F-V]. It is clear that Y(ℵ0)=F-V. Also it is easy to check that Y(ℵ0) is irreducible. We will apply Lemma A.4 in [27] to show that M is a uniform weak repeller. Since ℵ0 is a steady state solution, we can consider it to be a periodic orbit of period T=1. P(t,x), the fundamental matrix of the solutions for (20), is etY. Since the spectral radius of Y(ℵ0)=R0>1, the spectral radius of eY(ℵ0)>1. So condition 2 of Lemma A.4 is satisfied. Taking x=ℵ0, we get P(T,ℵ0)=eY(ℵ0) which is a primitive matrix, because Y(ℵ0) is irreducible, as mentioned in Theorem A.12(i) [28]. This satisfies the condition 1 of Lemma A.4. Thus, M is a uniform weak repeller and disease is weakly persistent. M is trivially closed and bounded relative to D and hence, compact. Therefore, by Theorem 1.3 [29], we have that M is a uniform strong repeller and disease is uniformly persistent.The epidemiological implication of Theorem5 is that the disease will persist in the population whenever R0>1. Numerical simulation results, depicted in Figure 3 using different initial conditions, shows convergence of solutions to the ℵ1 in-line with Theorem 5.Figure 3
Endemic equilibrium:R0=1.5545. Simulation shows the total infected population with different initial infected population sizes. The parameter values are given in the appendix.


## 2.4. Sensitivity Analysis

The asymptotic dynamics of the model are completely determined by the threshold quantityR0, which determines the prevalence of the disease. Since we have a deterministic model, the only uncertainty is generated by the input variation and parameters. Therefore, we present parameter-related global uncertainty and sensitivity analyses on R0. Parameter estimates can be uncertain because of many reasons including natural variation, error in measurements, or a lack of measuring techniques. Uncertainty analysis qualitatively decides which parameters are most influential in the model output (R0 in our case) and quantifies the degree of confidence in the existing data and parameter estimation. On the other hand, the sensitivity analysis identifies critical model parameters and quantifies the impact of each input parameter on the value of an output, in the presence of other input parameters.Ideally, uncertainty and sensitivity analyses should be performed simultaneously. Here we use the Latin-hypercube sampling based method to quantify uncertainty and sensitivity ofR0 as a function of 12 model parameters (μ, γ, ξ, α, κ, ψ, β, η, ζ, δa, δc and δq). For the sensitivity analysis, partial rank correlation coefficient (PRCC) measures the impact of the parameters on the output variable. PRCC provides a measure of monotonicity after the removal of the linear effects of all but one variable. PRCC method uses the rank transformation of the data (i.e., replacing the values with their ranks) to reduce the effects of nonlinearity. The Rank Correlation Coefficient (RCC) indicates the degree of monotonicity between the input and output variables. The resultant data are considered partially in some sense, that is, partial rank correlation coefficients (PRCC) are computed that take into account the correlations among other input variables. The basic reproduction number R0 is the output measure in the sensitivity and uncertainty analyses.The assumed distributions of the model parameters used in the two analyses are mentioned in the appendix. Our estimate ofR0 for Hep C from uncertainty analysis is 1.33 with 95% CI (1.1, 1.95) as shown in Figure 4. The probability that R0>1 is 90%. This suggests that Hep C will get endemic under the preset conditions. However, the time taken to reach that state could be large.Figure 4
Uncertainty analysis: the probability thatR0>1 is 90% with 95% confidence interval (1.1, 1.95). This suggests that hepatitis C will get endemic under the present conditions. However, the time taken to reach that state could be large. 10000 values were generated for each parameter according to their distributions and mean values. Values of parameters given in Appendix were used to calculate R0.The sensitivity analysis suggests that the most significant (PRCC values above 0.5 or below −0.5 in Figure5) sensitivity parameters to R0 are α, κ, β, and ζ. This suggests that these parameters need to be estimated with precision to capture the transmission dynamics of the Hep C. The analyses further suggest that the isolation strategy aimed to reduce the infected population yields the desired result.Figure 5
Sensitivity analysis: the proportion of chronically infected being quarantinedα, proportion of acute infections recovering spontaneously κ1, effective contact rate β, and modification parameter for infectiousness of quarantined ζ are the most significant parameters. This means that even a small error in the estimation of these parameters can greatly affect the value of R0 and hence the analysis of our model. Partial rank correlation coefficients (PRCC) are calculated with respect to R0. Parameters with modulus of PRCC values in excess of 0.5 are declared sensitive to R0.Since we are interested in the influence of critical model parameters on the basic reproductive number and hence the prevalence of chronic Hep C, we conduct numerical simulation to investigate it. In order to qualitatively measure the effect of isolation on the transmission dynamics of Hep C, a threshold analysis of the parameter associated with the isolation of chronically infected individuals is discussed (α). We computed the partial derivative of R0 with respect to this parameter.For the case of the isolation of chronically infected individuals, it is easy to see that(21)∂R0∂α=((βηk1+βζξ)(k1k2k3-αξk4)-β(ηk1k2+ξk3+ζαξ)(k1k3-ξk4))×((k1k2k3-αξk4)2)-1
which simplifies to
(22)∂R0∂α=βξk1k3[ζ(k2-α)+ηk4(k2-α)k3+ξk4k1-k3],
where
(23)k2-α=ψ+μ+δc>0.
Now we have
(24)∂R0∂α<0(>0)iffk3′<k3(>k3),
where
(25)k3′=ζ(k2-α)+ηk4(k2-α)k3+ξk4k1>0.
Thus, the isolation of chronically infected individuals will reduce R0 and, therefore, reduce disease burden (new infections, mortality, etc.) if k3′ does not exceed the threshold k3 (which is the rate of transfer of individuals out of the isolation state). This case is presented in Figure 6(a).Effect of isolation on the infected population.
(a)
Positive effect of isolation measures on the infected population withk3′=0.0414<k3=.2017
(b)
Negative effect of isolation measures on the infected population withk3′=0.0803>k3=.0517On the other hand, ifk3′>k3, then the use of isolation (of chronically infected individuals) will increase R0, and, consequently, increase disease burden (hence, the use of isolation is detrimental to the community in this case). This case is presented in Figure 6(b). This important result is summarized below.Lemma 6.
The use of isolation of the chronically infected individuals will have positive (negative) population-level impact ifk3′<(>)k3.Now we present the simulations of the critical parameters (as identified by sensitivity analysis) andR0. Figure 7 presents the dependence of the basic reproductive number on the parameters α and ξ, where α denotes the isolation rate of chronic and ξ denotes the progression rate to chronic from acute. From the contour plot, we see that if ξ is larger, then R0 is always greater than one, which implies that it is important to control the acute Hep C. Figure 7(b) shows that the basic reproductive number may be less than one if α and ξ can be restricted to a range, leading to the potential extinction of the disease.Plots of the basic reproductive numberR0 in terms of the parameters α and ξ, which show the estimated effects of α and ξ on R0. (a) A contour plot of the surface R0 for the values of R0=1,2,3. (b) Two surfaces, R0 and the constant 1, are plotted to show the curve on which R0=1.
(a)
(b)The dependence of basic reproductive numberR0 on the recovery rate κ, isolation rate α, and effective contact rate is explored in Figure 8. From Figure 8(a), it is clear that high isolation rate with low effective contact will result in smaller value of R0. Furthermore, R0 is very sensitive to β and basic reproductive number increases sharply when β is slightly increased. Therefore, keeping the effective contact rate low will result in disease extinction. In Figure 8(b), larger recovery rate of chronic individuals κ results in smaller values of R0. However, still the R0 increases as β increases but smoothly and not sharply as seen in Figure 8(a).Plots of the basic reproductive numberR0 in terms of the parameters β, α, and κ, which show the effects of β, α, and κ on R0. (a) Effect of isolation rate α and effective contact rate β on the R0. (b) Effect of recovery rate of chronic κ and effective contact rate β on the R0.
(a)
(b)


## 3. Optimal Control Strategies

Pontryagin and Boltyanskii [30] formulated the optimal control theory for models with underlying dynamics defined by a system of ordinary differential equations. The theory, application areas, and the numerical methods have progressed considerably. Pontryagin’s Maximum Principle allows us to adjust the control in a model to achieve the desired results. The control parameters are mostly functions of time appearing as coefficients in the model.Optimal control theory has been employed to make decisions involving epidemic and biological models. The desired results and performance of the control functions depend on the different situations. Lenhart’s HIV models [21, 22] used optimal control to design the treatment strategies. Jung et al. [23] provide a very good example of deciding how to divide the efforts between two treatment strategies (case holding and case finding) of the two-strain TB model. Yan et al. [24] used an optimal isolation strategy to fight the SARS epidemic. In [31] Joshi formulated two control functions as coefficients of the ODE system representing treatment effects in a two-drug regime in an HIV immunology model. The goal was to maximize the concentration of T cells while minimizing the toxic effects of the drug. The analytic and numerical results illustrated the level of two drugs to be used over the chosen time interval. The required balancing effect between two competing goals was well predicted by optimal control theory. Behncke [32] studied SIR models including vaccination, isolation, and health promotion campaign and obtained analytical results for optimal control. The optimal control intervention policies for stochastic epidemic models were treated by Clancy [33].Pontryagin’s Maximum Principle appends an adjoint system of differential equations with terminal boundary conditions to the original model (state system) of differential equations, in the attempt to characterize an optimal control. The optimality system, which characterizes the optimal controls, consists of the differential equations of the original model (state system) along with the adjoint differential equations (adjoint system). The adjoint system has the same number of equations as in the state system. The adjoint functions behave very similar to the Lagrange multipliers (appending constraints to the function of several variables to be maximized or minimized). The adjoint variables maximize or minimize the state variables with respect to the desired objective functional. The details of the necessary conditions for the adjoint and optimal controls are presented here [30, 34, 35]. For the application of these results see [21].Optimal control techniques are used to optimize the models given by system of differential equations. While formulating an optimal control problem, deciding how and where to introduce the control (through vaccination, drug treatment etc.) in the system of differential equation is very important. The formulation of the optimal control problem must be a reasonable and practical representation of the situation to be considered. The form of the optimal control depends heavily on the system being analyzed and the objective functional to be optimized. We will consider a quadratic dependence on the control in the objective functional.The questions of the existence and uniqueness of the optimal control in an optimality system can be dealt with using the Lipschitz properties of the differential equations [30, 34]. For a detailed example, see the work of Fister et al. [21] and Fleming and Rishel [34]. After establishing the existence and the uniqueness results, we can confidently continue to numerically solve the optimality system to get the desired optimal control.In this section, we consider two different time-dependent control strategies to prevent and control the spread of Hep C in the population. The first strategy is to introduce vaccination to our population and see the effects on the prevalence of the Hep C. We introduce a vaccination control, which is a function of time, in our model of differential equations to minimize the population of infected individual while keeping the cost of the vaccination to minimum through objective functional. Next strategy introduces control function for the isolation rate of chronically infected individuals to minimize the total infected population. The objective functional was designed carefully to minimize the infected population and the cost of the isolation facility.


### 3.1. Vaccination Control

The use of vaccination is an important disease control and prevention measure. Optimal vaccination control strategy for an SIR model has been devised using dynamic programming technique [36]. Also the optimal control strategies have been investigated for TB control [23, 37]. We have analyzed the dynamics of the model without vaccination in the last section which had a locally stable equilibrium. In this section we will explore the effects of a vaccination campaign on our deterministic model of Hep C by including vaccination. After including vaccination the model is given as follows:
(26)dSdt=Π+ωR-λ(1-v)S-(1-σ)λvS-μS-σvS,dAdt=λ(1-v)S+(1-σ)λvS+γ(1-f)Q-(ξ+κ+μ+δa)A,dCdt=ξA-(α+ψ1+μ+δc)C,dQdt=αC-(γ+μ+δq)Q,dRdt=γfQ+ψC+κA-(ω+μ)R,
where
(27)λ=β(ηA+C+ζQ)N.
In addition to the dynamics of the original model (1), now the susceptible population is being vaccinated at per capita rate v (vaccine coverage rate). The vaccination rate v includes both the medical inoculation and the educational campaigns to prevent Hep C. We also have to consider the partial efficiency of the vaccine, due to which only σ(0≤σ≤1) fraction of the vaccinated susceptible individuals go to the recovered class (σ is the vaccine efficacy). The remaining (1-σ) fraction of the vaccinated susceptible individuals goes to the acute class as a result of their contact with the infected individuals. When σ=0, it means that the vaccine has no effect at all and when σ=1, the vaccine is perfectly effective.It is well understood that in order to eradicate an epidemic we have to vaccinate a large fraction of the susceptible population. Take the example of smallpox, its eradication was achieved after an intensive worldwide vaccination campaign and a very high vaccination rate [38]. However, in some infectious diseases such as measles and TB, the disease persists despite the extensive usage of vaccine, mainly because of low vaccine efficacy [39] and vaccination campaigns that could not reach everyone. Different vaccination policies have been implemented in different parts of the world. Practically, the cost of vaccine is a very important factor influencing the policy. Hence, we need to find a balance between the rate of vaccinating susceptible individuals and the cost of the vaccination. Now we design an optimal vaccination strategy to minimize an objective functional that takes into account both the cost and the number of infectious individuals. Now let the vaccination rate be given as function of time v(t) in model (26). The control set V is
(28)V={v(t):0≤v(t)≤b,0≤t≤T,v(t)isLebesguemeasurable}.
The goal is to minimize the cost function defined as
(29)J[v]=∫0T‍P1C+P2A+P3Q+12Wv2(t).
This performance specification involves the numbers of individuals of acute, chronic, and isolated, respectively, as well as the cost of applying vaccination control (v). This cost also includes the cost for organization and management of vaccine and so forth. Hence, the cost function should be nonlinear. In this paper, a quadratic function is implemented for measuring the control cost [21–24]. The coefficients P1, P2, P3, and W are balancing cost factors due to scales and importance of the four parts of the objective function. We need to find an optimal control v*(t) such that
(30)J[v*]=min⁡v∈VJ[v].
The existence of a solution to the optimal control problem can be obtained by verifying sufficient conditions. We refer to the conditions in Theorem III.4.1 and its corresponding corollary in [34]. The boundedness of solutions to the system (26) for the finite time interval is needed to establish these conditions. Pontryagin’s Maximum Principle [30] provides the necessary conditions to be satisfied by the optimal vaccination v(t). This principle reduces (26), (28), and (29) into a problem of minimizing pointwise a Hamiltonian, H, with respect to v:
(31)H=P1C+P2A+P3Q+12Wv2(t)+∑i=1i=5‍λiki,
where ki represents the right hand side of model’s (26) ith differential equation. Using Pontryagin’s Maximum principle [30] and the optimal control existence result from [34], we have the following.Theorem 7.
There exists a unique optimalv*(t) which minimizes J over V. Also, there exists adjoint system of λi’s such that
(32)dλ1dt=(λ(1-v*)+(1-σ)v*λ+σv*+μ)λ1-((1-v*)λ+(1-σ)v*λ)λ2,dλ2dt=(βη(1-v*+(1-σ)v*)SN)λ1+(k1-βη(1-v*+(1-σ)v*)SN)λ2-ξλ3λ1-κλ5-P2,dλ3dt=(β(1-v*+(1-σ)v*)SN)λ1-β(1-v*+(1-σ)v*)SNλ2+k2λ3-αλ4-ψλ5-P1,dλ4dt=(βζ(1-v*+(1-σ)v*)SN)λ1-(k4+βζ(1-v*+(1-σ)v*)SN)λ2+k3λ4-γfλ5-P3,(33)dλ5dt=ωλ1+k5λ5,k1=(ξ+κ+μ+δa),k2=(α+ψ1+μ+δc),k3=(γ+μ+δq),k4=γ(1-f),k5=(ω+μ).
The transversality condition gives(34)λi(T)=0.
The vaccination control is characterized as
(35)v*(t)=min⁡[b,max⁡(0,SW(λλ1(1-σ)+σλ1+λλ2ssssssssssssssssssssssssssssssss-λλ1-λλ2(1-σ))SW)].Proof.
Clearly the integrand ofJ is convex with respect to v(t). Also the solutions of (26) are bounded as N(t)≤Π/μ for all time. Also it is easily verifiable that the state system (26) has the Lipschitz property with respect to the state variables. With these properties and using Corollary 4.1 of [34], we have the existence of the optimal control.
Since we have the existence of the optimal vaccination control, using Pontryagin’s Maximum Principle, we obtain(36)dλ1dt=-∂H∂S,λ1(T)=0,⋮dλ5dt=-∂H∂R,λ5(T)=0
evaluated at the optimal control, which results in the stated adjoint system (32). The optimality condition is
(37)∂H∂v|v*=0.
Therefore, on the set {t:0<v*(t)<b}, we obtain
(38)v*=SW(λλ1(1-σ)+σλ1+λλ2-λλ1-λλ2(1-σ)).
Considering the bounds on v*, we have the characterizations of the optimal control as in (35). Clearly the state and the adjoint functions are bounded. Also it is easily verifiable that state system and adjoint system have Lipschitz structure with respect to the corresponding variables; we obtain the uniqueness of the optimal control for sufficiently small time T [30, 34]. The uniqueness of the optimal control pair follows from the uniqueness of the optimality system, which consists of (26) and (32), with characterizations (35). There is a restriction on the length of the time interval in order to guarantee the uniqueness of the optimality system. This restriction on the length of the time interval is due to the opposite time orientations of (26) and (32); the state problem has initial values and the adjoint problem has final values. For example, see [21]. This restriction is very common in control problems [22].The following optimality system, consisting of 10 equations, characterizes the optimal vaccination control as defined in above theorem. It consists of five equations of (26) with initial condition {S(0)=S0,A(0)=A0,C(0)=C0,Q(0)=Q0,R(0)=R0} and five equations of (32) with transversality condition {λ1(T)=λ2(T)=λ3(T)=λ4(T)=λ5(T)=0}.Next, we discuss the numerical solutions of the optimality system and the corresponding optimal control pairs and the parameters. The optimal vaccination strategy is obtained by solving the optimality system (26) and (32), consisting of 10 ODEs from the state and adjoint equations. An iterative method is used for solving the optimality system. We start to solve the state equations with a guess for the control v(t) over the simulated time using a forward fourth-order Runge-Kutta scheme. The adjoint functions have final time conditions. Because of this transversality condition on the adjoint functions (32), the adjoint equations are solved by a backward fourth-order Runge-Kutta scheme using the current iteration solution of the state equations. Then, the controls are updated by using a convex combination of the previous control and the value from the characterizations in (35). This process is repeated and iteration is stopped if the values of unknowns at the previous iteration are very close to the ones at the present iteration.Numerical solutions to the optimal system (26) and (32) are carried out using MATLAB and are presented here. The parameter values and the initial conditions are given in the appendix. The parameter values used have R0>1 when the model without control is considered. Thus, the disease is not expected to die out without intervention strategies.Figure9(a) represents the control strategy to be employed for the optimal results. This control strategy minimizes both the cost and the infected population (A+C+Q). It is well understood that in order to eradicate an epidemic, a large fraction of susceptible population needs to be vaccinated. Therefore, an upper bound of b=0.7 was chosen for vaccination control v(t). The optimal control v is at its upper bound during the first 60 days and then v is steadily decreasing to 0. In fact, at the beginning of simulated time, the optimal control is staying at its upper bound in order to vaccinate as many susceptible individuals as possible to prevent the individuals from getting infected. The steady decrease of the v is determined by the balance between the cost of the infected individuals and the cost of the controls. Figure 9(b) shows the total infected population for the optimal vaccination control, constant vaccination, and without control. It is easy to see that the optimal control is much more effective for reducing the number of infected individuals and decreasing the time-span of the epidemic. As normally expected, in the early phase of the epidemic breakouts, keeping the vaccination controls at their upper bound will directly lead to the decrease of the number of the infected people.Simulations show the optimal vaccination control and its effectiveness. The left simulation presents the vaccination strategy to be followed to prevent the epidemic and disease spread. The right simulation presents comparison of the total chronically infected individuals under optimal and constant control. Clearly optimal strategy prevents the epidemic and retains the infected population to a minimum.
(a)
Optimal vaccination controlv(t)
(b)
Total infected population with different control strategiesFigure10 shows the cost associated with the optimal and constant control strategy. It is clear that the cost of optimal strategy is much less than the cost of constant strategy. In fact the costs differ by order of magnitude of ten. Figure 11 represents the population sizes of the distinct infected states (acute, chronic, and isolated) and with optimal control strategy, population sizes go to minimum in short period of time.Figure 10
Simulation of accumulated cost of different control strategies.Population sizes of different infected compartments under vaccination strategy. This simulation presents comparison of the total chronically infected, isolated, and acute individuals under optimal and constant control. Clearly optimal strategy prevents the epidemic and retains the infected population to a minimum.
(a)
Simulation of acute (Hep C) population with different control strategies
(b)
Simulation of chronic (Hep C) population with different control strategies
(c)
Simulation of isolated population with different control strategies


### 3.2. Isolation Control

Several approaches to Hep C vaccine development have now been studied and include synthetic peptides, DNA [40], recombinant E1 and E2 proteins [41], and prime-boost strategies [42]. Success of these approaches has been limited for a number of reasons including: the delivery of a limited number of protective viral epitopes, the inclusion of incorrectly folded recombinant proteins, the limited humoral and cell mediated responses that are associated with DNA vaccines, and the use of adjuvants with relatively poor potency. It is also now apparent that vaccine inducing strong T cell responses alone may not be sufficient to prevent hepatitis C infection [40]. An effective preventive vaccine against Hep C will, therefore, need to induce strong neutralizing and cellular immune responses [7].Since an effective vaccine is not available against all the genotypes of Hep C, we have to look for alternate strategies to control the spread of Hep C. The quarantine and isolation of those individuals with disease symptoms constitute what is probably the first infection control measure since the beginning of recorded human history [14]. In our model of Hep C, the isolation compartment was introduced to investigate the effect on the infected population size and results were discussed in the last section. Now we can attempt to control the isolation rate of the chronically infected individuals in order to control the Hep C. This section will explore the effects of isolation control rate, of chronically infected individuals, on the total size of the infected population. The model including the required control is as follows:
(39)dSdt=Π+ωR-λS-μS,dAdt=λS+γ(1-f)Q-(ξ+κ+μ+δa)A,dCdt=ξA-(ψ+u1+μ+δc)C,dQdt=u1C-(γ+μ+δq)Q,dRdt=γfQ+κA+ψC-(ω+μ)R,
where
(40)λ=β(ηA+C+ζQ)N.In addition to the dynamics of the original model (1), now we have a control parameter for the isolation rate of the chronically infected individuals labeled as u1. The control u1 represents the fraction of the chronically infected individuals that are being isolated in order to decrease the rate of the spread of infection.Now we design an optimal control strategy to minimize an objective functional that takes into account both the cost and the number of infectious individuals. Now let the isolation rate beu1(t) for model (39). The control set U is
(41)U={u1(t):0≤u1(t)≤b,0≤t≤T,u1(t)areLebesguemeasurable}.
The objective functional is defined as
(42)J[u1]=∫0T‍P1C+P2Q+12Wu12(t),
where we want to minimize the infectious individuals, while also keeping the cost of the isolation facilities low. W1, P1, P2, and P3 are the weight parameters. We need to find an optimal control u1*(t) such that
(43)J[u1*]=min⁡u1∈UJ[u1].
Pontryagin’s Maximum Principle [30] provides the necessary conditions to be satisfied by the optimal vaccination v(t). This principle reduces (39), (41), and (42) into a problem of minimizing pointwise a Hamiltonian H, with respect to u1, defined as
(44)H=P1C+P2Q+12Wu12(t)+∑i=1i=4‍λiki,
where ki represents the right hand side of model’s (39) ith differential equation. Using Pontryagin’s Maximum Principle [30] and the optimal control existence result from [34], we have the following.Theorem 8.
There exists a unique optimal controlu1*(t) which minimizes J over U. Also, there exists adjoint system of λi’s such that
(45)dλ1dt=(λ+μ)λ1-λλ2,dλ2dt=(βηSN)λ1+(k1-βηSN)λ2-ξλ3-κλ5,dλ3dt=(βηSN)λ1-βSNλ2+k2λ3-u1*λ4-ψλ5-P1,dλ4dt=(βζSN)λ1-(k4+βζSN)λ2+k3λ4-γfλ5-P2,dλ5dt=ωλ1+k5λ5,
where k2=(u1+ψ+μ+δc) and all other k’s are given above.
The transversality condition is(46)λi(T)=0.
The optimal treatment control pair is characterised as
(47)u1*(t)=min⁡[b,max⁡(0,C(λ3-λ4W)].The proof is identical to the proof of Theorem7.In this case the system of 10 differential equations, as stated above, characterizes the optimal treatment control pair, consisting of five equations of system (39) and five equations of system (45).Now, we discuss the numerical solutions of the optimality system and the corresponding optimal control pair and the parameters. The optimal treatment strategy is obtained by solving the optimality system (39) and (45), consisting of 10 ODEs from the state and adjoint system. The method is discussed in the last section.Figure12(a) represents the optimal isolation strategy to be employed to minimize the cost and the infected population. Considering the practical constraints, an upper bound of b=0.7 was chosen for the optimal isolation control u1(t). The optimal control u1 is at its upper bound during the first 70 days and then u1 is steadily decreasing to 0. In fact, at the beginning of simulated time, the optimal control is staying at its upper bound in order to isolate as many chronically infected individuals as possible to prevent the infected population from increasing. The steady decrease of u1 is determined by the balance between the cost of the infected individuals and the cost of the controls. Figure 12(b) shows the total infected population for the optimal control, constant control, and without control. It is clear that with the use of an optimal control strategy the disease can be eradicated in a short period of time. Figure 13 shows the cost associated with the optimal and constant control strategy. It is clear that the cost of optimal strategy is much less than the cost of constant strategy and in fact differs by order of magnitude of tens. Figure 14 represents the population sizes of the distinct infected states (acute, chronics and isolated) and with optimal control strategy, population sizes go to minimum in short period of time.Simulations show the optimal isolation control and its effectiveness. The left simulation presents the isolation strategy to be followed to prevent the epidemic and disease spread. The right simulation presents comparison of the total chronically infected individuals under optimal and constant control. Clearly optimal strategy prevents the epidemic and retains the infected population to a minimum.
(a)
Optimal isolation controlv(t)
(b)
Total infected population with different control strategiesFigure 13
Accumulated cost of different control strategies: simulation presents comparison of the cost incurred to implement optimal and different constant control strategies to control hepatitis C.Population Sizes of different infected compartments under optimal Isolation strategy. Simulation presents comparison of the total chronically infected and isolated individuals under optimal and constant control. Clearly optimal strategy prevents the epidemic and retains the infected population to a minimum.
(a)
Simulation of chronic (Hep C) population with different control strategies
(b)
Simulation of isolated population with different control strategies


## 3.1. Vaccination Control

The use of vaccination is an important disease control and prevention measure. Optimal vaccination control strategy for an SIR model has been devised using dynamic programming technique [36]. Also the optimal control strategies have been investigated for TB control [23, 37]. We have analyzed the dynamics of the model without vaccination in the last section which had a locally stable equilibrium. In this section we will explore the effects of a vaccination campaign on our deterministic model of Hep C by including vaccination. After including vaccination the model is given as follows:
(26)dSdt=Π+ωR-λ(1-v)S-(1-σ)λvS-μS-σvS,dAdt=λ(1-v)S+(1-σ)λvS+γ(1-f)Q-(ξ+κ+μ+δa)A,dCdt=ξA-(α+ψ1+μ+δc)C,dQdt=αC-(γ+μ+δq)Q,dRdt=γfQ+ψC+κA-(ω+μ)R,
where
(27)λ=β(ηA+C+ζQ)N.
In addition to the dynamics of the original model (1), now the susceptible population is being vaccinated at per capita rate v (vaccine coverage rate). The vaccination rate v includes both the medical inoculation and the educational campaigns to prevent Hep C. We also have to consider the partial efficiency of the vaccine, due to which only σ(0≤σ≤1) fraction of the vaccinated susceptible individuals go to the recovered class (σ is the vaccine efficacy). The remaining (1-σ) fraction of the vaccinated susceptible individuals goes to the acute class as a result of their contact with the infected individuals. When σ=0, it means that the vaccine has no effect at all and when σ=1, the vaccine is perfectly effective.It is well understood that in order to eradicate an epidemic we have to vaccinate a large fraction of the susceptible population. Take the example of smallpox, its eradication was achieved after an intensive worldwide vaccination campaign and a very high vaccination rate [38]. However, in some infectious diseases such as measles and TB, the disease persists despite the extensive usage of vaccine, mainly because of low vaccine efficacy [39] and vaccination campaigns that could not reach everyone. Different vaccination policies have been implemented in different parts of the world. Practically, the cost of vaccine is a very important factor influencing the policy. Hence, we need to find a balance between the rate of vaccinating susceptible individuals and the cost of the vaccination. Now we design an optimal vaccination strategy to minimize an objective functional that takes into account both the cost and the number of infectious individuals. Now let the vaccination rate be given as function of time v(t) in model (26). The control set V is
(28)V={v(t):0≤v(t)≤b,0≤t≤T,v(t)isLebesguemeasurable}.
The goal is to minimize the cost function defined as
(29)J[v]=∫0T‍P1C+P2A+P3Q+12Wv2(t).
This performance specification involves the numbers of individuals of acute, chronic, and isolated, respectively, as well as the cost of applying vaccination control (v). This cost also includes the cost for organization and management of vaccine and so forth. Hence, the cost function should be nonlinear. In this paper, a quadratic function is implemented for measuring the control cost [21–24]. The coefficients P1, P2, P3, and W are balancing cost factors due to scales and importance of the four parts of the objective function. We need to find an optimal control v*(t) such that
(30)J[v*]=min⁡v∈VJ[v].
The existence of a solution to the optimal control problem can be obtained by verifying sufficient conditions. We refer to the conditions in Theorem III.4.1 and its corresponding corollary in [34]. The boundedness of solutions to the system (26) for the finite time interval is needed to establish these conditions. Pontryagin’s Maximum Principle [30] provides the necessary conditions to be satisfied by the optimal vaccination v(t). This principle reduces (26), (28), and (29) into a problem of minimizing pointwise a Hamiltonian, H, with respect to v:
(31)H=P1C+P2A+P3Q+12Wv2(t)+∑i=1i=5‍λiki,
where ki represents the right hand side of model’s (26) ith differential equation. Using Pontryagin’s Maximum principle [30] and the optimal control existence result from [34], we have the following.Theorem 7.
There exists a unique optimalv*(t) which minimizes J over V. Also, there exists adjoint system of λi’s such that
(32)dλ1dt=(λ(1-v*)+(1-σ)v*λ+σv*+μ)λ1-((1-v*)λ+(1-σ)v*λ)λ2,dλ2dt=(βη(1-v*+(1-σ)v*)SN)λ1+(k1-βη(1-v*+(1-σ)v*)SN)λ2-ξλ3λ1-κλ5-P2,dλ3dt=(β(1-v*+(1-σ)v*)SN)λ1-β(1-v*+(1-σ)v*)SNλ2+k2λ3-αλ4-ψλ5-P1,dλ4dt=(βζ(1-v*+(1-σ)v*)SN)λ1-(k4+βζ(1-v*+(1-σ)v*)SN)λ2+k3λ4-γfλ5-P3,(33)dλ5dt=ωλ1+k5λ5,k1=(ξ+κ+μ+δa),k2=(α+ψ1+μ+δc),k3=(γ+μ+δq),k4=γ(1-f),k5=(ω+μ).
The transversality condition gives(34)λi(T)=0.
The vaccination control is characterized as
(35)v*(t)=min⁡[b,max⁡(0,SW(λλ1(1-σ)+σλ1+λλ2ssssssssssssssssssssssssssssssss-λλ1-λλ2(1-σ))SW)].Proof.
Clearly the integrand ofJ is convex with respect to v(t). Also the solutions of (26) are bounded as N(t)≤Π/μ for all time. Also it is easily verifiable that the state system (26) has the Lipschitz property with respect to the state variables. With these properties and using Corollary 4.1 of [34], we have the existence of the optimal control.
Since we have the existence of the optimal vaccination control, using Pontryagin’s Maximum Principle, we obtain(36)dλ1dt=-∂H∂S,λ1(T)=0,⋮dλ5dt=-∂H∂R,λ5(T)=0
evaluated at the optimal control, which results in the stated adjoint system (32). The optimality condition is
(37)∂H∂v|v*=0.
Therefore, on the set {t:0<v*(t)<b}, we obtain
(38)v*=SW(λλ1(1-σ)+σλ1+λλ2-λλ1-λλ2(1-σ)).
Considering the bounds on v*, we have the characterizations of the optimal control as in (35). Clearly the state and the adjoint functions are bounded. Also it is easily verifiable that state system and adjoint system have Lipschitz structure with respect to the corresponding variables; we obtain the uniqueness of the optimal control for sufficiently small time T [30, 34]. The uniqueness of the optimal control pair follows from the uniqueness of the optimality system, which consists of (26) and (32), with characterizations (35). There is a restriction on the length of the time interval in order to guarantee the uniqueness of the optimality system. This restriction on the length of the time interval is due to the opposite time orientations of (26) and (32); the state problem has initial values and the adjoint problem has final values. For example, see [21]. This restriction is very common in control problems [22].The following optimality system, consisting of 10 equations, characterizes the optimal vaccination control as defined in above theorem. It consists of five equations of (26) with initial condition {S(0)=S0,A(0)=A0,C(0)=C0,Q(0)=Q0,R(0)=R0} and five equations of (32) with transversality condition {λ1(T)=λ2(T)=λ3(T)=λ4(T)=λ5(T)=0}.Next, we discuss the numerical solutions of the optimality system and the corresponding optimal control pairs and the parameters. The optimal vaccination strategy is obtained by solving the optimality system (26) and (32), consisting of 10 ODEs from the state and adjoint equations. An iterative method is used for solving the optimality system. We start to solve the state equations with a guess for the control v(t) over the simulated time using a forward fourth-order Runge-Kutta scheme. The adjoint functions have final time conditions. Because of this transversality condition on the adjoint functions (32), the adjoint equations are solved by a backward fourth-order Runge-Kutta scheme using the current iteration solution of the state equations. Then, the controls are updated by using a convex combination of the previous control and the value from the characterizations in (35). This process is repeated and iteration is stopped if the values of unknowns at the previous iteration are very close to the ones at the present iteration.Numerical solutions to the optimal system (26) and (32) are carried out using MATLAB and are presented here. The parameter values and the initial conditions are given in the appendix. The parameter values used have R0>1 when the model without control is considered. Thus, the disease is not expected to die out without intervention strategies.Figure9(a) represents the control strategy to be employed for the optimal results. This control strategy minimizes both the cost and the infected population (A+C+Q). It is well understood that in order to eradicate an epidemic, a large fraction of susceptible population needs to be vaccinated. Therefore, an upper bound of b=0.7 was chosen for vaccination control v(t). The optimal control v is at its upper bound during the first 60 days and then v is steadily decreasing to 0. In fact, at the beginning of simulated time, the optimal control is staying at its upper bound in order to vaccinate as many susceptible individuals as possible to prevent the individuals from getting infected. The steady decrease of the v is determined by the balance between the cost of the infected individuals and the cost of the controls. Figure 9(b) shows the total infected population for the optimal vaccination control, constant vaccination, and without control. It is easy to see that the optimal control is much more effective for reducing the number of infected individuals and decreasing the time-span of the epidemic. As normally expected, in the early phase of the epidemic breakouts, keeping the vaccination controls at their upper bound will directly lead to the decrease of the number of the infected people.Simulations show the optimal vaccination control and its effectiveness. The left simulation presents the vaccination strategy to be followed to prevent the epidemic and disease spread. The right simulation presents comparison of the total chronically infected individuals under optimal and constant control. Clearly optimal strategy prevents the epidemic and retains the infected population to a minimum.
(a)
Optimal vaccination controlv(t)
(b)
Total infected population with different control strategiesFigure10 shows the cost associated with the optimal and constant control strategy. It is clear that the cost of optimal strategy is much less than the cost of constant strategy. In fact the costs differ by order of magnitude of ten. Figure 11 represents the population sizes of the distinct infected states (acute, chronic, and isolated) and with optimal control strategy, population sizes go to minimum in short period of time.Figure 10
Simulation of accumulated cost of different control strategies.Population sizes of different infected compartments under vaccination strategy. This simulation presents comparison of the total chronically infected, isolated, and acute individuals under optimal and constant control. Clearly optimal strategy prevents the epidemic and retains the infected population to a minimum.
(a)
Simulation of acute (Hep C) population with different control strategies
(b)
Simulation of chronic (Hep C) population with different control strategies
(c)
Simulation of isolated population with different control strategies


## 3.2. Isolation Control

Several approaches to Hep C vaccine development have now been studied and include synthetic peptides, DNA [40], recombinant E1 and E2 proteins [41], and prime-boost strategies [42]. Success of these approaches has been limited for a number of reasons including: the delivery of a limited number of protective viral epitopes, the inclusion of incorrectly folded recombinant proteins, the limited humoral and cell mediated responses that are associated with DNA vaccines, and the use of adjuvants with relatively poor potency. It is also now apparent that vaccine inducing strong T cell responses alone may not be sufficient to prevent hepatitis C infection [40]. An effective preventive vaccine against Hep C will, therefore, need to induce strong neutralizing and cellular immune responses [7].Since an effective vaccine is not available against all the genotypes of Hep C, we have to look for alternate strategies to control the spread of Hep C. The quarantine and isolation of those individuals with disease symptoms constitute what is probably the first infection control measure since the beginning of recorded human history [14]. In our model of Hep C, the isolation compartment was introduced to investigate the effect on the infected population size and results were discussed in the last section. Now we can attempt to control the isolation rate of the chronically infected individuals in order to control the Hep C. This section will explore the effects of isolation control rate, of chronically infected individuals, on the total size of the infected population. The model including the required control is as follows:
(39)dSdt=Π+ωR-λS-μS,dAdt=λS+γ(1-f)Q-(ξ+κ+μ+δa)A,dCdt=ξA-(ψ+u1+μ+δc)C,dQdt=u1C-(γ+μ+δq)Q,dRdt=γfQ+κA+ψC-(ω+μ)R,
where
(40)λ=β(ηA+C+ζQ)N.In addition to the dynamics of the original model (1), now we have a control parameter for the isolation rate of the chronically infected individuals labeled as u1. The control u1 represents the fraction of the chronically infected individuals that are being isolated in order to decrease the rate of the spread of infection.Now we design an optimal control strategy to minimize an objective functional that takes into account both the cost and the number of infectious individuals. Now let the isolation rate beu1(t) for model (39). The control set U is
(41)U={u1(t):0≤u1(t)≤b,0≤t≤T,u1(t)areLebesguemeasurable}.
The objective functional is defined as
(42)J[u1]=∫0T‍P1C+P2Q+12Wu12(t),
where we want to minimize the infectious individuals, while also keeping the cost of the isolation facilities low. W1, P1, P2, and P3 are the weight parameters. We need to find an optimal control u1*(t) such that
(43)J[u1*]=min⁡u1∈UJ[u1].
Pontryagin’s Maximum Principle [30] provides the necessary conditions to be satisfied by the optimal vaccination v(t). This principle reduces (39), (41), and (42) into a problem of minimizing pointwise a Hamiltonian H, with respect to u1, defined as
(44)H=P1C+P2Q+12Wu12(t)+∑i=1i=4‍λiki,
where ki represents the right hand side of model’s (39) ith differential equation. Using Pontryagin’s Maximum Principle [30] and the optimal control existence result from [34], we have the following.Theorem 8.
There exists a unique optimal controlu1*(t) which minimizes J over U. Also, there exists adjoint system of λi’s such that
(45)dλ1dt=(λ+μ)λ1-λλ2,dλ2dt=(βηSN)λ1+(k1-βηSN)λ2-ξλ3-κλ5,dλ3dt=(βηSN)λ1-βSNλ2+k2λ3-u1*λ4-ψλ5-P1,dλ4dt=(βζSN)λ1-(k4+βζSN)λ2+k3λ4-γfλ5-P2,dλ5dt=ωλ1+k5λ5,
where k2=(u1+ψ+μ+δc) and all other k’s are given above.
The transversality condition is(46)λi(T)=0.
The optimal treatment control pair is characterised as
(47)u1*(t)=min⁡[b,max⁡(0,C(λ3-λ4W)].The proof is identical to the proof of Theorem7.In this case the system of 10 differential equations, as stated above, characterizes the optimal treatment control pair, consisting of five equations of system (39) and five equations of system (45).Now, we discuss the numerical solutions of the optimality system and the corresponding optimal control pair and the parameters. The optimal treatment strategy is obtained by solving the optimality system (39) and (45), consisting of 10 ODEs from the state and adjoint system. The method is discussed in the last section.Figure12(a) represents the optimal isolation strategy to be employed to minimize the cost and the infected population. Considering the practical constraints, an upper bound of b=0.7 was chosen for the optimal isolation control u1(t). The optimal control u1 is at its upper bound during the first 70 days and then u1 is steadily decreasing to 0. In fact, at the beginning of simulated time, the optimal control is staying at its upper bound in order to isolate as many chronically infected individuals as possible to prevent the infected population from increasing. The steady decrease of u1 is determined by the balance between the cost of the infected individuals and the cost of the controls. Figure 12(b) shows the total infected population for the optimal control, constant control, and without control. It is clear that with the use of an optimal control strategy the disease can be eradicated in a short period of time. Figure 13 shows the cost associated with the optimal and constant control strategy. It is clear that the cost of optimal strategy is much less than the cost of constant strategy and in fact differs by order of magnitude of tens. Figure 14 represents the population sizes of the distinct infected states (acute, chronics and isolated) and with optimal control strategy, population sizes go to minimum in short period of time.Simulations show the optimal isolation control and its effectiveness. The left simulation presents the isolation strategy to be followed to prevent the epidemic and disease spread. The right simulation presents comparison of the total chronically infected individuals under optimal and constant control. Clearly optimal strategy prevents the epidemic and retains the infected population to a minimum.
(a)
Optimal isolation controlv(t)
(b)
Total infected population with different control strategiesFigure 13
Accumulated cost of different control strategies: simulation presents comparison of the cost incurred to implement optimal and different constant control strategies to control hepatitis C.Population Sizes of different infected compartments under optimal Isolation strategy. Simulation presents comparison of the total chronically infected and isolated individuals under optimal and constant control. Clearly optimal strategy prevents the epidemic and retains the infected population to a minimum.
(a)
Simulation of chronic (Hep C) population with different control strategies
(b)
Simulation of isolated population with different control strategies


## 4. Conclusions

This paper has discussed the transmission dynamics of the Hep C and eventually formulated an optimal control strategy to prevent disease spread. At first, a deterministic epidemic model for the spread of Hep C, which incorporates the possibility of an isolation state, is formulated. Global analysis of the equilibrium solution is performed. The existence of a disease-free equilibrium and an endemic equilibrium is shown. It is further demonstrated that the disease-free equilibrium is globally asymptotically stable for values of the basic reproduction number (R0<1). The disease uniformly persists when R0>1. In addition, it is shown that the endemic equilibrium is locally asymptotically stable for values of the basic reproduction number (R0>1) assuming constant total population. It is also shown that numerical simulations of the deterministic model agree with the theoretical results. Finally, the sensitivity and uncertainty analyses were performed along with numerical simulations to study the influence of vital parameters on the disease spread.The model is most sensitive to the control variablesα (proportion of chronic population being isolated) and β (effective contact rate); given the nature of the disease, controlling the isolation parameter (i.e., devising an effective isolation strategy) seems to be the most workable solution.Since we have a stable endemic equilibrium, we formulated (time-dependent) optimal control strategies to fight Hep C. We considered the two most effective strategies known to prevent disease spread, vaccination and isolation separately. Absence of an effective vaccine against all known genotypes of Hep C made us consider vaccination and isolation strategy separately. We assumed that the vaccination not only includes medical vaccination, but also the education and awareness schemes about Hep C.The numerical results show that the proportion that is isolated optimally with respect to time has a higher favorable impact (as compared to implementing a high but constant isolation rate) on keeping the cost of disease control low. However, it should be pointed out that the ideal time varying optimal strategy might not be applied easily. Still, it does provide a basis on which can be designed practical quasioptimal control strategies.Next we consider the possibility of an effective vaccine becoming available for Hep C. We formulate a possible optimal time varying vaccination strategy; our analysis shows that it is not required to vaccinate a constant proportion of the susceptible population over time, but rather the optimal policy (being cost effective) is to start with vaccinating a constant proportion of the susceptibles and then progressively reducing the fraction of the susceptibles vaccinated.Although the early diagnosis and treatment of Hep C virus infection are desirable to prevent spread of infection and to reduce the risk of progression of disease, the majority of acute individuals are asymptomatic and most infected persons are unaware of their exposure to the virus. More public awareness of Hep C may also increase the recovery from acute stage in the model. The goal of therapy is to prevent complications and death from the infection. Therapy along with an effective isolation strategy can greatly reduce the prevalence of Hep C.



---

*Source: 654050-2014-02-13.xml*