A theoretical model for the transmission dynamics of HIV/HSV-2 co-infection in the presence of poor HSV-2 treatment adherence

Herpes simplex virus type 2 (HSV-2) and human immunodeficiency virus (HIV) causes sexually transmitted infections (STIs) that are detrimental to human health [1]. HIV is a retrovirus that infects cells of the human system that leads to acquired immuno deficiency syndrome (AIDS). AIDS is a condition in humans in which progressive failure of the immune system allows life-threatening opportunistic infections and cancers to thrive. Since HIV was first identified in the United States in 1983, over 60 million people have been infected, and the WHO estimates that deaths due to AIDS exceeds 25 million people [1]. In 2015, an estimated 36.7 million people were living with HIV (including 1.8 million children), a global HIV prevalence of 0.8%, with the majority of this number living in low- and middle-income countries [2]. HSV-2 is a double strand DNA virus, with humans being the only natural hosts. There are 1.5 million new HSV-2 infections among 15–19 year olds in sub-Saharan Africa every year and an estimated 1.6 million infections annually in the USA. HSV-2 has been recognised as the most common cause of genital ulcer disease. HSV-2 is a significant factor for increased risk of acquisition and transmission of HIV [4]. In the United States of America (USA), annual health costs for STIs has reached US$17 Billion, with HSV-2 chewing up $541 million, making it the third most costly sexually transmitted infection (STI) after HIV and human papillomavirus (HPV) [5].

New strategies for the prevention of infection with HIV are needed, especially in sub-Saharan Africa [4]. The use of condoms remains low in the region, despite intensive educational campaigns, and sexually transmitted infections are highly prevalent, especially infection with HSV-2, with a prevalence of up to 80% [3, 4]. In HIV-uninfected individuals, HSV-2 causes genital ulceration and mucosal disruption, providing a direct portal for HIV entry [9]. A large number of CD4⁺ T-lymphocytes, HIV target cells, have been detected in herpetic lesions, and the presence of these cells could also increase susceptibility to HIV during sexual activity [31]. Observational studies suggest that HSV-2 infection doubles or triples the risk of acquiring HIV and may contribute to more than 50% of HIV infections in sub-Saharan Africa [8]. The epidemiological and biological association between HSV-2 and HIV has been the subject of numerous studies over the past two decades, with strong evidence supporting the hypothesis that HSV-2 increases the risk of acquisition among HSV-2-infected HIV-negative individuals [10, 13]. This biological synergy between these two viruses has led researchers to consider HSV-2-suppressive treatment as a biomedical prevention strategy to reduce the risk of HIV transmission [14]. Some researchers have seen poor HSV-2 treatment adherence as a possible factor that may explain the failure of some HSV-2 treatments to reduce systemic and genital HIV-1 levels [11]. Plummer et al. [12], conducted a qualitative study on HSV-2 treatment adherence in Tanzania and observed an 8% (7/86) adherence for HSV-2 patients who were on Acyclovir. Despite numerous research efforts that have been devoted to the study of HIV and HSV-2 co-infection, the aspect of poor HSV-2 treatment adherence and its impact on infection and spread of HIV has not yet been investigated. The global consensus currently is that, prevention of new infections is the key to reversing the HIV/AIDS epidemic, hence HSV-2 treatment with good treatment adherence can be another possible strategic method.

With HSV-2 being a pathway for the transmission of HIV, this has motivated many authors to look at HIV/HSV-2 co-infections. A number of mathematical models have looked into mathematical modelling of HIV and HSV-2 co-infections from a number of different issues [1, 17, 18, 19, 20, 25]. In contrast to other STIs such as chlamydia, gonorrhea, syphilis and trichomoniasis, infection with HSV-2 or HIV is lifelong, and, once established, there is currently no treatment to eliminate. To the best of our knowledge, no mathematical model has looked at the impact of poor HSV-2 treatment adherence on HIV/AIDS prevalence. The aim of this paper is to study an optimal control model of HIV and HSV-2 co-infection, by considering two control measures (namely monitoring and counselling individuals infected with HSV-2 only and monitoring and counselling individuals dually infected with HIV and HSV-2) that try to reduce the number of individuals who quit HSV-2 anti-viral treatment before completion. To achieve our aim, we shall apply the optimal control theory which has proven to be a successful tool in understanding ways to curtail the spread of infectious diseases by devising the optimal diseases intervention strategies. The method consists of minimising the cost of infection or the cost of implementing the controls, or both.

The paper is organised as follows. In the next section, the model is formulated. Section 3 is dedicated to the stability analysis of the model. The definition of an optimal control, model properties and proof of existence of optimal control, analysis and the numerical simulations are given in Section 4. The paper is concluded in Section 5.

The total sexually active population at time t, denoted by N, is sub-divided into mutually exclusive compartments, namely susceptibles (S), individuals infected with acute HSV-2 only and under antiviral treatment (I₁), individuals infected with acute HSV-2 only and not under anti-viral treatment (I₂), individuals infected with latent HSV-2 only after undergoing successful anti-viral treatment (Q₁), individuals infected with latent HSV-2 only after undergoing natural healing (Q₂), individuals infected with HIV only (H), individuals dually infected with HIV and acute HSV-2, also under HSV-2 anti-viral treatment (H_I1), individuals dually infected with HIV and acute HSV-2, also not under HSV-2 anti-viral treatment (H_I2), individuals dually infected with HIV and latent HSV-2 after undergoing successful HSV-2 anti-viral treatment (H_Q1), individuals dually infected with HIV and latent HSV-2 after undergoing natural healing (H_Q2) and the AIDS class (A), so that

N¯=S¯+I¯1+I¯2+Q¯1+Q¯2+H¯+H¯I1+H¯I2+H¯Q1+H¯Q2+A¯.$$\begin{array}{} \displaystyle \overline{N} = \overline{S} + \overline{I}_{1} + \overline{I}_{2} + \overline{Q}_{1} + \overline{Q}_{2} + \overline{H} + \overline{H}_{I1} + \overline{H}_{I2} + \overline{H}_{Q1} + \overline{H}_{Q2} + \overline{A}. \end{array} $$(1)

The recruitment of susceptibles is proportional to the population and is given by μN as applied in [1]. Both singly and dually infected individuals transmit either HSV-2 or HIV, not both at the same time. Susceptible individuals acquire infection following effective contact with HIV infected individuals at a rate λ_H, and acquire HSV-2 infection following effective contact with HSV-2 infectives at a rate λ_I. Natural mortality rate occurs in all classes at a rate μ and the AIDS class suffering an additional AIDS related death rate v. The force of infection associated with HSV-2 infection, denoted by λ_I, is given by

λ¯I=βI[(1−α1)(I¯2+(1−η1)I¯1)+(H¯I2+(1−η2)H¯I1)]N¯,$$\begin{array}{} \displaystyle \overline{\lambda}_{I} = \dfrac{\beta_{I} [(1-\alpha_{1})(\overline{I}_{2} + (1-\eta_{1})\overline{I}_{1}) + (\overline{H}_{I2} + (1-\eta_{2})\overline{H}_{I1})]}{\overline{N}}, \end{array} $$(2)

where, β_I is the effective contact rate for HSV-2 transmission, the modification parameter α₁ ∈ (0,1) accounts for the relative infectiousness of individuals dually infected with HIV and HSV-2, as compared to those infected with HIV and HSV-2 only. The modification parameters 0 < η₁,η₂ < 1 account for the assumed reduced likelihood for individuals in classes I₁ and H_I1 to pass on infection compared to classes I₂ and H_I2 respectively. This is due to the fact that individuals in treatment have reduced viral load compared to those who have failed to adhere to anti-viral treatment guidelines. Similarly, susceptibles acquire HIV infection following effective contact with HIV infected individuals at a rate λ_H, given by

λ¯H=βH[(1−α2)(H¯+H¯Q1+H¯Q2)+(H¯I2+(1−τ)H¯I1)]N¯.$$\begin{array}{} \displaystyle \overline{\lambda}_{H} = \dfrac{\beta_{H} [(1-\alpha_{2})(\overline{H} + \overline{H}_{Q1} + \overline{H}_{Q2}) + (\overline{H}_{I2} + (1-\tau)\overline{H}_{I1})]}{\overline{N}}. \end{array} $$(3)

In (3), β_H is the effective contact rate for HIV infection. α₂ ∈ (0,1) is a modification parameter accounting for increased likelihood of infectiousness by individuals in classes H_I1 and H_I2 relative to those in classes H, H_Q1 and H_Q2. τ ∈ (0,1) captures the assumption that individuals who are dually infected with HIV and HSV-2 and under treatment have reduced likelihood of passing on the infection as compared to the individuals who are dually infected with HIV and HSV-2 and also not under HSV-2 treatment. Upon being infected with HSV-2, the individuals infected with HSV-2 only will become latent at constant rate k₁, where as those dually infected with HIV and HSV-2 become latent at a rate k₂. Following an appropriate stimulus in individuals with latent HSV-2 and those dually infected with HIV and latent HSV-2, re-activation may occur [15]. The anti-viral treatment for individuals with acute HSV-2 only and individuals dually infected with HIV and HSV-2 is denoted by ψ. Since the anti-viral medication will also suppress re-activation of latent HSV-2, we assume that the re-activation rate of people with latent HSV-2 only and those dually infected with HIV and latent HSV-2 is at rate γ⁽¹⁾(ψ) and γ⁽²⁾(ψ), respectively. Both γ⁽¹⁾(ψ) and γ⁽²⁾(ψ) are decreasing functions of ψ [1]. Thus

γ(1)=γ(1)(0)w1w1+ψandγ(2)=γ(2)(0)w2w2+ψ,$$\begin{array}{} \displaystyle \gamma^{(1)} = \gamma^{(1)} (0) \dfrac{w_{1}}{w_{1} + \psi}\,\,\,\mbox{and}\,\,\,\gamma^{(2)} = \gamma^{(2)} (0) \dfrac{w_{2}}{w_{2} + \psi}, \end{array} $$(4)

in which the factor wiwi+ψ$\begin{array}{} \dfrac{w_{i}}{w_{i} + \psi} \end{array} $ for i = 1,2 represents reduced re-activation treatment. In the case where there is no treatment for the individuals infected with HSV-2 only, we have that ψ = 0; thus γ⁽¹⁾ (ψ₀) ≡ γ⁽¹⁾ (0) ≡ γ0(1).$\begin{array}{} \gamma^{(1)}_{0}. \end{array} $ Similarly, for the case where there is no treatment for the individuals dually infected with HIV and HSV-2, we also have that ψ = 0; thus γ⁽²⁾ (ψ₀) ≡ γ⁽²⁾ (0) ≡ γ0(2).$\begin{array}{} \gamma^{(2)}_{0}. \end{array} $ A proportion (1−θ₁) of individuals infected with HSV-2 only fail to adhere to treatment at rate δ₁ and then move to class I₂. Furthermore, a proportion (1−θ₂) of individuals dually infected with HIV and HSV-2 fail to adhere to HSV-2 anti-viral treatment at rate δ₂, and they move to class H_I2. σ ≥ 1 denotes the enhanced susceptibility to HIV infection for individuals infected with acute HSV-2. φ ≥ 1 denotes the enhanced susceptibility to HSV-2 infection for individuals infected with HIV only. It is worth noting that γ⁽²⁾ and γ⁽¹⁾ do not mean squared or to the power 1 respectively, but will be used as our notation through out the manuscript. The model flow diagram is depicted in Figure 1.

Fig. 1

From the above descriptions and assumptions on the dynamics of the epidemic above, the following are the model equations.

S′=μN¯−(λ¯H+λ¯I)S¯−μS¯,I1′=λ¯IS¯+γ(1)Q¯1−σλ¯HI¯1−(δ1(1−θ1)+κ1+ψ+μ)I¯1,I2′=δ1(1−θ1)I¯1−σλ¯HI¯2+γ0(1)Q¯2−(μ+κ1)I¯2,Q1′=(κ1+ψ)I¯1−λ¯HQ¯1−(γ(1)+μ)Q¯1,Q2′=κ1I¯2−λ¯HQ¯2−(μ+γ0(1))Q¯2,H′=λ¯HS¯−φλ¯IH¯−(μ+ϕ)H¯,HI1′=φλ¯IH¯−(δ(1−θ2)+ϕ+κ2+ψ+μ)H¯I1+γ(2)H¯Q1+σλHI1,HI2′=δ2(1−θ2)H¯I1+σλ¯HI¯2+γ0(2)H¯Q2−(μ+κ2+ϕ)H¯I2,HQ1′=(κ2+ψ)H¯I1−(μ+ψ+γ(2))H¯Q1+λ¯HQ¯1,HQ2′=λ¯HQ¯2+κ2H¯I2−(μ+ϕ+γ0(2))H¯Q2,A′=ϕ(H¯+H¯I1+H¯I2+H¯Q1+H¯Q2)−(μ+v)A¯.$$\begin{array}{} \displaystyle \left\{\begin{array}{llll} S' = \mu \overline{N} - (\overline{\lambda}_{H} + \overline{\lambda}_{I})\overline{S} - \mu \overline{S},\\\\ I_{1}' = \overline{\lambda}_{I} \overline{S} + \gamma^{(1)} \overline{Q}_{1} - \sigma \overline{\lambda}_{H} \overline{I}_{1} - (\delta_{1} (1-\theta_{1}) + \kappa_{1} + \psi + \mu)\overline{I}_{1}, \\\\ I_{2}' = \delta_{1} (1-\theta_{1}) \overline{I}_{1} - \sigma \overline{\lambda}_{H} \overline{I}_{2} + \gamma^{(1)}_{0} \overline{Q}_{2} - (\mu + \kappa_{1})\overline{I}_{2}, \\\\ Q_{1}' = (\kappa_{1} + \psi)\overline{I}_{1} - \overline{\lambda}_{H} \overline{Q}_{1} - (\gamma^{(1)} + \mu) \overline{Q}_{1}, \\\\ Q_{2}' = \kappa_{1} \overline{I}_{2} - \overline{\lambda}_{H} \overline{Q}_{2} - (\mu + \gamma^{(1)}_{0}) \overline{Q}_{2}, \\\\ H' = \overline{\lambda}_{H} \overline{S} - \varphi \overline{\lambda}_{I} \overline{H} - (\mu + \phi)\overline{H}, \\\\ H_{I1}' = \varphi \overline{\lambda}_{I} \overline{H} - (\delta (1-\theta_{2}) + \phi + \kappa_{2} + \psi + \mu)\overline{H}_{I1} + \gamma^{(2)} \overline{H}_{Q1} + \sigma \lambda_{H} I_{1}, \\\\ H_{I2}' = \delta_{2} (1-\theta_{2}) \overline{H}_{I1} + \sigma \overline{\lambda}_{H} \overline{I}_{2} + \gamma^{(2)}_{0} \overline{H}_{Q2} - (\mu + \kappa_{2} + \phi)\overline{H}_{I2}, \\\\ H_{Q1}' = (\kappa_{2} + \psi)\overline{H}_{I1} - (\mu + \psi + \gamma^{(2)})\overline{H}_{Q1} + \overline{\lambda}_{H} \overline{Q}_{1}, \\\\ H_{Q2}' = \overline{\lambda}_{H} \overline{Q}_{2} + \kappa_{2} \overline{H}_{I2} - (\mu + \phi + \gamma_{0}^{(2)})\overline{H}_{Q2}, \\\\ A' = \phi (\overline{H} + \overline{H}_{I1} + \overline{H}_{I2} + \overline{H}_{Q1} + \overline{H}_{Q2}) - (\mu + v)\overline{A}. \end{array} \right. \end{array} $$(5)

It is helpful to rescale model system (5) so that we have dimensionless variables. We let

S=S¯N¯,I1=I¯1N¯,I2=I¯2N¯,Q1=Q¯1N¯,Q2=Q¯2N¯,H=H¯N¯,HI1=H¯I1N¯,HI2=H¯I2N¯HQ1=H¯Q1N¯,HQ2=H¯Q2N¯,A=A¯N¯.$$\begin{array}{} \displaystyle \quad\, S = \dfrac{\overline{S}}{\overline{N}},\,I_{1} = \dfrac{\overline{I}_{1}}{\overline{N}},\,I_{2} = \dfrac{\overline{I}_{2}}{\overline{N}},\,Q_{1} = \dfrac{\overline{Q}_{1}}{\overline{N}},\,Q_{2} = \dfrac{\overline{Q}_{2}}{\overline{N}},\,H = \dfrac{\overline{H}}{\overline{N}},\,H_{I1} = \dfrac{\overline{H}_{I1}}{\overline{N}},\,H_{I2} = \dfrac{\overline{H}_{I2}}{\overline{N}} \\ \\\\\displaystyle H_{Q1} = \dfrac{\overline{H}_{Q1}}{\overline{N}},\,H_{Q2} = \dfrac{\overline{H}_{Q2}}{\overline{N}},\,A = \dfrac{\overline{A}}{\overline{N}}. \end{array} $$

We thus have the following re-scaled system

S′=μ−(λH+λI)S−μS,I1′=λIS+γ(1)Q1−σλHI1−(δ1(1−θ1)+κ1+ψ+μ)I1,I2′=δ1(1−θ1)I1−σλHI2+γ0(1)Q2−(μ+κ1)I2,Q1′=(κ1+ψ)I1−λHQ1−(γ(1)+μ)Q1,Q2′=κ1I2−λHQ2−(μ+γ0(1))Q2,H′=λHS−φλIH−(μ+ϕ)H,HI1′=φλIH−(δ2(1−θ2)+ϕ+κ2+ψ+μ)HI1+γ(2)HQ1+σλHI1,HI2′=δ2(1−θ2)HI1+σλHI2+γ0(2)HQ2−(μ+κ2+ϕ)HI2,HQ1′=(κ2+ψ)HI1−(μ+ϕ+γ(2))HQ1+λHQ1,HQ2′=λHQ2+κ2HI2−(μ+ϕ+γ0(2))HQ2,A′=ϕ(H+HI1+HI2+HQ1+HQ2)−(μ+v)A,$$\begin{array}{} \displaystyle \left\{\begin{array}{llll} S' = \mu - ({\lambda}_{H} + {\lambda}_{I}){S} - \mu {S},\\\\ I_{1}' = {\lambda}_{I} {S} + \gamma^{(1)} {Q}_{1} - \sigma {\lambda}_{H} {I}_{1} - (\delta_{1} (1-\theta_{1}) + \kappa_{1} + \psi + \mu){I}_{1}, \\\\ I_{2}' = \delta_{1} (1-\theta_{1}) {I}_{1} - \sigma {\lambda}_{H} {I}_{2} + \gamma^{(1)}_{0} {Q}_{2} - (\mu + \kappa_{1}){I}_{2}, \\\\ Q_{1}' = (\kappa_{1} + \psi){I}_{1} - {\lambda}_{H} {Q}_{1} - (\gamma^{(1)} + \mu) {Q}_{1}, \\\\ Q_{2}' = \kappa_{1} {I}_{2} - {\lambda}_{H} Q_{2} - (\mu + \gamma^{(1)}_{0}) {Q}_{2}, \\\\ H' = {\lambda}_{H} S - \varphi {\lambda}_{I} H - (\mu + \phi){H}, \\\\ H_{I1}' = \varphi {\lambda}_{I} H - (\delta_{2} (1-\theta_{2}) + \phi + \kappa_{2} + \psi + \mu){H}_{I1} + \gamma^{(2)} {H}_{Q1} + \sigma \lambda_{H} I_{1}, \\\\ H_{I2}' = \delta_{2} (1-\theta_{2}) {H}_{I1} + \sigma {\lambda}_{H} {I}_{2} + \gamma^{(2)}_{0} {H}_{Q2} - (\mu + \kappa_{2} + \phi){H}_{I2}, \\\\ H_{Q1}' = (\kappa_{2} + \psi){H}_{I1} - (\mu + \phi + \gamma^{(2)}){H}_{Q1} + {\lambda}_{H} {Q}_{1}, \\\\ H_{Q2}' = {\lambda}_{H} {Q}_{2} + \kappa_{2} {H}_{I2} - (\mu + \phi + \gamma_{0}^{(2)}){H}_{Q2}, \\\\ A' = \phi ({H} + {H}_{I1} + {H}_{I2} + {H}_{Q1} + {H}_{Q2}) - (\mu + v){A}, \end{array} \right. \end{array} $$(6)

where λ_I = β_I [(1−α₁)(I₂ + (1−η₁) I₁) + (H_I2 + (1−η₂)H_I1)] and

λ_H = β_H [(1−α₂)(H + H_Q1 + H_Q2) + (H_I2 + (1−τ) H_I1)].

Since model system (6) monitors human populations, all the state variables and parameters of the model are non-negative. Consider the biologically-feasible region

Ω={(S,I1,I2,Q1,Q2,H,HI1,HI2,HQ1,HQ2,A)∈R+11:N≤1}.$$\begin{array}{} \displaystyle \Omega = \lbrace ( S,I_{1},I_{2},Q_{1},Q_{2},H,H_{I1},H_{I2},H_{Q1},H_{Q2},A) \in \mathbb{R}_{+}^{11} : N \leq 1 \rbrace. \end{array} $$(7)

The following steps are followed to establish the positive invariance of Ω. The rate of change of the total population, obtained by adding all the equations in model system (6), is given by

dNdt=μ−μN−vA≤μ−μN.$$\begin{array}{} \displaystyle \dfrac{dN}{dt} = \mu - \mu N - v A \leq \mu - \mu N. \end{array} $$(8)

It is easy to see that whenever N > 1, then dNdt$\begin{array}{} \dfrac{dN}{dt} \end{array} $ < 0. Since dNdt$\begin{array}{} \dfrac{dN}{dt} \end{array} $ is bounded by μ − μN, a standard comparison theorem [6] can be used to show that

N(t)≤N(0)exp−μt+(1−exp−μt).$$\begin{array}{} \displaystyle N(t) \leq N(0) \exp^{-\mu t} + (1 - \exp^{-\mu t}). \end{array} $$(9)

In particular, N(t) ≤ 1 if N(0) ≤ 1. Thus, every solution of the model system (6) with initial conditions in Ω remain there for t > 0. Thus Ω is positive invariant and attracting. Hence it is sufficient to consider the dynamics of the flow generated by model system (6) in Ω. In this region, the model can be considered as been epidemiologically and mathematically well-posed [7].

Before analysing the full model system (6), it is instructive to gain insights into the dynamics of the models for HIV only (HIV-only model) [1, 21, 38] and HSV-2 only (HSV-2 only model) [22].

In this section, we present an optimal control problem, describing our goal and the restrictions of the epidemics, in order to investigate an effective campaign on reducing HIV/AIDS and HIV/HSV-2 co-infections. HSV-2 anti-viral treatment adherence is low in some parts of the world especially in sub-Saharan Africa where HIV/AIDS is at its highest. Monitoring and counselling of individuals infected with HSV-2 only and those dually infected with HIV/HSV-2, to HSV-2 anti-viral treatment can be a possible aid in the reduction of HIV/AIDS. We will reconsider model system (6) and introduce two time-dependent control variables. We seek to reduce the number of individuals who quit HSV-2 treatment before completion. There will be a lot of costs generated during the control process. So it is advisable to balance between the costs and the non-adherence effects.

Assuming the control variables in the set

U={(u1,u2):[0,T]→R2|ui(t)is a Lebesgue measure on[0,Ui],i=1,2.},$$\begin{array}{} \displaystyle U = \lbrace (u_{1},u_{2}) : [0,T] \rightarrow \mathbb{R}^{2} \vert u_{i}(t) \,\,\mbox{is a Lebesgue measure on}\,\,[0,U_{i}], i =1,2. \rbrace, \end{array} $$(15)

in which all control variables are bounded and Lebesgue measurable and U_i, i = 1,2 is denoted to be the upper bound of the control variables. In our controls, u₁ denotes a time-dependent control effort on counselling and monitoring individuals with HSV-2 only, while u₂ is the same but for individuals dually infected with HIV and HSV-2.

In view of this, our optimal control problem is to minimise the objective functional given by

J(u1,u2)=∫0T[z1I2+z2HI2+b1u1I2+b2u2HI2+12(c1u12+c2u22)]dt$$\begin{array}{} \displaystyle J(u_{1},u_{2}) = \int_{0}^{T} \bigg [z_{1}I_{2} + z_{2}H_{I2} + b_{1}u_{1}I_{2} + b_{2}u_{2}H_{I2} + \dfrac{1}{2} (c_{1}u_{1}^{2} + c_{2}u_{2}^{2}) \bigg] dt \end{array} $$(16)

subject to the system

S′=μ−(λH+λI)S−μS,I1′=λIS+γ(1)Q1−σλHI1−(δ1(1−u1)+κ1+ψ+μ)I1,I2′=δ1(1−u1)I1−σλHI2+γ0(1)Q2−(μ+κ1)I2,Q1′=(κ1+ψ)I1−λHQ1−(γ(1)+μ)Q1,Q2′=κ1I2−λHQ2−(μ+γ0(1))Q2,H′=λHS−φλIH−(μ+ϕ)H,HI1′=φλIH−(δ2(1−u2)+ϕ+κ2+ψ+μ)HI1+γ(2)HQ1+σλHI1,HI2′=δ2(1−u2)HI1+σλHI2+γ0(2)HQ2−(μ+κ2+ϕ)HI2,HQ1′=(κ2+ψ)HI1−(μ+ϕ+γ(2))HQ1+λHQ1,HQ2′=λHQI2+κ2HI2−(μ+ϕ+γ02)HQ2,A′=ϕ(H+HI1+HI2+HQ1+HQ2)−(μ+v)A,$$\begin{array}{} \displaystyle \left\{\begin{array}{llll} S' = \mu - ({\lambda}_{H} + {\lambda}_{I}){S} - \mu {S},\\\\ I_{1}' = {\lambda}_{I} {S} + \gamma^{(1)} {Q}_{1} - \sigma {\lambda}_{H} {I}_{1} - (\delta_{1} (1- u_{1}) + \kappa_{1} + \psi + \mu){I}_{1}, \\\\ I_{2}' = \delta_{1} (1- u_{1}) {I}_{1} - \sigma {\lambda}_{H} {I}_{2} + \gamma^{(1)}_{0} {Q}_{2} - (\mu + \kappa_{1}){I}_{2}, \\\\ Q_{1}' = (\kappa_{1} + \psi){I}_{1} - {\lambda}_{H} {Q}_{1} - (\gamma^{(1)} + \mu) {Q}_{1}, \\\\ Q_{2}' = \kappa_{1} {I}_{2} - {\lambda}_{H} Q_{2} - (\mu + \gamma^{(1)}_{0}) {Q}_{2}, \\\\ H' = {\lambda}_{H} S - \varphi {\lambda}_{I} H - (\mu + \phi){H}, \\\\ H_{I1}' = \varphi {\lambda}_{I} H - (\delta_{2} (1- u_{2}) + \phi + \kappa_{2} + \psi + \mu){H}_{I1} + \gamma^{(2)} {H}_{Q1} + \sigma \lambda_{H} I_{1}, \\\\ H_{I2}' = \delta_{2} (1- u_{2}) {H}_{I1} + \sigma {\lambda}_{H} {I}_{2} + \gamma^{(2)}_{0} {H}_{Q2} - (\mu + \kappa_{2} + \phi){H}_{I2}, \\\\ H_{Q1}' = (\kappa_{2} + \psi){H}_{I1} - (\mu + \phi + \gamma^{(2)}){H}_{Q1} + {\lambda}_{H} {Q}_{1}, \\\\ H_{Q2}' = {\lambda}_{H} {Q}_{I2} + \kappa_{2} {H}_{I2} - (\mu + \phi + \gamma_{0}^{2}){H}_{Q2}, \\\\ A' = \phi ({H} + {H}_{I1} + {H}_{I2} + {H}_{Q1} + {H}_{Q2}) - (\mu + v){A}, \end{array} \right. \end{array} $$(17)

where z_i, b_i, and c_i (i = 1,2) are (positive) balancing coefficients transferring the integrals into a monetary quantity over a finite period of T months. The z_i values represent the weights of the individuals infected with acute HSV-2 only and not under anti-viral treatment and the individuals dually infected with HIV and HSV-2, also not under HSV-2 anti-viral treatment, respectively. The terms with b₁ and b₂ represent the costs associated with intervention strategies in monitoring and counselling. We assume that the costs are proportional to the quadratic form of their corresponding control functions. The objective functional in (16) also includes some quadratic terms with coefficients c₁ and c₂ to indicate potential nonlinearities in the costs.

Our objective is to find an optimal control pair (u1∗(t),u2∗(t))$\begin{array}{} (u_{1}^{*}(t),u_{2}^{*}(t)) \end{array} $ in order to seek the minimum value of the objective functional J(u1∗(t),u2∗(t)),$\begin{array}{} J(u_{1}^{*}(t),u_{2}^{*}(t)), \end{array} $ such that

J(u1∗(t),u2∗(t))=min{J(u1(t),u2(t))|u1(t),u2(t)∈U}$$\begin{array}{} \displaystyle J(u_{1}^{*}(t),u_{2}^{*}(t)) = \mbox{min} \lbrace J (u_{1}(t),u_{2}(t)) \vert u_{1}(t),u_{2}(t) \in U \rbrace \end{array} $$(18)

subject to the system given by (17).

We now derive necessary conditions that the control pair and corresponding states must satisfy. By using the same method as in [33], the existence of the optimal control can be proved. In the above minimising problem, we can easily verify that the objective functional is convex on the closed, convex control set U. The optimal system is bounded, which determines the compactness needed for the existence of the optimal control.

In order to find an optimal solution of model system (17), first let us define the Hamiltonian functions H for the optimal control system (17) as

H(t,X,U,λ)=L+λ1dSdt+λ2dI1dt+λ3dI2dt+λ4dQ1dt+λ5dQ2dt+λ6dHdtλ7dHI1dt+λ8dHI2dt+λ9dHQ1dt+λ10dHQ2dt+λ11dAdt,$$\begin{array}{} \displaystyle H(t,X,U,\lambda) = L + \lambda_{1}\dfrac{dS}{dt} + \lambda_{2}\dfrac{dI_{1}}{dt} + \lambda_{3}\dfrac{dI_{2}}{dt} + \lambda_{4}\dfrac{dQ_{1}}{dt} + \lambda_{5}\dfrac{dQ_{2}}{dt} + \lambda_{6}\dfrac{dH}{dt} \\\displaystyle\qquad\qquad\qquad\,\,\,\,\, \lambda_{7}\dfrac{dH_{I1}}{dt} + \lambda_{8}\dfrac{dH_{I2}}{dt} + \lambda_{9}\dfrac{dH_{Q1}}{dt} + \lambda_{10}\dfrac{dH_{Q2}}{dt} + \lambda_{11}\dfrac{dA}{dt}, \end{array} $$(19)

where L is the Lagrangian function L = z₁I₁ + z₂H_I1 + b₁u₁I₂ + b₂u₂H_I2 + 12(c1u12+c2u22).$\begin{array}{} \dfrac{1}{2} (c_{1}u_{1}^{2} + c_{2}u_{2}^{2}). \end{array} $ With the existence of the optimal control system, we now present and discuss the adjoint system and the characterisations of the optimal control system.

For simplicity, we denote

X(t)=(S(t),I1(t),I2(t),Q1(t),Q2(t),HH1(t),HI2(t),QI1(t),QI2(t),A(t))$$\begin{array}{} \displaystyle X(t) = (S(t),I_{1}(t),I_{2}(t),Q_{1}(t),Q_{2}(t),H_{H1}(t),H_{I2}(t),Q_{I1}(t),Q_{I2}(t),A(t)) \end{array} $$(20)

and λ = (λ₁, λ₂,λ₃, λ₄,λ₅, λ₆,λ₇, λ₈,λ₉, λ₁₀,λ₁₁).

4.1

Numerical Simulations

This section is dedicated to the numerical study of the control strategies against HIV and HSV-2 co-infection. Numerical results in this section were generated using the forward and backward sweep methods (FBSM). For detailed explanation on FBSM we refer the reader to [36].

In order to illustrate the results of the foregoing analysis, we have simulated model system (6) using the parameters in Table 1. Unfortunately, the scarcity of the data on HSV-2 and HIV/AIDS correlation with a focus on HSV-2 treatment adherence limits our ability to calibrate; nevertheless, we assume some of the parameters in the realistic range for illustrative purposes. We use parameters obtained from previous research by Feng et al. [1], Foss et al. [25] and Abu-Raddad et al. [26]. In this work, the authors report that the model fits well with some of the parameters/data from the previous research by Abu-Raddad et al. [26] on HIV/HSV-2 interaction model. Reliable data on HSV-2 and HIV/AIDS correlation with a focus on HSV-2 treatment adherence would enhance our understanding and aid in the possible intervention strategies to be implemented.

Table 1

Model parameters and their baseline values. The time unit is in months

Definition	Symbol	Baseline values (Range)	Source
Effective contact rate for HSV-2 infection	βI	0.2	[1, 25]
Effective contact rate for HIV infection	βH	0.3 (0.11-0.95)	[30]
Rate of acute HSV-2 becoming latent	κ1	2.3805(2.083-2.678)	[26]
Rate of acute HSV-2 and HIV infected becoming latent	κ2	1.683(1.458-1.875)	[26]
Treatment quitting rate	δ1	0.3	[26, 27]
Treatment quitting rate (dually infected)	δ2	0.3	Assumed
Proportion of HSV-2 patients who successfully complete treatment	θ1	0.7(0-1)	Assumed
Proportion of patients who successfully complete treatment (dually infected)	θ2	0.8(0-1)	Assumed
Reactivation rate with an effect of treatment (HSV-2 only)	γ(1)(ψ)	Varies	[1]
Reactivation rate with an effect of treatment (dually infected)	γ(2)(ψ)	Varies	[1]
Baseline reactivation rate of latent HSV-2	γ0(1)$\begin{array}{} \gamma_0^{(1)} \end{array} $(ψ)	0.339-0.436	[26]
Baseline reactivation rate of latent HSV-2 (dually infected)	γ0(2)$\begin{array}{} \gamma_0^{(2)} \end{array} $(ψ)	0.365-0.469	[26]
Enhanced susceptibility of people with acute HSV-2 to HIV infection	σ	2.7-3.1	[29]
Enhanced susceptibility to HSV-2 infection by HIV infectives	φ	≥ 1	[1]
Treatment rate of acute HSV-2	ψ	Varied	[1]
Rate of progression from HIV to AIDS	ϕ	0.0104	[26]
Average sexual lifespan	μ	0.004(0.003-0.005)	[25, 26]
AIDS-related death rate	v	0.03	[32]
Modification parameter	α1	(0−1)	Assumed
Modification parameter	α2	(0−1)	Assumed
Modification parameter	η1	(0−1)	Assumed
Modification parameter	η2	(0−1)	Assumed
Modification parameter	τ	(0−1)	Assumed

The estimated parameters, among which are the weights on costs c₁ and c₂, the weights on I₂ (individuals infected with acute HSV-2 only and not under antiviral treatment) and H_I2 (individuals dually infected with HIV and acute HSV-2, also not under HSV-2 anti-viral treatment); have been chosen for illustrative purposes. Thus z₁ = 5000, z₂ = 5000, b₁ = 1, b₂ = 1, c₁ = 0.5, c₂ = 0.5, and the following initial population levels are used in all the simulations S(0) = 0.6, I₁(0) = 0.25, I₂(0) = 0.0, Q₁(0) = 0.0, Q₂(0) = 0.0, H(0) = 0.01, H_I1(0) = 0.0, H_I2(0) = 0.0, Q_I1(0) = 0.0, Q_I2(0) = 0.0, A(0) = 0,0 over a period of 36 months (3 years).

Figure 2(a) illustrates the impact of controls on the population of the individuals infected with with HIV-only. The population of the HIV-only have moderate increases for the whole period under study for both cases in the presence and absence of controls. It is worth noting that the controls are effective for a period from 6 to 36 months. Further, we note that after 36 months, the controls have an effect of reducing the HIV cases by approximately 35%.

Fig. 2

Numerical results in figure 2(b) illustrates the impact of optimal counselling and monitoring on the individuals infected with HSV-2 only and those dually infected with HIV/HSV-2 on the AIDS cases. Overall we observe that the population of the AIDS cases increases for both cases ( with and without the controls ) for the whole period under review. For the period 4–36 months, we observe an increase of AIDS cases in the absence of controls as compared to the presence of controls. Thus, the controls start showing some effectiveness on the AIDS cases after 4 months.

Figure 2(a) illustrates the proportion of individuals who are dually infected with HIV and acute HSV-2, also under HSV-2 anti-viral treatment, in the presence and absence of controls. For both cases, in the absence and presence of controls, there is a sharp increase for the period 0 – 15 months. It is worth noting that both cases reach the maximum value at 15 months. The maximum proportion reached in the presence of controls is approximately 0.004 and the maximum proportion in the absence of controls is approximately 0.01. For the period 15–36 months for both cases, we note that they both drop uniformly with the H_I1 proportion higher in the absence of controls as compared to their presence.

Figure 2(b) illustrates the proportion of individuals who are dually infected with HIV and acute HSV-2, also not under HSV-2 anti-viral treatment, in the presence and absence of controls. In the presence of controls we observe that the population is lower, than in the absence of controls for the whole period under study. It is worth noting that, the controls start being effective after 2 months.

Figure 4 shows the evolution of the control profiles over time. An interesting result is that the second control, u₂, starts at the lower bound, while u₁ starts at the upper bound. This can be explained by our assumption that, in a given community, an individual can be infected by either HIV or HSV-2 but not both at the same time. The results suggest that more effort should be devoted to optimal counselling and monitoring of individuals dually infected with HIV and acute HSV-2 control u₂ which is feasible for up to 36 months, as compared to individuals infected with acute HSV-2 only control u₁ which is feasible for about 27 months. However, if optimal counselling and monitoring is implemented for 27 months or less, then both controls can be feasible.

Fig. 3

Fig. 4

HSV-2 is the most prevalent disease in most parts of the world especially in sub-Saharan Africa. HSV-2 is a significant factor for increased risk of acquisition and transmission of HIV and also being the leading disease in causing genital ulcers. Thus, proper HSV-2 treatment adherence, might be beneficial in the reduction of HIV/AIDS. In this paper, we formulated and analysed a deterministic compartmental model for the transmission dynamics of HIV and HSV-2. We considered the epidemiological synergy between sexually transmitted HIV and HSV-2. We calculate the basic reproduction number of the model. Applying the comparison theorem by Lakshmikantham, we managed to prove the global stability of the disease-free equilibrium. The Centre Manifold theory has been used to show that the endemic equilibrium point is locally asymptotically stable when the associated reproduction number is greater than unity. Introducing two time dependent controls to our model, we then formulated an optimal control problem for model system (6) to investigate optimal monitoring and counselling which seek to minimise the number of individuals infected with HIV/AIDS and those dually infected with HIV and HSV-2. Optimal monitoring and counselling is applied so as to try and reduce the number of individuals who quit HSV-2 anti-viral treatment before completion. We proved the existence of the optimal control and characterised the controls using Pontryagin’s Maximum Principle. From the illustrations in this study it is clear that time dependent intervention strategies can lead to the reduction of HIV/AIDS and those dually infected by HIV and HSV-2. It is worth noting that, as much as the HIV cases and the HIV/HSV-2 dual cases are reduced, they do not converge to zero since there are other factors influencing their prevalence within the population. Numerical analysis also suggests that the implementation of the controls should be denoted to both since they are feasible for a longer period considered in this study. Further, the optimal control results suggests that more effort should be devoted to monitoring and counselling of the individuals dually infected with HIV and HSV-2 compared to those infected with HSV-2 only. It would be vital to implement these controls at minimum costs, hence we suggest the monitoring and counselling can be disseminated at various health care centres. Various methods that target HIV and AIDS patients, can now involve some information on the dangers of poor HSV-2 treatment adherence.

Our model has several limitations, which should be acknowledged. Limited data exist on the co-infection of HSV-2 and HIV/AIDS, particularly mathematical modelling study on HSV-2 is low. Therefore, some of our numerical estimates remain uncertain, particularly those that are most sensitive to HSV-2 related parameters, which we took from other previous research. More data sets and experimental studies are needed to include more realistic biological processes in the models. However just like any other model, we cannot say the model is complete it can be extended to include resource limited or resource given communities, with HIV treatment compartments included.