# Mathematical Modeling and Analysis of an Alcohol Drinking Model with the Influence of Alcohol Treatment Centers

Authors: Bouchaib Khajji; Abderrahim Labzai; Omar Balatif; Mostafa Rachik
Journal: International Journal of Mathematics and Mathematical Sciences
 (2020)

---

## Abstract

In this paper, we present a continuous mathematical modelPMHTrTpQ of alcohol drinking with the influence of private and public addiction treatment centers. We study the dynamical behavior of this model and we discuss the basic properties of the system and determine its basic reproduction number R0. We also study the sensitivity analysis of model parameters to know the parameters that have a high impact on the reproduction number R0. The stability analysis of the model shows that the system is locally as well as globally asymptotically stable at drinking-free equilibrium E0 when R0≤1. When R0>1, drinking present equilibrium E∗ exists and the system is locally as well as globally asymptotically stable at alcohol present equilibrium E∗.


---

## Body

## 1. Introduction

Alcoholism is a social phenomenon that affects all social classes and people with different educational levels and age groups and can have an impact on many aspects of life. It is a chronic disorder that causes the person to drink uncontrollably; as a result, it can destroy people’s relationships and make them less productive as well as harm their physical and mental health. Long-term abuse can also cause severe damage to the brain and liver. Alcohol has several social, economic, and health effects on the individual and society as a whole including heart disease, waste of money, poverty, crime, family disintegration, and liver disease and also causes ulcers, diabetes complications, sexual problems, birth defects, bone loss, vision problems, increased risk of cancer, and suppressed immune function which need treatment in hospitals and private or public addiction treatment centers. Alcoholism remains one of the most frequent and important topics discussed by the world community due to its dire consequences [1].According to the World Health Organization (WHO), in its report “World Status Report on Alcohol and Health 2018,” alcohol addiction causes the death of three million people worldwide every year, which represents one in every 20 deaths. For example, alcohol consumption kills more than AIDS, tuberculosis, and violence combined [1].In 2016, alcohol has killed around 3 million people worldwide, compared with 3.3 million in 2012. Three quarters of the deaths are in men. According to these figures, there is a clear reduction in the number of deaths caused by alcohol addiction. In the European region and in the Americas region, the harmful consequences including illness and injury resulted from the use of alcohol are detrimental in comparison with other regions. Several diseases, more than 200, are related to alcohol consumption.28% of the 3 million deaths ascribed to alcohol were linked to traffic accidents, violence, suicides, and other violent acts, 21% to digestive disorders, and 19% to cardiovascular diseases. The other cases of deaths are resulted from infectious diseases, cancers, mental disorders, and other health problems [1].Many studies and research in social, medical, and political sciences have focused on this topic and other related topics ([2–4] and the references cited therein). But the mathematical studies and research on this topic are still limited and most of them have focused on the statistical aspect of the phenomenon ([5–8]).Mathematical models can be used to analyze the spread of infectious diseases or the social behavior of individuals [9–20]. As regard to drinking, several different mathematical models have been formulated and studied to help in reducing the number of drinkers [17, 21–24]. Agrawal et al. [22] developed a nonlinear SHTR mathematical model of alcohol abuse with a nonlinear incidence rate. The stability analysis of the model that they proposed shows that the system is locally asymptotically stable at alcohol-free equilibrium E0 when R0≤1. When R0>1, alcohol equilibrium E∗ exists and the system becomes locally asymptotically stable at E∗ and E0 becomes unstable. Manthey et al. [23] proposed a mathematical model to study the dynamics of campus drinking as an epidemiological model. According to their results, the reproductive numbers were not sufficient to predict whether drinking behavior would persist on campus and that the pattern of recruiting new members played a significant role in the reduction of campus alcohol problems. Sharma et al. [17] developed a mathematical model of alcohol abuse and discussed the existence, local, global stability of drinking-free, endemic equilibria, and sensitivity analysis of R0. They demonstrated that backward bifurcation can occur when R0=1. Sharma et al. [24] developed a two-stage (four compartments) model for youths with serious drinking problems and their treatment having admitted and having not admitted. They also analyzed the stability of all the equilibria. Adu [21] used a nonlinear SHTR mathematical model to study the dynamics of drinking epidemic. He divided the population into four compartments: nondrinkers S, heavy drinkers H, drinkers in treatment T, and recovered drinkers R. He discussed the existence and stability of drinking-free and endemic equilibria.Besides these works, we will study the analysis of a mathematical alcohol modelPMHTrTpQ with the following additions:(i)
CompartmentTp represents the number of the poor heavy drinkers who join public addiction treatment centers, which may not have sophisticated equipment and do not provide high quality treatment services in comparison with the private ones, especially in the developing countries, and those individuals who do not have the financial capacity to go to private centers(ii)
CompartmentTr represents the number of the rich heavy drinkers who join private addiction treatment centers that encompass special facilities and provide advanced treatment for such individuals with sufficient financial capacity(iii)
The death rate induced by the heavy drinkersδThroughout this research, we present a nonlinear mathematical alcohol modelPMHTrTpQ with private and public addiction treatment centers that describes the dynamics of drinkers’ classes with the influence of the private and public addiction treatment centers. Hence, the awareness programs and the addiction treatment centers of alcohol drinking will alert the susceptible individuals and heavy drinkers so they decline drinking. We assume that the heavy drinkers can recover from heavy alcohol drinking due to media, health reasons, joining addiction treatment centers, prohibition, and hike tax on alcohol beverages. A fraction of heavy drinkers’ population will join the private addiction treatment centers, whereas the rest will join the public addiction treatment centers. The goal is to analyze the impact of the addiction treatment centers and awareness programs on the drinkers. We also analyze the effect of addiction treatment centers on system stability. Our results show that addiction treatment centers and awareness programs are an effective measure in controlling heavy drinking. By using Routh–Hurwitz criteria and constructing Lyapunov functions, the local and the global stability of alcohol-free equilibrium and alcohol equilibrium are obtained. We also study the sensitivity analysis of model parameters to know the parameters that have a high impact on the reproduction number R0. The drinker classes of this model are divided into six compartments: potential drinkers P, moderate drinkers M, heavy drinkers H, the rich heavy drinkers who join private addiction treatment centers Tr, the poor heavy drinkers who join public addiction treatment centers Tp, and quitters of drinking Q.The rest of the paper is organized as follows. In section2, the formulation of the model and some basic properties are derived. In section 3, equilibria of the proposed model are obtained and their stability discussed. In section 4, the global stability of equilibrium point is discussed. The problem of sensitivity’s parameters is discussed in section 5. Some numerical simulations and discussions are given in section 6. Lastly we give the conclusion of the whole paper.


## 2. Model Formulation

### 2.1. Model Description

We propose a continuous modelPMHTrTpQ to describe the interaction between drinkers classes. The population is divided into six compartments: potential drinkers Pt, moderate drinkers Mt, heavy drinkers Ht, rich heavy drinkers Trt, poor heavy drinkers Tpt, and quitters of drinking Qt.The graphical representation of the proposed model is shown in Figure1.Figure 1
Diagram of drinkers’ class.We consider the following system of six nonlinear differential equations:(1)dPtdt=b−β1PMN−μP,dMtdt=β1PMN−μ+β2M,dHtdt=β2M−μ+δ+α1+α2+α3H,dTrtdt=α1H−γ1+μTr,dTptdt=α2H−γ2+μTp,dQtdt=γ1Tr+γ2Tp+α3H−μQ,where P0≥0,M0≥0,H0≥0,Tr0≥0,Tp0≥0, and Q0≥0 are the given initial states.


#### 2.1.1. Potential DrinkersP

(2)
d
P
t
d
t
=
b
−
β
1
P
M
N
−
μ
P
.The potential drinkersPt represents individuals whose age is over adolescence and adulthood and may become drinkers. This compartment is increased by the recruitment rate denoted by b and decreased by an effective contact with the moderate drinkers at β1 rate and natural death μ. It is assumed that potential drinkers can acquire drinking behavior and can become moderate drinkers through effective contact with moderate drinkers in some social occasions such as weddings, celebrating graduation ceremonies, week-end parties, and end of the year celebration. In other words, it is assumed that the acquisition of a drinking behavior is analogous to acquiring disease infection.


#### 2.1.2. Moderate DrinkersM

(3)
d
M
t
d
t
=
β
1
P
M
N
−
μ
+
β
2
M
.The compartmentM is composed of moderate drinkers who can control their consumption during some events and occasions or in a way that is unapparent to their social environment. This category of drinkers does not face any problems or negative consequences. Friends or family do not complain about their intake of alcohol. A moderate drinker does not think about drinking very often or often feel a need to drink. They may go out to have a few drinks and are able to handle their alcohol consumption without experiencing a loss of control. Alcohol does not dominate their thoughts and they do not need to set limits when they drink. They are not prone to extreme mood swings, fighting, or being violent. It is increased by potential drinkers who turn to be moderate drinkers at β1 rate. This compartment is decreased when moderate drinkers become heavy drinkers at a rate β2 and also by natural death at rate μ.


#### 2.1.3. Heavy DrinkersH

(4)
d
H
t
d
t
=
β
2
M
−
μ
+
δ
+
α
1
+
α
2
+
α
3
H
.The compartmentH is composed of heavy drinkers suffering from addiction to alcohol. When an individual becomes an alcoholic, they face a great difficulty to control or set limits for their consumption. The majority of alcoholics begin as potential drinkers and then turn to moderate drinkers. Alcohol seems to exert a control on the alcoholic’s life. Their job, their family, social circle, and health are all endangered. Despite these negative consequences, the alcoholic is unable to quit drinking. The alcoholics may begin to disclaim that they have a problem; this disclaim can make it even more difficult for the person to get help. Alcohol addiction is considered to be a disease; it changes chemicals in the addict’s brain and has made alcohol the most important thing in their life. At the time a person is an alcoholic, they will usually need to get help at a rehab to overcome their addiction. This compartment becomes larger as the number of heavy drinkers increases by the rate β2 and decreases when some of them give up drinking at a rate α3 as well as when they join treatment centers whether private or public at rates α1 and α2. In addition, this compartment decreases by natural death μ and due to deaths caused by diseases resulted from excessive alcohol intake at a rate δ.


#### 2.1.4. Rich Heavy DrinkersTr

(5)
d
T
r
t
d
t
=
α
1
H
−
γ
1
+
μ
T
r
.The compartmentTr contains the number of heavy drinkers who take advantage of their financial potentials to join private treatment centers of alcohol addiction that are very often well equipped and provide good and quality services. This compartment is increased by the rate α1 and decreased by the rates γ1 and μ.


#### 2.1.5. Poor Heavy DrinkersTp

(6)
d
T
p
t
d
t
=
α
2
H
−
γ
2
+
μ
T
p
.The compartmentTp represents the number of heavy drinkers who join public treatment centers of alcohol addiction which may not provide advanced treatment and that are marked by a shortage of equipment and low-quality services, especially in the developing countries, and it also contains individuals who do not have the financial capacity to join the private centers. This compartment is increased by the rate α2 and decreased by the rates γ2 and μ.


#### 2.1.6. Quitters of DrinkingQ

(7)
d
Q
t
d
t
=
γ
1
T
r
+
γ
2
T
p
+
α
3
H
−
μ
Q
.The compartmentQ encompasses the individuals who quit drinking. It is increased with the recruitment of individuals who have been treated in treatment centers of alcohol addiction at rates α1 and α2. It is also increases at the rate α3 of those who quit alcohol without resort to treatment centers and decreases at the rate μ due to natural deaths.The total population size at timet is denoted by Nt with Nt=Pt+Mt+Ht+Trt+Tpt+Qt.


### 2.2. Basic Properties

#### 2.2.1. Invariant Region

It is necessary to prove that all solutions of system (1) with positive initial data will remain positive for all times t>0. This will be established by the following lemma.Lemma 1.
All feasible solutionPt,Mt,Ht,Trt,Tpt, and Qt of system equation (1) are bounded by the region:(8)Ω=P,M,H,Tr,Tp,Q∈IR+6:P+M+H+Tr+Tp+Q≤bμ.Proof.
From system equation (1),(9)dNtdt=dPtdt+dMtdt+dHtdt+dTrtdt+dTptdt+dQtdt,dNtdt=b−μNt−δH,implies that(10)dNtdt≤b−μNt,and it follows that(11)Nt≤bμ+N0e−μt,where N0 is the initial value of total number of people; thus,(12)limt⟶+∞supNt≤bμ,then(13)Pt+Mt+Ht+Trt+Tpt+Qt≤bμ.
Hence, for the analysis of model (1), we get the region which is given by the set(14)Ω=P,M,H,Tr,Tp,Q∈IR+6:P+M+H+Tr+Tp+Q≤bμ,which is a positively invariant set for (1), so we only need to consider the dynamics of system (1) on the set Ω nonnegative of solutions.


#### 2.2.2. Positivity of Solutions of the Model

Theorem 1.
IfP0≥0,M0≥0,H0≥0,Tr0≥0,Tp0≥0, and Q0≥0, then the solution of system equation (1) Pt,Mt,Ht,Trt,Tpt, and Qt are positive for all t>0.
Proof. From the first equation of system (1), we have(15)dPtdt=b−β1MtN+μPt.It is assumed that(16)At=β1MtN+μ.We multiply equation (15) by exp∫0tAsds and we find(17)dPtdt∗exp∫0tAsds=b−AtPt∗exp∫0tAsds,which implies that(18)dPtdt∗exp∫0tAsds+AtPt∗exp∫0tAsds=b∗exp∫0tAsds.Therefore,(19)ddtPt∗exp∫0tAsds=b∗exp∫0tAsds.Taking integral with respect tos from 0 to t, we obtain(20)Pt∗exp∫0tAsds−P0=b∗∫0texp∫0wAsdsdw.Multiplying equation (20) by exp−∫0tAsds, we obtain(21)Pt−P0∗exp−∫0tAsds=b∗exp−∫0tAsds∗∫0texp∫0wAsdsdw.Then, we can obtain(22)Pt=P0∗exp−∫0tAsds+b∗exp−∫0tAsds∗∫0texp∫0wAsdsdw≥0.So, the solutionPt is positive.Similarly, from the second equation of system (1), we have(23)Mt≥M0exp−∫0tBsds≥0,where(24)Bt=β1PtN−β2+μ.Similarly, from the third, fourth, and fifth equations of system (1), we have(25)Ht≥H0exp−μ+δ+α1+α2+α3t≥0,Trt≥Tr0exp−μ+γ1t≥0,Tpt≥Tp0exp−μ+γ2t≥0,Qt≥Q0exp−μt≥0.Therefore, we can see thatPt>0, Mt>0, Ht>0, Trt>0, Tpt, and Qt>0 ∀t≥0, and this completes the proof.□The first three equations in system (1) are independent of the variables Tr,Tp, and Q. Hence, the dynamics of equation system (1) is equivalent to the dynamics of the equation system:(26)dPtdt=b−β1PMN−μP,dMtdt=β1PMN−β2M−μM,dHtdt=β2M−μ+δ+α1+α2+α3H.


## 2.1. Model Description

We propose a continuous modelPMHTrTpQ to describe the interaction between drinkers classes. The population is divided into six compartments: potential drinkers Pt, moderate drinkers Mt, heavy drinkers Ht, rich heavy drinkers Trt, poor heavy drinkers Tpt, and quitters of drinking Qt.The graphical representation of the proposed model is shown in Figure1.Figure 1
Diagram of drinkers’ class.We consider the following system of six nonlinear differential equations:(1)dPtdt=b−β1PMN−μP,dMtdt=β1PMN−μ+β2M,dHtdt=β2M−μ+δ+α1+α2+α3H,dTrtdt=α1H−γ1+μTr,dTptdt=α2H−γ2+μTp,dQtdt=γ1Tr+γ2Tp+α3H−μQ,where P0≥0,M0≥0,H0≥0,Tr0≥0,Tp0≥0, and Q0≥0 are the given initial states.


### 2.1.1. Potential DrinkersP

(2)
d
P
t
d
t
=
b
−
β
1
P
M
N
−
μ
P
.The potential drinkersPt represents individuals whose age is over adolescence and adulthood and may become drinkers. This compartment is increased by the recruitment rate denoted by b and decreased by an effective contact with the moderate drinkers at β1 rate and natural death μ. It is assumed that potential drinkers can acquire drinking behavior and can become moderate drinkers through effective contact with moderate drinkers in some social occasions such as weddings, celebrating graduation ceremonies, week-end parties, and end of the year celebration. In other words, it is assumed that the acquisition of a drinking behavior is analogous to acquiring disease infection.


### 2.1.2. Moderate DrinkersM

(3)
d
M
t
d
t
=
β
1
P
M
N
−
μ
+
β
2
M
.The compartmentM is composed of moderate drinkers who can control their consumption during some events and occasions or in a way that is unapparent to their social environment. This category of drinkers does not face any problems or negative consequences. Friends or family do not complain about their intake of alcohol. A moderate drinker does not think about drinking very often or often feel a need to drink. They may go out to have a few drinks and are able to handle their alcohol consumption without experiencing a loss of control. Alcohol does not dominate their thoughts and they do not need to set limits when they drink. They are not prone to extreme mood swings, fighting, or being violent. It is increased by potential drinkers who turn to be moderate drinkers at β1 rate. This compartment is decreased when moderate drinkers become heavy drinkers at a rate β2 and also by natural death at rate μ.


### 2.1.3. Heavy DrinkersH

(4)
d
H
t
d
t
=
β
2
M
−
μ
+
δ
+
α
1
+
α
2
+
α
3
H
.The compartmentH is composed of heavy drinkers suffering from addiction to alcohol. When an individual becomes an alcoholic, they face a great difficulty to control or set limits for their consumption. The majority of alcoholics begin as potential drinkers and then turn to moderate drinkers. Alcohol seems to exert a control on the alcoholic’s life. Their job, their family, social circle, and health are all endangered. Despite these negative consequences, the alcoholic is unable to quit drinking. The alcoholics may begin to disclaim that they have a problem; this disclaim can make it even more difficult for the person to get help. Alcohol addiction is considered to be a disease; it changes chemicals in the addict’s brain and has made alcohol the most important thing in their life. At the time a person is an alcoholic, they will usually need to get help at a rehab to overcome their addiction. This compartment becomes larger as the number of heavy drinkers increases by the rate β2 and decreases when some of them give up drinking at a rate α3 as well as when they join treatment centers whether private or public at rates α1 and α2. In addition, this compartment decreases by natural death μ and due to deaths caused by diseases resulted from excessive alcohol intake at a rate δ.


### 2.1.4. Rich Heavy DrinkersTr

(5)
d
T
r
t
d
t
=
α
1
H
−
γ
1
+
μ
T
r
.The compartmentTr contains the number of heavy drinkers who take advantage of their financial potentials to join private treatment centers of alcohol addiction that are very often well equipped and provide good and quality services. This compartment is increased by the rate α1 and decreased by the rates γ1 and μ.


### 2.1.5. Poor Heavy DrinkersTp

(6)
d
T
p
t
d
t
=
α
2
H
−
γ
2
+
μ
T
p
.The compartmentTp represents the number of heavy drinkers who join public treatment centers of alcohol addiction which may not provide advanced treatment and that are marked by a shortage of equipment and low-quality services, especially in the developing countries, and it also contains individuals who do not have the financial capacity to join the private centers. This compartment is increased by the rate α2 and decreased by the rates γ2 and μ.


### 2.1.6. Quitters of DrinkingQ

(7)
d
Q
t
d
t
=
γ
1
T
r
+
γ
2
T
p
+
α
3
H
−
μ
Q
.The compartmentQ encompasses the individuals who quit drinking. It is increased with the recruitment of individuals who have been treated in treatment centers of alcohol addiction at rates α1 and α2. It is also increases at the rate α3 of those who quit alcohol without resort to treatment centers and decreases at the rate μ due to natural deaths.The total population size at timet is denoted by Nt with Nt=Pt+Mt+Ht+Trt+Tpt+Qt.


## 2.1.1. Potential DrinkersP

(2)
d
P
t
d
t
=
b
−
β
1
P
M
N
−
μ
P
.The potential drinkersPt represents individuals whose age is over adolescence and adulthood and may become drinkers. This compartment is increased by the recruitment rate denoted by b and decreased by an effective contact with the moderate drinkers at β1 rate and natural death μ. It is assumed that potential drinkers can acquire drinking behavior and can become moderate drinkers through effective contact with moderate drinkers in some social occasions such as weddings, celebrating graduation ceremonies, week-end parties, and end of the year celebration. In other words, it is assumed that the acquisition of a drinking behavior is analogous to acquiring disease infection.


## 2.1.2. Moderate DrinkersM

(3)
d
M
t
d
t
=
β
1
P
M
N
−
μ
+
β
2
M
.The compartmentM is composed of moderate drinkers who can control their consumption during some events and occasions or in a way that is unapparent to their social environment. This category of drinkers does not face any problems or negative consequences. Friends or family do not complain about their intake of alcohol. A moderate drinker does not think about drinking very often or often feel a need to drink. They may go out to have a few drinks and are able to handle their alcohol consumption without experiencing a loss of control. Alcohol does not dominate their thoughts and they do not need to set limits when they drink. They are not prone to extreme mood swings, fighting, or being violent. It is increased by potential drinkers who turn to be moderate drinkers at β1 rate. This compartment is decreased when moderate drinkers become heavy drinkers at a rate β2 and also by natural death at rate μ.


## 2.1.3. Heavy DrinkersH

(4)
d
H
t
d
t
=
β
2
M
−
μ
+
δ
+
α
1
+
α
2
+
α
3
H
.The compartmentH is composed of heavy drinkers suffering from addiction to alcohol. When an individual becomes an alcoholic, they face a great difficulty to control or set limits for their consumption. The majority of alcoholics begin as potential drinkers and then turn to moderate drinkers. Alcohol seems to exert a control on the alcoholic’s life. Their job, their family, social circle, and health are all endangered. Despite these negative consequences, the alcoholic is unable to quit drinking. The alcoholics may begin to disclaim that they have a problem; this disclaim can make it even more difficult for the person to get help. Alcohol addiction is considered to be a disease; it changes chemicals in the addict’s brain and has made alcohol the most important thing in their life. At the time a person is an alcoholic, they will usually need to get help at a rehab to overcome their addiction. This compartment becomes larger as the number of heavy drinkers increases by the rate β2 and decreases when some of them give up drinking at a rate α3 as well as when they join treatment centers whether private or public at rates α1 and α2. In addition, this compartment decreases by natural death μ and due to deaths caused by diseases resulted from excessive alcohol intake at a rate δ.


## 2.1.4. Rich Heavy DrinkersTr

(5)
d
T
r
t
d
t
=
α
1
H
−
γ
1
+
μ
T
r
.The compartmentTr contains the number of heavy drinkers who take advantage of their financial potentials to join private treatment centers of alcohol addiction that are very often well equipped and provide good and quality services. This compartment is increased by the rate α1 and decreased by the rates γ1 and μ.


## 2.1.5. Poor Heavy DrinkersTp

(6)
d
T
p
t
d
t
=
α
2
H
−
γ
2
+
μ
T
p
.The compartmentTp represents the number of heavy drinkers who join public treatment centers of alcohol addiction which may not provide advanced treatment and that are marked by a shortage of equipment and low-quality services, especially in the developing countries, and it also contains individuals who do not have the financial capacity to join the private centers. This compartment is increased by the rate α2 and decreased by the rates γ2 and μ.


## 2.1.6. Quitters of DrinkingQ

(7)
d
Q
t
d
t
=
γ
1
T
r
+
γ
2
T
p
+
α
3
H
−
μ
Q
.The compartmentQ encompasses the individuals who quit drinking. It is increased with the recruitment of individuals who have been treated in treatment centers of alcohol addiction at rates α1 and α2. It is also increases at the rate α3 of those who quit alcohol without resort to treatment centers and decreases at the rate μ due to natural deaths.The total population size at timet is denoted by Nt with Nt=Pt+Mt+Ht+Trt+Tpt+Qt.


## 2.2. Basic Properties

### 2.2.1. Invariant Region

It is necessary to prove that all solutions of system (1) with positive initial data will remain positive for all times t>0. This will be established by the following lemma.Lemma 1.
All feasible solutionPt,Mt,Ht,Trt,Tpt, and Qt of system equation (1) are bounded by the region:(8)Ω=P,M,H,Tr,Tp,Q∈IR+6:P+M+H+Tr+Tp+Q≤bμ.Proof.
From system equation (1),(9)dNtdt=dPtdt+dMtdt+dHtdt+dTrtdt+dTptdt+dQtdt,dNtdt=b−μNt−δH,implies that(10)dNtdt≤b−μNt,and it follows that(11)Nt≤bμ+N0e−μt,where N0 is the initial value of total number of people; thus,(12)limt⟶+∞supNt≤bμ,then(13)Pt+Mt+Ht+Trt+Tpt+Qt≤bμ.
Hence, for the analysis of model (1), we get the region which is given by the set(14)Ω=P,M,H,Tr,Tp,Q∈IR+6:P+M+H+Tr+Tp+Q≤bμ,which is a positively invariant set for (1), so we only need to consider the dynamics of system (1) on the set Ω nonnegative of solutions.


### 2.2.2. Positivity of Solutions of the Model

Theorem 1.
IfP0≥0,M0≥0,H0≥0,Tr0≥0,Tp0≥0, and Q0≥0, then the solution of system equation (1) Pt,Mt,Ht,Trt,Tpt, and Qt are positive for all t>0.
Proof. From the first equation of system (1), we have(15)dPtdt=b−β1MtN+μPt.It is assumed that(16)At=β1MtN+μ.We multiply equation (15) by exp∫0tAsds and we find(17)dPtdt∗exp∫0tAsds=b−AtPt∗exp∫0tAsds,which implies that(18)dPtdt∗exp∫0tAsds+AtPt∗exp∫0tAsds=b∗exp∫0tAsds.Therefore,(19)ddtPt∗exp∫0tAsds=b∗exp∫0tAsds.Taking integral with respect tos from 0 to t, we obtain(20)Pt∗exp∫0tAsds−P0=b∗∫0texp∫0wAsdsdw.Multiplying equation (20) by exp−∫0tAsds, we obtain(21)Pt−P0∗exp−∫0tAsds=b∗exp−∫0tAsds∗∫0texp∫0wAsdsdw.Then, we can obtain(22)Pt=P0∗exp−∫0tAsds+b∗exp−∫0tAsds∗∫0texp∫0wAsdsdw≥0.So, the solutionPt is positive.Similarly, from the second equation of system (1), we have(23)Mt≥M0exp−∫0tBsds≥0,where(24)Bt=β1PtN−β2+μ.Similarly, from the third, fourth, and fifth equations of system (1), we have(25)Ht≥H0exp−μ+δ+α1+α2+α3t≥0,Trt≥Tr0exp−μ+γ1t≥0,Tpt≥Tp0exp−μ+γ2t≥0,Qt≥Q0exp−μt≥0.Therefore, we can see thatPt>0, Mt>0, Ht>0, Trt>0, Tpt, and Qt>0 ∀t≥0, and this completes the proof.□The first three equations in system (1) are independent of the variables Tr,Tp, and Q. Hence, the dynamics of equation system (1) is equivalent to the dynamics of the equation system:(26)dPtdt=b−β1PMN−μP,dMtdt=β1PMN−β2M−μM,dHtdt=β2M−μ+δ+α1+α2+α3H.


## 2.2.1. Invariant Region

It is necessary to prove that all solutions of system (1) with positive initial data will remain positive for all times t>0. This will be established by the following lemma.Lemma 1.
All feasible solutionPt,Mt,Ht,Trt,Tpt, and Qt of system equation (1) are bounded by the region:(8)Ω=P,M,H,Tr,Tp,Q∈IR+6:P+M+H+Tr+Tp+Q≤bμ.Proof.
From system equation (1),(9)dNtdt=dPtdt+dMtdt+dHtdt+dTrtdt+dTptdt+dQtdt,dNtdt=b−μNt−δH,implies that(10)dNtdt≤b−μNt,and it follows that(11)Nt≤bμ+N0e−μt,where N0 is the initial value of total number of people; thus,(12)limt⟶+∞supNt≤bμ,then(13)Pt+Mt+Ht+Trt+Tpt+Qt≤bμ.
Hence, for the analysis of model (1), we get the region which is given by the set(14)Ω=P,M,H,Tr,Tp,Q∈IR+6:P+M+H+Tr+Tp+Q≤bμ,which is a positively invariant set for (1), so we only need to consider the dynamics of system (1) on the set Ω nonnegative of solutions.


## 2.2.2. Positivity of Solutions of the Model

Theorem 1.
IfP0≥0,M0≥0,H0≥0,Tr0≥0,Tp0≥0, and Q0≥0, then the solution of system equation (1) Pt,Mt,Ht,Trt,Tpt, and Qt are positive for all t>0.
Proof. From the first equation of system (1), we have(15)dPtdt=b−β1MtN+μPt.It is assumed that(16)At=β1MtN+μ.We multiply equation (15) by exp∫0tAsds and we find(17)dPtdt∗exp∫0tAsds=b−AtPt∗exp∫0tAsds,which implies that(18)dPtdt∗exp∫0tAsds+AtPt∗exp∫0tAsds=b∗exp∫0tAsds.Therefore,(19)ddtPt∗exp∫0tAsds=b∗exp∫0tAsds.Taking integral with respect tos from 0 to t, we obtain(20)Pt∗exp∫0tAsds−P0=b∗∫0texp∫0wAsdsdw.Multiplying equation (20) by exp−∫0tAsds, we obtain(21)Pt−P0∗exp−∫0tAsds=b∗exp−∫0tAsds∗∫0texp∫0wAsdsdw.Then, we can obtain(22)Pt=P0∗exp−∫0tAsds+b∗exp−∫0tAsds∗∫0texp∫0wAsdsdw≥0.So, the solutionPt is positive.Similarly, from the second equation of system (1), we have(23)Mt≥M0exp−∫0tBsds≥0,where(24)Bt=β1PtN−β2+μ.Similarly, from the third, fourth, and fifth equations of system (1), we have(25)Ht≥H0exp−μ+δ+α1+α2+α3t≥0,Trt≥Tr0exp−μ+γ1t≥0,Tpt≥Tp0exp−μ+γ2t≥0,Qt≥Q0exp−μt≥0.Therefore, we can see thatPt>0, Mt>0, Ht>0, Trt>0, Tpt, and Qt>0 ∀t≥0, and this completes the proof.□The first three equations in system (1) are independent of the variables Tr,Tp, and Q. Hence, the dynamics of equation system (1) is equivalent to the dynamics of the equation system:(26)dPtdt=b−β1PMN−μP,dMtdt=β1PMN−β2M−μM,dHtdt=β2M−μ+δ+α1+α2+α3H.


## 3. Equilibria and Their Stability Analysis

### 3.1. Equilibrium Point

The standard method is used to analyze model (26). In this model, there are two equilibrium points, that is, drinking-free equilibrium point and drinking present equilibrium point. The equilibrium points are found by setting the right-hand side of equations (1)–(3) equal to zero.The drinking-free equilibriumE0b/μ,0,0 is achieved in the absence of drinking M=H=0.The drinking present equilibriumE∗P∗,M∗,H∗ is achieved when drinkers exist M≠0andH≠0, where(27)P∗=bμR0,M∗=bR0−1β1,H∗=bβ2R0−1β1μ+α1+α2+α3+δ,R0=β1μ+β2.R
0 is the basic reproduction number that measures the average number of new drinkers generated by single drinker in a population of potential drinkers. The value of R0 will indicate whether the epidemic could occur or not. The reproduction basic number can be determined by using the next generation matrix method formulated in [25].


### 3.2. Local Stability Analysis

Now we proceed to study the stability behavior of equilibriaE0 and E∗.


#### 3.2.1. The Drinking-Free Equilibrium

In this section, we analyze the local stability of the drinking-free equilibrium.Theorem 2. The drinking-free equilibriumE0b/μ,0,0 of system (26) is asymptotically stable if R0<1 and unstable if R0>1.Proof.
The Jacobian matrix atE is given by(28)JE=−β1MN−μ−β1PN0β1MNβ1PN−β2−μ00β2−μ−α1−α2−α3−δ.
The Jacobian matrix for the drinking-free equilibrium is given by(29)JE0=−μ−β100β1−β2−μ00β2−μ−α1−α2−α3−δ,where P0=b/μ=N.
The characteristic equation of this matrix is given bydetJE0−λI3=0, where I3 is a square identity matrix of order 3.
Therefore, eigenvalues of the characteristic equation ofJE0 are(30)λ1=−μ,λ2=−β2+μ−β1=−μ+β21−β1μ+β2,λ3=−μ+α1+α2+α3+δ,R0=β1μ+β2.
Therefore, all the eigenvalues of the characteristic equationJE0 are clearly real and negative if R0<1.We conclude that the drinking-free equilibrium is locally asymptotically stable ifR0<1 and unstable if R0>1.


#### 3.2.2. Drinking Present Equilibrium

In this section, we analyze the local stability of drinking present equilibrium.To find the drinking present equilibrium of the system of equation (26) we consider dPt/dt=0,dMt/dt=0, and dHt/dt=0 provided that at least one of the infected compartments is nonzero. We evaluate the equilibrium of system (22) by setting the right-hand side of equation of system (26) equal to zero and then solve for P∗,M∗, and H∗.From the second equation in system (26), we have(31)P∗=Nμ+β2β1,P∗=bμR0.From the first equation in system (22), we have(32)M∗=bR0−1β1.Also, equation (3) in system (26) gives(33)H∗=bβ2R0−1β1μ+α1+α2+α3+δ.Let the following theorem analyzes the local stability of the drinking present equilibrium whenR0>1.Theorem 3.
The drinking present equilibriumE∗ is locally asymptotically stable if R0>1 and unstable otherwise.Proof.
We presentE∗P∗,M∗,H∗ as drinking present equilibrium of system (26) and P∗≠0,M∗≠0, and H∗≠0.The Jacobian matrix is(34)JE∗=−β1M∗N−μ−β1P∗N0β1M∗Nβ1P∗N−β2−μ00β2−μ−α1−α2−α3−δ,where(35)P∗=bμR0,M∗=bR0−1β1,H∗=bβ2R0−1β1μ+α1+α2+α3+δ.We see that the characteristic equationPλ of JE∗ has an eigenvalue λ1=−μ+α1+α2+α3+δ whose real part is negative. So, in order to determine the stability of the drinking present equilibrium of model (26), we discuss the roots of the following equation φλ:(36)φλ=λ2+a1λ+a2,where(37)a1=β1M∗N+μ>0,a2=β12P∗NM∗N>0.By Routh–Hurwitz criterion, system (26) is locally asymptotically stable if a1>0 and a2>0.Thus, the drinking present equilibriumE∗ of system (26) is locally asymptotically stable.


## 3.1. Equilibrium Point

The standard method is used to analyze model (26). In this model, there are two equilibrium points, that is, drinking-free equilibrium point and drinking present equilibrium point. The equilibrium points are found by setting the right-hand side of equations (1)–(3) equal to zero.The drinking-free equilibriumE0b/μ,0,0 is achieved in the absence of drinking M=H=0.The drinking present equilibriumE∗P∗,M∗,H∗ is achieved when drinkers exist M≠0andH≠0, where(27)P∗=bμR0,M∗=bR0−1β1,H∗=bβ2R0−1β1μ+α1+α2+α3+δ,R0=β1μ+β2.R
0 is the basic reproduction number that measures the average number of new drinkers generated by single drinker in a population of potential drinkers. The value of R0 will indicate whether the epidemic could occur or not. The reproduction basic number can be determined by using the next generation matrix method formulated in [25].


## 3.2. Local Stability Analysis

Now we proceed to study the stability behavior of equilibriaE0 and E∗.


### 3.2.1. The Drinking-Free Equilibrium

In this section, we analyze the local stability of the drinking-free equilibrium.Theorem 2. The drinking-free equilibriumE0b/μ,0,0 of system (26) is asymptotically stable if R0<1 and unstable if R0>1.Proof.
The Jacobian matrix atE is given by(28)JE=−β1MN−μ−β1PN0β1MNβ1PN−β2−μ00β2−μ−α1−α2−α3−δ.
The Jacobian matrix for the drinking-free equilibrium is given by(29)JE0=−μ−β100β1−β2−μ00β2−μ−α1−α2−α3−δ,where P0=b/μ=N.
The characteristic equation of this matrix is given bydetJE0−λI3=0, where I3 is a square identity matrix of order 3.
Therefore, eigenvalues of the characteristic equation ofJE0 are(30)λ1=−μ,λ2=−β2+μ−β1=−μ+β21−β1μ+β2,λ3=−μ+α1+α2+α3+δ,R0=β1μ+β2.
Therefore, all the eigenvalues of the characteristic equationJE0 are clearly real and negative if R0<1.We conclude that the drinking-free equilibrium is locally asymptotically stable ifR0<1 and unstable if R0>1.


### 3.2.2. Drinking Present Equilibrium

In this section, we analyze the local stability of drinking present equilibrium.To find the drinking present equilibrium of the system of equation (26) we consider dPt/dt=0,dMt/dt=0, and dHt/dt=0 provided that at least one of the infected compartments is nonzero. We evaluate the equilibrium of system (22) by setting the right-hand side of equation of system (26) equal to zero and then solve for P∗,M∗, and H∗.From the second equation in system (26), we have(31)P∗=Nμ+β2β1,P∗=bμR0.From the first equation in system (22), we have(32)M∗=bR0−1β1.Also, equation (3) in system (26) gives(33)H∗=bβ2R0−1β1μ+α1+α2+α3+δ.Let the following theorem analyzes the local stability of the drinking present equilibrium whenR0>1.Theorem 3.
The drinking present equilibriumE∗ is locally asymptotically stable if R0>1 and unstable otherwise.Proof.
We presentE∗P∗,M∗,H∗ as drinking present equilibrium of system (26) and P∗≠0,M∗≠0, and H∗≠0.The Jacobian matrix is(34)JE∗=−β1M∗N−μ−β1P∗N0β1M∗Nβ1P∗N−β2−μ00β2−μ−α1−α2−α3−δ,where(35)P∗=bμR0,M∗=bR0−1β1,H∗=bβ2R0−1β1μ+α1+α2+α3+δ.We see that the characteristic equationPλ of JE∗ has an eigenvalue λ1=−μ+α1+α2+α3+δ whose real part is negative. So, in order to determine the stability of the drinking present equilibrium of model (26), we discuss the roots of the following equation φλ:(36)φλ=λ2+a1λ+a2,where(37)a1=β1M∗N+μ>0,a2=β12P∗NM∗N>0.By Routh–Hurwitz criterion, system (26) is locally asymptotically stable if a1>0 and a2>0.Thus, the drinking present equilibriumE∗ of system (26) is locally asymptotically stable.


## 3.2.1. The Drinking-Free Equilibrium

In this section, we analyze the local stability of the drinking-free equilibrium.Theorem 2. The drinking-free equilibriumE0b/μ,0,0 of system (26) is asymptotically stable if R0<1 and unstable if R0>1.Proof.
The Jacobian matrix atE is given by(28)JE=−β1MN−μ−β1PN0β1MNβ1PN−β2−μ00β2−μ−α1−α2−α3−δ.
The Jacobian matrix for the drinking-free equilibrium is given by(29)JE0=−μ−β100β1−β2−μ00β2−μ−α1−α2−α3−δ,where P0=b/μ=N.
The characteristic equation of this matrix is given bydetJE0−λI3=0, where I3 is a square identity matrix of order 3.
Therefore, eigenvalues of the characteristic equation ofJE0 are(30)λ1=−μ,λ2=−β2+μ−β1=−μ+β21−β1μ+β2,λ3=−μ+α1+α2+α3+δ,R0=β1μ+β2.
Therefore, all the eigenvalues of the characteristic equationJE0 are clearly real and negative if R0<1.We conclude that the drinking-free equilibrium is locally asymptotically stable ifR0<1 and unstable if R0>1.


## 3.2.2. Drinking Present Equilibrium

In this section, we analyze the local stability of drinking present equilibrium.To find the drinking present equilibrium of the system of equation (26) we consider dPt/dt=0,dMt/dt=0, and dHt/dt=0 provided that at least one of the infected compartments is nonzero. We evaluate the equilibrium of system (22) by setting the right-hand side of equation of system (26) equal to zero and then solve for P∗,M∗, and H∗.From the second equation in system (26), we have(31)P∗=Nμ+β2β1,P∗=bμR0.From the first equation in system (22), we have(32)M∗=bR0−1β1.Also, equation (3) in system (26) gives(33)H∗=bβ2R0−1β1μ+α1+α2+α3+δ.Let the following theorem analyzes the local stability of the drinking present equilibrium whenR0>1.Theorem 3.
The drinking present equilibriumE∗ is locally asymptotically stable if R0>1 and unstable otherwise.Proof.
We presentE∗P∗,M∗,H∗ as drinking present equilibrium of system (26) and P∗≠0,M∗≠0, and H∗≠0.The Jacobian matrix is(34)JE∗=−β1M∗N−μ−β1P∗N0β1M∗Nβ1P∗N−β2−μ00β2−μ−α1−α2−α3−δ,where(35)P∗=bμR0,M∗=bR0−1β1,H∗=bβ2R0−1β1μ+α1+α2+α3+δ.We see that the characteristic equationPλ of JE∗ has an eigenvalue λ1=−μ+α1+α2+α3+δ whose real part is negative. So, in order to determine the stability of the drinking present equilibrium of model (26), we discuss the roots of the following equation φλ:(36)φλ=λ2+a1λ+a2,where(37)a1=β1M∗N+μ>0,a2=β12P∗NM∗N>0.By Routh–Hurwitz criterion, system (26) is locally asymptotically stable if a1>0 and a2>0.Thus, the drinking present equilibriumE∗ of system (26) is locally asymptotically stable.


## 4. Global Stability

### 4.1. Global Stabilty of the Alcohol-Free Equilibrium

To show that system (26) is globally asymptotically stable, we use the Lyapunov function theory for both the drinking-free equilibrium and the drinking present equilibrium. First, we present the global stability of the drinking-free equilibrium E0 when R0≤1.Theorem 4.
The drinking-free equilibruimE0 is globally asymptotically stable if R0≤1 and unstable otherwise.Proof.
Consider the following Lyapunov function:(38)V=cM,where c positive constant.
The derivative ofVP,M,Hwith respect to t gives(39)dVdt=cdMdt=cβ1PN−μ+β2M=cμ+β2β1Pμ+β2N−1M≤R0−1M,where c=1/μ+β2 and R0=β1/μ+β2.
So,dV/dt≤0 if R0≤1.
Furthermore,dV/dt=0 if and only if M=0. Hence, by LaSalle’s invariance principle [26], E0 is globally asymptotically stable.


### 4.2. Global Stability of the Drinking-Present Equilibrium

The final result of the global stability ofE∗ in this section is as follows.Theorem 5.
The drinking-present equilibrium pointE∗ is globally asymptotically stable if R0>1.Proof.
Consider the Lyapunov functionV:(40)V:Γ⟶IR,VP,M=c1P−P∗lnPP∗+c2M−M∗lnMM∗,where c1 and c2 are positive constants to be chosen later and Γ=P,M∈Γ/P>0,M>0. Then, the time derivative of the Lyapunov function is given by(41)dVP,Mdt=−bc1P−P∗2PP∗+β1Nc2−c1P−P∗M−M∗.
Then, the time derivative of the Lyapunov function is given by(42)dVP,Mdt=−bc1P−P∗2PP∗+β1Nc2−c1P−P∗M−M∗.
Forc1=c2=1, we have(43)dVP,Mdt=−bP−P∗2PP∗≤0.
Also, we obtain(44)dVP,Mdt=0⇒P=P∗.
Hence, by LaSalle’s invariance principle [26], the drinking present equilibrium point E∗ is globally asymptotically stable on Γ.


## 4.1. Global Stabilty of the Alcohol-Free Equilibrium

To show that system (26) is globally asymptotically stable, we use the Lyapunov function theory for both the drinking-free equilibrium and the drinking present equilibrium. First, we present the global stability of the drinking-free equilibrium E0 when R0≤1.Theorem 4.
The drinking-free equilibruimE0 is globally asymptotically stable if R0≤1 and unstable otherwise.Proof.
Consider the following Lyapunov function:(38)V=cM,where c positive constant.
The derivative ofVP,M,Hwith respect to t gives(39)dVdt=cdMdt=cβ1PN−μ+β2M=cμ+β2β1Pμ+β2N−1M≤R0−1M,where c=1/μ+β2 and R0=β1/μ+β2.
So,dV/dt≤0 if R0≤1.
Furthermore,dV/dt=0 if and only if M=0. Hence, by LaSalle’s invariance principle [26], E0 is globally asymptotically stable.


## 4.2. Global Stability of the Drinking-Present Equilibrium

The final result of the global stability ofE∗ in this section is as follows.Theorem 5.
The drinking-present equilibrium pointE∗ is globally asymptotically stable if R0>1.Proof.
Consider the Lyapunov functionV:(40)V:Γ⟶IR,VP,M=c1P−P∗lnPP∗+c2M−M∗lnMM∗,where c1 and c2 are positive constants to be chosen later and Γ=P,M∈Γ/P>0,M>0. Then, the time derivative of the Lyapunov function is given by(41)dVP,Mdt=−bc1P−P∗2PP∗+β1Nc2−c1P−P∗M−M∗.
Then, the time derivative of the Lyapunov function is given by(42)dVP,Mdt=−bc1P−P∗2PP∗+β1Nc2−c1P−P∗M−M∗.
Forc1=c2=1, we have(43)dVP,Mdt=−bP−P∗2PP∗≤0.
Also, we obtain(44)dVP,Mdt=0⇒P=P∗.
Hence, by LaSalle’s invariance principle [26], the drinking present equilibrium point E∗ is globally asymptotically stable on Γ.


## 5. Sensitivity Analysis ofR0

Sensitivity analysis is commonly used to determine the model robustness to parameter values, that is, to help us know the parameters that have a high impact on the reproduction numberR0.Using the approach in Chitnis et al. [27], we calculate the normalized forward sensitivity indices of R0. Let(45)ϒmℛ0=mR0∗∂ℛ0∂m,denote the sensitivity index of R0 with respect to the parameter m. We obtain(46)ℛ0=β1μ+β2,ϒβ1ℛ0=1,ϒβ2ℛ0=−β2μ+β2,ϒμℛ0=−μμ+β2.From the abovementioned discussion, we observe that the basic reproduction numberℛ0 is most sensitive to changes in β1. If β1 increases ℛ0 will also increase with the same proportion, and if β1 decreases in the same proportion, μ and β2 will have an inversely proportional relationship with ℛ0. So, an increase in any of them will bring about a decrease in ℛ0; however, the size of the decrease will be proportionally smaller. Recall that μ is the natural death rate of the population. Given ℛ0’s sensitivity to β1, it seems sensible to focus efforts on the reduction of β1. In other words, this sensitivity analysis tells us that prevention is better than cure. Efforts to increase prevention are more effective in controlling the spread of habitual drinking than efforts to increase the numbers of individuals accessing treatment (Table 1).Table 1
The academic parameters used for model (1).
Parameter
Description
Sensitivity index
β
1
The effective contact rate
+
1
β
2
Coefficient of transmission theM at H
−
0.32
μa
The natural death rate
−
0.68


## 6. Numerical Simulations

In this section, we illustrate some numerical solutions of model (1) for different values of the parameters. The resolution of system (1) was created using the Gauss-Seidel-like implicit finite-difference method developed by Gumel et al. [28], presented in [29] and denoted the GSS1 method. We use the following different initial values such that P+M+H+Tr+Tp+Q=1000.We use the following parameters in Table2.Table 2
The academic parameters used for model (1).
P
1
0
P
2
0
P
3
0
M
1
0
M
2
0
M
3
0
H
1
0
H
2
0
H
3
0
b
300
600
800
100
300
500
50
150
200
65
T
1
r
0
T
2
r
0
T
3
r
0
T
1
p
0
T
2
p
0
T
3
p
0
Q
1
0
Q
2
0
Q
3
0
γ
1
60
30
10
35
20
10
18
10
6
0.001
μ
β
1
β
2
N
δ
γ
2
α
2
α
3
α
1
0.065
0.02
0.14
1000
0.002
0.002
0.001
0.001
0.001We begin by a graphic representation of the drinking-free equilibriumE0, and we use the same parameters and different initial values given in Table 2, R0=0.37 and R0<1.From these figures, using the different values of initial variablesP0, M0, H0, T0r, T0p, and Q0, we obtained the following remarks (Figure 2):(i)
The number of potential drinkers increases and approaches the numberP0=1000 (see Figure 2(a))(ii)
The number of moderate drinkers decreases and approaches zero (see Figure2(b))(iii)
The number of the heavy drinkers increases at first, after that it decreases and approaches zero (see Figure2(c))(iv)
The number of the poor heavy drinkers decreases and approaches zero (see Figure2(d))(v)
The number of the rich heavy drinkers decreases and approaches zero (see Figure2(e)),(vi)
The number of the quitters of drinking decreases and approaches zero (see Figure2(f))Figure 2
WhenR0 < 1, the drinking-free equilibrium E0 is globally asymptotically stable.
(a)
(b)
(c)
(d)
(e)
(f)Therefore, the solution curves to the equilibriumE0P0,0,0,0,0,0 when R0<1. Hence, model (1) is globally asymptotically stable.We use the following parameters in Table3.Table 3
The academic parameters used for model (1).
P
1
0
P
2
0
P
3
0
M
1
0
M
2
0
M
3
0
H
1
0
H
2
0
H
3
0
b
200
300
600
300
400
500
100
300
500
65
T
1
r
0
T
2
r
0
T
3
r
0
T
1
p
0
T
2
p
0
T
3
p
0
Q
1
0
Q
2
0
Q
3
0
γ
1
60
30
5
35
20
6
18
10
6
0.001
μ
β
1
β
2
N
δ
γ
2
α
2
α
3
α
1
0.065
0.75
0.14
1000
0.002
0.002
0.001
0.001
0.001Also, we begin by a graphic representation of the drinking present equilibriumE∗ and we use the same parameters and different initial values given in Table 3, R0=3.66 and R0>1.From these figures, using the different values of initial variablesP0,M0,H0,T0r,T0p, and Q0, we obtained the following remarks (Figure 3):(i)
The number of potential drinkers increases at first, then it decreases slightly and approaches the valueP∗=273 (see Figure 3(a))(ii)
The number of moderate drinkers decreases rapidly at first, then increases slightly and approaches the valueM∗=230 (see Figure 3(b))(iii)
The number of heavy drinkers increases and approaches the valueH∗=461 (see Figure 3(c))(iv)
The number of poor heavy drinkers decreases and approaches the valueTp∗=6.89 (see Figure 3(d))(v)
The number of rich heavy drinkers decreases and approaches the valueTr∗=7.02 (see Figure 3(e))(vi)
The number of quitters of drinking decreases and approaches the valueQ∗=7.03 (see Figure 3(f))Figure 3
WhenR0 < 1, the drinking-present equilibrium E∗ is globally asymptotically stable.
(a)
(b)
(c)
(d)
(e)
(f)Therefore, the solution curves to the equilibriumE∗P∗,M∗,H∗,Tp∗,Tr∗,Q∗ when R0>1. Hence, model (1) is globally asymptotically stable.


## 7. Conclusion

In this paper, we have presented a continuous mathematical modelPMHTrTpQ of alcohol drinking with the influence of private and public addiction treatment centers and the dynamical behavior of the model is studied. We have found R0=β1/μ+β2 as basic reproduction number of system (1), which helps us to determine the dynamical behavior of the system. We also studied the sensitivity analysis of model parameters to know the parameters that have a high impact on the reproduction number R0. We used the stability analysis theory for nonlinear systems to analyze the mathematical drinking model and to study both the local and global behavior of drinking dynamics. Local asymptotic stability for the drinking-free equilibrium E0 can be obtained, if the threshold quantity R0≤1. On the other hand, if R0>1, then the alcohol present equilibrium E∗ is locally asymptotically stable. A Lyapunov function was used to show global stability of E0. E0 is globally asymptotically stable if R0≤1. Also a Lyapunov function was used to show global stability of E∗. E∗ is globally asymptotically stable if R0>1.



---

*Source: 4903168-2020-02-10.xml*