Estimating the dynamics of the drinking epidemic model with control interventions: A sensitivity analysis

Yasir Nadeem Anjam; Kashaf Ishfaq; Salman Arif Cheema; Sufian Munawar; Najma Saleem; Muhammad Farman

doi:10.1515/nleng-2024-0010

Artikel Open Access

Estimating the dynamics of the drinking epidemic model with control interventions: A sensitivity analysis

Yasir Nadeem Anjam , Kashaf Ishfaq , Salman Arif Cheema , Sufian Munawar , Najma Saleem und Muhammad Farman

Veröffentlicht/Copyright: 12. Juli 2024

Veröffentlicht von

Veröffentlichen auch Sie bei De Gruyter Brill

Manuskript einreichen Informationen für Autor*innen

Aus der Zeitschrift Nonlinear Engineering Band 13 Heft 1

Abstract

This article presents a non-linear mathematical model that captures the dynamics of drinking prevalence within a population. The model is analyzed under an optimal control framework, dividing the total population into four compartments: susceptible, heavy drinker, drinker in treatment, and recovered classes. The model’s validity is affirmed through considerations of positivity, boundedness, reproduction number, stability, and sensitivity analysis. Stability theory is employed to explore both local and global stabilities. Sensitivity analysis identifies parameters with a significant impact on the reproduction number ( R 0 ), with maximum sensitivity observed in parameters related to drinking transmission and transitions from heavy drinking to treatment stages. These parameters exhibit sensitivity indices of ( 0.538 , 1 ) , indicating that a 10% increase in these parameters would result in a ( 5.38 , 1 ) increase in the threshold quantity. The study introduces an optimal control strategy that involves awareness campaigns and treatment as control variables. These controls aim to minimize the number of heavy drinkers while maximizing the number of recovered individuals. Pontryagin’s maximum principle is used to solve optimal control problems. Additionally, the research explores various parametric settings for each compartment, enriching the study environment. The effectiveness of the proposed control scheme is evaluated through rigorous numerical simulations, highlighting its competitive edge. The results, validated using MATLAB simulations, are detailed throughout the article.

Keywords: drinking model; stability; sensitivity; optimal control theory; numerical simulations

1 Introduction

The ever-evolving research literature nominates the association of several physical, psychological, and socio-economic issues with the access use of alcohol. For example, Lees et al. [1] provided an enlightened account documenting the linkages of excessive alcohol and brain functioning ranging from compromised decision-making to memory loss. Similarly, World Health Organization [2,3] endorsed an intriguing association between chronic alcohol use and the development of medical complications, including a higher risk of liver disease, increased odds of cardiovascular issues, heart problems, high blood pressure, weakened immune system, and an increased risk of certain types of cancer. Furthermore, a thought-provoking investigation of Ferrari et al. [4] found vivid links connecting alcohol use with the streaming of negative emotions, such as depression, self-harm, suicidal attitude, anxiety, and alcohol dependence or addiction. Heavy drinking can also impair cognitive function, leading to difficulties with memory, attention, and decision-making. The main factors of alcohol drinking include income levels, education, employment status, and social norms. Research has found that individuals with lower incomes and lower levels of education are more likely to engage in heavy drinking [5]. Unemployment and job instability have also been associated with increased alcohol consumption [6]. Social norms and peer pressure can also influence drinking behaviors, as alcohol consumption is often seen as a way to fit in or cope with stress. The negative impulses of alcohol consumption are known to have drastic effects on more fragile groups of society, i.e., children. Excess drinking in mothers during pregnancy is found to be related to fetal alcohol spectrum disorders in children. Moreover, in a wider perspective, WHO Global Status Report [7] reported an astonishing estimate of 30% of absenteeism associated with alcohol, in Costa Rica. However, in the UK, the share of alcohol usage in fatal accidents is found to be 60%. The degree of burden of alcohol consumption on the health care systems can be witnessed from the documentation of Johnson et al. [8] reporting a 47% increase in the rate of all alcohol-related emergency visits from 2006 to 2014. From a global perspective, the unchecked use of alcohol is found to be related to the instigation of more than 200 medical complexities, including liver cirrhosis, malignancies, and cardiovascular disease [9,10]. The resulting number of deaths is counted as higher as 3.3 million fatalities per annum, contributing to a notable share of 6% of all deaths, globally [3]. More alarmingly, excessive alcohol consumption is becoming an emerging trend in adolescents aged 15–19 years [11]. Almost 46% of the world’s adolescents aged 15–19 years reported having overused alcohol. In 2016, the National Survey on Drug Use and Health reported that 85.6% of people aged 18 and older drank alcohol, in the USA.

Dynamic modeling, employing mathematical frameworks such as differential equations, unveils the temporal evolution of systems, articulating variable interactions for precise predictions and analysis across diverse fields. Mathematical representations unravel patterns, identify parameters, and enhance our understanding of complex systems. The impact of dynamic modeling spans scientific progress and technological innovations, underscoring its significance across disciplines [12]. Chen has made significant contributions to dynamic modeling, particularly in integrated energy systems and natural gas dynamics [13]. His work involves developing robust state estimators and Kalman filter-based approaches to enhance the accuracy and reliability of state predictions [14]. Focusing on complex energy networks and pipeline systems, Chen’s research addresses critical challenges and offers insights for advancing dynamic modeling methodologies in the energy sector.

Mathematical studies serve as invaluable tools for evaluating, testing, and implementing strategies on both short- and long-term scales, particularly in addressing chronic relapsing diseases such as alcoholism. There has been limited research on mathematical modeling methods for social problems such as alcohol and drug use, despite these issues often being referred to as epidemics, in contrast to research on other epidemic problems. However, in recent years, an increasing number of mathematicians have dedicated their efforts to developing mathematical models aimed at understanding and addressing alcohol-related concerns. Khan and Khan [15] proposed a mathematical model for drinking using numerical analysis of fractional order. Chinnadurai et al. [16] focused on employing mathematical modeling to analyze strategies for controlling alcohol consumption, particularly examining their impact on socioeconomically disadvantaged populations. Through the optimal control approach, Imken and Fatmi [17] introduced a novel mathematical model to investigate alcohol consumption within the diabetic population, especially under antidiabetic drug treatments. ur Rahman et al. [18] proposed a fractional mathematical model for drinking under Atangana–Baleanu Caputo derivatives. Bunonyo et al. [19] research focused on using the eigenvalue method to develop a mathematical model for the time-dependent concentration of alcohol in the human bloodstream, contributing to a deeper understanding of alcohol metabolism dynamics and potentially informing strategies for alcohol-related risk assessment and management. Huo and Wang’s nonlinear mathematical model [20] illustrated the impact of awareness campaigns on binge drinking, demonstrating their effectiveness in mitigating alcohol-related issues. Ma et al. [21] tested a mathematical model of alcoholism as a communicable disease, incorporating awareness campaigns and a time delay, using optimal control techniques. Wang et al. [22] proposed and investigated a nonlinear model of alcoholism with optimal control to prevent interactions between susceptible and infected individuals. Sharma et al. [23] established a mathematical model of alcohol consumption, examining the stability and existence of an endemic equilibrium without drinking, and analyzing the sensitivity of R 0 [24,25]. Manthey et al. [26] proposed a mathematical model to study the dynamics of campus drinking as an epidemiological model. Huo et al. [27] proposed a new social epidemic model of alcoholism with media coverage aimed at encouraging people to abstain from drinking. Huo and Song [28] distinguished between heavy drinkers who admit to drinking and those who do not, developing a two-stage model for problematic drinking that accounts for the shift in drinkers’ status from susceptible individuals to admitting drinkers. This article discusses optimal control analysis of disease transmission within a community, as presented in the study by Anjam et al. [29].

In recent times, there has been a growing emphasis on employing optimal control strategies to improve viral transmission management. These applications extend beyond medical contexts and encompass various areas such as policy development, engineering risk assessment, emergency preparedness, and control program refinement [30,31]. Optimal control strategies are designed to identify effective protocols for addressing complex problems, utilizing mathematical techniques such as dynamic programming, Pontryagin’s principle, or model predictive control [32]. Utilizing optimal control techniques in mathematical models offers benefits by enabling the identification of efficient control strategies to achieve desired outcomes, including improved system performance, maximized benefits, or reduced costs [33]. Verma et al. proposed a mathematical model to investigate coronavirus dynamics, considering the impact of lockdown as an epidemiological measure [34]. They also studied a nonlinear smoking model to analyze interactions between smokers and smoking quitters [35]. Moreover, optimal control strategies have been introduced for coronavirus disease 2019, human immunodeficiency virus, and Dengue through mathematical models that capture relationships among these infectious diseases [36]. For further insights, references such as Omame et al. [37] sightsee backward bifurcation and optimal control in a co-infection model for severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) and zika virus, while Omame et al. [38] delve into SARS-CoV-2 and HBV co-dynamics modeling.

Motivated by the importance of the aforementioned issue, this study delves into the application of optimal control theory to devise a new mathematical framework capable of illustrating drinking prevalence trends in a community. The goal of the optimal control scheme is to establish feasible control laws that optimize specific performance metrics. This article introduces a proposed mathematical model for a drinking epidemic, categorizing the entire population into four classes: susceptible drinkers, individuals in treatment, heavy drinkers, and recovered drinkers. Initially, we validate the proposed model by examining fundamental properties such as boundedness, positivity, reproduction number, equilibrium points, stability, and conducting sensitivity analyses. We conduct comprehensive stability analyses at equilibrium points to assess both local and global stabilities, outlining the legitimacy of the model through detailed sensitivity analysis and documenting local and global features. Subsequently, we develop the drinking epidemic model using control variables derived from qualitative optimal control principles, aiming to maximize susceptible and recovered drinkers while minimizing heavy drinkers and optimizing the aforementioned compartments. The optimal control involves two control variables: a social awareness campaign u 1 ( t ) and treatment assistance for reducing the drinking population u 2 ( t ) . Employing Pontryagin’s maximum principle (Pontryagin et al. [39]), we establish a proof for the existence of the developed optimal control problem and analyze it to identify the necessary optimality conditions. The parametric setting covers a range of scenarios.

This article is organized into the following sections. Section 2 concentrates on model analysis, while Section 3 delves into sensitivity analysis. Section 4 discusses the stability aspects of the model under consideration. Section 5 elaborates on the integration of optimal control, emphasizing the implementation of more effective health interventions. Finally, Section 6 summarizes the key findings of this study.

1.1 Model formulation

Understanding the transmission patterns of diseases and creating effective control strategies heavily relies on mathematical models. Therefore, it is essential to focus on delineating the epidemiological facets of the disease and pinpointing crucial, adjustable parameters to bolster disease control measures.

This study delves into exploring a mathematical model aimed at understanding the drinking epidemic within the human population [40]. It highlights several key limitations in the model analysis, including simplification of real-world complexities, potential oversights in parameter selection, reliance on assumed relationships, and the absence of dynamic external factors. These constraints may limit the model’s ability to fully capture the complexities of human behavior and societal dynamics, potentially impacting the generalizability of the findings. Therefore, caution is advised in interpreting the predictions, as real-world applicability may be constrained by the inherent simplifications in the modeling approach. Additionally, certain assumptions were made in developing the proposed drinking epidemic model. The epidemic occurs within a closed environment, where factors such as sex, race, and social status do not influence the likelihood of developing heavy drinking habits. Members interact homogeneously, meaning they have an equal degree of interaction with one another. Heavy drinking behavior can be transmitted to non-drinkers when they come into contact with individuals who regularly consume large amounts of alcohol. Individuals undergoing treatment for their drinking habits may transition through stages of recovery and susceptibility, potentially leading to a return to heavy drinking. Conversely, those who have successfully stopped drinking enter the recovery compartment.

At any time t ≥ 0 , the mathematical model divides the entire population N ( t ) into four compartments, enhancing accuracy and applicability by systematically categorizing individuals based on shared characteristics related to drinking behavior. This structured approach allows for a more nuanced analysis of diverse factors influencing alcohol consumption within specific subgroups. By recognizing and delineating distinct patterns and tendencies, the model can provide targeted insights into the dynamics of drinking behavior, offering a comprehensive understanding of the population’s heterogeneous nature. This segmentation facilitates precision in research, intervention strategies, and policy development, contributing to a more tailored and effective approach to addressing alcohol-related issues at a population level. This division not only improves the model’s accuracy in representing the heterogeneous nature of drinking behavior but also enhances its applicability by providing a framework for understanding how different factors contribute to the overall dynamics of alcohol consumption in a population.

Therefore, the model classes include the susceptible compartment S ( t ) , encompassing individuals who either abstain from alcohol or engage in moderate drinking without compromising their physical health. The heavy drinker compartment H ( t ) represents individuals who engage in binge drinking, resulting in severe adverse effects on their physical health. The drinker in treatment compartment T ( t ) represents individuals undergoing therapeutic interventions, such as medication, following a bout of alcoholism, and the recovered compartment R ( t ) , signifying individuals who successfully recover from alcoholism through treatment and subsequently maintain sustained abstinence from alcohol. This structured approach enables us to account for the varying levels of susceptibility, severity of drinking habits, treatment interventions, and recovery rates within the population. It allows us to model the transitions between these compartments over time, reflecting the real-world complexities of drinking behavior and recovery processes. By incorporating specific parameters such as recruitment rates, transmission rates, death rates, and recovery rates into the model, we can simulate and analyze the interactions and flow patterns between these compartments more accurately. Moreover, the model formulation incorporates specific parameters such as the recruitment rate of the susceptible class b , the transmission rate from the susceptible class to the heavy drinker class α , the transmission rate from the recovered class to the susceptible class η , the natural death rate μ , the drinking-induced death rate in the class of heavy drinkers δ 1 , the drinking-induced death rate in the treatment class δ 2 , the proportion of individuals entering the drinker in treatment class ϕ , and the recovery rate γ . Considering the aforementioned factors, the schematic diagram in Figure 1 depicts the flow pattern of interactions between different states related to drinking behavior.

Figure 1

The schematic diagram elucidates the dynamics of transmission in drinking behavior.

Hence, the total population at any given time t is expressed as:

N ( t ) = S ( t ) + H ( t ) + T ( t ) + R ( t ) .

Therefore, the fundamental model for the drinking epidemic is governed by a system of nonlinear differential equations, which include:

(1) d S d t = b − α S H − μ S + η R , d H d t = α S H − ( μ + δ 1 + ϕ ) H , d T d t = ϕ H − ( μ + δ 2 + γ ) T , d R d t = γ T − ( μ + η ) R ,

subject to the initial conditions,

S ( 0 ) ≥ 0 , H ( 0 ) ≥ 0 , T ( 0 ) ≥ 0 , and R ( 0 ) ≥ 0 .

The specifics of each term that impact the system, as outlined in Model (1), are identified in Table 1.

Table 1

Explanation for each term in the equations of Model (1) is presented

Terms of each equation	Explanation
b	Recruitment rate of the susceptible people
α S H	The rate of transmission from the susceptible class to individuals with heavy drinking habits
η R	Transmission rate of recovered people
μ S	The rate of natural mortality among individuals in the susceptible population.
δ 1 H	The mortality rate caused by drinking among individuals in the heavy drinker population
ϕ H	Proportion of drinking people of heavy drinker class
δ 2 T	The mortality rate attributed to drinking among individuals in the drinker in treatment population
γ T	Recovered rate of drinker in treatment people
μ H	Mortality rate among heavy drinkers
μ T	Natural death rate of drinker in treatment people
μ R	Natural death rate of recovered people

2 Qualitative analysis of the proposed model

In the forthcoming section, our analysis will focus on the mathematical formulation encapsulated in Model (1). Our examination will delve into the aspects of boundedness and positivity within the model. Additionally, we will determine the equilibrium points and compute the basic reproduction number.

2.1 Positivity

Positivity constraints are imposed on both the initial conditions and parameters in the devised model, ensuring that their values remain non-negative or greater than zero throughout the modeling process. This imposition enhances the model’s reliability by minimizing the occurrence of unrealistic outcomes, thereby improving its applicability and making it more representative of real-world scenarios and behaviors. Consequently, we rearrange Model (1) in accordance with these constraints.

ψ ˙ ( t ) = G ( ψ ( t ) ) ,

where ψ ( t ) = ( ψ 1 , ψ 2 , ψ 3 , ψ 4 ) T ≔ ( S , H , T , R ) T , ψ ( 0 ) = ( S ( 0 ) , H ( 0 ) , T ( 0 ) , R ( 0 ) ) T ε ℛ + 4 , and

G ( ψ ) = G 1 ( ψ ) G 2 ( ψ ) G 3 ( ψ ) G 4 ( ψ ) = b − α S H − μ S + η R α S H − ( μ + δ 1 + ϕ ) H ϕ H − ( μ + δ 2 + γ ) T γ T − ( μ + η ) R .

In a stress-free scenario, G i ( ψ ) ∣ ψ i = 0 ≥ 0 for i = 1 , … , 4 . Consistent with the seminal discovery by Nagumo [41], the solution of Model (1) with an initial condition ψ 0 ∈ R + 4 , denoted as ψ ( t ) = ψ ( t ; ψ 0 ) , ensures that ψ ( t ) ∈ R + 4 for all t > 0 .

2.2 Boundedness

Boundedness in the drinking model refers to inherent constraints within the study, such as the exclusive focus on specific demographics or contexts. These limitations risk compromising the model’s reliability by narrowing its perspective and potentially overlooking broader behavioral patterns. As a result, reduced generalizability restricts the model’s applicability to specific circumstances and populations, which may compromise its ability to capture the nuanced complexity of real-world drinking behaviors. To understand the aforementioned model, we need the bounds on the dependent variables involved. To achieve this, we must determine the feasible region of attraction as outlined in the following theorem.

Theorem 2.1

There is a positive ψ for nonzero such that all solutions meet ψ > ( S ( t ) , H ( t ) , T ( t ) , R ( t ) ) for a long time t.

Proof

The solutions of System (1) are greater than zero, now in the first compartment of System (1) as

d S ( t ) d t = b − α S H − μ S + η R ≤ b − μ S .

Thus, d S ( t ) d t < 1 + b μ for undue time t , let us say t > t 0 for an excessive time t . Define R 1 ( t ) = S ( t ) + H ( t ) + T ( t ) . Now, differentiating of R 1 with the respect to the solution of Model (1) yields

d R 1 ( t ) d t = b − μ S − ( μ + δ 1 ) H − ( μ + δ 2 ) T , ≤ − h R 1 ( t ) + b ,

where h = min ( μ , ( μ + δ 1 ) , ( μ + δ 2 ) ) . Remember S ( t ) ≤ 1 + b μ for any t > t 0 . ψ 1 exists, reliant merely on the parameters of the considered model, such that R 1 ( t ) ≤ ψ 1 for eventually t > t 0 and H ( t ) and T ( t ) are bounded above. Subsequently, the second and third compartments of Model (1), R are eventually bounded above, and let ψ be the maximum. This displays that Model (1) is destructive.

Ω = ( S , H , T , R ) 1 + b μ ≥ S ≥ 0 , H ≥ 0 , ψ ≥ R .

Clearly, Ω is convex.□

2.3 Equilibrium points

Equilibrium points denote the stable-state prevalence of individuals engaged in drinking behavior, wherein the influx of new drinkers is counterbalanced by the rates of recoveries and mortality. There exist two equilibrium points: the drinking-free equilibrium point, and the drinking-present equilibrium point. For the computation of the drinking model, the process involves setting the left-hand side of Model (1) equations equal to zero.

(2) 0 = b − α S H − μ S + η R , 0 = α S H − ( μ + δ 1 + ϕ ) H , 0 = ϕ H − ( μ + δ 2 + γ ) T , 0 = γ T − ( μ + η ) R .

2.3.1 Drinking-free equilibrium point ( E 0 )

In the modeling of drinking behaviors, a drinking-free equilibrium point ( E 0 ) represents a stable state with sustained absence of alcohol consumption, offering insights into effective prevention or mitigation measures. Understanding the dynamics around this equilibrium enhances comprehension of influential factors governing the spread of drinking behaviors. Analyzing these points is crucial for devising precise strategies in public health and intervention efforts. Therefore, the drinking-free equilibrium point ( E 0 ) is derived by setting H = 0 , T = 0 , and R = 0 in Model (1) as follows:

E 0 = ( S 0 , H 0 , T 0 , R 0 ) = b μ , 0 , 0 , 0 .

2.3.2 Drinking endemic equilibrium point ( E 1 )

An endemic equilibrium point ( E 1 ) in the modeling of drinking behaviors is a stable state characterized by persistent alcohol consumption prevalence, providing insights into long-term influencing factors. Examining the dynamics around this point enhances understanding of the complex interplay governing the sustained spread of drinking behaviors. Analysis of endemic equilibrium points informs strategic interventions for addressing prolonged patterns of alcohol consumption in public health initiatives. The drinking endemic equilibrium point is written as E 1 = ( S * , H * , T * , R * ) ≠ 0 . Upon solving the equations from Model (1), we arrive at the following outcomes:

S * = μ + δ 1 + ϕ α , H * = μ ( μ + δ 1 + ϕ ) − α b α [ η − ( μ + δ 1 + ϕ ) ] , T * = ϕ [ μ ( μ + δ 1 + ϕ ) − α b ] α ( μ + δ 2 + γ ) [ η − ( μ + δ 2 + ϕ ) ] , R * = γ ϕ [ μ ( μ + δ 1 + ϕ ) − α b ] α ( μ + η ) ( μ + δ 2 + γ ) [ η − ( μ + δ 1 + ϕ ) ] ,

and ϕ = μ + δ 2 + γ , γ = μ + η .

2.4 Reproduction number

The basic reproduction number R 0 is a significant metric for assessing the transmission potential of heavy drinkers. It quantifies the total number of secondary heavy drinkers generated from the introduction of a single drinker individual into a susceptible population. In biological terms, it indicates the infectiousness and transmissibility of a pathogen.

We utilize the next-generation matrix method [42,43] to calculate the reproduction number. The matrices F and V are employed to represent the recruitment of new drinkers and the internal and external transmission terms associated with compartments of individuals engaged in drinking behavior, one obtains,

F = α S H 0 ,

V = − ( μ + δ 1 + ϕ ) H ϕ H − ( μ + δ 2 + γ ) T .

The Jacobian of F and V matrix are given as

(3) F = α 0 0 0 and V = − ( μ + δ 1 + ϕ ) 0 ϕ − ( μ + δ 2 + γ ) .

The basic reproduction number R 0 is calculated as the spectral radius of the next-generation matrix F V − 1 , defined as follows:

F V − 1 = α μ + δ 1 + ϕ 0 0 0 .

Hence, the basic reproduction number R 0 is obtained through the following derivation:

R 0 = α μ + δ 1 + ϕ .

If R 0 < 1 , it indicates that the disease is unlikely to propagate extensively in the population and might eventually diminish. Conversely, if R 0 > 1 , the disease can spread, potentially causing an outbreak.

3 Sensitivity analysis

Sensitivity analysis plays a crucial role in assessing the impact of input parameters on the dynamics of the epidemic drinking model. The parameters selected for sensitivity analysis in this study align closely with the model’s objectives, which aim to understand the dynamics of drinking behavior and its transmission within a population. These carefully chosen parameters represent key factors influencing the spread of drinking behaviors, including the transmission rate of heavy drinkers, natural death rate, drinking-induced death rates, and the proportion of individuals entering treatment. The sensitivity indices used in the drinking model to calculate the reproduction number ( R 0 ) are derived using the methodology proposed by Chitnis et al. [44]. We specifically compute the normalized forward sensitivity indices of the parameters m in relation to R 0 . Let us delve into this analysis further

Δ m R 0 = ∂ R 0 ∂ m m R 0 .

The sensitivity estimation for the reproduction number concerning various parameters is presented as follows:

Δ α R 0 = ∂ R 0 ∂ α α R 0 = 1 , Δ μ R 0 = ∂ R 0 ∂ μ μ R 0 = − 0.192 , Δ δ 1 R 0 = ∂ R 0 ∂ δ 1 δ 1 R 0 = − 0.270 , Δ ϕ R 0 = ∂ R 0 ∂ ϕ ϕ R 0 = − 0.538 .

Sensitivity analysis is commonly used to identify parameters significantly impacting the reproduction number R 0 . A positive value indicates a direct proportional relationship with the reproduction number, while a negative value indicates an inversely proportional relationship. In Table 2, the parameter ( α ) has a positive sign, indicating a directly proportional relationship with the reproduction number. Thus, an escalation in the reproduction number corresponds to an increase in parameters like the transmission rate of heavy drinkers ( α ). Reducing the parameter ( α ) can effectively control drinking in the population. Conversely, the reproduction number is negatively associated with parameters such as the natural death rate ( μ ), the proportion of drinkers entering treatment at the heavy drinker stage ( ϕ ), and the drinking-induced death rate of heavy drinker individuals ( δ 1 ). Increasing either the natural death rate of the population ( μ ), the drinking-induced death rate of heavy drinker individuals ( δ 1 ), or the proportion of drinkers entering treatment ( ϕ ) is neither ethical nor practical. Considering the sensitivity of R 0 to the drinking transmission rate of heavy drinkers ( α ), it is sensible to prioritize efforts on reducing this rate due to its direct impact on reproduction. In essence, this sensitivity analysis highlights the efficacy of prevention over cure. Enhancing preventive measures proves more effective in curbing the proliferation of habitual drinking compared to expanding access to treatment for individuals.

Table 2

Indices of sensitivity and parameters affecting the reproduction number R 0

Parameters	Sensitivity indices
α	+
μ	‒
δ 1	‒
ϕ	‒

A notable revelation from this sensitivity analysis is the substantial influence of the transmission rate of heavy drinkers on the reproduction number ( R 0 ). This discovery underscored the pivotal role that heavy drinkers have in propagating drinking behaviors within the population. Moreover, the sensitivity analysis underscored the significance of treatment interventions and their ability to lower heavy drinking rates, highlighting the importance of focused intervention strategies in addressing alcohol consumption. In simpler terms, this analysis suggests that prevention is more effective than trying to fix the issue later. Prioritizing prevention efforts is better for controlling the spread of drinking than solely focusing on increasing treatment numbers. The initial conditions for the model compartments, such as S ( 0 ) = 0.50 , T ( 0 ) = 0.15 , H ( 0 ) = 0.25 , and R ( 0 ) = 0.1 , as well as the values of compartmental parameters, are provided in Table 3 for numerical simulations and graphical representation of the model’s dynamic behavior.

Table 3

Parameters and its value of Model (1)

Parameters	Description	Value	Sources
b	Recruitment rate of the susceptible people	0.25	Estimated
α	Transmission rate from susceptible class to heavy drinker people	0.7	[40]
η	Transmission rate from recovered people to susceptible people	0.1	[40]
μ	Natural death rate	0.25	[45]
δ 1	Drinking-induced death rate of heavy drinker people	0.35	[45]
ϕ	Proportion of drinking entering of drinker in treatment class	0.7	[40]
δ 2	Drinking-induced death rate of drinker in treatment people	0.3	[45]
γ	Recovered rate of drinker in treatment people	0.09	[40]

The effect of different significant model parameters on the value of R 0 is shown graphically in Figure 2. Figure 2(a) describes the reproduction number R 0 sensitivity versus drinking transmission rate ( α ) and drinking including death rate ( δ 1 ). The influence of the drinking transmission rate ( α ) on the reproduction number is positive, whereas the impact of the drinking-related death rate ( δ 1 ) is negative. Although the effect of the drinking-related death rate is almost negligible, the influence of the drinking transmission rate ( α ) is quite significant, so it increased the value of R 0 by almost 0.8. A similar interpretation can be seen in Figure 2(f) for ( α ) and ( μ ) on R 0 . Figure 2(b) depicts the sensitivity of the reproduction number ( R 0 ) concerning the drinking-related death rate ( δ 1 ) and the natural death rate ( μ ). Both parameters are decreased since they negatively impact the reproduction number ( R 0 ). Therefore, both ( δ 1 ) and μ have a similar effect on R 0 , reaching a value close to 0.25. Figure 2(c) illustrates the impact of the parameters ( α ) and ( ϕ ) on R 0 . It is observed that R 0 increases significantly with these parameters, reaching a maximum value close to 0.8 or less than 1. This implies that the disease is unlikely to occur. Figure 2(d) illustrates how variations in the drinking-related death rate ( δ 1 ) and the proportion of heavy drinkers in the treatment group ( ϕ ) influence the reproduction number ( R 0 ). Both values decrease because they have a negative impact on R 0 . So, both ( δ 1 ) and ( ϕ ) have a similar effect on R 0 , making it approach a value near 0.5. Similarly, interpretation can be seen in Figure 2(e) for the proportional of heavy drinkers in the treatment compartment ( α ) and natural death rate ( μ ) on ( R 0 ).

$Figure 2 Fluctuation of different compartmental parameters and their effect on the basic reproduction number R 0 {R}_{0} : (a) R 0 {R}_{0} vs α \alpha and δ 1 {\delta }_{1} , (b) R 0 {R}_{0} vs δ 1 {\delta }_{1} and μ \mu , (c) R 0 {R}_{0} vs α \alpha and ϕ \phi , (d) R 0 {R}_{0} vs ϕ \phi and δ 1 {\delta }_{1} , (e) R 0 {R}_{0} vs ϕ \phi and μ \mu , and (f) R 0 {R}_{0} vs α \alpha and μ \mu .$

Figure 2

Fluctuation of different compartmental parameters and their effect on the basic reproduction number R 0 : (a) R 0 vs α and δ 1 , (b) R 0 vs δ 1 and μ , (c) R 0 vs α and ϕ , (d) R 0 vs ϕ and δ 1 , (e) R 0 vs ϕ and μ , and (f) R 0 vs α and μ .

4 Analysis of stability

In this segment, we explore the local stability and global stability of the developed model at equilibrium points. The local stability is assessed by analyzing the eigenvalues of the Jacobian matrix at the equilibrium point. An equilibrium point is deemed locally stable if all eigenvalues have negative real parts. Global stability entails studying the system’s behavior across its whole range, frequently necessitating progressive mathematical techniques such as Lyapunov analysis.

4.1 Local stability

A locally stable system tends to revert to its equilibrium state after experiencing minor disruptions, showcasing resilience and a tendency to maintain its original configuration near the equilibrium point. Conducting a local stability analysis within a drinking model aims to assess how minor disturbances can affect the persistence or reduction of drinking behavior. In the context of drinking models, the concept of local stability at equilibrium points becomes a crucial framework for understanding the complex dynamics of alcohol consumption and its spread within a population. In this discussion, we explore the locally stable drinking-free equilibrium point and the endemic equilibrium point in the devised model, utilizing the following theorem [29].

Theorem 4.1

The drinking-free equilibrium point E 0 is locally asymptotically stable with the condition R 0 < 1 , whereas unstable with the condition R 0 > 1 .

Proof

The Jacobian matrix of drinking-free equilibrium point E 0 = ( S 0 , H 0 , T 0 , R 0 ) = ( b μ , 0 , 0 , 0 ) for the purposed model is given as

J ( E 0 ) = − ( α H + μ ) − α S 0 η α H α S − ( μ + δ 1 + ϕ ) 0 0 0 ϕ − ( μ + δ 2 + γ ) 0 0 0 γ − ( μ + η ) .

Evaluating the Jacobian matrix at the drinking-free equilibrium yields

J ( E 0 ) = − μ − α 0 η 0 α − ( μ + δ 1 + ϕ ) 0 0 0 ϕ − ( μ + δ 2 + γ ) 0 0 0 γ − ( μ + η ) .

The resulting eigenvalues are λ 1 = − μ < 0 and λ 2 = − ( μ + δ 2 + ϕ ) < 0 , whereas

J ( E 0 ) − λ I = α − ( μ + δ 1 + ϕ ) − λ 0 ϕ − ( μ + δ 2 + γ ) − λ = 0 ,

λ 4 = − ( μ + δ 2 + γ ) < 0 . λ 1 , λ 2 , and λ 4 are real and negative. Thus, λ 3 = α − ( μ + δ 1 + ϕ ) < 0 if and only if R 0 < 1 . Using the Routh–Hurwitz criterion [46,47], we confirm that each eigenvalue of the polynomial equation possesses a non-positive real part when R 0 < 1 . Consequently, E 0 is locally asymptotically stable.□

Theorem 4.2

The endemic equilibrium point E 1 in Model (1)exhibits local asymptotic stability under the condition R 0 > 1 ; otherwise, it is unstable.

Proof

The Jacobin matrix is computed as

(4) J ( E 1 ) = − α μ ( μ + δ 1 + ϕ ) − α b α [ η − ( μ + δ 1 + ϕ ) ] + μ − α μ + δ 1 + ϕ α 0 η α μ ( μ + δ 1 + ϕ ) − α b α [ η − ( μ + δ 1 + ϕ ) ] α μ + δ 1 + ϕ α − ( μ + δ 1 + ϕ ) 0 0 0 ϕ − ( μ + δ 2 + γ ) 0 0 0 γ − ( μ + η ) ,

where

C 11 = − α μ ( μ + δ 1 + ϕ ) − α b α [ η − ( μ + δ 1 + ϕ ) ] + μ , C 12 = − α μ + δ 1 + ϕ α , C 13 = 0 , C 14 = η , C 21 = α μ ( μ + δ 1 + ϕ ) − α b α [ η − ( μ + δ 1 + ϕ ) ] , C 22 = α μ + δ 1 + ϕ α − ( μ + δ 1 + ϕ ) , C 23 = 0 , C 24 = 0 , C 31 = 0 , C 32 = ϕ , C 33 = − ( μ + δ 2 + γ ) , C 34 = 0 , C 41 = 0 , C 42 = 0 , C 43 = γ , C 44 = − ( μ + η ) .

Substituting C i j into equation (4), we obtain

(5) J ( E 1 ) = C 11 C 12 0 C 14 C 21 C 22 0 0 0 C 32 C 33 0 0 0 C 43 C 44 .

The auxiliary equation of (5) can be calculated as follows:

(6) ∣ ( J ( E 1 ) − λ I ) ∣ = C 11 − λ C 12 0 C 14 C 21 C 22 − λ 0 0 0 C 32 C 33 − λ 0 0 0 C 43 C 44 − λ = 0 .

λ 4 + ( C 11 + C 22 + C 33 + C 44 ) λ 3 + ( C 12 C 21 − C 33 C 44 − C 11 C 44 − C 22 C 44 − C 11 C 33 − C 22 C 33 − C 11 C 22 ) λ 2 + ( C 11 C 33 C 44 + C 22 C 33 C 44 + C 11 C 22 C 44 − C 12 C 21 C 33 − C 12 C 21 C 44 ) λ 1 + ( C 12 C 21 C 33 C 44 + C 14 C 21 C 32 C 43 − C 11 C 22 C 33 C 44 ) = 0 .

The aforementioned equation can be written as

λ 4 + D 1 λ 3 + D 2 λ 2 + D 3 λ 1 + D 4 = 0 ,

where

(7) D 1 = ( C 11 + C 22 + C 33 + C 44 ) , D 2 = ( C 12 C 21 − C 33 C 44 − C 11 C 44 − C 22 C 44 − C 11 C 33 − C 22 C 33 − C 11 C 22 ) , D 3 = ( C 11 C 33 C 44 + C 22 C 33 C 44 + C 11 C 22 C 44 − C 12 C 21 C 33 − C 12 C 21 C 44 ) , D 4 = ( C 12 C 21 C 33 C 44 + C 14 C 21 C 32 C 43 − C 11 C 22 C 33 C 44 ) .

Applying the Routh–Hurwitz criterion to fourth-order polynomials [46], we ensure that D 1 > 0 , D 4 > 0 , D 1 D 2 − D 3 > 0 , and ( D 1 D 2 − D 3 ) D 3 − D 1 2 D 4 > 0 hold true only if R 0 > 1 . The non-positive eigenvalue of the auxiliary equation confirms the local asymptotic stability of E 1 as per the Routh–Hurwitz criterion.□

4.2 Global stability

Global stability in a drinking model investigates whether, under various conditions, the system converges to a stable state regarding drinking behavior across its entire range. Global stability analysis in the drinking model is essential for understanding the sustained, long-term dynamics of population drinking behavior. By examining the system’s equilibrium conditions and considering the influence of broader societal factors, this analysis provides a comprehensive understanding of how the population’s drinking patterns evolve over an extended period. Such insights are crucial for developing enduring public health strategies and policies tailored to address persistent trends and dynamics in alcohol consumption within the population. Global stability analyses offer a holistic framework for effective long-term planning and management of drinking dynamics. In this discussion, we investigate the globally stable equilibrium points for both drinking-free and endemic states in the proposed model, employing the methodology outlined in the study by Anjam et al. [29].

Theorem 4.3

The drinking-free equilibrium point E 0 is globally asymptomatically stable with the condition R 0 < 1 and unstable with the condition R 0 > 1 .

Proof

The global stability of the devised model at point E 0 is instigated through the construct of the Lyapunov function such as

(8) W ( t ) = ( S 0 − S ) + H + T .

Upon computing the time derivative of equation (8) and applying it to the system of equations in (1), the following result is obtained:

(9) d W d t = d S d t + d H d t + d T d t , d W d t = b − μ S − ( μ + δ 1 ) H − ( μ + δ 2 ) T , d W d t = − [ μ ( S − S 0 ) + ( μ + δ 1 ) H + ( μ + δ 2 ) T ] < 0 .

It is verifiable that d W d t < 0 if R 0 < 1 . Also, d W d t = 0 if S = S 0 and H = T = 0 . The Lyapunov LaSalle’s invariant principle [48,49] indicates that E 0 achieves global asymptotic stability.□

Theorem 4.4

The endemic equilibrium point E 1 of Model (1)exhibits global asymptotic stability when R 0 > 1 ; otherwise, it is unstable.

Proof

The global stability is established by formulating the Lyapunov function at the endemic equilibrium points E 1 = ( S * , H * , T * , R * ) in the following manner:

(10) Z = 1 2 [ ( S − S * ) + ( H − H * ) + ( T − T * ) ] 2 .

Upon calculating the time derivatives of the aforementioned equation and utilizing Model (1), the resulting expression is as follows:

d Z d t = [ ( S − S * ) + ( H − H * ) + ( T − T * ) ] × [ b − μ S − μ + δ 1 H − ( μ + δ 2 ) T ] , d Z d t = [ ( S − S * ) + ( H − H * ) + ( T − T * ) ] × [ μ S * R 0 − μ S − ( μ + δ 1 ) H − ( μ + δ 2 ) T ] ,

(11) d Z d t = − [ ( S − S * ) + ( H − H * ) + ( T − T * ) ] × [ μ ( S − S * R 0 ) + ( μ + δ 1 ) H + ( μ + δ 2 ) T ] .

As a consequence, d Z d t ≤ 0 for each value ( S * , H * , T * , R * ) , whereas d Z d t = 0 holds only for S = S * , H = H * , and T = T * . Consequently, the endemic equilibrium E 1 is the only positively invariant set contained in [ ( S , H , T , R ) , S = S * , H = H * , T = T * , R = R * ] . As a result, the positive E 1 exhibits global asymptotic stability.□

4.3 Results from numerical stability analysis

Runge–Kutta (RK) method is highly favored in epidemiological modeling because of its accuracy and robustness in solving the differential equations that govern disease dynamics. These models typically involve a system of coupled differential equations that describe interactions among different population compartments. The fourth-order RK method is particularly noteworthy for its high accuracy in approximating solutions to these equations, guaranteeing reliable and precise numerical results. Its robust convergence properties enhance both stability and efficiency, making it highly effective in solving diverse mathematical problems, especially those involving intricate and dynamic infectious disease dynamics. In this context, we employ the RK method to solve our proposed deterministic model, leveraging its ability to compute solutions through a series of approximations based on weighted function evaluations. Notably, the fourth-order RK method delivers accurate results with minimal model evaluations, making it a preferred choice for obtaining numerical solutions in Model (1). Hence, we obtain

(12) S i + 1 − S i l = b − α S i + 1 H i − μ S i + 1 + η R i , H i + 1 − H i l = α S i H i + 1 − ( μ + δ 1 + ϕ ) H i + 1 , T i + 1 − T i l = ϕ H i − ( μ + δ 2 + γ ) T i + 1 , R i + 1 − R i l = γ T i − ( μ + η ) R i + 1 .

To numerically solve the proposed model, we initially focus on the equation for the susceptible compartment S in Model (12). This process results in

(13) S i + 1 − S i = l b − l α S i + 1 H i − l μ S i + 1 + l η R i , S i + 1 + l α S i + 1 H i + l μ S i + 1 = S i + l b + l η R i , S i + 1 = S i 1 + l α H i + l μ + l b 1 + l α H i + l μ + l η R i 1 + l α H i + l μ .

Similarly, the remaining compartments of the proposed model are solved, respectively, and give,

H i + 1 − H i l = α S i H i + 1 − ( μ + δ 1 + ϕ ) H i + 1 , H i + 1 − H i = l α S i H i + 1 − l ( μ + δ 1 + ϕ ) H i + 1 , ⋮ H i + 1 = H i 1 − l α S i + l ( μ + δ 1 + ϕ ) . T i + 1 − T i l = ϕ H i − ( μ + δ 2 + γ ) T i + 1 , T i + 1 − T i = l ϕ H i − l ( μ + δ 2 + γ ) T i + 1 , ⋮ T i + 1 = T i 1 + l ( μ + δ 2 + γ ) + l ϕ H i 1 + l ( μ + δ 2 + γ ) , R i + 1 − R i l = γ T i − ( μ + η ) R i + 1 , R i + 1 − R i = l γ T i − l ( μ + η ) R i + 1 , ⋮ R i + 1 = R i 1 + l ( μ + η ) + γ T i 1 + l ( μ + η ) .

4.3.1 Algorithm

Step 1: S ( 0 ) = H ( 0 ) = T ( 0 ) = R ( 0 ) = 0 .

Step 2: for i = 1 , 2 , 3 , … , n − 1

S i + 1 = S i 1 + l α H i + l μ + l b 1 + l α H i + l μ + l η R i 1 + l α H i + l μ , H i + 1 = H i 1 − l α S i + l ( μ + δ 1 + ϕ ) , T i + 1 = T i 1 + l ( μ + δ 2 + γ ) + l ϕ H i 1 + l ( μ + δ 2 + γ ) , R i + 1 = R i 1 + l ( μ + η ) + γ T i 1 + l ( μ + η ) .

Step 3: for i = 1 , 2 , 3 , … , n − 1 , write S * ( t i ) = S * , H * ( t i ) = H * , T * ( t i ) = T * , and R * ( t i ) = R * .

For the graphical depiction of the stability analysis, MATLAB software is utilized to run the aforementioned findings. The graphs presented in the figures below are generated using compartmental initial conditions for the devised model: S ( 0 ) = 0.50 , T ( 0 ) = 0.15 , H ( 0 ) = 0.25 , and R ( 0 ) = 0.1 , along with the parametric values specified in Table 3. These graphs serve to visually characterize the dynamics and stability characteristics of the model under varying initial conditions and parameter settings. Moreover, the time is taken between 0 and 30 units in terms of days with the initial population of susceptible individuals S ( t ) , drinkers in treatment individuals T ( t ) , heavy drinkers individuals H ( t ) , and recovered individuals R ( t ) .

Figures 3, 4, 5, 6 depict the stability dynamics of all population classes under various initial conditions. The susceptible class varies in magnitude with values of 0.5, 0.35, 0.2, and 0.05, showing a direct relationship among the size of the susceptible population and the time required to achieve stability (Figure 3). Similarly, Figure 4 illustrates stability achievements across different sizes of the heavy drinker class, with values of 0.25, 0.20, 0.15, and 0.10. The heavy drinker population decreases over time, stabilizing after 10 days. The population in the drinkers in treatment class is shown for various prevalence levels, represented by values such as 0.15, 0.13, 0.11, and 0.09 in Figure 5. The number of individuals in treatment for drinking initially rises over time, followed by a rapid decline as time progresses. Figure 6 illustrates the correlation between stability and the varying extent of the recovered population, considering values such as 0.1, 0.09, 0.08, and 0.07. The population experiences an initial rapid decrease, stabilizing with constant numbers after 15 days.

$Figure 3 Simulation results illustrate how changes in the susceptible class S ( t ) {\mathbb{S}}\left(t) manifest their impacts.$

Figure 3

Simulation results illustrate how changes in the susceptible class S ( t ) manifest their impacts.

$Figure 4 Simulation results depict how variations in the heavy drinker class H ( t ) {\mathbb{H}}\left(t) manifest their effects.$

Figure 4

Simulation results depict how variations in the heavy drinker class H ( t ) manifest their effects.

$Figure 5 Simulation results depict how variations in the drinker in treatment class T ( t ) {\mathbb{T}}\left(t) manifest their effects.$

Figure 5

Simulation results depict how variations in the drinker in treatment class T ( t ) manifest their effects.

$Figure 6 Simulation results illustrate how changes in the recovered class R ( t ) {\mathbb{R}}\left(t) manifest their impacts.$

Figure 6

Simulation results illustrate how changes in the recovered class R ( t ) manifest their impacts.

5 Optimal control theory

This section delves into the application of the drinking model with regard to optimal control strategies. Optimal control is utilized to intercede in the transmission of drinking within a specific community. The integration of an optimal control strategy, focusing on awareness campaigns ( u 1 ( t ) ) and treatment ( u 2 ( t ) ), enhances the model’s effectiveness in reducing heavy drinking rates by orchestrating a coordinated and targeted approach. Awareness campaigns play a pivotal role in educating the population about the risks and consequences of heavy drinking, promoting healthier behaviors, and reducing the social acceptability of excessive alcohol consumption. This can lead to a shift in social norms and attitudes towards drinking, thereby influencing individual behaviors and reducing heavy drinking rates over time. On the other hand, treatment interventions target individuals who are already engaged in heavy drinking behaviors. By providing access to effective treatment options, such as counseling, therapy, or medication, the model can help individuals reduce their alcohol consumption, manage cravings, and address underlying issues contributing to their drinking habits. This approach not only aids in reducing heavy drinking rates but also promotes overall health and well-being among affected individuals. Furthermore, the combination of awareness campaigns and treatment interventions in the optimal control strategy allows for a comprehensive approach to addressing heavy drinking rates. Awareness campaigns create a supportive environment for individuals seeking help, while treatment interventions offer practical solutions for those struggling with heavy drinking. Together, these strategies work synergistically to reduce heavy drinking rates and promote healthier behaviors within the population. Moreover, the control variable u 1 ( t ) pertains to the implementation of a public awareness program targeting susceptible and treated individuals. This program aims to safeguard them against excessive alcohol consumption and facilitate successful avoidance. On the other hand, the u 2 ( t ) control variable represents the treatment applied to heavy drinker individuals. Therefore, heavy drinkers represent the proportion of individuals transitioning from regular alcohol consumption. Upon applying the u 2 ( t ) variable, they transition to the group of those permanently abstaining from alcohol. This control variable determines the optimal values for decision variables, representing the adjustable action variable to achieve the desired result. By selecting the optimal value for the control variable, the best course of action can be determined. For emerging a control strategy, we bring into usage the optimal control theory [50–52]. The adaptation of the control scheme interning susceptible drinkers class S ( t ) , drinkers in treatment class T ( t ) , heavy drinkers class H ( t ) , and recovered class R ( t ) is discussed in next lines. The differential formulation of the model’s equation (1), incorporating the control variable, is now expressed as follows:

(14) d S d t = b − ( 1 − u 1 ) α S H − μ S + η R , d H d t = ( 1 − u 1 ) α S H − ( μ + δ 1 + ϕ + u 2 ) H , d T d t = ϕ H − ( μ + δ 2 + γ + u 1 ) T , d R d t = γ T − ( μ + η ) R + u 2 H + u 1 T ,

whereas the initial conditions are,

S ( 0 ) ≥ 0 , H ( 0 ) ≥ 0 , T ( 0 ) ≥ 0 , and R ( 0 ) ≥ 0 .

Therefore, we consider an optimal control problem, wherein the objective function is specified as

(15) J ( u 1 ( t ) , u 2 ( t ) ) = ∫ 0 t A 1 * S ( t ) + A 2 * H ( t ) + A 3 * T ( t ) + A 4 * R ( t ) + 1 2 A 5 * u 1 2 ( t ) + 1 2 A 6 * u 2 2 ( t ) d t ,

where A 1 * , A 2 * , A 3 * , A 4 * , A 5 * , and A 6 * represent the weight constant allied with each class. The objective function is minimized with respect to control pairs such as

(16) J { u 1 * , u 2 * } = min { J ( u 1 ( t ) , u 2 ( t ) ) , u 1 ( t ) , u 2 ( t ) ∈ U } .

Here, U = { ( u 1 ( t ) , u 2 ( t ) ) \ u i ( t ) } is a control set with Lebesgue measurable onto [ 0 , 1 ] and permissible range such as 0 ≤ u i ( t ) ≤ 1 , where i = { 1 , 2 } .

Furthermore, Pontryagin’s maximum principle [53,54] is used to find the optimal control strategy u that minimizes (or maximizes) the objective functional J by considering the behavior of the Hamiltonian function H at each point in the state variables.

Theorem 5.1

For the optimal control Problem (14), there exists a u * ( t ) = ( u 1 * ( t ) , u 2 * ( t ) ε U ) such that

min ( u 1 ( t ) , u 2 ( t ) ∈ U ) J ( u 1 ( t ) , u 2 ( t ) ) = J ( u 1 ( t ) , u 2 ( t ) ) .

Proof

A variety of techniques are considered to quantify the effectiveness of the control scheme [55,56], and consequently, nonnegative values are witnessed for control variables along with state variables. The convexity and the closeness of control variables describe the degree of solidity necessary to validate the control system. Moreover, it is noteworthy that the integrand exhibits a convex nature, representing an objective function in the form of A 1 * S ( t ) + A 2 * H ( t ) + A 3 * T ( t ) + A 4 * R ( t ) + 1 2 A 5 * u 1 2 ( t ) + 1 2 A 6 * u 2 2 ( t ) as certification of the proof, where A i * , i = 1 , 2 , 3 , … , 6 represents the weight constant. We employ the Pontryagin maximum principle [39,53,54] to address our proposed problem. The Hamiltonian function is provided as

(17) H = L ( S ( t ) , H ( t ) , T ( t ) , R ( t ) , u 1 , u 2 ) + λ 1 d S ( t ) d t + λ 2 d H ( t ) d t + λ 3 d T ( t ) d t + λ 4 d R ( t ) d t .

Moreover, the presence of non-trivial vector functions like λ ( t ) = ( λ 1 ( t ) , λ 2 ( t ) , λ 3 ( t ) , λ 4 ( t ) ) becomes apparent, particularly when ( y * , u * ) is regarded as the optimal procedure for the implemented control problem:

(18) d y d t = ∂ H ( t , y , u , λ ) ∂ u , 0 = ∂ H ( t , y , u , λ ) ∂ u , λ ′ ( t ) = − ∂ H ( t , y , u , λ ) ∂ y .

The outcomes are achieved by the employment of the required condition to the Hamiltonian function.□

5.1 Existence of the optimal control problem

In this section, we establish the existence of optimal control by ensuring that all conditions of the control system are met at the initial time t = 0 . To achieve this, we maintain a bounded Lebesgue measurable control [56,57] with initial conditions and an upward bounded solution of the system. The problem of optimal control is tackled through an exploration of Lagrangian and Hamiltonian methods. The optimal control problem in Lagrangian formulation is represented by the following equation:

(19) L { S ( t ) , H ( t ) , T ( t ) , R ( t ) , u 1 ( t ) , u 2 ( t ) } = A 1 * S ( t ) + A 2 * H ( t ) + A 3 * T ( t ) + A 4 * R ( t ) + 1 2 A 5 * u 1 2 ( t ) + 1 2 A 6 * u 2 2 ( t ) .

The minimum optimal control value is determined by formulating the Hamiltonian function H in the following manner:

H = L ( S ( t ) , H ( t ) , T ( t ) , R ( t ) , u 1 , u 2 ) + λ 1 d S ( t ) d t + λ 2 d H ( t ) d t + λ 3 d T ( t ) d t + λ 4 d R ( t ) d t .

The resultant system, in terms of adjoint variable λ 1 , λ 2 , λ 3 , λ 4 , and optimal control variable u 1 and u 2 , is written as

(20) λ 1 ′ = − { A 1 * − ( 1 − u 1 ) α H ( λ 1 − λ 2 ) − λ 1 μ } , λ 2 ′ = − { A 2 * − ( 1 − u 1 ) α S ( λ 1 − λ 2 ) − λ 2 ( μ + δ 1 + ϕ + u 2 ) + λ 4 u 2 + λ 3 ϕ } , λ 3 ′ = − { A 3 * − λ 3 ( μ + δ 2 + γ + u 1 ) + λ 4 ( γ + u 1 ) } , λ 4 ′ = − { A 4 * + λ 1 η − λ 4 ( η + μ ) } , u 1 ( t ) = ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * , u 2 ( t ) = H ( λ 2 − λ 4 ) A 6 * .

Theorem 5.2

Given an optimal control u 1 * , u 2 * and solution S * ( t ) , H * ( t ) , T * ( t ) , and R * ( t ) of the equivalent state system (14), there exist adjoint variable λ m ( t ) , m = 1 , … , 4 such as,

(21) λ 1 ′ = − { A 1 * − ( 1 − u 1 ) α H ( λ 1 − λ 2 ) − λ 1 μ } , λ 2 ′ = − { A 2 * − ( 1 − u 1 ) α S ( λ 1 − λ 2 ) − λ 2 ( μ + δ 1 + ϕ + u 2 ) + λ 4 u 2 + λ 3 ϕ } , λ 3 ′ = − { A 3 * − λ 3 ( μ + δ 2 + γ + u 1 ) + λ 4 ( γ + u 1 ) } , λ 4 ′ = − { A 4 * + λ 1 η − λ 4 ( η + μ ) } ,

using the transversality criterion λ m ( t ) , m = 1 , … , 4 .

Proof

Let us assume that S ( t ) = S * ( t ) , H ( t ) = H * ( t ) , T ( t ) = T * ( t ) , and R ( t ) = R * ( t ) , and let also describe Hamiltonian respecting state variables such as S ( t ) , H ( t ) , T ( t ) , and R ( t ) . The altered adjoint system under the condition of transversality is now given as

(22)□ λ 1 ′ = − { A 1 * − ( 1 − u 1 ) α H ( λ 1 − λ 2 ) − λ 1 μ } , λ 2 ′ = − { A 2 * − ( 1 − u 1 ) α S ( λ 1 − λ 2 ) − λ 2 ( μ + δ 1 + ϕ + u 2 ) + λ 4 u 2 + λ 3 ϕ } , λ 3 ′ = − { A 3 * − λ 3 ( μ + δ 2 + γ + u 1 ) + λ 4 ( γ + u 1 ) } , λ 4 ′ = − { A 4 * + λ 1 η − λ 4 ( η + μ ) } .

Theorem 5.3

The optimal control pair ( u 1 * ( t ) , u 2 * ( t ) ) over the region U is given by

(23) u 1 * ( t ) = max min ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * , 1 , 0 , u 2 * ( t ) = max min H ( λ 2 − λ 4 ) A 6 * , 1 , 0 .

Proof

The application of optimal conditions provided outcomes such as

(24) ∂ H ∂ u 1 = A 5 * u 1 ( t ) + λ 1 α H S − λ 2 α H S + λ 4 T − λ 3 T , ∂ H ∂ u 2 = A 6 * u 2 ( t ) − λ 2 H + λ 4 H .

The control variable is also solved such as

(25) u 1 * ( t ) = ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * , u 2 * ( t ) = H ( λ 2 − λ 4 ) A 6 * .

One may note that the condition of control space remains expressible as under

u 1 * ( t ) = 0 , if ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * ≤ 0 , ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * , if 0 < ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * < 1 , 1 , if ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) , T A 5 * ≥ 1 . u 2 * ( t ) = 0 , if H ( λ 2 − λ 4 ) A 6 * ≤ 0 , H ( λ 2 − λ 4 ) A 6 * , if 0 < H ( λ 2 − λ 4 ) A 6 * < 1 , 1 , if H ( λ 2 − λ 4 ) A 6 * ≥ 1 .

In simplification, the control variables are documented as follows:

(26) u 1 * ( t ) = max min ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * , 1 , 0 , u 2 * ( t ) = max min H ( λ 2 − λ 4 ) A 6 * , 1 , 0 .

The optimal system is now given as:

(27) d S * d t = b − 1 − max min ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * , 1 , 0 α S * H * − μ S * + η R * , d H * d t = 1 − max min ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * , 1 , 0 α S * H * − μ + δ 1 + ϕ − max min H ( λ 2 − λ 4 ) A 6 * , 1 , 0 H * , d T * d t = ϕ H * − μ + δ 2 + γ − max min ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * , 1 , 0 T * , d R * d t = γ T * − ( μ + η ) R * + max min H ( λ 2 − λ 4 ) A 6 * , 1 , 0 H * + max min ( λ 2 − λ 1 ) α S H + ( λ 3 − λ 4 ) T A 5 * , 1 , 0 T * .

The state variable and control variable in the optimal formation are achieved by the employment of an adjoint variable, the optimal system along with the initial conditions. The positivity of the second-order derivative of the objective function with respect to control variables is obvious. Therefore, the resultant Hamiltonian is documented as

□ H * = A 1 * S * ( t ) + A 2 * H * ( t ) + A 3 * T * ( t ) + A 4 * R * ( t ) + 1 2 A 5 * u 1 2 ( t ) + 1 2 A 6 * u 2 2 ( t ) + ∑ m = 1 4 λ m g m ( S * , H * , T * , R * ) .

5.2 Computational simulation and discussion

In this section, we aim to validate the proposed control strategy through numerical evaluation. We enrich the strategic framework by integrating a wide range of parametric settings. The numerical analysis commences with initializing the compartmental initial conditions from the study by Adu et al. [40] for our proposed model, comprising S ( 0 ) = 0.50 , T ( 0 ) = 0.15 , H ( 0 ) = 0.25 , and R ( 0 ) = 0.1 . Furthermore, we ensure consistency by adopting literature-based parameter values as outlined in Table 3. We address the optimization problem by employing the RK fourth-order technique and applying the transversality criterion within the time interval of [0,25]. Additionally, we assume the weight constants to be A 1 = 0.6610000 , A 2 = 0.54450 , A 3 = 0.0090030 , A 4 = 0.44440 , A 5 = 0.3550 , and A 6 0.5560 . Figures 7, 8, 9, 10 illustrate the dynamic effectiveness of our model under both controlled and uncontrolled scenarios across all four compartments over time in days. The advantages of the optimal control strategy are evident across all population compartments. The steeper rise of the susceptible population with respect to control is obvious in Figure 7. Similarly, the class of heavy drinkers is noted to be reduced more sharply with control as compared to without control in Figure 8. Moreover, the effectiveness of employed control can be seen for the compartment of drinkers in treatment from Figure 9. Finally, Figure 10 depicts the utility of the control scheme for the class of recovered drinkers. Moreover, Figures 11 and 12 provide insights into the effectiveness of control variables, namely, the social awareness campaign u 1 ( t ) and treatment for reducing the drinking population u 2 ( t ) , both with and without control interventions over time in days.

$Figure 7 Dynamics behavior of susceptible individual S ( t ) {\mathbb{S}}\left(t) with and without control are depicted in the graph.$

Figure 7

Dynamics behavior of susceptible individual S ( t ) with and without control are depicted in the graph.

$Figure 8 Dynamics behavior of heavy drinker individual H ( t ) {\mathbb{H}}\left(t) with and without control are depicted in the graph.$

Figure 8

Dynamics behavior of heavy drinker individual H ( t ) with and without control are depicted in the graph.

$Figure 9 Dynamics behavior of drinker in treatment individual T ( t ) {\mathbb{T}}\left(t) with and without control are depicted in the graph.$

Figure 9

Dynamics behavior of drinker in treatment individual T ( t ) with and without control are depicted in the graph.

$Figure 10 Dynamics behavior of recovered individual R ( t ) {\mathbb{R}}\left(t) with and without control are depicted in the graph.$

Figure 10

Dynamics behavior of recovered individual R ( t ) with and without control are depicted in the graph.

$Figure 11 Dynamics of control variable u 1 {u}_{1} are depicted in the graphic.$

Figure 11

Dynamics of control variable u 1 are depicted in the graphic.

$Figure 12 Dynamics of control variable u 2 {u}_{2} are depicted in the graphic.$

Figure 12

Dynamics of control variable u 2 are depicted in the graphic.

6 Conclusion

This study is related to the development of a novel deterministic mathematical model elaborating the dynamics of the drinking epidemic model under the launch of an optimal control scheme. The investigation focuses on stratifying the population into four compartments: the susceptible class, heavy drinker class, drinkers in treatment class, and recovered drinker class. To ensure generalizability, a wide range of parametric settings are employed, encompassing the recruitment rate of the susceptible class, transmission rates from the susceptible class to the heavy drinker class, and from the recovered class to the susceptible class, natural death rate, drinking-induced death rates in the heavy drinker and treatment classes, the proportion of individuals entering the treatment class, and recovery rate. Throughout the investigation, a comprehensive analysis is provided, covering aspects such as positivity, boundedness, equilibrium points, basic reproduction number, stability, and sensitivity. The model ensures boundedness and positivity by incorporating constraints and selecting appropriate functional forms. This precautionary approach prevents unrealistic growth of variables, aligning the model with real-world constraints and enhancing its reliability. The careful management of these aspects improves the model’s applicability, offering policymakers an accurate representation that supports well-informed decision-making within practical limitations. The analysis of stability in the model indicates that when R 0 < 1 , the system is both locally and globally asymptotically stable at the drinking-free equilibrium E 0 . Conversely, if R 0 > 1 , the system exhibits an endemic equilibrium E 1 , at which point it is locally and globally asymptotically stable. Furthermore, the sensitivity analysis to see the impact of different compartmental parameters on the reproduction number is discussed. The sensitivity analysis showed that the drinking transmission rates of heavy drinkers ( α ) are highly sensitive to R 0 , indicating their significant influence on the spread of drinking behaviors. This suggests that these parameters could be considered appropriate for inclusion as control variables in optimal control analysis. Furthermore, this finding underscores the critical role of heavy drinkers in driving the spread of drinking behaviors within the population. Additionally, the sensitivity analysis underscores the importance of treatment interventions and their potential to reduce heavy drinking rates, emphasizing the value of targeted intervention strategies in curbing alcohol consumption.

We initiate the optimal control process by integrating two time-dependent variables: a social awareness campaign and treatment, both acting as control variables. An optimal control problem is formulated by defining an objective function aimed at determining the optimal values for these control variables to minimize overall costs. By applying Pontryagin’s maximum principle, we derive significant conditions for the optimal solution. Our research investigates the viability of two different optimal control methodologies. Numerical simulations conducted using MATLAB validate the effectiveness of the proposed control strategies. Graphical results reveal that the simultaneous use of both control variables is more effective in reducing the transmission flow of drinking behavior. These findings are extensively discussed, considering scenarios with and without optimal control. Through rigorous analytical procedures, we establish the utility of the control strategy in achieving optimal conditions across all four compartments. The susceptible population is estimated to show a gradual increase with control implementation, while a notable decrease in the proportions of heavy drinkers and individuals in treatment is observed with optimal control. Additionally, the recovery rate is significantly faster with the use of optimal control. Furthermore, our research advocates that simulations incorporating time-dependent controls are supplementary cost-effective equated to those utilizing time-independent controls.

In future endeavors, extending the model to include fractional-order dynamics alongside optimal control remains an intriguing avenue. This approach is expected to provide deeper insights into the dynamics of drinking behavior in society by leveraging available information more effectively. Additionally, there will be a strong emphasis on assessing the proposed scheme’s applicability across diverse geographical regions, socio-economic strata, and varying healthcare access levels. Understanding how these strategies interact or counteract each other will help formulate more effective plans for disease control on a global scale. Moreover, future research may also focus on integrating real-time data and advanced machine-learning techniques to enhance the accuracy and predictive capabilities of the model.

Funding information: This research was supported by the Deanship of Research Development, Prince Mohammad Bin Fahd, University, Alkhubar 31952, Saudi Arabia.
Author contributions: All authors have accepted responsibility for the entire content of this manuscript and have approved its submission. YNA and KI designed the concept, while KI and SAC carried out the methodology. MF and SM developed the software model code and conducted the simulations. YNA, SM, and NS prepared the manuscript with contributions from all coauthors. The authors applied the SDC approach to determine the sequence of authors.
Conflict of interest: The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this article.
Data availability statement: Data sharing is not applicable to this article as no datasets were generated or analysed during the current study.

References

[1] Lees B, Meredith LR, Kirkland AE, Bryant BE, Squeglia LM. Effect of alcohol use on the adolescent brain and behavior. Pharmacol Biochem Behav. 2020;192:172906. 10.1016/j.pbb.2020.172906Suche in Google Scholar PubMed PubMed Central

[2] World Health Organization. Global status report on alcohol and health; 2014. p. 1–68.Suche in Google Scholar

[3] James LM, Van Kampen E, Miller RD, Engdahl BE. Risk and protective factors associated with symptoms of post-traumatic stress, depression, and alcohol misuse in OEF/OIF veterans. Military Med. 2013;178(2):159–65. 10.7205/MILMED-D-12-00282Suche in Google Scholar PubMed

[4] Ferrari AJ, Charlson FJ, Norman RE, Patten SB, Freedman G, Murray CJ, et al. Burden of depressive disorders by country, sex, age, and year: findings from the global burden of disease study 2010. PLoS Med. 2013;10(11):1001547. 10.1371/journal.pmed.1001547Suche in Google Scholar PubMed PubMed Central

[5] Keyes KM, Hasin DS. Socio-economic status and problem alcohol use: the positive relationship between income and the DSM-IV alcohol abuse diagnosis. Addiction. 2008;103(7):1120–30. 10.1111/j.1360-0443.2008.02218.xSuche in Google Scholar PubMed PubMed Central

[6] Catalano R, Dooley D, Wilson G, Hough R. Job loss and alcohol abuse: a test using data from the epidemiologic Catchment area project. J Health Soc Behav. 1993;34:215–25. 10.2307/2137203Suche in Google Scholar

[7] WHO Global Status Report. Social Problems associated with alcohol use, Alcohol consumption and the workplace. 2004. p. 59–60.Suche in Google Scholar

[8] Johnson RC, Schwitters SY, Wilson JR, Nagoshi CT, McClearn GE. A cross-ethnic comparison of reasons given for using alcohol, not using alcohol or ceasing to use alcohol. J Stud Alcohol. 1985;46(4):283–8. 10.15288/jsa.1985.46.283Suche in Google Scholar PubMed

[9] World Health Organization. Global status report on alcohol and health 2018. World Health Organization; 2019. Suche in Google Scholar

[10] Anderson P. The UK government’s alcohol strategy. Drugs Educat Prevent Policy. 2012;19(5):360–1. 10.3109/09687637.2012.698425Suche in Google Scholar

[11] WHO global status report on alcohol and health Geneva. 2014.Suche in Google Scholar

[12] Smith JC. Basic interdiction models. Wiley encyclopedia of operations research and management science. Hoboken: John Wiley and Sons; 2010. 10.1002/9780470400531.eorms0089Suche in Google Scholar

[13] Chen L, Li Y, Huang M, Hui X, Gu S. Robust dynamic state estimator of integrated energy systems based on natural gas partial differential equations. IEEE Trans Industry Appl. 2022;58(3):3303–12. 10.1109/TIA.2022.3161607Suche in Google Scholar

[14] Ma P, Yao N, Yang X. Service quality evaluation of terminal express delivery based on an integrated SERVQUAL-AHP-TOPSIS approach. Math Problems Eng. 2021;2021:1–10. 10.1155/2021/8883370Suche in Google Scholar

[15] Khan FM, Khan ZU. Numerical analysis of fractional order drinking mathematical model. J Math Tech Model. 2024;1(1):11–24. Suche in Google Scholar

[16] Chinnadurai K, Athithan S, Kareem MG. Mathematical modelling on alcohol consumption control and its effect on poor population. IAENG Int J Appl Math. 2024;54(1):1–9. Suche in Google Scholar

[17] Imken I, Fatmi NI. A new mathematical model of drinking alcohol among diabetes population taking anti-diabetic drugs: an optimal control approach. Commun Math Biol Neurosci. 2024;2024:25. Suche in Google Scholar

[18] ur Rahman M, Arfan M, Shah Z, Alzahrani E. Evolution of fractional mathematical model for drinking under Atangana-Baleanu Caputo derivatives. Phys Scr. 2021;96(11):115203. 10.1088/1402-4896/ac1218Suche in Google Scholar

[19] Bunonyo KW, Ebiwareme L, Iworiso PB. Mathematical modeling of time-dependent concentration of alcohol in the human bloodstream using the eigenvalue method. TWIST. 2024;19(1):58–64. 10.37745/ijmss.13/vol11n4119Suche in Google Scholar

[20] Huo HF, Wang Q. Modelling the influence of awareness programs by media on the drinking dynamics. Abstr Appl Anal. 2014;2014:1–8. 10.1155/2014/938080Suche in Google Scholar

[21] Ma SH, Huo HF, Meng XY. Modelling alcoholism as a contagious disease: a mathematical model with awareness programs and time delay. Discrete Dynams Nature Soc. 2015;2015(2):1–13. 10.1155/2015/260195Suche in Google Scholar

[22] Wang XY, Huo HF, Kong QK, Shi WX. Optimal control strategies in an alcoholism model. Abstr Appl Anal. 2014;2014:1–18. 10.1155/2014/954069Suche in Google Scholar

[23] Sharma S, Samanta GP. Drinking as an epidemic: a mathematical model with dynamic behaviour. J Appl Math Inform. 2013;31(1–2):1–25. 10.14317/jami.2013.001Suche in Google Scholar

[24] Michel AN, Hou L, Liu D. Stability of dynamical systems. Birkhaüser, Boston, MA. 2008. Suche in Google Scholar

[25] Huo HF, Chen YL, Xiang H. Stability of a binge drinking model with delay. J Biol Dynam. 2017;11(1):210–25. 10.1080/17513758.2017.1301579Suche in Google Scholar PubMed

[26] Manthey JL, Aidoo AY, Ward KY. Campus drinking: an epidemiological model. J Biol Dynam. 2008;2(3):346–56. 10.1080/17513750801911169Suche in Google Scholar PubMed

[27] Huo HF, Huang SR, Wang XY, Xiang H. Optimal control of a social epidemic model with media coverage. J Biol Dynam. 2017;11(1):226–43. 10.1080/17513758.2017.1321792Suche in Google Scholar PubMed

[28] Huo HF, Song NN. Global stability for a binge drinking model with two stages. Discr Dynam Nat Soc. 2012;2012:1–15. 10.1155/2012/829386Suche in Google Scholar

[29] Anjam YN, Shahid I, Emadifar H, Arif Cheema S, ur Rahman M. Dynamics of the optimality control of transmission of infectious disease: a sensitivity analysis. Scientif Reports. 2024;14(1):1041. 10.1038/s41598-024-51540-7Suche in Google Scholar PubMed PubMed Central

[30] Biegler LT, Cervantes AM, Wächter A. Advances in simultaneous strategies for dynamic process optimization. Chem Eng Sci. 2002;57(4):575–93. 10.1016/S0009-2509(01)00376-1Suche in Google Scholar

[31] Feng Z, Yang Y, Xu D, Zhang P, McCauley MM, Glasser JW. Timely identification of optimal control strategies for emerging infectious diseases. J Theoretic Biol. 2009;259(1):165–71. 10.1016/j.jtbi.2009.03.006Suche in Google Scholar PubMed PubMed Central

[32] Lewis FL, Vrabie D, Syrmos VL. Optimal control. Hoboken: John Wiley and Sons; 2012. 10.1002/9781118122631Suche in Google Scholar

[33] Bryson AE. Applied optimal control: optimization, estimation and control. New York: Routledge; 2018. 10.1201/9781315137667Suche in Google Scholar

[34] Verma V, Agarwal M, Verma A. A mathematical model for the novel coronavirus with effect of lockdown. Int J Model Simulat Scientif Comput. 2023;14(03):2350005. 10.1142/S1793962323500058Suche in Google Scholar

[35] Verma V. Optimal control analysis of a mathematical model on smoking. Model Earth Syst Environ. 2020;6(4):2535–42. 10.1007/s40808-020-00847-1Suche in Google Scholar

[36] Omame A, Raezah AA, Diala UH, Onuoha C. The optimal strategies to be adopted in controlling the co-circulation of COVID-19, Dengue and HIV: Insight from a mathematical model. Axioms. 2023;12(8):773. 10.3390/axioms12080773Suche in Google Scholar

[37] Omame A, Abbas M, Onyenegecha CP. Backward bifurcation and optimal control in a co-infection model for SARS-CoV-2 and ZIKV. Results Phys. 2022;37:105481. 10.1016/j.rinp.2022.105481Suche in Google Scholar PubMed PubMed Central

[38] Omame A, Abbas M. Modeling SARS-CoV-2 and HBV co-dynamics with optimal control. Phys A Stat Mech Appl. 2023;615:128607. 10.1016/j.physa.2023.128607Suche in Google Scholar PubMed PubMed Central

[39] Pontryagin LS. Mathematical theory of optimal processes. New York: Routledge; 2018. 10.1201/9780203749319Suche in Google Scholar

[40] Adu IK, Osman MA, Yang C. Mathematical model of drinking epidemic. Br J Math Comput Sci. 2017;22(5):1–10. 10.9734/BJMCS/2017/33659Suche in Google Scholar

[41] Nagumo M. Über die lage der integralkurven gewöhnlicher differentialgleichungen. Proceedings of the Physico-Mathematical Society of Japan. 3rd Series. 1942; vol. 24. p. 551–9. Suche in Google Scholar

[42] Van den Driessche P, Watmough J. Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission. Math Biosci. 2002;180(1–2):29–48. 10.1016/S0025-5564(02)00108-6Suche in Google Scholar

[43] Van den Driessche P. Reproduction numbers of infectious disease models. Infect Disease Model. 2017;2(3):288–303. 10.1016/j.idm.2017.06.002Suche in Google Scholar PubMed PubMed Central

[44] Chitnis N, Hyman JM, Cushing JM. Determining important parameters in the spread of malaria through the sensitivity analysis of a mathematical model. Bulletin Math Biol. 2008;70:1272–96. 10.1007/s11538-008-9299-0Suche in Google Scholar PubMed

[45] Satana TS, Kassaye MT. Mathematical modeling and analysis of alcoholism epidemics: a case study in Ethiopia. J Pure Appl Math. 2023;7:70–80. 10.21203/rs.3.rs-1450151/v1Suche in Google Scholar

[46] Merkin DR. Introduction to the theory of stability. New York (NJ): Springer Science and Business Media; 2012. p. 24. Suche in Google Scholar

[47] Zhang DC, Shi B. Oscillation and global asymptotic stability in a discrete epidemic model. J Math Anal Appl. 2003;278(1):194–202. 10.1016/S0022-247X(02)00717-5Suche in Google Scholar

[48] Ma X, Zhou Y, Cao H. Global stability of the endemic equilibrium of a discrete SIR epidemic model. Adv Differ Equ. 2013;2013:1–9. 10.1186/1687-1847-2013-42Suche in Google Scholar PubMed PubMed Central

[49] La Salle JP. The stability of dynamical systems. Philadelphia (PA): Society for Industrial and Applied Mathematics; 1976. p. 25. Suche in Google Scholar

[50] Khajji B, Kouidere A, Balatif O, Rachik M. Mathematical modeling, analysis and optimal control of an alcohol drinking model with liver complication. Commun Math Biol Neurosci. 2020;2020:1–29. Suche in Google Scholar

[51] Khajji B, Labzai A, Kouidere A, Balatif O, Rachik M. A discrete mathematical modeling of the influence of alcohol treatment centers on the drinking dynamics using optimal control. J Appl Math. 2020;2020:1–3. 10.1155/2020/9284698Suche in Google Scholar

[52] Lee S, Jung E. Optimal control intervention strategies in low-and high-risk problem drinking populations. Socio-Econom Plan Sci. 2010;44(4):258–65. 10.1016/j.seps.2010.07.006Suche in Google Scholar

[53] Hwang CL, Fan LT. A discrete version of Pontryagin’s maximum principle. Operat Res. 1967;15(1):139–46. 10.1287/opre.15.1.139Suche in Google Scholar

[54] Guibout V, Bloch A. A discrete maximum principle for solving optimal control problems. In 2004 43rd IEEE Conference on Decision and Control; 2004. vol. 2. p. 1806–11. 10.1109/CDC.2004.1430309Suche in Google Scholar

[55] Zaman G, Kang YH, Jung IH. Optimal treatment of an SIR epidemic model with time delay. BioSystems. 2009;98(1):43–50. 10.1016/j.biosystems.2009.05.006Suche in Google Scholar PubMed

[56] Zakary O, Rachik M, Elmouki I. On the analysis of a multi-regions discrete SIR epidemic model: an optimal control approach. Int J Dynam Control. 2017;5:917–30. 10.1007/s40435-016-0233-2Suche in Google Scholar PubMed PubMed Central

[57] Todorov E, Jordan MI. Optimal feedback control as a theory of motor coordination. Nature Neurosci. 2002;5(11):1226–35. 10.1038/nn963Suche in Google Scholar PubMed

Received: 2023-12-15

Revised: 2024-03-26

Accepted: 2024-04-15

Published Online: 2024-07-12

This work is licensed under the Creative Commons Attribution 4.0 International License.

Artikel in diesem Heft

https://doi.org/10.1515/nleng-2024-0010

Schlagwörter für diesen Artikel

drinking model; stability; sensitivity; optimal control theory; numerical simulations

Creative Commons

BY 4.0