Mathematical Model of HIV/AIDS Considering Sexual Preferences Under Antiretroviral Therapy, a Case Study in San Juan de Pasto, Colombia.

While several studies on human immunodeficiency virus (HIV)/acquired immunodeficiency syndrome (AIDS) in the homosexual and heterosexual population have demonstrated substantial advantages in controlling HIV transmission in these groups, the overall benefits of the models with a bisexual population and initiation of antiretroviral therapy have not had enough attention in dynamic modeling. Thus, we used a mathematical model based on studying the impacts of bisexual behavior in a global community developed in the PhD thesis work of Espitia (2021). The model is governed by a nonlinear ordinary differential equation system, the parameters of which are calibrated with data from the cumulative cases of HIV infection and AIDS reported in San Juan de Pasto in 2019. Our model estimations show which parameters are the most influential and how to modulate them to decrease the HIV infection.

INTRODUCTION

In the beginning of the 80's, the Center for Disease Control and Prevention (CDC) caught the world's attention because of its presentation of the acquired immunodeficiency syndrome (AIDS) caused by the human immunodeficiency virus (HIV). The scientific community optimistically believed that an effective vaccine would soon be developed. However, almost 40 years later, neither a vaccine nor a successful treatment could stop the transmission of this virus.

An alternative treatment called “antiretroviral therapy” (ART) is currently the best option for a long-lasting viral suppression and subsequently for a reduction of mortality. Nowadays, available drugs do not completely eradicate the HIV infection, but they decrease virus replication and consequently reduce the morbidity and mortality, which are generally associated with AIDS.

According to the United Nations Program on HIV/AIDS (UNAIDS, 2021), since the beginning of the epidemic and until 2019, 75.7 million people contracted HIV infection, and 32.7 million people died from AIDS. At the end of 2019, 38 million people were living with HIV worldwide, and 1.7 million people contracted HIV infection in that year; besides, more than 690,000 people died from AIDS-related illnesses.

In 2014, UNAIDS proposed the Fast-Track Treatment target called “90-90-90 intervention,” in which by 2020, it was required that 90% of individuals living with HIV infection worldwide should be diagnosed, 90% of these should have enrolled in ART, and 90% of individuals receiving ART should have a suppressed viral load. In Colombia in 2019, ART coverage was 85.63%, and 72.09% had undetectable viral load, <50 copies/mm3. The department of Nariño had a coverage of 90.26% of ART with 72.97% of undetectability (CAC 2019).

One of the first investigations in HIV/AIDS mathematical modeling was carried out by Anderson et al. (1986) who presented two epidemiological deterministic mathematical models: the first one is sexual transmission in susceptible infectious people living with AIDS and noninfectious zeropositive; the second one is sexual transmission between susceptible individuals, two classes of infectious, people living with AIDS and, noninfectious zeropositive. Since then, a large number of variations and improvements have been made using current mathematical and medicinal tools. These models can help predict the future behavior of the disease and offer a better understanding of its epidemiological patterns.

The following review is a list of the main articles about HIV/AIDS from 2012 to 2019: Kaur et al. (2012), explored the spread of HIV in children, adults, and seniors taking into consideration horizontal and vertical transmissions. Sun et al. (2013) investigated the epidemic among men who have sex with men enrolled in public health programs such as condom use and ART, dividing the infected population according to CD4+ cell counts in the blood. Bhunu et al. (2014) studied the dynamics of prostitution in the transmission of the disease based on an epidemiological mathematical model in susceptible infected individuals living with AIDS, in prostitutes and nonprostitutes receiving ART and considering viral load in the force of infection. Bhunu (2015) claimed that the care for individuals living with HIV/AIDS is more than the provision of ART, and the homelessness on HIV/AIDS transmission is described through a deterministic mathematical model.

The PhD thesis of Afassinou (2016) introduced a compartmental mathematical model incorporating pre-exposure prophylaxis in the susceptible population and infected individuals. A research about injectable drug users was presented by Yang et al. (2017); they studied female sex workers and senior male clients with heterosexual transmission using an epidemiological deterministic model. Omondi et al. (2018) also approached the dynamics in the sexually mature age group aged 15 years in a population of susceptible and infected individuals considering CD4+ T cell counts and individuals enrolled in ART. Finally, Omondi et al. (2019) studied the HIV transmission dynamics using a compartmental model. They considered heterosexual transmission in two age groups in Kenya; young adults aged 15–24 and adults aged 25 and older, each population was subdivided into susceptible and infected people receiving ART treatment.

The previous literature review shows a pattern of HIV/AIDS mathematical models, which is often considered in a heterosexual population and a few times in heterosexuals with homosexuals separately. However, populations consisting of heterosexuals, bisexuals, and homosexuals at the same time are not included; it means that bisexual behavior is not studied in this mathematical modeling. Thus, we investigate the influence of sexual intercourse among these groups on the transmission of HIV and subsequent AIDS always under ART.

The proposed model is an original result of an extensive literary review of researches about the following: first, epidemiological HIV/AIDS mathematical models, second, biological researches about HIV behavior and later AIDS, and third, psychological aspects regarding homosexual and bisexual contacts. Thus, our model considers one of its main contributions the homosexual, heterosexual, and bisexual contacts with different sexual preferences. Sensibility of parameters in the model was analyzed according to the Latin Hypercube Sampling Method with the aim of knowing which ones have greater sensitiveness to the variations on the initial conditions.

BACKGROUND

A contextualization of sexual contact and its transmission mechanisms is provided in this section, it is provided a contextualization about sexual contact and its transmission mechanisms. Heterosexual and homosexual contact is contemplated as the main route of HIV/AIDS contagion, and bisexual contact in heterosexuals and homosexual individuals as another way to get infected. In this regard, the transmission between exclusively homosexual women will not be taken into account due to the absence of the exchange of body fluids, which is the principal mechanism of sexual transmission. However, there are some cases of transmission mainly due to the exchange of sex toys, for more information, see works such as Lloyd and May (1996), Mastro and De Vincenzi (1996), Kwakwa and Ghobrial (2003), and Chan et al. (2014).

Sexual contact between homosexual and heterosexual men is justified by Rosario et al. (2006) and Thompson and Morgan (2008), where the authors state that: “there exist individuals that change their sexual behavior depending on the situation or at different stages in their life. A possibly common and transient example of situational sexuality is the person who self identifies as heterosexual, but will sexually interact with a member of the same sex when lacking other opportunities. Less transient but also possibly common, a person who self identifies as homosexual or lesbian (either at the time, or later) may sexually interact with a member of the opposite sex if a same-sex relationship seems unfeasible.” As a result, in this model, sexual interaction between homosexual men and heterosexual men will be assumed.

The CDC estimates that HIV rates in men having sex with men (MSM) are higher than heterosexual contacts. In part, these differences reflect the fact that an individual MSM can engage in both insertive and receptive sexual roles (versatility), while exclusively heterosexual men and women each engage in only one of these roles (Beloqui, 2008; Simon et al., 2008; Glick et al., 2012).

METHODS

A compartmental model of HIV/AIDS is represented by a system of ordinary differential equations (ODE) formed by observing the flow of individuals from one compartment to the others. As an example, the target population is divided into three compartments, where the number of individuals changes with respect to time; the susceptible class comprises individuals who were not exposed to infection, the infected class consists of infected individuals who infect others, and the class of individuals living with AIDS involves individuals who have already developed the syndrome as a result of HIV infection. The infection occurs due to the interaction of susceptible individuals with infected ones, and depends on the initial hypothesis of the model.

The presented model was made under the advice of the HIV/AIDS infectious disease specialist Dr. Alexandre Naime Barbosa, who evaluated the dynamics hypotheses; this model is part of the doctoral thesis work in applied mathematics at the Institute of Mathematics, Statistics and Scientific Computing in State University of Campinas by Espitia (2021). The model assumed that the only way to transmit the HIV virus is through sexual relationships between heterosexuals, homosexuals, and bisexual individuals; moreover, the effect of ART among infected people is studied making use of a deterministic mathematical model of nonlinear ODE.

The total population is divided in homosexual men and heterosexual men and women; therefore, into eight classes as described in Table 1.

Table 1.

Description of State Variables in the Model

Population	Description
Sh(t)	Susceptible homosexual men
Ih(t)	Untreated infected homosexual men
Sw(t)	Susceptible women
Iw(t)	Untreated infected women
Sm(t)	Susceptible heterosexual men
Im(t)	Untreated infected heterosexual men
T(t)	Treated with antiretrovirals
A(t)	People living with AIDS

AIDS, acquired immunodeficiency syndrome.

Figure 1 represents the transmission dynamics between the three studied sexual preferences. Each vertex of the triangle represents one population, and the sides of the triangle denote the different forms of transmission between the populations involved. To begin, the exclusive transmission among homosexual men is illustrated by the upper circular arrow in purple labeled . Then, the transmission between homosexual and heterosexual men and the transmission between homosexual men and women are represented by blue lines labeled and , respectively. Finally, heterosexual transmission between men and women is in the red line represented by . The direction of the arrows represents the sense of the analyzed contagion; nonetheless, contagions can biologically occur in all directions, in this work, we assume only those directions shown in Figure 1, because with the acceptance of all the directions, besides demanding a cumbersome mathematical treatment also causing the model to bring about qualitatively different equilibrium dynamics.

FIG. 1.

Triangle transmission between homosexual men, heterosexual men, and women.

Consequently, it is assumed that the only form of contagion among homosexuals is among themselves, and that heterosexual people become infected due to contact with homosexual men or heterosexuals of the opposite sex. Thus, the blue lines have only one direction, while the red line between heterosexual men and women has two.

Mathematical model

To be able to model these contacts, the infection forces of which are shown in Figure 1, these mathematical expressions must be able to represent the way in which susceptible individuals are becoming infected. represents the product between the probability of contagion and the rate of the number of sexual partners, . In these infection forces, it was assumed that individuals under treatment and individuals living with AIDS under treatment are the least infectious because, due to ART, their viral load is undetectable, and therefore, the contagion is almost null.

Hypotheses assumed for modeling

H1 Constant recruitment in all susceptible classes is assumed. H2 Sexual transmission in discordant couples is considered. H3 Homosexual individuals get infected among themselves. The HIV transmission in a female population occurs through sexual relationships with infected heterosexual men or through sexual contact with infected homosexual men. Susceptible heterosexual men can get infected by infected women or infected homosexual men. H4 No gender difference is taken into account for treated individuals as well as for individuals living with AIDS. H5 Individuals living with AIDS could be treated or untreated, noting that an individual who developed AIDS during a hospital treatment will be diagnosed and enrolled in ART. H6 Both natural mortality in all classes and induced mortality for individuals living with AIDS are considered.

The parameters used in this model are as follows: represents constant recruitment in susceptible individuals. is the proportion of men, while is that of heterosexual men. Thus, recruitments in susceptible populations of homosexual men, women, and heterosexual men are , , and , respectively. p is the proportion of initially treated individuals. denotes the natural mortality rate. d is the parameter for induced disease mortality rate. represents the AIDS development rate in treated individuals due to ART resistance or treatment failure. is the departure rate from all compartments of infected individuals. means probability of sexual transmission and denotes sexual partner rate. The diagram and ODE system is shown in Eq. (2), and its dynamics is governed by a nonlinear and autonomous system of ODE represented in Eq. (2).

The basic reproduction number, denoted by R0, represents the expected number of secondary cases produced by a typically infected individual in a completely susceptible population. If , then one infected individual can generate less than one new infection over the course of the infectious period; therefore, the disease cannot grow. However, if , then one infected individual produces more than one new infection, and thus, the disease invades the population and it is considered an epidemic. For important features of R0, we recommend Holland (2007). The next-generation method exposed in Van den Driessche (2017) is used, the results of which indicate that this basic reproduction number for the ODE model presented in Figure 2 is given by the following:

FIG. 2.

Diagram and ODE system for triangle transmission model. ODE, ordinary differential equations.

A forward bifurcation is present in this model as shown in the diagram of Figure 3. Two dotted vertical lines have been illustrated; the black one represents , while the green one marks .

FIG. 3.

Bifurcation diagram of untreated infected heterosexual men, in function of transmission probability in homosexual contact.

Local and global stability of stationary points was demonstrated in the PhD thesis work (Espitia, 2021). Thus, we claim that stability does not depend on initial conditions.

DATA AND PARAMETERS FOR HIV/AIDS MODEL IN SAN JUAN DE PASTO

Data adjustment for HIV infection in San Juan de Pasto was made using government bibliographical sources such as data from the Departmental Health Institute of Nariño (IDSN, 2019) and the National Administrative Department of Statistics (DANE, 2019), as well as scientific references. Table 2 summarizes active cases in San Juan de Pasto from 1989 to 2019.

Table 2.

Active Cases Reported in San Juan de Pasto from 1989 to 2019

Year	1989	1990	1991	1992	1993	1994	1995	1996	1997	1998	1999
Individuals	1	2	3	6	3	16	11	8	8	21	21
Year	2000	2001	2002	2003	2004	2005	2006	2007	2008	2009	2010
Individuals	24	25	41	44	50	54	75	67	35	36	35
Year	2011	2012	2013	2014	2015	2016	2017	2018	2019		Total
Individuals	30	39	61	54	73	67	101	101	121		1233

Source: Departmental Institute of Health of Nariño.

With the aim of knowing the dynamical behavior and the effects of varying the parameters in the ODE system, we carry out numerical simulation following Runge-Kutta's mathematical method of order 4 determining the baseline values of parameters known from the literature that correspond to available experimental data and biological facts. The unknown parameters were estimated intuitively based on the general behavior of the country and department data. The parameter values are approximated as follows: , which is defined as the recruitment in the susceptible population, is estimated using data from the National Census of Population and Housing (CNPV 2019); with regard to San Juan de Pasto, this value is given by 333 individuals per year. The natural death rate, , is estimated by DANE (2019) and based on life expectancy, which is reported as 77.3 years of age, and thus, . The AIDS-related death rate, d, based on the life expectancy of individuals living with AIDS, as stated by Nettleman (2019), is on average 3 years, and therefore, .

Sources such as Fernández (2019) and Choi et al. (2020) were useful in estimating parameters such as sexual partners, . ART parameters were adapted from CAC (2019). All parameter values are listed in Table 3, and the initial conditions were assumed by DANE (2019) and UNAIDS (2021); those are shown in equation (2).

Table 3.

Parameters of the Human Immunodeficiency Virus/Acquired Immunodeficiency Syndrome Model for Case Study in San Juan de Pasto

Parameters	Value	Unit	Source
Ψ	333	Individuals/year	CNPV (2019)
θ	0.48	Nondimension	DANE (2019)
γ	0.92	Nondimension	DANE (2019), Williams (2019)
p	0.90	Nondimension	CAC (2019), Cornejo et al. (2020)
μ	0.0129	1/year	DANE (2019)
d	0.3333	1/year	Nettleman (2019)
δ	0.018	1/year	Apenteng and Ismail (2019)
α	0.3333	1/year	Akudibillah et al. (2018)
βs	0.02	Nondimension	Ostadzad et al. (2016)
βh	0.44	Nondimension	UNAIDS (2021)
βhw	0.018	Nondimension	Assumed
βhm	0.25	Nondimension	Assumed
cs	4	1/year	Fernández (2019)
ch	7	1/year	Fernández (2019)
Chw	2	1/year	Fernández (2019)
Chm	1	1/year	Fernández (2019)

RESULTS

We present one main scenario and three possible additional scenarios, which, in addition to allowing a better understanding of the dynamics of HIV infection, are a fundamental tool in the design and analysis of public health policies. A common behavior in all scenarios is the decline of susceptible populations due to infected people, and thus, we only present changes in infected populations.

Scenarios

Main scenario

The parameter values in Table 3 correspond to a basic reproduction number in homosexuals and a basic reproduction number in heterosexuals . Figure 4 is a numerical simulation of populations for a period of 10 years, noticing that susceptible populations decline over time.

FIG. 4.

On the left, the main scenario using parameters in table description of parameters. On the right, disease-free equilibrium when the basic reproduction number in homosexuals is less than one, both in logarithm scale.

The disease-free equilibrium is defined as the point at which no disease is present in the population, which is generally represented in epidemiological models when the basic reproduction number is less than one. Getting to this point means controlling the epidemic and eradicating the infection in the long term, which is the main goal of every epidemiological model. As a result, the right side of Figure 4 shows the infected population when the parameter of homosexual partners is , and thus, ; it follows that the population reaches a disease-free equilibrium, and, as a consequence, the infected population is eradicated.

Scenario 1: modifying the basic reproduction number in heterosexuals,

We simulated the effect of considering several heterosexual partners, c (keeping the other parameters as in Table 3), and thus, is modified. The first modification is when , in this way and are assumed, noting that the infected heterosexual men population declines in approximately 4.5 years, and the infected women population declines in about 6 years, and thus heterosexual men decrease faster than infected women. The second modification is when or 24; in this case we have , verifying that the populations of infected heterosexual men, women, and treated individuals grow over time. Figure 5 depicts this behavior.

FIG. 5.

Scenario 1: modifying heterosexual partners’ rate.

Scenario 2: modifying the departure rate of infected individuals

The departure rate of infected individuals has the greatest variability of populations; the parameter is the most sensitive* one for the model. Susceptible populations decline, while infected individuals, treated populations, and people living with AIDS have changes, as shown in Figure 6.

FIG. 6.

Scenario 2: modifying departure rate from infected individuals.

Scenario 3: modifying sexual partners between homosexual men with women, and homosexual men with heterosexual men

Figure 7 shows the effect of modifying sexual partners in a bisexual contact for infected heterosexual men and women.

FIG. 7.

Scenario 3: modifying bisexual partner.

DISCUSSIONS

The dynamics of the HIV/AIDS epidemic, to a large extent, depends on changes of the basic reproduction number for homosexuals, . This was also evidenced by modifying several parameters in the scenarios above. The basic reproduction number in heterosexuals, , is less influential because in the Background section, it was mentioned that the probability of HIV infection in homosexuals is higher than in heterosexuals. In addition, investigations such as Cohen et al. (2011) and Del Romero et al. (2016) permit to conclude that the number of sexual partners of homosexuals is larger than the heterosexuals, thus, . It suggests that HIV infection can be controlled or eliminated from the community if control programs are directed toward reducing to values less than one. The model showed the persistence of the disease when .

The dynamics of HIV/AIDS is, in general, too complex to allow intuitive predictions, and requires both the support of mathematical modeling for quantitative assessments and the understanding of the system functioning; furthermore, one of the most difficult tasks of mathematical modeling was the model parameters. Moreover, the proposed HIV/AIDS model tries to be as approximate as possible using real parameter values for the behavior in San Juan de Pasto. The emphasis is not on the accuracy of the scenarios, but on the actions that can be taken now as a result of seeing the state of the epidemic in the future; these actions can involve, among other things, the prevention of new infections, provision of ART, and educational campaigns such as the reduction of the number of sexual partners or use of condom for self-protection.

The model presented in this work should be treated with circumspection due to the assumptions made surrounding the estimation of the model parameters. Nonetheless, the model provides useful insights about the dynamics of the epidemic.

CONCLUSIONS

A deterministic mathematical model of nonlinear differential equations was performed to model the dynamics of HIV/AIDS in an adult population considering sexual preferences such as the following: exclusive homosexual men contact, contact between homosexual men and women, contact between homosexual and heterosexual men, and finally, heterosexual contact, according to scientific references that validate these interactions. Thus, we consider three susceptible populations: homosexual men, women and heterosexual men, and consequently infected individuals under ART and people living with AIDS.

HIV/AIDS epidemic data were collected from the Departmental Health Institute of Nariño (IDSN, 2019) and the National Administrative Department of Statistics (DANE, 2019) to model HIV infection in San Juan de Pasto, Colombia. Then, we performed numerical simulations according to Runge-Kutta's mathematical method of order 4, where it was possible to know the effects of manipulating the parameters corresponding to the model, and in this way, it is expected to contribute to the design of public health policies in favor of the region.

Based on numerical simulations with city-specific parameters and considering the three different scenarios, it can be given three referent conclusions related to the dynamical behavior about parameters involved in the model. First, the number of sexual partners and the departure rate in infected individuals are the most sensitive parameters in the dynamics of infection because infected populations have big changes when these values are slightly modified. Second, scenario 2 shows that the best way to decrease heterosexual contagion in San Juan de Pasto is to increase the departure rate in infected individuals . Third, scenarios 1 and 3 show that increasing the number of sexual partners in homosexual or heterosexual contact implies an increase in the number of infections.

Thus, the best way to reduce contagion and consequently to reach a disease-free equilibrium is mainly decreasing the number of homosexual partners since they are the most affected population by the virus and the most likely to get infected and to spread the infection. Increasing the departure rate of infected individuals leads to a decrease in infected heterosexual men and infected women, but it is not enough to prevent and curb the rate of contagion.

With the population parameters of San Juan de Pasto, we performed several numerical simulations modifying parameters that return the basic reproduction number greater than or less than one, it suggests that: when and , there is a general decline in the HIV infection over a period of few years, but the infection persists. Thus, it can be concluded that the most important observation from our findings is that in the population of the city, there is a short-term rise of HIV infection, in which there is a significant increase of new HIV infections followed by a substantial decline in the generation of new infections.