Model analysis for an HIV infectious disease using elasticity and sensitivity techniques

1.   Introduction

The human immunodeficiency virus (HIV) weakens an individual's immune system by infecting the CD4⁺ T cells. On the other hand, acquired immune deficiency syndrome (AIDS) describes the symptoms that follow HIV infection, which results in a compromised immune system. The United Nations Program on HIV/AIDS (UNAIDS) and the World Health Organization (WHO) reported that 36.7 million people worldwide were living with HIV/AIDS in 2016. Despite significant advancements in science and effective health intervention techniques, HIV remains one of the most devastating diseases in human history. HIV is mainly transmitted from mother to child during pregnancy, birth, or breastfeeding, as well as through unprotected sexual contact, sharing injectors and other injection tools with HIV-positive individuals, and tainted blood transfusions. The virus causes a serious infection once it enters the body, which frequently presents flu-like symptoms. The immune system becomes increasingly compromised as the disease grows more serious, thus increasing the likelihood that the person with the illness will get additional illnesses, known as opportunistic infections, that are uncommon for healthy people.

There are four stages to HIV disease:

(1) Individuals who are living with HIV (PLHIV) and have progressed to stage one might experience mild symptoms such as fever, diarrhea, and flu-like symptoms.

(2) Individuals who have progressed to stage two may exhibit signs such as skin disorders, swollen lymph nodes, and Tuberculosis (TB).

(3) Stage three symptoms might include lymph node-related TB and other symptoms affecting the mucous membranes.

(4) In the last stage, systemic meningoencephalitis may strike a person. A frequently term for stage four is AIDS.

Although there is currently no vaccine or cure for HIV, antiretroviral therapies can reduce the disease's progression as well as improve the survival rates to levels that are close to a normal life [1]. HIV infection-related morbidity and death have been successfully decreased by highly active antiretroviral therapy. HIV is one of the most dangerous and devastating infectious diseases. The most important immune system cells, CD4⁺ cells, are harmed by HIV infection. The virus gradually damages the human immune system, making the infected individual more susceptible to disease [2].

One important technique for performing time and cost-cutting analyses of real-world problems is mathematical modeling. Additionally, it is an excellent tool to characterize, predict, analyze, and suggest potential preventative measures for many infectious diseases [3]. In recent years, researchers from various scientific fields have become interested in the topic of modeling in infection disease such as [4]–[7].

In the field of mathematical modeling, analyzing the model dynamics of HIV compartments became an interesting topic locally and globally. A dynamic system is represented by a mathematical model derived from mathematical concepts, which is essential for the dynamical system of HIV transmission prediction, examination, and control. While a model can be modified using controlling functions, developing a model demands several assumptions and parameters [8]. A sensitivity analysis makes it possible to investigate which input variables typically cause output variation, as well as how uncertainty in the input variables affects the model outputs [9]. The most important aspect of the biological procedure is mathematical modeling, which involves expressing hypotheses as either mathematical statements or models that are subsequently utilized to predict situations and/or provide decisions. Models of disease transmission, statistical models, dynamical systems, and other types of mathematical models are among them. There are various steps in converting an idea into a theoretical model, which is followed by a quantitative model. It is clear that our idea is represented by arrows and boxes in a theoretical model, and that chemical kinetics modeling is used to convert physical and biological realities into mathematical representations. Additionally, there are some other real world problems that can be expressed in terms of mathematical equations [10]–[14]. The field of mathematical modeling of HIV infections has been proposed in recent decades such as [15], [16]. The authors showed that the suggested models could fit well with some clinical HIV infection data sets by examining the local and global stability of the endemic states. In [17], the authors investigated the effects of HIV therapy, screening, and educational campaigns on the dynamics of HIV transmission.

To discuss the progression of HIV and predict the peak of the pandemic, we propose a five-compartmental mathematical model (HEICV model), shown in Figure 1, that includes uninfected CD4⁺ T cells (H), HIV-exposed CD4⁺ T cells (E), an HIV-infected CD4⁺ T cells (I), cytotoxic T lymphocyte (CTL) cells (C), and free HIV particles in the blood V(t).

This research aims to investigate the transmission dynamics of HIV/AIDS infection by numerical simulation and a sensitivity analysis for different model initial states and parameters. Moreover, the basic reproduction number (R₀), is calculated to identify the critical model transmissions.

The following are the main contributions of the present investigation:

(1) Obtaining the basic reproduction number (R₀) and evaluating the stability of the HIV/AIDS-free equilibrium point.

(2) A sensitivity analysis aids in determining which variable is very sensitive regarding the parameters. By identifying the parameters that impact the most of spread and control of HIV, we can concentrate on the most critical aspects of the epidemic.

(3) The concept of elasticity is important in describing the sensitivity between R₀ and the model parameters in the infectious disease models.

The article is structured as follows. In Section 2, we introduce the modified model of HIV infectious transmission provided in [18]. This model is defined by nonlinear ordinary differential equations with constant rates based on the mass action law. In Section 3, the basic reproduction numbers are calculated for both the disease-free equilibrium point, E₀, and the endemic equilibrium point, E₁. In Section 5, we discuss the local stability of the disease-free equilibrium point. Accordingly, the methods of local sensitivity and elasticity are used to identify the critical model parameters; see Sections 6 and 7. Moreover, some numerical results are computed for the different initial states and parameters in Section 8. Finally, we conclude the main results in Section 9.

2.   Mathematical formulation of the problem

In this section, we extend an HIV infectious disease model given in [18]; the model describes the interaction between the HIV viruses, CD4⁺ T cells, infected cells, and the CTL immune response. In this model, we have included the natural death rates to the compartments that represent exposed people (E), infected people (I), free virus particles (V), and CTL in order to create a more accurate model of HIV dynamics. This makes it possible for the model to represent the natural mortality that happens separate from HIV-related variables. The model includes six variables: the total number of healthy CD4⁺ T cells (H), exposed CD4⁺ T cells (E), infected CD4⁺ T cells (I) free virus (V), and CTL cells (C). Additionally, the model includes thirteen parameters that describe the model transmissions. According to the model transmission rates shown in Figure 1, healthy CD4⁺ T cells are produced at a rate of λ, die with a natural death rate of µ₁H, and become infected by the virus at a transmission rate of k1HVH+V. Exposed CD4⁺ T cells die by HIV viruses at a rate of d₁E, die with a natural death at a rate of µ₁E, increase at a rate of k1HVH+V, and become infected at a rate of k₂E. Infected CD4⁺ T cells increase at a rate of k₂E, die by HIV viruses at a rate of d₂I, and die with a natural death rate µ₁I. The infected cells are killed by the CTL response at a rate d₃CI. In addition, CTL cells expand in response to a viral antigen derived from infected cells at a rate of γCI, decrease in the absence of antigenic stimulation at a rate of d₄C and die with a natural death rate of µ₂C. Finally, free HIV viruses are produced by infected cells at a rate of αI, decay at a rate of d₅V, and die with a natural death rate of µ₃C.

The dynamics of an HIV infection disease with a CTL response are expressed by the following nonlinear system of differential equations:

with the initial model states H(0) = H₀, E(0) = E₀, I(0) = I₀, D(0) = D₀, C(0) = C₀ and V(0) = V₀. The model compartments and transmission rates with their biological meaning and estimate values are given in Table 1.

3.   Next generation matrix for infection diseases

An essential epidemiological concept used to measure the potential for an HIV infectious disease transmission, is the basic reproduction number, denoted by R₀ [1], [19]. In this work, we use the next-generation matrix approach to calculate the basic reproduction number, R₀, for the developed compartmental model. Recently, this method has been the most widely used to calculate this threshold parameter. This approach can be used when a model compartment includes both asymptomatic and symptomatic cells. This section explains the compartmental model for infection transmission. Based on the next-generation matrix approach, a compartment is referred to as being infected if there are infected individuals within it. Individuals in this class demonstrate both asymptomatic and symptomatic behaviors. Assume that there are s components in an infectious disease model (x₁, x₂, ..., x_s).

Consider a system of ordinary differential equations for an infectious disease model as follows:

Assume that x=(x1,x2,⋯,xs)T are compartments; these compartments can be divided into two classes: m infected compartments (x1,x2,⋯,xm) and n non-infected (healthy) compartments (xm+1,xm+2,⋯,xs), where s = m + n. This means that

Let $\mathcal{F}_i(x)$ represent the new infection rate in the i^th compartment, and assume the following:

Where Vi− is the rate at which people are transferred out of the i^th compartment, and Vi+ is the rate at which people are transferred into the i^th compartment (see Figure 2).

Furthermore, we denote the infected compartments by xI∈Rm and non-infected compartments by xN∈Rn. Then, the model Eq 3.1 takes the following form:

Now, we determine the matrices $\mathcal{H}$ and $\mathcal{W}$ with components as follows:

Where $\mathcal{W}$ is a non-singular matrix, $\mathcal{H}$ is non-negative, and the free equilibrium of Eq 3.1 is E₀. The spectral radius, or dominant eigenvalue, of the next generation matrix $\mathcal{H}\mathcal{W}^{-1}$ yields the basic reproduction number R₀, as demonstrated by the following expression:

Where ρ represents the spectral radius. The readers can see more details about the idea of the next generation method in [20]–[22].

Understanding the dynamic behavior of infectious disease models and determining the model key parameters are greatly aided by the next generation approach of computing R₀.

4.   Disease equilibrium E: The basic reproduction number R₀

Two kinds of equilibrium points exist for system 2.1 when we take the right-hand side and equate it to zero. The free disease equilibrium (FDE) point is provided below:

additionally, the endemic equilibrium point is provided below:

where

The basic reproduction number (R₀) of system 2.1 for E₀ is calculated as follows:

Let

and

Then,

At the infection-free equilibrium E₀ and using Eq 3.5, we have the following:

where η=µ1+d1+k2, σ=d2+µ1, β=µ2+d4 and ρ=µ3+d5.

Therefore, the next generation matrix can be given as follows:

Thus, $\mathcal{H} \mathcal{W}^{-1}$ has three eigenvalues:

The basic reproduction number is given by the spectral radius of $\mathcal{H} \mathcal{W}^{-1}$, which takes the following form:

In addition, for system 2.1 and using the endemic equilibrium point E1=(H*,E*,I*,C*,V*), the matrices $\mathcal{H}$ and $\mathcal{W}$ are defined as follows

Thus, the next generation matrix can be given as follows:

Therefore, HW⁻¹ for E₁ has three eigenvalues:

The basic reproduction number of the endemic equilibrium point is given as follows:

5.   Local stability of the disease-free equilibrium point

In this section, we study the local stability of the disease-free equilibrium point E0=(λµ1,0,0,0,0) by calculating the Jacobian matrix of system 2.1, which is provided as follows:

Theorem 5.1. The disease-free-equilibrium, E₀ is locally asymptotically stable if the basic reproduction number R₀ < 1 and unstable if R₀ > 1. The value of R₀ is determined using Eq 4.1.

Proof. The local stability of the disease-free-equilibrium can be determined from the Jacobian matrix. In addition, the Jacobian matrix for the disease-free equilibrium is given by the following:

Therefore, the eigenvalues are diagonal elements λ1=−µ1<0, λ2=−η<0, λ3=−σ<0, λ4=−β<0, and λ5=−ρ+k1k1αησ<0 if ρ>k1k1αησ when R₀ < 1 and λ₅ > 0 when R₀ > 1.

Thus, λ₅ < 0 if R₀ < 1 and λ₅ > 0 if R₀ > 1. Clearly, it can be seen that when R₀ < 1, all the eigenvalues have a negative real part; therefore, the disease-free equilibrium, E₀, is locally asymptotically stable and unstable when R₀ > 1.

6.   Elasticity of R₀

In the context of epidemiological models, elasticity refers to how sensitive the model outcomes are to changes in particular parameters. Understanding the relative change of the model parameter responsible for transmission might help select the most effective control techniques. Accordingly, R₀ is associated with disease transmissions, and the sensitivity indicates which parameters have a significant influence on R₀. Let us consider a transmissible disease model with p parameters (ω_j) for j = 1, 2, ..., p. According to the model parameters, ω_j, the sensitivity index of R₀ is ∂R0∂ωj. The elasticity index, also known as the normalized sensitivity index, is an additional measure that calculates the relative change of R₀ with regard to ω_j, represented by ΥωjR0, which is defined as follows:

Although the value of the elasticity index establishes the proportional importance of the parameter, the sign indicates whether R₀ grows (positive sign) or decreases (negative sign) with the parameter.

The elasticity index for each parameters of Model 2.1 is as follows:

We can use the baseline parameter values given in Table 1 to compute the model elasticity.

The sensitivity index of R₀ (elasticity) with regard to the model parameters is determined using MATLAB codes, and the findings are shown in Table 2 and Figure 3. Such indices may have a positive sign such as (k₁, k₂ and α), or negative sign such as (µ₁, µ₃, d₁, d₂ and d₅). These indicate either a direct or indirect relationship between R₀ and the model parameters. For example, positive values indicate that the value of R₀ increases when these parameters increase, and the virus spreads more quickly. Furthermore, a strong positive association is found for the parameters k₁ and α, indicating that R₀ is extremely sensitive to the transmission rate between healthy CD4⁺ T cells and infected CD4⁺ cells. Additionally, a negative sign indicates that raising these parameters can reduce the value of the basic reproduction number. Eventually, we may conclude that the parameters λ, µ₂, d₃, d₄, and γ have a smaller impact on the spreading virus.

7.   Sensitivity analysis

Sensitivity analyses have been extensively used in systems biology models during the few last decades. This method can be applied into infectious disease models to determine the most sensitive parameters and variables. Local and global sensitivities are the two most common types of sensitivities used to study biological systems [23]. Let us consider a transmissible disease model with p parameters (ω_j) for j = 1, 2, ..., p, and s compartments (x_i) for i = 1, 2, ..., s. The differential equation system that represents the model balancing equations appears as follows:

where $x \in \mathbb{R}^s$ and $\omega \in \mathbb{R}^p$.

A local sensitivity is the changes in the state variables x_i regarding the parameters ω_j. From a mathematical perspective, the first order derivatives represent the time-dependent sensitivities of x_i with regard to each parameter value:

A direct sensitivity analysis, commonly known as a “forward sensitivity analysis,” is an additional technique to compute the derivatives. It is possible to solve ordinary differential equations for the sensitivity coefficients using this method. The differential equations represent the state variables x_i's first order derivatives over time with regard to the model parameters ω_j.

The Jacobian matrix for performing a sensitivity analysis is described as follows, using the following chain rule of the system:

Where $W_{\omega_j}, \mathcal{J}$, and $\mathcal{S}$ are defined by the following:

The input parameter ω_j and the starting state of the output variables x_i determine the system's initial circumstances 7.5.

Three well-known methods (non-normalization, half-normalization, and full-normalization) can be used to calculate the local sensitivity. We provide ςωjxi as a measured sensitivity coefficient for each x_i with respect to each parameter k_j in order to compute the sensitivity analysis.

The full-normalization sensitivities can be found as follows:

The half-normalization sensitivities are defined as follows:

The non-normalization sensitivities are given as follows:

This approach can be applied in simulations to calculate the local sensitivities for full-, half-, and non-normalization in order to demonstrate how the input parameters affect the output variables. The model-sensitive analysis is a significant aspect that can be investigated in the context of HIV. Additionally, MATLAB's SimBiology Toolbox can be used to calculate this approach. Additional information, improvements, and sensitivity analysis applications in the field of systems biology can be found in [20], [23]–[25].

As mentioned before, our aim in this section is to determine the crucial model parameters based on the sensitivity analysis; more thorough and accurate work is required for the described HIV model presented in Eq 2.1. We use the estimated values of the model parameters and variables listed in Table 1 for the computational cases. Computational findings are produced by the use of SimBiology Toolbox for MATLAB in three different types: full-normalizations, half-normalizations, and non-normalizations, see Figures 4–6. Interestingly, the results assist us to understand the model better and allow us identify which model parameters are crucial.

8.   Numerical results and discussions

Using the Matlab software, we can calculate the elasticity of R₀ with respect to the model parameters. The results are shown in the third column of Table 2 and are represented as a bar graph (see Figure 3). Table 2 shows that all parameters have an impact on the basic reproduction number R₀, which can be positive such as (k₁, k₂ and α), or negative such as (µ₁, µ₃, d₁, d₂ and d₅). These show whether R₀ and the model parameters are related either directly or indirectly. For instance, positive results show that the virus spreads more quickly when these parameters are increased and the value of R₀ grows. In addition, R₀ is highly sensitive to the rate of transmission between healthy CD4⁺ T cells and infected CD4⁺ cells; A negative sign suggests that increasing these parameters can lower the fundamental reproduction number. Finally, we may draw the conclusion that the parameters λ, µ₂, d₃, d₄, and γ have little effect on the virus's propagation. The simulated findings shown indicate that the majority of the model classes demonstrate a sensitivity with respect to the critical parameters. The suggested techniques of local sensitivity used here indicate how the input parameters affect the output variables. For example, through the first full-normalization technique of the sensitivity analysis, the variable V (those have free HIV particles in the blood) is more sensitive to d₅ compared to the other relevant parameters; see Figure 4). Alternatively, in Figure 5, variables E (HIV-exposed CD4⁺ T cells) and I (HIV-infected CD4⁺ T cells) are more sensitive to the parameter k₁ as compared to others parameters. Furthermore, Figure 6 shows that variable H (uninfected CD4⁺ T cells) is sensitive to the parameter µ₁, and variable C (CTL cells) is sensitive to the model parameters µ₂ and d₅, as they have a smaller sensitivity compared to the other parameters.

Accordingly, there are some interesting points based on the influence of each included parameter over the model states using all methods of the local sensitivity analysis. First, even when all model parameters are included, the non-normalization approach is unable to accurately determine the effect of the model parameters on the model variables; see Figure 6. In order to understand the effect of all parameters, we can exclude some parameters such as k₁, µ₁, µ₂, and d₄ from the model computational simulations; see Figure 6b. Consequently, as shown in Figure 5, the computational results based on the local sensitivity for a half-normalization simulation indicates that most of the parameters significantly influence the model variable sensitivities.

Consequently, the same set of parameter values given in Table 1 is used to show the dynamics of the model Equation 2.1 for the disease-free and endemic steady states. The basic reproduction number R₀ < 1 if λ = 10, k₁ = 1.1, d₅ = 0.9, and α = 45. The model dynamics for the given disease-free equilibrium converge to the disease-free equilibrium point, as illustrated graphically in Figure 8. However, R₀ > 1 when λ = 10, k₁ = 0.05, d₅ = 2, and α = 25, and Figure 7 shows the model dynamics to an endemic steady state.

9.   Conclusions

This study presented a mathematical model that captured the dynamics of HIV transmission and the impact of antiretroviral therapy (ART). The model divides the population into five compartments: healthy individuals (H), exposed individuals (E), infected individuals (I), individuals receiving ART (C), and free virus particles (V). The model equations and assumptions were biologically grounded and reflected the current understanding of HIV epidemiology. The key findings of this study are as follows:

(1) Disease equilibria and basic reproduction number (R₀): The analysis of the disease-free and endemic equilibria showed that the stability of these equilibria is determined by the basic reproduction number R₀. When R₀ < 1, the disease-free equilibrium is stable, and the infection dies out. When R₀ > 1, the endemic equilibrium is stable, and the infection persists in the population.

(2) Sensitivity analysis: Using various normalization techniques, the sensitivity analysis identified the key parameters that have the greatest impact on the value of R₀ and the spread of HIV. Parameters such as the transmission rate (k₁), the rate of progression from exposed to infected (k₂), and the infectivity of treated individuals (α) were found to be the most influential.

These results align with the current understanding of HIV epidemiology and the factors that drive disease transmission. The sensitivity analysis provided valuable insights into the relative importance of different model parameters, which can inform the development of targeted interventions and control strategies. Overall, the manuscript presents a comprehensive mathematical modeling approach to study the dynamics of HIV transmission and the effects of ART. The findings contribute to the existing knowledge on HIV epidemiology and can aid in the design of more effective strategies for managing the HIV epidemic.

Based on the findings, we recommend the following future research directions:

(1) Incorporating spatial heterogeneity: The current model assumes a well-mixed population. Extending the model to incorporate spatial heterogeneity, such as the distribution of high-risk and low-risk regions, could provide insights into the localized dynamics of HIV transmission and the impact of targeted interventions.

(2) Exploring optimal control strategies: Building upon the sensitivity analysis, future studies could focus on developing optimal control strategies that target the most influential parameters. This could involve the use of an optimal control theory to determine the best combination of interventions, such as increased testing, an improved treatment adherence, and behavioral change campaigns.

(3) Validating model with real-world data: The model parameters used in this study were obtained from the literature. Validating the model against real-world epidemiological data, such as the incidence and prevalence rates, would strengthen the model's predictive capabilities and its applicability to specific geographic regions or populations.

(4) Investigating the impact of comorbidities: The current model does not consider the potential impact of comorbidities, such as TB or hepatitis, on the dynamics of HIV transmission and disease progression. Incorporating these factors could lead to a more comprehensive understanding of the challenges faced in HIV management.

(5) Exploring the role of social determinants: Future research could investigate the influence of social determinants of health, such as socioeconomic status, education, and access to healthcare, on the spread of HIV and the effectiveness of intervention strategies.

By addressing these future research directions, the scientific community can further enhance the understanding of HIV epidemiology and develop more targeted and effective strategies to control the HIV pandemic.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.