Analysis of an HIV-HCV simultaneous infection model with time delay

A novel mathematical delay model for simultaneous infection of HIV and hepatitis C virus is formulated and dynamically analyzed. Basic properties of the model are established and proved. Basic reproductive threshold is systematically calculated as the maximum of three subthreshold parameters. A disease free equilibrium is determined to be globally asymptotically stable for all values of the delay when the threshold is less than unity. However, when the threshold is greater than one, endemic equilibrium emerged which is shown to be locally asymptotically stable for any length of delay. Although the delay has no effect on stabilities of equilibria points, however, it is found to reduce the infectivity of the viruses as the length of the delay is increased. Epidemiological interpretations of the results and numerical simulations illustrating them are given. DOI:10.46481/jnsps.2021.109


Introduction
Around 2.75 million people who have human immunodeficiency virus (HIV) are coinfected with hepatitis C virus (HCV) globally and and on average HIV-infected individuals are six time more likely to get HCV infection than HIV-uninfected [1]. Worldwide, HIV and HCV are public health challenges which affected several populations. Across the continents, there are about 37 million and 115 people infected with HIV and chronic HCV infections, respectively.
In Dhutta and Gupta [25], a mathematical model is proposed for the dynamics of HIV/AIDS with incorporation of weak CD4 + T cells. They computed the local stability of the infection-free and infection equilibria for the model when the 1 Hassan & Hussaini / J. Nig. Soc. Phys. Sci. 3 (2021) 1-11 2 valuable reproduction number is less than and greater than one respectively. Furthermore, using Lyapunov's second method and the geometric approach, they define novel conditions for the global stability of the equilibria. Further, Maimuna [27] formulated a mathematical model for the spread of HIV with an Anti-Retroviral Therapy (ART) intervention. The author established that dynamics of the model with no ART depends on the basic reproduction number. However, sensitivity analysis revealed that the basic reproduction number decrease with increasing number of infected humans who follow the ART treatment. Jia et al [26], constructed a dynamic model for HCV transmission and prevalence based on the reported data from China and determined the most influential parameters to evaluate the effectiveness of control measures. Finally, in [28], the authors applied ideas from the ecology of infectious diseases to model the transmission of HCV in a population of injection drug users. The authors suggested that modelling HCV as an indirectly transmitted infection facilitates a more nuanced understanding of disease dynamics. In addition, sensitivity analysis of parameters on the value of R 0 was carried out and parameters related to an interaction with the environmental reservoir are found to be more influential.
Time delay in epidemiology is mostly incorporated to represent developmental stages, incubation period or wining of immunity [29,30]. It is an important component that cannot be neglected as it affects the dynamics of models causing instabilities of equilibria resulting in Hopf bifurcations [32,31,33]. In addition, an important threshold value (the basic reproduction number), as can be seen later in this article, is expressed in terms of time delay [34,30,35]. Many delay models of HCV/hepatitis B virus (HBV) have been presented in the literature to explore the effects of delay as incubation period. Banerjee et al. [36] presented an intracellular time delay model for hepatitis C virus which has shown that the delay doesn't affect the stability of steady-state, however, destabilized the infected equilibrium with resulting in Hopf bifurcation. Gourley et al. [37] considered the dynamical properties of HBV -delayed model with standard incidence formulation. They realized that the infected steady-state, expressed in terms of delay, is globally stable regardless of the time delay length. In [33], Zhao and Xu presented a delay HCV model using Beddigton-DeAngels formulation and obtained global stabilities of equilibria when the threshold parameters satisfy certain conditions.
In this study, we begin by developing a model using delay differential equation for HIV and HCV co-infection which includes simultaneous transmission from person-to-person. Unlike in [5,14,15,24], we introduce delay in order to capture the incubation periods of both HIV and HCV infections.
In this research, the model is presented in Section 2 while the analysis is given in Section 3 . In Section 4, we present illustrative graphs for the dynamical properties of the model. Lastly, conclusion and summary of our findings are reported in Section 5.

Model Formulation
The total human population at time t, N(t), is divided into four different classes of healthy (susceptible, (S (t))), HIV-only infected individuals (I H (t)), HCV-only infected individuals (I C (t)), individuals simultaneously infected with both HIV and HCV (I HC (t)), so that The susceptible population (S (t)) is increased by the recruitment of people into the community at a constant rate Π (all newly-recruited individuals are assumed to be susceptible to both infections) and by recovery from HCV (at a rate ψ). The population is decreased by infection with HIV only (at a rate λ H ) or HCV only (at a rate λ C ) or simultaneous infection with HIV and HCV (at a rate λ HC ) (see, for instance, [19,17] and references therein for the biology). However, the infections with HIV and HCV do not occur as soon as the infectives get into contact with the susceptible individuals. Rather, there are time lapses or delays, representing the latent periods, in which the infectives progress to fully infectious people. Thus we have I H (t − τ 1 )e −µτ 1 representing the HIV only infectives that can only infect after the elapse of τ 1 , (the latent period for HIV) and e −µτ 1 , the survival probability over the latent period. Similarly, I C (t − τ 2 )e −µτ 2 represent the HCV only infectives that can infect after the elapse of τ 2 , (the latent period for HCV), e −µτ 2 is the survival probability for natural death over the latent period. Lastly, I HC (t − τ 3 )e −µτ 3 stand for the HIV and HCV infectives that can infect after the elapse of τ 3 , (the latent period for simultaneous infection of HIV and HCV), while e −µτ 3 is the survival probability over the latent period. Thus, In (1), β H , β C and β HC represent the effective contact rates for HIV, HCV and simultaneous HIV/HCV infections, respectively. The parameter η > 1 accounts for the assumed increase in infectiousness of dually-infected individuals due to high immunosuppression caused by the simultaneous infection of HIV and HCV, compared to singly-infected with HIV or HCV. For mathematical simplification we assume that τ 1 = τ 2 = τ 3 = τ.
The susceptible population further decreased by natural death (at a rate µ; natural death occurs in all compartments at the same rate).
The population of individuals infected with HIV-only (I H ) is generated by HIV infection following the effective contact (at the rate λ H ) and by dually-infected individuals (I HC ) recovered from HCV (at a rate κψ; where the modification parameter 0 < κ < 1 models the reduced likelihood of dually-infected 2 person to recover from HCV, in comparison to those infected with HCV only). It is decreased by HCV infection at rate λ C , HIV-related death at rate δ C and natural death.
The population of individuals infected with HCV-only (I C ) is generated by HCV infection (at the rate λ C ). It is decreased by HIV infection at rate λ H ), recovery from HCV at rate ψ, HCV-related death at rate δ H and natural death.
The class of individuals simultaneously infected (I HC ) is generated by HIV and HCV simultaneous infections at rate λ HC , by HIV-only infected individuals (I H (t)) who are HCV infected at rate θ 1 λ C and by HCV-only infected individuals (I C (t)) who are HIV infectious at rate θ 2 λ H , where the parameter θ i ≥ 1 i = 1, 2 captures the fact that singly-infected individuals have weak immune-system (due to illness) compared to wholly-susceptible individuals. The population suffers death due to the co-infection of both diseases (at a rate δ HC ).
The following system of equations describes the model based on the aforementioned assumptions and notations while the parameters are presented in Table 1: with initial data where φ 1 (t), · · · , φ 4 (t) are continuous functions on the interval [−τ, 0] and

Basic properties
In this part, we present results for basic qualitative properties of the model as follows.
Theorem 1. The solution (S (t), I H (t), I C (t), I HC (t)) of the model (6), with initial data (7), exists for all time, t > 0 and is unique. Furthermore, the solution is nonnegative for all t > 0.
Proof. For existence and uniqueness of solution, we start as follows: Let X(t) = (S (t), I H (t), I C (t), I HC (t)). The system (6), can be represented as where X t (θ) = X(t + θ), f is Lipschitz and continuous in X. It follows from Theorems 2.1, 2.3 in [38] that the model (6) has a unique solution (S (t), I H (t), I C (t), I HC (t)) satisfying the initial data (7). The positivity of solution is proved by the method of contradiction as shown below: Suppose the conclusion is not true. Under the given initial function, there exists a time, t 1 ∈ [0, +∞) where S (t) will changes sign, at least once, so that S (t) > 0 for all t ∈ (0, t 1 ), For the first case, consider the first equation in (6), thus This is a contradiction to the earlier assumption that S (t 1 ) = 0. Therefore, t 1 doesn't exists, hence S (t) > 0 for all t > 0.
Similarly for the second argument, consider the second equation in (6): This also contradicts the earlier assumption that I H (t 2 ) = 0. Therefore, there is no such time t 2 exists. Hence I H (t) > 0 for all t > 0. Similar arguments can be applied to other two equations to show that I C > 0 and I HC > 0 for all t > 0.
Theorem 2. The biologically-feasible region Ω is positivelyinvariant for the model (6), where Assumed increase in infectiousness of dually infected to single infection 0 < κ < 1 Modification parameter for reduced dually infected to single infection µ Natural death rate in all classes δ H , δ C , δ HC HIV, HCV, HIV/HCV induced death rates Modification parameters for singly-infected to wholly susceptibles Proof. Adding the equations in (6) gives (8) so that, with positivity of I H , I C , I HC . It follows from (9), using Gronwall lemma [39], that Therefore, Ω is positively invariant.

Existence and stability of equilibria
At equilibrium, equations in (6) are equated to zero. In the absence of diseases, the DFE is obtained to be It should be noted that at equilibrium, HC . Hence the following gives the implicit values for the endemic equilibrium: where The local asymptotic stability of a given equilibrium of the model (6), is determined by linearizing the system about the equilibrium point [40,41], and showing that all the roots of the transcendental polynomial have negative real parts. Thus linearizing the system (6) about the following variables: a 1 a 2 a 2 a 3 a 4 a 5 a 6 a 7 a 8 a 9 a 10 a 11 a 12 a 13 a 14 a 15 ...d 4 ), be constants and z is a complex number. Equation (13) can be express as: However, before we determine the stabilities of equilibria, it is important to calculate the basic reproduction number of the model, denoted by R 0 . It can be recall that, susceptible humans acquire HIV, HCV infections and HIV-HCV simultaneously following effective contact with infected individuals with these diseases. It follows that, the number of HIV infectives generated by an infected human (near the DFE) is given by the product of the infection rate ( β H N * ), the probability that an infected human survives natural death (e −µτ ) and the average duration of stay in the infected HIV class ( 1 K 1 ). Hence, the average number of HIV infections generated by infected individual is given by Similarly, the number of HCV infections generated by an infected individual with HCV (near the DFE), is given by the product of the infection rate ( β C N * ), the probability that an infected human survived natural death (e −µτ ) and the average duration of stay in the infected HCV class ( 1 K 2 ). Therefore, the average number of HCV infected generated by infectious individual is given by Furthermore, the number of HIV-HCV simultaneous infections generated by an infected individual with the two diseases (near the DFE), is given by the product of the infection rate of infected HIV-HCV ( β HC N * ), survival probability (e −µτ ) and the average duration of stay in the class I HC ( 1 K 3 ). Thus, the average number of HIV-HCV infections generated by infected individual is given by Finally, the maximum of the three sub thresholds in (15), (16) and (17) gives the value of R 0 , i.e. the average number of new infections generated by infected individuals (with HIV, HCV and HIV-HCV) is given by noting that at DFE, S * = N * = Π µ . It is worth stating here that, the formulation of R 0 in (18) is derived from the procedure in [12] used in a model that have more than one strains/diseases, to be the maximum of various sub thresholds for the infections. Moreover, similar approaches are also applied in [46,42].

Local asymptotic stability of DFE
At DFE, we claim the stability properties as follows: Theorem 3. The DFE E 0 , is locally asymptotically stable (LAS) in Ω whenever R 0 < 1 for all τ ≥ 0 and unstable if R 0 > 1.
Proof. At DFE, the transcendental equation (14) is reduced to Expanding equation (19), we have where First we determine the stability of E 0 when τ = 0. From (19), when τ = 0, it is obvious that the eigenvalue −µ is negative, Further, we determine the stability of E 0 when τ > 0. Here we investigate whether there is a root z = iy of (19), y ∈ R + , that may cross the imaginary axis and cause stability switches. Substituting z = iy in (19) or equivalently in (22), we find y that must satisfy where P 1 (iy) = (y 4 − P 1 iy 3 − P 2 y 2 + P 3 iy + P 4 ), Q 1 (iy) = (−Q 1 iy 3 − Q 2 y 2 + Q 3 iy + Q 4 ). However, (22) is implicit in e −zτ , therefore the process of separating real, imaginary parts and squaring (23) may not lead to explicit expression in y. Without loss of generality, one can see from (19) that if R 0 < 1 and G(y) = +∞ as y → ∞. Therefore, when τ > 0, equation (19) has no positive root on (0, +∞), hence the DFE is LAS for all τ ≥ 0 if R 0 < 1 and unstable if R 0 > 1. The biological implication of this result is that the diseases can be completely eliminated whenever the initial population is close to the basin of attraction of the DFE.

Global asymptotic stability of DFE
To ensure complete eradication of the diseases in the population irrespective of the population size started with, we present and prove the following result.
Proof. The global properties can be established using a Comparison Theorem [23] and the approach in [43]. The equations of the three infected components in model (6) can be expressed using matrix-vector form as where F =  are the matrices for new infection terms and transitions respectively.
Since S ≤ N at every time t in Ω and by letting V = V 1 + V 2 , it follows from (24) that It can be shown that the spectral radius of F V −1 gives the value of R 0 in (18). Furthermore, since the eigenvalues of the matrix (F − V) have negative real parts if R 0 < 1 [43,44], consequently, the linearized differential system (25) is stable whenever R 0 < 1. As a result, lim t→∞ (I H (t), I C (t), I HC (t)) → (0, 0, 0).
Therefore, from the first equation in (6) and substituting I H (t) = I C (t) = I HC (t) = 0 from (26), for any t, it follows that so that Thus, By LaSalles invariance principle [45] E 0 is the largest invariance set in Ω, hence is GAS whenever R 0 < 1.

Stability analysis of endemic equilibrium
Here, due to the complexity of transcendental equation, we give conditions for stability of any unique endemic equilibrium (EE), whenever it exists.
From (12) the transcendental equation (14) is simplified to β C e −µτ S N a 11 a 12 −a 13 where P 2 (z) and Q 2 (z) are polynomials of degree 4 and 3 respectively, with real coefficients and hence have no common imaginary roots whenever R 0 > 1.
Furthermore, substituting z = iy, as purely imaginary root, we define which is a polynomial of degree 8 whose leading coefficient is positive. Suppose (31) has at least zero positive root, then the following stability result for E * * can therefore be stated Theorem 5. If the polynomial G 2 (y) has (a) no positive root, then E * * is locally asymptotically stable for all delays τ ≥ 0 whenever R 0 > 1.
(b) at least one simple positive root, then as delay increases there will be n ∈ Z + number of stability switches for fixed parameter values and the endemic equilibrium E * * is locally asymptotically stable for 0 ≤ τ < τ * if R 0 > 1, here, y and the subscripts represents the real and imaginary parts of P 2 and Q 2 .

Interior endemic equilibrium
Here, we consider one of the interior equilibria with simultaneous infection of HIV and HCV only. In the absence of HIV and HCV only infections, equation (12) will be reduced to provided R HC e −µτ β HC µ+e −µτ β HC . Similarly, the transcendental equation (30) is simplified to be with When τ = 0, (33) yields so that after simplification, we have (p 1 p 4 − p 2 p 3 ) + (p 1 p 5 + p 3 p 5 ) = δ HC K 2 3 (R HC − 1) 2 + µK 3 (R HC − 1). (35) From (35), it can be seen that all the coefficients of G 2 (z) are positive whenever R HC > 1. Therefore, G 2 (z) is stable (all the roots have negative real parts).

Analysis for impact of delay
The impact of the time delay is analyzed using threshold approach to show whether increase (decrease) of the delay will have effect on the infectivity of the two diseases in the community. We use the threshold quantity, R 0 , basic reproduction number, as a function of the delay τ. Since R 0 = 1 is a threshold value for the infection in the population, we determine the critical delay value τ crit above (below) which the diseases have effect as follows:   The rate of change of R 0 with respect to τ is given by Since the derivative in (41) is negative, it indicates that R 0 decreases with increase of the time delay τ and vice varsa. This implies that, increasing the delay beyond the threshold value, with fix parameter values, will results in decreasing the number of infected individuals with the two diseases and those that are simultaneously infected. This result can be summarized as Proposition 1. The increase in time delay in the model (1) will reduces the infectivity of the two diseases whenever the delay is greater than the threshold value τ crit .

Numerical simulations
Numerical experiments are conducted to illustrate the dynamics as described in theoretical results and display the impact of delay in the model. Figure 1 (a) depicts the global asymptotic stability for the DFE using parameter values from Table 2, except for β H = β C = 0.002, β HC = 0.0025, τ = 3 and different initial conditions so that R 0 = 0.16. It can be seen that the susceptible and total populations of infectives converges to Π µ and zero asymptotically, respectively. In Figures 1 (b)-(d), the local stability of endemic equilibrium is displayed using values from Table 2 except for in (b) τ = 1, β H = β C = 0.02, β HC = 0.025, (c) τ = 10 β H = β C = 0.002, β HC = 0.0025, and (d) τ = 30 , β H = β C = 0.04, β HC = 0.045, as proved in Theorem 5(a), so that R 0 = 1.60, R 0 = 1.42 and R 0 = 1.45 respectively. In Figure 2 (a), it is shown that the length of delay (τ = 30) doesn't cause any oscillation in the stability of the endemic equilibria as observed in some delay models. Figure 2 (b), illustrates the stability if one of the interior equilibrium, E * * 3 using parameter values in Table 2 while β H = β C = 0.0002, β HC = 0.025, so that R HC = 1.62. In this case, one or simultaneous infection will occur. The effect of delay is illustrated in Figure 2 (c) to support Proposition 1 in which increasing the delay, τ from 1, 15, 25 to 35, caused remarkable decrease in the total number of infected individuals with HIV and HCV.

Concluding remarks
In this work, a delay model of simultaneous infection of HIV and HCV is formulated and dynamically analyzed. The novelty (results) of the research is summarized and discussed below: (i) Basic properties (existence, boundedness and positivity of solution)of the model are stated and proved as Theorems 1 and 2.
(ii) The basic reproduction threshold is systematically obtained as the maximum of subthreshold values of the individual viruses; a threshold parameter above which the viruses will persist in the population.
(iii) A disease free equilibrium is found to be globally asymptotically stable when the basic reproduction threshold is less than unity for any length of time delay. Under this case, the diseases will die out in the population whenever the basic reproduction number is brought and maintained below one irrespective of the initial populations started with. This is shown in Figure 1(a).
(iv) However, if the basic reproduction threshold is greater than unity, endemic equilibrium exists which is shown to be locally asymptotically stable for all values of the incubation period (delay) as shown in Figures 1(b -d). This implies that whenever the initial population is within the basin of attraction of the endemic equilibrium and the basic reproduction number is greater than one, the diseases will persist in the population no matter the length of delay.
(v) Increasing the length of time delay is observed to be decreasing the number of cumulative infectious people when the delay is above a critical value. It is worth remarking here that although the time delay has no effect on the stability of equilibria however, it do affect the infectivity of the viruses as shown in the formulation of basic reproduction number. The implication of this result is, if the incubation period can be extended by say intervention, the number of infected people will be reduced.
(vi) Numerical simulations, using data from the literature are used to illustrate our results.